DLT.rar_DLT_DLT C_visual c

#include <iostream.h> #include <math.h> #include <stdio.h> #define DataType long double #define N 12 int GaussEliminant(DataType a[N][N],DataType b[N],DataType c[N]);//利用高斯消元法解多元一次线性方程组 int main(int argc,char *argv[]) { DataType control[7][4]= // 6个控制点的(XYZ)坐标 { {0,0,0,0}, {0,1.311,1.751,1.951}, {0,5.121,5.531,5.941}, {0,10.511,10.751,10.191}, {0,15.311,15.551,15.971}, {0,20.711,20.551,20.981}, {0,25.311,25.541,25.951}, }; DataType photo_control[7][3] = //6个控制点在图像内的坐标(XY) { {0,0,0}, {0,1.211,1.531}, {0,5.141,5.551}, {0,10.161,10.751}, {0,15.181,15.591}, {0,20.121,20.541}, {0,25.151,25.561}, }; DataType b[N]; DataType L[N]; for(int ii=0; ii<N; ii++) { b[ii] = 0;//b 是方程右边 L[ii] = 0;//L 是L系数 } DataType a[N][N];//相当于方程Al = -B 左侧的A a[0]行与列不参与运算 for (ii=1; ii<=11; ii++)//初始化为0 for(int jj=1;jj<=11;jj++) a[ii][jj] =0; //以下是给a[11][11]赋值过程 for ( ii=1; ii<=11; ii++)//ii 行下标 for( int jj=1;jj<=11;jj++)//jj 列下标 { if (ii%2)//奇数行 { if (jj<=3)//奇数行的前3项与实际控制点坐标相同 { a[ii][jj] = control[ii/2+1][jj]; } else if (jj==4)//奇数行的第4项为1 { a[ii][jj] = 1; } else if (jj>4&&jj<8)//第4到第8项为0 { a[ii][jj] = 0; } else if (jj>8&&jj<N)//第8到第11项为一乘积 { a[ii][jj] = control[ii/2+1][jj-8]*photo_control[ii/2+1][1]; } } else if (ii%2 == 0) { if (jj<=4) { a[ii][jj] = 0; } else if (jj>4&&jj<8) { a[ii][jj] = control[ii/2][jj-4]; } else if (jj==8) { a[ii][jj] = 1; } else if (jj>8&&jj<N) { a[ii][jj] = control[ii/2][jj-8]*photo_control[ii/2][2]; } } }//N维数组A赋值过程结束 //给B数组赋值 for ( ii=1;ii<N;ii++) { if (ii%2) { b[ii] = photo_control[(ii+1)/2][1]; } else if(ii%2==0) { b[ii] = photo_control[ii/2][2]; } } char szBuff[50];//临时数组用于格式化输出DOUBLE cout<<"方程左边矩阵A为："; for (ii=1;ii<N;ii++) { cout<<endl; for (int jj=1;jj<N;jj++) { sprintf(szBuff,"%f",a[ii][jj]); cout<<szBuff<<" "; } } cout<<endl<<"方程右边矩阵B为："<<endl; for (ii=0;ii<N;ii++) { cout<<b[ii]<<endl; } GaussEliminant(a,b,L); return 0; } int GaussEliminant(DataType a[N][N],DataType b[N],DataType c[N]) { DataType d=0,t=0,sum=0; int k=1,l=0,i=0,j=0,m,n; n=N-1;//n代表方程的组数 while(k<=n) { //*****选主元素***** d=a[k][k]; l=k; for(i=k+1;i<=n;i++) { if(fabs(a[i][k])>fabs(d)) { d=a[i][k];//现在d中存放的是：a[k]列中最大的系数 就是对角的那一列 l=i;// l 中存放：最大数的行标 } } if(l!=k)//如果对角元素在这一列中不是最大的 { for(j=k;j<=n;j++)//则将对角元素所在式 与 最大元素所在式进行 系数 交换 //两式对调 { //目的是：使最大元素在对角线上 t=a[l][j]; a[l][j]=a[k][j]; a[k][j]=t; } //下面三步是交换方程 等号 右面的系数 t=b[k]; b[k]=b[l]; b[l]=t; } //*****消元***** for(j=k+1;j<=n;j++) { a[k][j]=a[k][j]/a[k][k]; } b[k]=b[k]/a[k][k]; for(i=k+1;i<=n;i++) { for(j=k+1;j<=n;j++) { a[i][j]=a[i][j]-a[i][k]*a[k][j]; } j=1; b[i]=b[i]-a[i][k]*b[k]; } k++; } //*****回代***** for(i=n-1;i>=1;i--) { sum=0; for(j=i+1;j<=n;j++) { sum=sum+a[i][j]*b[j]; } b[i]=b[i]-sum; } cout<<endl<<endl<<"最后结果:"<<endl; for(int loop=1;loop<=n;loop++) { c[loop] = b[loop]; cout<<"x"<<loop<<"="<<c[loop]<<endl; } return 0; }