关于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)表示图像横向和纵向的分辨率。图像的分辨率就是图像每个像素所能表示的面积，一般都是正方形的格网，所以也就是没两个相邻像元坐标的差值。基于这个原理，使用两个坐标进行验证。假设当前点行列号坐标为A（i，j），相邻的右侧点坐标为B（i+1，j）。分别计算A和B的横向地理坐标，并计算差值，即：

dX = XgeoB - XgeoA = [GT(0) + (i+1)*GT(1) + j*GT(2)] - [GT(0) + i*GT(1) + j*GT(2)] = GT(0) -GT(0) + (i+1)*GT(1) - i*GT(1) + j*GT(2) - j*GT(2) = (i+1)*GT(1) - i*GT(1) = GT(1)

同理可以得到 dY= GT(5)。

对于一个普通的标准图像来说（这里的标准图像是指GT(2)和GT(4)都为0），如图1所示，图像的行列号坐标为XOY，每个网格代表一个图像像素区域，i表示列号，j表示行号，淡蓝色右下角的行列坐标为（i，j），图中红色方块纵向长度为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<int>(dCol);iRow = static_cast<int>(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等分带投影，那么六参数的单位就是米。

关于GDAL计算图像坐标的几个问题相关推荐

（一）GDAL计算图像坐标与像素坐标之间的关系
1. 原理使用GDAL处理地理图像时,不可避免的会遇到一个问题,图像的地理坐标问题,因为有了这个地理坐标,地理图像才和普通图像有了最本质的区别,那么在使用GDAL时,如果处理与地理坐标相关的信息呢? ...
python dataset[trans_python gdal根据图像坐标生成矢量框（含图像坐标转地理坐标）...
要生成矢量框需要将图像坐标转换为地理坐标或者投影坐标,以下代码是生成了满足条件的1000*1000区域对应的矢量框,关键在于红色字体部分. # -*- coding: utf-8 -*- import ...
学习笔记——【python】GetGeoTransform()使用，gdal截取图像，使用GDAL进行影像投影坐标、地理坐标、图上坐标的转换
1. GetGeoTransform()使用.gdal截取图像 GetGeoTransform() GeoTransform[0],左上角横坐标(应该是投影坐标) GeoTransform[2],行旋 ...
相机计算坐标公式_相机标定后，从图像坐标（像素）向物理坐标的转换计算，...
各位大侠,本人用1组35张图片,标定相机参数.分别用MATLAB和OPENCV分别标定,标定结果基本一致. 现在,我希望从标定的一副图像中(使用这幅图像标定的相机外参),计算图像上几个红点之间的在3维 ...
python dataset[trans_科学网—Python GDAL 图像坐标，投影坐标，经纬度坐标三者映射及运行错误解决 - 吴妍潼的博文...
题记: 写该博客是因为自己经常遇到这个问题,而我发现网络上关于这方面浏览量高的一些代码竟然都有误,每次照搬都被虐得很惨.有一些同志在某些博客下方留言说代码有问题,而博主没有回应,也没有更改错误.为了自 ...
python运行报错无法映射_科学网-Python GDAL 图像坐标，投影坐标，经纬度坐标三者映射及运行错误解决-吴妍潼的博文...
题记: 写该博客是因为自己经常遇到这个问题,而我发现网络上关于这方面浏览量高的一些代码竟然都有误,每次照搬都被虐得很惨.有一些同志在某些博客下方留言说代码有问题,而博主没有回应,也没有更改错误.为了自 ...
gdal读写图像分块处理
转自赵文原文 gdal读写图像分块处理(精华版) Review: 用gdal,感觉还不如直接用C++底层函数对遥感数据进行处理.因为gdal进行太多封装,如果你仅仅只是Geotif等格式进行处理,IO ...
使用OpenCV计算图像的轮廓矩的代码
图像处理开发需求.图像处理接私活挣零花钱,请加微信/QQ 2487872782 图像处理开发资料.图像处理技术交流请加QQ群,群号 271891601 轮廓矩用来干嘛?说实话,我是没有找到相关文章专门 ...
图像坐标：我想和世界坐标谈谈(A) 【计算机视觉学习笔记--双目视觉几何框架系列】
玉米竭力用轻松具体的描述来讲述双目三维重建中的一些数学问题.希望这样的方式让大家以一个轻松的心态阅读玉米的<计算机视觉学习笔记>双目视觉数学架构系列博客.这个系列博客旨在捋顺一下已标定的双 ...

关于GDAL计算图像坐标的几个问题

关于GDAL计算图像坐标的几个问题相关推荐

最新文章

热门文章