pointMat.write(i,j,point[k]);
}
}
cal.transpos(&POINTMat,&POINTTMat);//矩阵求逆
//最小二乘法求矩阵系数,结果存在了solveMat中,为6行二列的矩阵 _Matrix mul1Mat(6,6); mul1Mat.init_matrix();
cal.multiply(&POINTTMat,&POINTMat,&mul1Mat); _Matrix mul2Mat(6,2); mul2Mat.init_matrix();
cal.multiply(&POINTTMat,&pointMat,&mul2Mat); _Matrix invMat(6,6); invMat.init_matrix();
cal.inverse(&mul1Mat,&invMat);
cal.multiply(&invMat,&mul2Mat,&solveMat); 2. 遥感图像的几何纠正变换
首先确定输出地图坐标的边界坐标,通过给出的输入图像左上角端点的坐标值求解出原始图像的各个端点坐标,通过多项式进行转换得到输入地图的端点坐标。
Mat image;
image = imread(\ if(image.empty()) {return;}
getArea(X1,Y1,X2,Y2);//得到输出图像的地理坐标范围 int row,col; row = image.rows; col = image.cols;
//得到输出图像的行列号,已知该图像分辨率为30m int M = (Y2 - Y1)/ 30 + 1; int N = (X2 - X1) / 30 + 1;
unsigned char * im = new unsigned char[M*N];//存储输出图像像素 for(i = 0 ; i < M;i++) {
for(j = 0 ; j k++; } im[i*N + j] = 0; }//先赋值为0 3. 数字图像亮度值重采样 数字图像亮度重采样的方法有很多,这里采用了两种方法。 (1)最邻近像元法,采取距离被采样点最近的像素值作为采样亮度。 for( i = 1 ; i <= M ;i++)//重采样 { for( j = 1 ; j<= N;j++) { double tempX = X1 + (j - 1)*30; double tempY = Y2 - (i - 1)*30; double geo[6] = {1,tempX,tempY,tempX*tempY,tempX*tempX,tempY*tempY}; _Matrix geoMat(1,6); geoMat.init_matrix(); int what; int whatk = 0; for (what = 0 ; what < 6;what++) { geoMat.write(0,what,geo[whatk]); whatk ++; } _Matrix imgMat(1,2); imgMat.init_matrix(); cal.multiply(&geoMat,&solveMat,&imgMat); double imgX = imgMat.arr[0];//原始图像坐标 double imgY = imgMat.arr[1]; int t1 = (int(imgX+0.5) - 800); int t2 = (-(int(imgY+0.5 ) + 3600));//最邻近插值 if(t1 >= 0 &&t1 < col && t2 >= 0 && t2 < row) { im[(i - 1) * N + (j - 1)] = image.data[t1 + t2* col] ; } geoMat.free_matrix(); imgMat.free_matrix(); } } (2)双线性内插法。对被采样点周围的四个像素的亮度值进行双线性内插,具体原理可以用图1所示的三角形线性函数来表示。 for(i=0;i for(j=0;j calculate(&x1,&y1,Xp,Yp,a,b); if(x1>=xa&&x1<=xc && y1>=ya&&y1<=yb){ double r,c; c = x1 - xa; r = yb - y1; double deltax,deltay; deltax = int(r)+1-r; deltay = int(c)+1-c; double Wx1,Wy1,Wx2,Wy2; Wx1 = 1-deltax; Wy1 = 1-deltay; Wx2 = deltax; Wy2 = deltay; int I11,I12,I21,I22; I11 = lpData[nWidth*int(r) + int(c)]; I12 = lpData[nWidth*int(r) + int(c)+1]; I21 = lpData[nWidth*(int(r)+1) + int(c)]; I22 = lpData[nWidth*(int(r)+1) + int(c)+1]; lpRData[N*i + j] = Wx1*Wy1*I11 + Wx2*Wy1*I21 + Wx1*Wy2*I12 + Wx2*Wy2*I22;} } else lpRData[N*i + j] = 0;} 3 软件使用方法 1. 打开Microsoft Visual Studio 2010 软件,新建一个Win32控制台应用程序。 2. 配置相应的OpenCV环境。 首先将OpenCV2.2 文件夹放在与工程根目录平行的目录之下。 回到工程,在“解决方案资源管理器”中,右键单击工程名,弹出菜单中点选“属性”。 在弹出的属性页对话框中,选择“VC++目录”-点击“包含目录”右边的倒三角-点击“编辑”。 在弹出的包含目录对话框中,分别加入“...\\OpenCV2.2\\include”,“...\\OpenCV2.2\\include\\opencv”,“...\\OpenCV2.2\\include\\opencv2”三个文件夹。最后点击“确定”。 点击“库录”右边的倒三角-点击“编辑”。在弹出的库目录对话框中,加入“...\\OpenCV2.2\\lib”文件夹。最后点击“确定” 在属性页对话框中,选择“链接器”-“输入”。点击“附加依赖项”右边的倒三角-点击“编辑”。在弹出的“附加依赖项”对话框中,分别填入如下文件名: opencv_core220d.lib opencv_highgui220d.lib opencv_video220d.lib opencv_ml220d.lib opencv_legacy220d.lib opencv_imgproc220d.lib 点击“确定”。最后点击“应用”-“确定”。 通过上述方法的操作,即完成了OpenCV环境的配置。 3. 程序的编写。基于多项式的遥感图像几何纠正处理可以大致分为三个步骤:通过地面已知控制点解算多项式系数、遥感图像的纠正变换、数字图像亮度值的重采样。 4.实验结果与分析 基于多项式的遥感图像几何纠正处理前后的效果图如下: 图2:几何纠正效果图 从这张结果图像可以看出,原先具有几何变形的图像得到了一定程度的纠正,图像信息得到了良好的还原。 在亮度重采样的方法选择上,我们采用了最邻近像元法和双线性插值法两种方法进行对比。效果图如下: 图3:最邻近像元法 图4:双线性插值法 从对比上可以看出,图像的几何变形都得到了较为完善的纠正,目视判读的纠正结果基本一致,但是仔细对比可以发现后者相对于前者在图像细节上要模糊一些,但是根据理论后者的几何精度更高。因此纠正方法的选择还是应该综合考虑实际因素。 5 讨论和总结 5.1 小组分工 本实习小组共有贾丹、刘爽、石东博和吴京航四人。 贾丹主要负责答辩PPT 的主要制作。 刘爽主要负责答辩PPT的参与制作,实习报告的参与编写。 石东博主要负责实习报告的主要编写,答辩PPT的参与制作。 吴京航主要负责用于答辩展示的代码编写。 此外组内四人都已各自独立完成彩色图像的几何变形、基于多项式的遥感图像几何纠正的代码的编写。 5.2 实习总结 本次遥感原理课程设计的编程部分为期五天,共六个课时,我们的收获有很多。