QQ登录

只需一步,快速开始

上位机MFC实现单像空间后方交会源代码

[ 复制链接 ]
一、程序流程图
2019-10-11_101425.jpg
二、部分代码
两个结构体,分别用来保存外方位元素和原始数据。
  1. struct EOEO                // elements of exterior orientation
  2. {
  3.         double Xs;
  4.         double Ys;
  5.         double Zs;
  6.         double omega;
  7.         double phi;
  8.         double kappa;
  9. };
复制代码
  1. 几个主要函数

  2. //////////////////////////////////////////////////////////////////////////
  3. // 函数功能:初始化坐标数据
  4. // 参数说明:sd:保存原始数据的结构体数组,PhotographicScale:摄影比例尺,focus:摄影机主距
  5. //                         filename:保存坐标数据的文件名
  6. //////////////////////////////////////////////////////////////////////////
  7. void InitData(SourceData* sd, char* filename)
  8. {
  9.         //打开文件以备读取数据
  10.         cout <<"开始读取数据..."<<endl;
  11.         ifstream datafile(filename,ios::in | ios::nocreate);
  12.         if( !datafile)
  13.         {
  14.                 cout<<"打开文件失败, 文件可能不存在"<<endl;
  15.                 system("pause");
  16.                 _exit(1);
  17.         }

  18.         memset(sd, 0, N*sizeof(SourceData));
  19.        
  20.         int i;
  21.         for (i = 0; i != N; i ++)
  22.         {
  23.                 datafile >> sd[i].x >> sd[i].y >> sd[i].X >> sd[i].Y >> sd[i].Z;
  24.         }
  25.         datafile.close();
  26.         cout <<"读取数据完毕..."<<endl;
  27. }


  28. //////////////////////////////////////////////////////////////////////////
  29. // 函数功能:检查改正数是否已达精度要求
  30. // 参数说明:deta:保存改正数的数组
  31. //////////////////////////////////////////////////////////////////////////
  32. bool CheckPrecision(double* deta)
  33. {
  34.         //2.9088820866572159615394846141477e-5 即 0.1′的弧度表示
  35.         bool ret;
  36.         ret = (fabs(deta[0])<0.000001 && fabs(deta[1])<0.000001        && fabs(deta[2])<0.000001 && \
  37.                 fabs(deta[3])<2.90888208656e-5 && fabs(deta[4])<2.90888208656e-5 && \
  38.                 fabs(deta[5])<2.90888208656e-5);
  39.         return ret;
  40. }

  41. //////////////////////////////////////////////////////////////////////////
  42. // 函数功能:迭代器,计算的主体部分
  43. // 参数说明:sd:保存原始数据的结构体数组,PhotographicScale:摄影比例尺,focus:摄影机主距
  44. //////////////////////////////////////////////////////////////////////////
  45. void Iterator(SourceData sd[N], double PhotographicScale, double Focus)
  46. {
  47.         double phi, omega, kappa, Xs, Ys, Zs;
  48.         phi = omega = kappa =0.0;
  49.         Xs = Ys = Zs = 0.0;

  50.         for (int k=0;k<N;k++)
  51.         {
  52.                 sd[k].x /= 1000.0;
  53.                 sd[k].y /= 1000.0;
  54.                 Xs += sd[k].X;
  55.                 Ys += sd[k].Y;
  56.         }
  57.         Xs /= N;
  58.         Ys /= N;
  59.        
  60.         double f = Focus / 1000.0;
  61.         double m = PhotographicScale;
  62.         Zs =m*f;
  63.        
  64.         cout<<endl<<"六元素初始化值:"<<endl;
  65.         cout <<"Xs: "<<Xs<<'\t'<<"Ys: "<<Ys<<'\t'<<"Zs: "<<Zs<<endl;
  66.         cout <<"phi: "<<phi<<"\t\t"<<"omega: "<<omega<<"\t"<<"kappa: "<<kappa<<endl<<endl;

  67.         //声明并初始化六元素改正数矩阵
  68.         double deta[6] = {1,1,1,1,1,1};

  69.         double x0(0);
  70.         double y0(0);//内方位元素
  71.         double X0[N] = {0.0};
  72.         double Y0[N] = {0.0};
  73.         double Z0[N] = {0.0};

  74.         //声明旋转矩阵
  75.         double R[9];

  76.         double A[2*6] = {0.0};
  77.         double AT[6*2] = {0.0};

  78.         double Buf1[36] = {0.0};
  79.         double Buf2[36] = {0.0};//ATA累加

  80.         double Buf3[6] = {0.0};
  81.         double Buf4[6] = {0.0};//ATL累加
  82.        
  83.         double Buf5[8*6] = {0.0};//存储8×6的A矩阵,没办法
  84.         double Buf6[8*1] = {0.0};//存储8×1的L矩阵,同上
  85.         double V[8*1] = {0.0};

  86.         double ATA[36] = {0.0};
  87.         double ATL[6] = {0.0};
  88.         double L[2] = {0.0};

  89.         int iCount = 0;//迭代次数
  90.         cout << "开始迭代计算..."<<endl;
  91.         while(!CheckPrecision(deta))
  92.         {
  93.                 cout <<endl<<"第 "<< ++iCount<<"次迭代:"<<endl;
  94.                 if (iCount == MAXITERATION)
  95.                 {
  96.                         cout << ">>>迭代次数超限,可能不收敛"<<endl;
  97.                         break;
  98.                 }

  99.                 //每次迭代之前必须清空两个保存累加值的矩阵ATA与ATL
  100.                 for (int i = 0; i != 36; i ++)
  101.                 {
  102.                         ATA[i] = 0.0;
  103.                         if (i < 6)
  104.                         {
  105.                                 ATL[i] = 0.0;
  106.                         }
  107.                 }

  108.                 //计算旋转矩阵
  109.                 R[0] = cos(phi)*cos(kappa) - sin(phi)*sin(omega)*sin(kappa);
  110.                 R[1] = -cos(phi)*sin(kappa) - sin(phi)*sin(omega)*cos(kappa);
  111.                 R[2] = -sin(phi)*cos(omega);
  112.                 R[3] = cos(omega)*sin(kappa);
  113.                 R[4] = cos(omega)*cos(kappa);
  114.                 R[5] = -sin(omega);
  115.                 R[6] = sin(phi)*cos(kappa) + cos(phi)*sin(omega)*sin(kappa);
  116.                 R[7] = -sin(phi)*sin(kappa) + cos(phi)*sin(omega)*cos(kappa);
  117.                 R[8] = cos(phi)*cos(omega);

  118.                 for (i = 0; i != N; i ++)
  119.                 {
  120.                         Z0[i] = R[2]*(sd[i].X-Xs) + R[5]*(sd[i].Y-Ys) + R[8]*(sd[i].Z-Zs);
  121.                         X0[i] = x0 - f*(R[0]*(sd[i].X-Xs) + R[3]*(sd[i].Y-Ys) + R[6]*(sd[i].Z-Zs)) / Z0[i];
  122.                         Y0[i] = y0 - f*(R[1]*(sd[i].X-Xs) + R[4]*(sd[i].Y-Ys) + R[7]*(sd[i].Z-Zs)) / Z0[i];

  123.                         A[0] = ((R[0]*f + R[2]*(sd[i].x-x0)) / Z0[i]);
  124.                         A[1] = ((R[3]*f + R[5]*(sd[i].x-x0)) / Z0[i]);
  125.                         A[2] = ((R[6]*f + R[8]*(sd[i].x-x0)) / Z0[i]);
  126.                         A[3] = ((sd[i].y-y0)*sin(omega) - ((sd[i].x-x0)*((sd[i].x-x0)*cos(kappa)\
  127.                                 - (sd[i].y-y0)*sin(kappa))/f + f*cos(kappa))*cos(omega));
  128.                         A[4] = (-f*sin(kappa) - (sd[i].x-x0)*((sd[i].x-x0)*sin(kappa) + (sd[i].y-y0)*cos(kappa))/f);
  129.                         A[5] = (sd[i].y - y0);
  130.                         A[6] = ((R[1]*f + R[2]*(sd[i].y-y0)) /Z0[i]);
  131.                         A[7] = ((R[4]*f + R[5]*(sd[i].y-y0)) /Z0[i]);
  132.                         A[8] = ((R[7]*f + R[8]*(sd[i].y-y0)) /Z0[i]);
  133.                         A[9] = (-(sd[i].x-x0)*sin(omega) - ((sd[i].y-y0)*((sd[i].x-x0)*cos(kappa)\
  134.                                 - (sd[i].y-y0)*sin(kappa))/f - f*sin(kappa))*cos(omega));
  135.                         A[10] = (-f*cos(kappa) - (sd[i].y-y0)*((sd[i].x-x0)*sin(kappa) + (sd[i].y-y0)*cos(kappa))/f);
  136.                         A[11] = (-(sd[i].x-x0));

  137.                         //该循环保存A矩阵,最后评定精度用
  138.                         for (int l=0;l<12;l++)
  139.                         {
  140.                                 Buf5[12*i+l] = A[l];
  141.                         }
  142.                        
  143.                         //所谓的逐步法化,即要在循环内部就将ATA计算出来累加,下面的L矩阵类似
  144.                         MatrixTranspose(A,AT,2,6);
  145.                         MatrixMultiply(AT,A,Buf1,6,2,6);
  146.                         MatrixCopy(ATA, Buf2, 36);
  147.                         MatrixAdd(Buf1,Buf2,ATA,36);        // 为逐步法化后的ATA矩阵累加

  148.                         L[0] = (sd[i].x - X0[i]);
  149.                         L[1] = (sd[i].y - Y0[i]);

  150.                         //保存L矩阵,最后评定精度用
  151.                         for (l=0;l<2;l++)
  152.                         {
  153.                                 Buf6[2*i+l] = L[l];
  154.                         }

  155.                         MatrixMultiply(AT,L,Buf3,6,2,1);
  156.                         MatrixCopy(ATL, Buf4, 6);
  157.                         MatrixAdd(Buf3,Buf4,ATL,6);                       
  158.                 }//for

  159.                 //“逐步法化”的另一处不同,出循环即可直接计算ATA逆乘ATL
  160.                 MatrixInversion(ATA,6);
  161.                 MatrixMultiply(ATA,ATL,deta,6,6,1);

  162.                 //deta即为改正数
  163.                 Xs += deta[0];
  164.                 Ys += deta[1];
  165.                 Zs += deta[2];
  166.                 phi += deta[3];
  167.                 omega += deta[4];
  168.                 kappa += deta[5];

  169.                 cout<<"改正数值:"<<endl;
  170.                 for (i=0; i != 6; i ++)
  171.                 {
  172.                         cout <<"deta["<<i<<"]: "<<deta[i]<<endl;
  173.                 }
  174.                 cout << endl<<"六元素值:"<<endl;
  175.                 cout <<"Xs: "<<Xs<<endl<<"Ys: "<<Ys<<endl<<"Zs: "<<Zs<<endl;
  176.                 cout <<"kappa: "<<kappa<<endl<<"omega: "<<omega<<endl<<"phi: "<<phi<<endl<<endl;
  177.         }//while

  178.         EOEO eoeo;
  179.         eoeo.kappa = kappa;
  180.         eoeo.omega = omega;
  181.         eoeo.phi = phi;
  182.         eoeo.Xs = Xs;
  183.         eoeo.Ys = Ys;
  184.         eoeo.Zs = Zs;

  185.         cout << ">>>正常退出迭代"<<endl<<endl<<endl;

  186.         //精度评定
  187.         double Q[6] = {0.0};
  188.         for (int h=0;h<6;h++)
  189.         {
  190.                 Q[h] = ATA[h*6+h];
  191.         }

  192.         MatrixMultiply(Buf5,deta,V,8,6,1);//V=Ax-L
  193.         MatrixMinus(V,Buf6,V,8);

  194.         double m0(0);//单位权中误差
  195.         double VSum(0);//[vv],即平方和

  196.         int i;
  197.         for (i=0;i<8;i++)
  198.         {
  199.                 VSum += V[i]*V[i];
  200.         }
  201.         m0=sqrt(VSum / (2*N - 6));//中误差m0

  202.         double M[6] = {0.0};//保存六个值的中误差
  203.         for (i = 0; i != 6; i ++)
  204.         {
  205.                 M[i] = m0 * sqrt(Q[i]);
  206.                 if (i >= 3)
  207.                 {
  208.                         M[i] = M[i]*180*3600/PI;
  209.                 }
  210.         }

  211.         OutputResult(&eoeo, R, M, m0);
  212.         cout<<endl<<"解算全部完成"<<endl<<endl;
  213. }

