GDAL读取图像与投影变换(Python)

安装

先安装好Anaconda和配置好的Vscode,然后安装GADL的轮子,(https://www.lfd.uci.edu/~gohlke/pythonlibs/#gdal) 在网站上面下载的轮子要匹配自己的python版本和操作系统,比如我的python是3.8版本,64位操作系统,我应该下载的轮子就是这个: 在Vscode的命令行中用下面命令切换到要安装的位置,这个位置里面放好刚下的轮子

cd path
cd E:

然后直接用pip命令即可安装 pip install GDAL‑3.3.2‑cp38‑cp38‑win_amd64.whl

一、读取图像

下面是一个类,可以放在一个新的文件中作为头文件,每次读取图像可以直接从头文件中引入类

class Grid:
    # 读取图像文件
    def read_img(self, filename):
        data = gdal.Open(filename)  # 打开文件
        im_width = data.RasterXSize  # 读取图像行数
        im_height = data.RasterYSize  # 读取图像列数
        im_geotrans = data.GetGeoTransform()
        # 仿射矩阵,左上角像素的大地坐标和像素分辨率。
        # 共有六个参数,分表代表左上角x坐标;东西方向上图像的分辨率;如果北边朝上,地图的旋转角度,0表示图像的行与x轴平行;左上角y坐标;
        # 如果北边朝上,地图的旋转角度,0表示图像的列与y轴平行;南北方向上地图的分辨率。
        im_proj = data.GetProjection()  # 地图投影信息
        im_data = data.ReadAsArray(0, 0, im_width, im_height)  # 此处读取整张图像
        # ReadAsArray(<xoff>, <yoff>, <xsize>, <ysize>)
        # 读出从(xoff,yoff)开始,大小为(xsize,ysize)的矩阵。
        del data  # 释放内存,如果不释放,在arcgis,envi中打开该图像时会显示文件被占用
        return im_proj, im_geotrans, im_data

调用上面类的语句如下:

a = Grid() #初始化类
prj, geotrans, data = a.read_img(r'path')  # 读取数据,path表示文件路径,是一个字符串,读取遥感影像的时候需要与影像头文件放在一处
print(prj, geotrans, data, data.shape, sep='\n')

我在调用了上述代码之后的运行结果如下:

我这里读的数据是一个2232行,2228列,共44个波段的数据,而上述函数将图像信息传递到了一个数组中,数组有44个元素,每个元素是一个列表,这个列表有2232个元素,表示图像的每一行,每一行也是一个数组,长度是2228,表示这一行有2228个像素。 如果想要得到第三个波段第100行第300列的像素值,可以直接对数组进行索引,在上面的语句中就是: data[2][99][299] 注意下标从0开始 上述部分参考自这篇博客:https://blog.csdn.net/weixin_44839513/article/details/103168166

下面对代码进行拆解:

  1. data = gdal.Open(filename) 打开文件,这里的data是一个数据集,包含了图像的一些地理投影信息,不能用访问数组的方式访问data
  2. im_width = data.RasterXSize # 读取图像行数 im_height = data.RasterYSize # 读取图像列数 读取图像的行列数,从这两句可以就看出来图像的行列数信息是在data里面的
  3. im_data = data.ReadAsArray(0, 0, im_width, im_height) 这句将图像转换为一个数组,上面数组先出现的变量是行数,然后是列数。

至此,不考虑投影信息,读取一张图像可以用下面的简化代码:

DataSet = gdal.Open(filename)  # 打开文件
im_width = DataSet.RasterXSize  # 读取图像行数
im_height = DataSet.RasterYSize  # 读取图像列数
im_data = DataSet.ReadAsArray(0, 0, im_width, im_height) #读取整张图像

二、投影变换

且听下回分解


本文章使用limfx的vscode插件快速发布