Python与开源GIS:访问数据集的数据

通过访问前面介绍的GDAL的一些函数,可以概括地知道影像的获取时间、处理时间、空间分辨率、 影像大小等一些信息。

但是为了对遥感影像进行处理,需要进一步能够访问遥感影像中的数据, 即影像中像元的灰度值。

GDAL提供了下面两个函数来访问是影像的数值。

  • ReadRaster() 读取图像数据(以二进制的形式)
  • ReadAsArray() 读取图像数据(以数组的形式)

下面看一下具体的用法。

    >>> help(dataset.ReadRaster)
    ReadRaster(self, xoff, yoff, xsize, ysize, buf_xsize=None, buf_ysize=None, buf_type=None, band_list=None)
    >>> help(dataset.ReadAsArray)
    ReadAsArray(self, xoff=0, yoff=0, xsize=None, ysize=None)

这是两个非常重要的函数,它们直接读取图像的数据, 从而可以对栅格数据进行分析。 可以看到两个函数的帮助中有许多的参数。 解释一下:

  • xoff,yoff:指定想要读取的部分原点位置在整张图像中距离全图原点的位置(以像元为单位)。
  • xsize,ysize: 指定要读取部分图像的矩形的长和宽(以像元为单位)。
  • buf_xsize,buf_ysize:可以在读取出一部分图像后进行缩放。那么就用这两个参数来定义缩放后图像最终的宽和高,gdal将帮你缩放到这个大小。
  • buf_type:可以对读出的数据的类型进行转换(比如原图数据类型是short,你要把它们缩小成byte)。
  • band_list:适应多波段的情况。可以指定要读取的哪几个波段。

这里简单看一下如何获取GeoTIFF文件中的数据。

    >>> dataset2 = gdal.Open("/bk/gdata/img_landsat_subset.tif")
    >>> dataset2.ReadAsArray(50,50,3,3)
    array([[[ 45,  40,  43],
            [ 51,  45,  45],
            [ 54,  50,  44]],
           [[148, 151, 168],
            [157, 153, 160],
            [170, 163, 156]],
           [[121, 119, 132],
            [130, 127, 127],
            [143, 135, 127]]], dtype=uint8)
    >>> dataset2.ReadRaster(50,50,3,3)
        '-(+3--62,\x94\x97\xa8\x9d\x99\xa0\xaa\xa3\x9cyw\x84\x82\x7f\x7f\x8f\x87\x7f'

我们就把图像左上角位于(3340, 3780),宽高都为10个像元的数据读取出来了。

这两个函数返回的结果不一样,其中ReadAsArray()读出的是numpy的数组,类型为int16; 而ReadData()读出二进制型。