中小型网站建设怎么样,湖北今天发生的重大新闻,百度大数据分析,网站建设及相关流程图关于GDAL计算图像坐标的几个问题_gdal读取菱形四角点坐标-CSDN博客
这篇文章写的很好#xff0c;讲清楚了图像行列号与图像点坐标#xff08;x,y#xff09;对应关系#xff0c;以及图像行列号如何转为地理坐标的#xff0c;转载一下做个备份。
1.关于GDAL计算图像坐标的…关于GDAL计算图像坐标的几个问题_gdal读取菱形四角点坐标-CSDN博客
这篇文章写的很好讲清楚了图像行列号与图像点坐标x,y对应关系以及图像行列号如何转为地理坐标的转载一下做个备份。
1.关于GDAL计算图像坐标的几个问题
使用GDAL处理地理图像时不可避免的会遇到一个问题图像的地理坐标问题因为有了这个地理坐标地理图像才和普通图像有了最本质的区别那么在使用GDAL时如何处理与地理坐标相关的信息呢下面进行简单的说明。
1如何使用行列号计算图像的地理坐标或者如何通过地理坐标来定位在图像的某个位置
2如何获取图像的四至范围或者如果通过指定的地理范围计算图像的所在区域
要解决上面三个问题首先需要知道和了解GDAL的数据模型其中里面有个非常重要的就是投影和六参数。这两个可以使用GDALDataset类中的GeoTransform()函数和GetProjectionRef()函数来进行获取。第一个参数获取的是图像的六参数我自己起的名字是一个仿射变化的参数第二个是图像的投影也就是空间参考系统。下面先说说第一个六参数六参数其实是图像行列号坐标和地理坐标转换的一组转换系数。下面是用GT来表示六参数图像行列号与图像的地理坐标之间的数学关系式如下 Xgeo GT(0) Xpixel*GT(1) Yline*GT(2)Ygeo GT(3) Xpixel*GT(4) Yline*GT(5) 上式中Xgeo和Ygeo表示的图像的地理坐标Xpixel表示图像的列号Yline表示图像的行号GT(i)就是上面所说的六参数一共是六个值。这六个值大致可以分为三组GT(0)和GT(3)是第一组表示图像左上角的地理坐标GT(1)和GT(5)是第二组表示图像横向和纵向的分辨率一般这两者的值相等符号相反横向分辨率为正数纵向分辨率为负数GT(2)和GT(4)是第三组表示图像旋转系数对于一般图像来说这两个值都为0。
为什么说图像的GT(0)和GT(3)表示图像左上角的坐标对于图像行列号坐标系统来说坐标的原点在左上角所以左上角的行列号是(0,0)将坐标带入上式可以得到 Xgeo GT(0)Ygeo GT(3) 所以说GT(0)和GT(3)表示图像左上角的坐标。
GT(1)和GT(5)表示图像横向和纵向的分辨率。图像的分辨率就是图像每个像素所能表示的面积一般都是正方形的格网所以也就是没两个相邻像元坐标的差值。基于这个原理使用两个坐标进行验证。假设当前点行列号坐标为Aij相邻的右侧点坐标为Bi1j。分别计算A和B的横向地理坐标并计算差值即 dX XgeoB - XgeoA [GT(0) (i1)*GT(1) j*GT(2)] - [GT(0) i*GT(1) j*GT(2)] GT(0) -GT(0) (i1)*GT(1) - i*GT(1) j*GT(2) - j*GT(2) (i1)*GT(1) - i*GT(1) GT(1) 同理可以得到 dY GT(5)。
对于一个普通的标准图像来说这里的标准图像是指GT(2)和GT(4)都为0如图1所示图像的行列号坐标为XOY每个网格代表一个图像像素区域i表示列号j表示行号淡蓝色右下角的行列坐标为ij图中红色方块纵向长度为dy横向长度为dx分别为图像的分辨率图中O点的地理坐标就是GT(0),GT(3)。 图1 一个标准的图像行列号坐标及其地理坐标说明
有了上面的说明那么就可以很简单的来进行图像的行列号与地理坐标进行相互转换具体的代码如下共有两个一个正算一个反算。
bool Projection2ImageRowCol(double *adfGeoTransform, double dProjX, double dProjY, int iCol, int iRow)
{try{double dTemp adfGeoTransform[1]*adfGeoTransform[5] - adfGeoTransform[2]*adfGeoTransform[4];double dCol 0.0, dRow 0.0;dCol (adfGeoTransform[5]*(dProjX - adfGeoTransform[0]) - adfGeoTransform[2]*(dProjY - adfGeoTransform[3])) / dTemp 0.5;dRow (adfGeoTransform[1]*(dProjY - adfGeoTransform[3]) - adfGeoTransform[4]*(dProjX - adfGeoTransform[0])) / dTemp 0.5;iCol static_cast(dCol);iRow static_cast(dRow);return true;}catch(...){return false;}
}bool ImageRowCol2Projection(double *adfGeoTransform, int iCol, int iRow, double dProjX, double dProjY)
{//adfGeoTransform[6] 数组adfGeoTransform保存的是仿射变换中的一些参数分别含义见下//adfGeoTransform[0] 左上角x坐标 //adfGeoTransform[1] 东西方向分辨率//adfGeoTransform[2] 旋转角度, 0表示图像 北方朝上//adfGeoTransform[3] 左上角y坐标 //adfGeoTransform[4] 旋转角度, 0表示图像 北方朝上//adfGeoTransform[5] 南北方向分辨率try{dProjX adfGeoTransform[0] adfGeoTransform[1] * iCol adfGeoTransform[2] * iRow;dProjY adfGeoTransform[3] adfGeoTransform[4] * iCol adfGeoTransform[5] * iRow;return true;}catch(...){return false;}
}
现在我们再回到之前开头的两个问题。对于第一个问题实际就是图像行列号坐标与地理坐标的相互转换上面的代码就可以用来解决。对于第二个问题可以转换为第一个问题第二个问题其实就是两个点的坐标转换分别是左上角点和右下角点。比如第一个如何计算图像的四至范围图像的四至范围从图1中可以看出图像的四至其实就是图像左上角坐标和右下角的坐标为起来的矩形区域那么就分别将左上角和右下角的行列号按照上面的公式进行转换即可得到四至范围对与第二问就是将这个坐标进行反算得到。 关于地理图像的坐标问题就说到这里关于上面的投影信息空间参考信息需要说明一下这个六参数里面的坐标范围是和空间参考是一一对应的比如空间参考是一个WGS84椭球体那么这个六参数一般的单位就是度如果是北京54等分带投影那么六参数的单位就是米。