复制代码
三、计算结果
1.png
四、本程序的两种拓展形式
1、  将空间后方交会计算封装成类
为了实现空间后交计算的移植性,这里我将该过程封装成一个类CResection,并定义了相应的接口。
该类接受原始数据的接口为带参构造函数,其原型为:
CResection ();
输出接口为两个类成员函数:
       OutputResult();     //该成员适用于windows console application 程序调用时调用iostream库输出结果。
       SaveResult();        //该成员适用于任何情况下将计算结果保存为文件,用户指定文件路径与文件名。
       以上三个函数类型声明为public以满足从外部环境声明类对象调用的需要。另包含3个私有成员函数Iterator (),
CheckPrecision(),InitData() 进行内部计算,编码与运行时对用户不可见。
       详细请参看工程CResection代码。
用户调用该类时,只需将CResection.h与CResection.cpp添加进目标工程并在文件适当位置包含该类头文件即可。
使用实例:
#include “CResection.h”
//…………………
CResection res (50000,154.23, “SourceData.txt”);
res.SaveResult (“Results.txt”);
res.OutputResult();       //仅在控制台程序可用
2、  将该功能封装成动态链接库
封装成动态链接库,只要将类CResection稍加修改即可,对于代码能实现可靠的隐蔽性。
调用该动态链接库,只需在工程设置中link标签页下添加CResection.lib,
并在需要调用该库的文件中包含CResection.h,而后即可声明该类对象进行空间后交计算。
详细请参看工程CResection_dll代码
请点击此处下载

请先注册会员后在进行下载

已注册会员,请先登录后下载

文件名称:上位机MFC实现单像空间后方交会源代码.rar 
文件大小:1.79 MB  售价:3金币
下载权限: 不限 以上或 VIP会员   [购买捐助会员]   [充值积分]   有问题联系我

  

halcon从自学到接项目视频教程,另外再赠送全网最全资源  

  

欢迎围观我录制的一套halcon自学视频教程(进入)

  

上位机VC MFC程序开发精典实例大全源码与视频讲解配套下载408例

  

经历1年的编程与录制点击进入查看

  

如果您认可,可联系功能定制!

  

如果您着急,充值会员可直接联系发您资料!

  

QQ联系我

微信扫扫联系我

  

回复

使用道具 举报

快速回复 返回列表 客服中心 搜索