单像空间后方交会c++程序代码

武汉大学摄影测量原理课后练习空间后方交会作业，C++编程实现单像空间后方交会算法，并根据实验数据输出结果，并进行精度评定

//单张相片的空间后方交会程序
//说明：程序未作底片变形系统误差改正，像点坐标直接使用且内方位元素X0=Y0=0像点坐标为等精度不相关观测值
#include 
#include //标准输出输入流
#include //用于数学计算
#include  //调用函数exit（）需包含此文件
#include //控制输出格式需包含此文件
#include “matrixfunction.h“//矩阵运算头文件
int N=4;//控制点总数
/////////////////////////////////读入数据函数/////////////////////////////////////////

void getdata（char *filenamedouble *xdouble *ydouble *x1double *y1double *z1）
//读取数据文件（文件名，像x，像y，物x1，物y1，物z1）
{
	ifstream filein（filenameios::in|ios::nocreate）;  
    // 以输入方式打开一个文件，ios::nocreate表示原文件不存在时，不创建
	if（!filein）
	{  cout<<“打开文件失败!文件可能不存在!“<	   exit（1）;
	}
	for（int i=0;i		filein>>x[i]>>y[i]>>x1[i]>>y1[i]>>z1[i];//读入每行（每个控制点）的像方和物方坐标
	filein.close（）;
}

int main（）
{	
	int ij;//循环所用辅助参数
	double a1a2a3b1b2b3c1c2c3;//旋转矩阵九个元素
	double f=153.24;//摄影机主距（mm）
	double X0=0;	double Y0=0;//内方位元素为0
	double M=50000;//摄影比例尺分母
	double xsyszsfaiomgkaf;//6个外方位元素
	double *A=new double [2*N*6];//存放外方位元素系数A的矩阵2i*6
	double det[6];//6个外方位元素改正值矩阵
	double *L=new double [2*N];
	double *gpx=new double [N];double *gpy=new double [N];double *gpz=new double [N];//存储控制点物方坐标
	double *ipx=new double [N];double *ipy=new double [N];//存储控制点像方坐标
	double *ipx0=new double [N];double *ipy0=new double [N];//存储由共线条件解算的近似像点坐标
	double *AT=new double[2*N*6];//存储外方位元素系数转置矩阵AT6*2i
	double *TEP1=new double[36];double *TEP2=new double[6];double *TEP3=new double [36];//中间临时矩阵
	double Q[6];//存储6个权倒数
	double *V=new double [2*N];//存储像点观测值改正数
	double Vsum=0m0;//单位权中误差
	double m[6];//存储六个外方位元素的中误差
	////////////////////////读入观测数据/////////////////////////////
	getdata（“sourcefile.txt“ipxipygpxgpygpz）;
	////////////////////////确定未知数的初始值////////////////////////
	f/=1000;//将主距f单位从mm划成m
	fai=omg=kaf=0.0;
	xs=ys=0.0;zs=M*f;
	for（i=0;i	{
		ipx[i]/=1000;
		ipy[i]/=1000;
		xs+=gpx[i];
		ys+=gpy[i];
	}
	xs/=N;
	ys/=N;
	//////////////////////迭代///////////////////////
	while （fabs（det[0]）>0.000001||fabs（det[1]）>0.000001||fabs（det[2]）>0.000001||fabs（det[3]）>2.90888208656e-5
		||fabs（det[4]）>2.90888208656e-5||fabs（det[5]）>2.90888208656e-5）	//2.90888208656e-5为0.1分换算成的弧度值
	{
		/////////////////计算旋转矩阵
		a1=cos（fai）*cos（kaf）-sin（fai）*sin（omg）*sin（kaf）;
		a2=-cos（fai）*sin（kaf）-sin（fai）*sin（omg）*cos（kaf）;
		a3=-sin（fai）*cos（omg）;
		b1=cos（omg）*sin（kaf）;
		b2=cos（omg）*cos（kaf）;
		b3=-sin（omg）;
		c1=sin（fai）*cos（kaf）+cos（fai）*sin（omg）*sin（kaf）;
		c2=-sin（fai）*sin（kaf）+cos（fai）*sin（omg）*cos（kaf）;
		c3=cos（fai）*cos（omg）;

		///////////////////////////逐点计算///////////////////////////

		for （i=0;i		{
			////////////////根据共线方程由初始外方位元素计算像点近似值//////////////////
			ipx0[i]=X0-f*（a1*（gpx[i]-xs）+b1*（gpy[i]-ys）+c1*（gpz[i]-zs））/（a3*（gpx[i]-xs）+b3*（gpy[i]-ys）+c3*（gpz[i]-zs））;
			ipy0[i]=Y0-f*（a2*（gpx[i]-xs）+b2*（gpy[i]-ys）+c2*（gpz[i]-zs））/（a3*（gpx[i]-xs）+b3*（gpy[i]-ys）+c3*（gpz[i]-zs））;
			//////////////////////逐点计算误差方程的系数矩阵A///////////////

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----

     文件     262239  2009-10-25 21:11  Space resection\Debug\Spaceresection.exe

     文件     304600  2009-10-25 21:11  Space resection\Debug\Spaceresection.ilk

     文件      41070  2009-10-25 21:11  Space resection\Debug\Spaceresection.obj

     文件     343320  2009-10-25 21:11  Space resection\Debug\Spaceresection.pch

     文件     484352  2009-10-25 21:11  Space resection\Debug\Spaceresection.pdb

     文件      50176  2009-10-25 21:11  Space resection\Debug\vc60.idb

     文件      69632  2009-10-25 21:11  Space resection\Debug\vc60.pdb

     文件       2949  2009-10-20 11:51  Space resection\matrixfunction.h

     文件        433  2009-10-25 21:11  Space resection\result.txt

     文件        244  2009-10-04 15:36  Space resection\sourcefile.txt

     文件       8103  2009-10-22 14:17  Space resection\Spaceresection.cpp

     文件       3497  2009-10-20 21:54  Space resection\Spaceresection.dsp

     文件        553  2009-10-20 22:04  Space resection\Spaceresection.dsw

     文件      33792  2009-10-25 21:11  Space resection\Spaceresection.ncb

     文件      53760  2009-10-25 21:11  Space resection\Spaceresection.opt

     文件       1227  2009-10-25 21:11  Space resection\Spaceresection.plg

     目录          0  2009-10-25 21:11  Space resection\Debug

     目录          0  2009-10-25 21:11  Space resection

     文件        244  2009-10-04 15:36  后方交会数据.txt

----------- ---------  ---------- -----  ----

              1660191                    19