了解GDAL的图像处理/Python

GDAL是一个操做各类栅格地理数据格式的库。包括读取、写入、转换、处理各类栅格数据格式（有些特定的格式对一些操做如写入等不支持）。它使用了一个单一的抽象数据模型就支持了大多数的栅格数据（GIS对栅格，矢量，3D数据模型的抽象能力实在使人叹服）。固然除了栅格操做，这个库还同时包括了操做矢量数据的另外一个有名的库ogr（ogr这个库另外介绍），这样这个库就同时具有了操做栅格和矢量数据的能力，买一送一，这么合算的买卖为何不作。最最最重要的是这个库是跨平台的，开源的！现在这个库对各类数据格式的支持强大到使人啧啧的地步了。若是你对他的强大有什么怀疑的话，看看这里一大串的GDAL所支持格式清单，吓到了吧！再看看它的主页最后那些使用了它做为底层数据处理的软件列表吧！其中你能够不知道GRASS，你也能够不知道Quantum GIS (QGIS)，可是你总该知道Google Earth吧！
你即便不玩GIS，这个库也是满有用的。首先，哪一个库支持这么多栅格（图片）格式，哪一个库在C/C++/python/ruby/VB/java/C#(这个暂时不彻底支持)下都能用，并且都同样用？退一步讲，3S软件又不必定要用在3S下（不少医学影像就是用PCI软件来处理的）。再退一步，你的生活即便和3S一点关系都没有，栅格数据又不仅仅只有GIS下才用到。你大可用这个库来读取jpg，gif，tif，xpm等格式。并且对各类格式支持得不是通常的好，很大一部分非标准格式照样支持得很是好。我曾经在java下玩过jai，以及一系列jai的扩展库，一些图像格式在不少图片浏览器中均可以正确读取（有的甚至不是非标准格式），用jai死活就读不出来！
这个库的python版和其余的python库结合的很好。最直接、明显的支持是使用Numeric库来进行数据读取和操做。各类矩阵魔术能够发挥得淋漓尽致（图像其实就是矩阵）。并且按个人观点，python对矩阵的操做比其余的语言有明显的优点。写出来的东西比其余语言写出来的短小的多，并且好看得多。而且python的弱类型在处理栅格数据格式类型的时候代码量比强类型的语言少了数倍（不用double，byte，short等等分开处理，这简直就是先天上的优点）。因此我就喜欢用python作图像的处理。因此就连GIS界的微软ESRI也直接在ARCGIS9中用python来做栅格数据的导入导出。一句话，真是太方便啦！

2. 安装
2.1. windows下的安装
官方安装文档在这里。下面是我本身的实践步骤：
先去http://www.gdal.org/dl/下一个版本，解压。 打开控制台，输入： “D:\Program Files\Microsoft Visual Studio .NET 2003\Vc7\bin\vcvars32.bat" 注册vc的编译环境。
打 开gdal文件夹下的nmake.opt修改GDAL_HOME = "C:\warmerda\bld"把路径改到须要把gdal安装的地方。不改也能够。这里须要添加python支持，因此修改PY_INST_DIR = (GDAL_HOME)\pymod把路径改为python下的Lib\site-packages文件夹下。PYDIR = "C:\Software\Python24" 改为python的安装路径。 下面的参数爱改什么就把前面的#删除（要看您有没有那些库的源码），注意一下路径就能够了。我是都没改。 后面就依次运行       Toggle line numbers       1 nmake /f makefile.vc       2 nmake /f makefile.vc install       3 nmake /f makefile.vc devinstall       最后最后，还要去GDAL_HOME目录下的bin文件夹下把gdal13.dll （也有多是gdal12.dll）copy到PY_INST_DIR路径下       到此处就完成安装gdal（python）的工做。       最后须要注意一下，gdal在.net2005下只能顺利编译1.2，1.3以上版本不能顺利编译，有一个地方指针转换出错。多是2005的编译器比以往的严厉一点吧。另外,安装了QGIS,对编译也有一些影响,主要是proj库的冲突,致使一个找不到"d:/program.obj"文件的错误,若是你有静态编译过proj,那么你能够打开nmake.opt修改有关proj的设置,若是搞不定,就卸载QGIS,而后编译,编译后再安装QGIS.呵呵,还好QGIS的体积没有ArcGIS那么可怕.       2.2. linux下的安装       linux下的安装就更简单了。直接       Toggle line numbers       1 ./configure       2 make       3 su       4 make install       5 ldconfig       就ok(默认就已经支持python)。固然在第一步的时候须要看看是否依赖的库都安装了。若是缺乏，就去安装一个。若是对configure的条件不理解，就用./configure --help看看具体状况。       2.3. 安装其余驱动       这里讲一个安装hdf4的驱动的例子（默认状况下gdal是不安装hdf4的），其余驱动应该和这个也差不了多少吧，能够做为其余的参考。完整步骤以下：       在windows下的安装：       从ftp://ftp.ncsa.uiuc.edu/HDF/HDF/HDF%5FCurrent/bin/windows/下载42r1-win.ZIP，解压。       编辑gdal根目录下的nmake.opt，找到“# Uncomment the following and update to enable NCSA HDF Release 4 support.”这一行       把下面两行前面的#去掉，而后改为：       HDF4_DIR = D:\warmerda\42r1-win\release       #HDF4_LIB = /LIBPATH:<?XML:NAMESPACE PREFIX = "[default] http://www.w3.org/1998/Math/MathML" NS = "http://www.w3.org/1998/Math/MathML" />(GDAL_HOME)\pymod把路径改为python下的Lib\site-packages文件夹下。PYDIR = "C:\Software\Python24" 改为python的安装路径。 下面的参数爱改什么就把前面的#删除（要看您有没有那些库的源码），注意一下路径就能够了。我是都没改。 后面就依次运行       Toggle line numbers       1 nmake /f makefile.vc       2 nmake /f makefile.vc install       3 nmake /f makefile.vc devinstall       最后最后，还要去GDAL_HOME目录下的bin文件夹下把gdal13.dll （也有多是gdal12.dll）copy到PY_INST_DIR路径下       到此处就完成安装gdal（python）的工做。       最后须要注意一下，gdal在.net2005下只能顺利编译1.2，1.3以上版本不能顺利编译，有一个地方指针转换出错。多是2005的编译器比以往的严厉一点吧。另外,安装了QGIS,对编译也有一些影响,主要是proj库的冲突,致使一个找不到"d:/program.obj"文件的错误,若是你有静态编译过proj,那么你能够打开nmake.opt修改有关proj的设置,若是搞不定,就卸载QGIS,而后编译,编译后再安装QGIS.呵呵,还好QGIS的体积没有ArcGIS那么可怕.       2.2. linux下的安装       linux下的安装就更简单了。直接       Toggle line numbers       1 ./configure       2 make       3 su       4 make install       5 ldconfig       就ok(默认就已经支持python)。固然在第一步的时候须要看看是否依赖的库都安装了。若是缺乏，就去安装一个。若是对configure的条件不理解，就用./configure --help看看具体状况。       2.3. 安装其余驱动       这里讲一个安装hdf4的驱动的例子（默认状况下gdal是不安装hdf4的），其余驱动应该和这个也差不了多少吧，能够做为其余的参考。完整步骤以下：       在windows下的安装：       从ftp://ftp.ncsa.uiuc.edu/HDF/HDF/HDF%5FCurrent/bin/windows/下载42r1-win.ZIP，解压。       编辑gdal根目录下的nmake.opt，找到“# Uncomment the following and update to enable NCSA HDF Release 4 support.”这一行       把下面两行前面的#去掉，而后改为：       HDF4_DIR = D:\warmerda\42r1-win\release       #HDF4_LIB = /LIBPATH:(HDF4_DIR)\lib hd421m.lib
HDF4_LIB = (HDF4DIR)\dll\hd421m.lib(HDF4DIR)\dll\hd421m.lib(HDF4_DIR)\dll\hm421m.lib \
(HDF4DIR)\lib\hd421.lib(HDF4DIR)\lib\hd421.lib(HDF4_DIR)\lib\hm421.lib
用HDF4_LIB=/LIBPATH：这种形式彷佛能够创建gdal的库，可是往下编译会出错。并且要把(HDF4DIR)\dll和(HDF4DIR)\dll和(HDF4_DIR)\lib拷贝到同一个目录下，否则会提示找不到库
你也能够试一试在D:\Program Files\Microsoft Visual Studio .NET 2003\Common7\Tools\vsvars32.bat文件中添加HDF4_LIB路径到“@set LIB=”这行的末尾（不要忘记;的分割符）。
而后找一下"INC="这行,把 -I(HDF4DIR)\include加到下一行的末尾（应该也能够在vsvars32.bat中添加路径，不过要重启命令行）。而后编译吧！祝你好运。注意：上面的HDF4DIR是我本机的路径，你要根据你本身的路径进行设置（想起个人一个老师说过的话：“抄人家的做业能够，不要连名字也一块儿抄走啊”），下面的(HDF4DIR)\include加到下一行的末尾（应该也能够在vsvars32.bat中添加路径，不过要重启命令行）。而后编译吧！祝你好运。注意：上面的HDF4DIR是我本机的路径，你要根据你本身的路径进行设置（想起个人一个老师说过的话：“抄人家的做业能够，不要连名字也一块儿抄走啊”），下面的(HDF4_DIR)能够不用改，那个是变量，会自动替代HDF4_DIR 路径。
编译成功后，要HDF4能运行，还须要两个库，一个是zlib,一个是szip，能够到下面两个连接去下载一个
ftp://hdf.ncsa.uiuc.edu/lib-external/zlib/1.2/bin
ftp://hdf.ncsa.uiuc.edu/lib-external/szip/2.0/bin/windows
把这两个库下载后解压，而后设置PATH系统变量，使得它们在默认状态下也能够被动态连接成功 。
在Linux下比在Windows下简单：
只要用./configure --help察看一下打开HDF4的编译开关（包括库路径，头文件路径等，看清楚），而后在./configure 后面加上那个开关以及hdf4的安装路径后就能够了。在configure后gdal会提示是否支持HDF4。 编译后也要把zlib和szip的动态连接库设置好 。
到此你已经能够用C/C++来操做gdal读写hdf4的格式了！
最后，为了让Python可以支持hdf的读写，别忘了把从新生成安装后的pymod目录下的内容，还有gdal13，还有那两个hdf的库，还有zlib，szip的库拷贝到Python的Lib\site-packages目录下 。
2.4. 下载
若是你实在玩不转，能够在这里下载已经编译好的gdal1.3.2程序库 以及其依赖的其余库，其中包括hdf4，hdf5支持，以及proj，geos插件。注意，这里的geos是静态连接的，注意版权（geos是LGPL的license）。hdf4和hdf5用的是release版本。这里是个人nmake配置文件，你能够对照你的实际状况参考一下。
3. 快速开始
其实在主站的教程里已经有python的示例了。可是咱们仍是按照本身的思路来开始吧。
第一步就是打开一个数据集。对于“数据集”这个名词你们可能不会太习惯，可是对于通常的格式来讲，一个“数据集”就是一个文件，好比一个gif文件就是一个以gif为扩展名的文件。可是对于众多RS数据来讲，一个数据集包含的绝对不只仅是一个文件。对于不少RS数据，他们把一张图像分红数个图像文件，而后放在一个文件夹中，用一些额外的文件来组织它们之间的关系，造成一个“数据集”。若是你不理解，那么就算了，当成jpg或者gif文件好了。
下面咱们打开一个tiff文件(GeoTIFF)。这个文件是我从GRASS的示例数据spearfish中导出的一个同名影像数据。
Toggle line numbers
1 >>> import gdal
2 >>> dataset = gdal.Open("j:/gisdata/gtif/spot.tif")
3 >>> dir(dataset)
4 ['AddBand', 'AdviseRead', 'BuildOverviews', 'FlushCache', 'GetDescription', 'Get
5 Driver', 'GetGCPCount', 'GetGCPProjection', 'GetGCPs', 'GetGeoTransform', 'GetMe
6 tadata', 'GetProjection', 'GetProjectionRef', 'GetRasterBand', 'GetSubDatasets',
7 'RasterCount', 'RasterXSize', 'RasterYSize', 'ReadAsArray', 'ReadRaster', 'Refr
8 eshBandInfo', 'SetDescription', 'SetGCPs', 'SetGeoTransform', 'SetMetadata', 'Se
9 tProjection', 'WriteRaster', '__del__', '__doc__', '__init__', '__module__', '_b
10 and', '_o']
11 >>>
这样咱们就打开了这个文件。而且咱们能够看到能够供咱们调用的函数们(更具体的API列表能够看这里)。如今咱们不作修改，不作添加，因此只要带有Set开头的函数以及有Write开头的函数咱们暂时都无论。由于RS影像必然要和地理上的位置挂上钩，才能把图像正确铺展到一个坐标系中。其中的信息和对应关系有点复杂，不适合在快速开始中介绍，咱们暂时也先无论。这里须要注意的就是几个函数。
GetDescription 得到栅格的描述信息。
Toggle line numbers
1 >>> dataset.GetDescription()
2 'j:/gisdata/gtif/spot.tif'
3 >>>
看来这里的图像描述是图像的路径名，可是这是和各类不一样数据集相关的，不一样数据集可能有不一样的描述。这要看读取驱动的实现做者的高兴了。
RasterCount 得到栅格数据集的波段数。
GetRasterBand 得到栅格数据集的波段。
Toggle line numbers
1 >>> dataset.RasterCount
2 1
3 >>> band = dataset.GetRasterBand(1)
4 >>>
Band这个词能够翻译成“波段。
这里咱们获取了第一个波段(红色值组成的表)。注意！这里的波段获取和一般的C数组获取不同，开始是1不是0。获取了波段，咱们就能够在下面的操做中读取这个波段的全部数值。
RasterXSize 图像的宽度(X方向上的像素个数)
RasterYSize 图像的高度(Y方向上的像素个数)
Toggle line numbers
1 >>> dataset.RasterXSize
2 950
3 >>> dataset.RasterYSize
4 700
5 >>>
能够看出咱们的图像大小是950*700。仍是很小的一张图。
ReadRaster 读取图像数据(以二进制的形式)
ReadAsArray 读取图像数据(以数组的形式)
Toggle line numbers
1 >>> help(dataset.ReadRaster)
2 Help on method ReadRaster in module gdal:
3 ReadRaster(self, xoff, yoff, xsize, ysize, buf_xsize=None, buf_ysize=None, buf_t
4 ype=None, band_list=None) method of gdal.Dataset instance
5 >>> help(dataset.ReadAsArray)
6 Help on method ReadAsArray in module gdal:
7 ReadAsArray(self, xoff=0, yoff=0, xsize=None, ysize=None) method of gdal.Dataset
8 instance
9 >>>
这两个函数很重要，它们直接读取图像的数据，能够看到两个函数的帮助中有一大溜的参数。解释一下： xoff,yoff,xsize,ysize 你可能不想读取整张图像。只想读取其中的一部分。那么就用xoff，yoff指定想要读取的部分原点位置在整张图像中距离全图原点的位置。用xsize和ysize指定要读取部分图像的矩形大小。
buf_xsize buf_ysize 你能够在读取出一部分图像后进行缩放。那么就用这两个参数来定义缩放后图像最终的宽和高，gdal将帮你缩放到这个大小。
buf_type 若是你要读取的图像的数据类型不是你想要的(好比原图数据类型是short，你要把它们缩小成byte)，就能够设置它。
band_list 这就适应上面多波段的状况。你能够指定读取的波段序列。要哪几个波段，不要哪几个波段，你说了算。
举个例子吧：
Toggle line numbers
1 >>> dataset.ReadAsArray(230,270,10,10)
2 array([[255, 255, 255, 232, 232, 255, 255, 255, 255, 222],
3 [255, 255, 255, 255, 255, 255, 210, 110, 11, 122],
4 [255, 255, 255, 255, 255, 255, 210, 255, 11, 243],
5 [201, 255, 255, 255, 255, 200, 200, 110, 122, 243],
6 [111, 211, 255, 201, 255, 255, 100, 11, 132, 243],
7 [255, 100, 100, 100, 110, 100, 110, 111, 122, 243],
8 [255, 255, 255, 255, 255, 255, 122, 222, 255, 255],
9 [255, 255, 255, 255, 255, 255, 243, 243, 255, 255],
10 [255, 255, 255, 255, 255, 255, 255, 255, 255, 255],
11 [255, 255, 255, 255, 255, 255, 255, 255, 255, 255]],'b')
12 >>> dataset.ReadRaster(230,270,10,10)
13 '\xff\xff\xff\xe8\xe8\xff\xff\xff\xff\xde\xff\xff\xff\xff\xff\xff\xd2n\x0bz\xff\
14 xff\xff\xff\xff\xff\xd2\xff\x0b\xf3\xc9\xff\xff\xff\xff\xc8\xc8nz\xf3o\xd3\xff\x
15 c9\xff\xffd\x0b\x84\xf3\xffdddndnoz\xf3\xff\xff\xff\xff\xff\xffz\xde\xff\xff\xff
16 \xff\xff\xff\xff\xff\xf3\xf3\xff\xff\xff\xff\xff\xff\xff\xff\xff\xff\xff\xff\xff
17 \xff\xff\xff\xff\xff\xff\xff\xff\xff'
18 >>>
咱们就把图像中位于230,270，宽度10高度10的数据读取出来了。
咱们看完了数据集的主要函数。彷佛已经够用了。的确，若是只是为了显示图像，这些的确已经够了。可是若是须要更多信息，咱们就不得不进入波段操做数据(实际上咱们大多数时候都须要进入band获取信息)。下面咱们如今来看看刚才读取出来的那个band有些什么东西能够供咱们操做的(具体的API列表看这里)。
Toggle line numbers
1 >>> dir(band)
2 ['AdviseRead', 'Checksum', 'ComputeBandStats', 'ComputeRasterMinMax', 'DataType'
3 , 'Fill', 'FlushCache', 'GetDefaultHistogram', 'GetDescription', 'GetHistogram',
4 'GetMaximum', 'GetMetadata', 'GetMinimum', 'GetNoDataValue', 'GetOffset', 'GetO
5 verview', 'GetOverviewCount', 'GetRasterColorInterpretation', 'GetRasterColorTab
6 le', 'GetScale', 'GetStatistics', 'ReadAsArray', 'ReadRaster', 'SetDefaultHistog
7 ram', 'SetDescription', 'SetMetadata', 'SetNoDataValue', 'SetRasterColorInterpre
8 tation', 'SetRasterColorTable', 'WriteArray', 'WriteRaster', 'XSize', 'YSize', '
9 __doc__', '__init__', '__module__', '_o']
10 >>>
挑几个有用的吧。
Toggle line numbers
1 >>> band.XSize
2 950
3 >>> band.YSize
4 700
5 >>> band.DataType
6 1
7 >>>
不用解释了吧，波段图像的宽和高(象元为单位)。DataType，图像中实际数值的数据类型。具体数据类型定义在gdalconst模块里。使用的时候用import gdalconst引入。
Toggle line numbers
1 >>> import gdalconst
2 >>> dir(gdalconst)
3 ['CE_Debug', 'CE_Failure', 'CE_Fatal', 'CE_None', 'CE_Warning', 'CPLES_Backslash
4 Quotable', 'CPLES_CSV', 'CPLES_SQL', 'CPLES_URL', 'CPLES_XML', 'CPLE_AppDefined'
5 , 'CPLE_AssertionFailed', 'CPLE_FileIO', 'CPLE_IllegalArg', 'CPLE_NoWriteAccess'
6 , 'CPLE_None', 'CPLE_NotSupported', 'CPLE_OpenFailed', 'CPLE_OutOfMemory', 'CPLE
7 _UserInterrupt', 'CXT_Attribute', 'CXT_Comment', 'CXT_Element', 'CXT_Literal', '
8 CXT_Text', 'DCAP_CREATE', 'DCAP_CREATECOPY', 'DMD_CREATIONDATATYPES', 'DMD_CREAT
9 IONOPTIONLIST', 'DMD_EXTENSION', 'DMD_HELPTOPIC', 'DMD_LONGNAME', 'DMD_MIMETYPE'
10 , 'GA_ReadOnly', 'GA_Update', 'GCI_AlphaBand', 'GCI_BlackBand', 'GCI_BlueBand',
11 'GCI_CyanBand', 'GCI_GrayIndex', 'GCI_GreenBand', 'GCI_HueBand', 'GCI_LightnessB
12 and', 'GCI_MagentaBand', 'GCI_PaletteIndex', 'GCI_RedBand', 'GCI_SaturationBand'
13 , 'GCI_Undefined', 'GCI_YellowBand', 'GDT_Byte', 'GDT_CFloat32', 'GDT_CFloat64',
14 'GDT_CInt16', 'GDT_CInt32', 'GDT_Float32', 'GDT_Float64', 'GDT_Int16', 'GDT_Int
15 32', 'GDT_TypeCount', 'GDT_UInt16', 'GDT_UInt32', 'GDT_Unknown', 'GF_Read', 'GF_
16 Write', 'GPI_CMYK', 'GPI_Gray', 'GPI_HLS', 'GPI_RGB', 'GRA_Bilinear', 'GRA_Cubic
17 ', 'GRA_CubicSpline', 'GRA_NearestNeighbour', '__builtins__', '__doc__', '__file
18 __', '__name__']
19 >>>
那些GDT开头的就是数值数据类型。
Toggle line numbers
1 >>> band.GetNoDataValue()
2 65535.0
3 >>> band.GetMaximum()
4 >>> band.GetMinimum()
5 >>> band.ComputeRasterMinMax()
6 (1.0, 255.0)
7 >>>
Maximum 是表示在本波段数值中最大的值，Minimum固然就是表示本波段中最小的值啦。咱们能够看到在一开始这两个都没有值。由于对于文件格式不会有固有的最大最小值。因此咱们经过函数ComputeRasterMinMax计算获得了。注意！这里的最大最小值不包括“无心义值”！也就是上面显示的NoDataValue。须要解释一下“无心义值”。不要觉得0或者255在任何状况下都无心义。在不少状况下0，255须要和其余值同样表示一个实际意义。虽然可能它最终会被显示得和黑色同样。而一些位置上的点要表示的意思是“什么也不是”，它在那个位置上只是为了占一个位置，使得总体图像看起来像个矩形而已。在作实际应用的时候两种值的处理将会彻底不同。因此须要设置无心义值，来和其余的值区别开来。而用ComputeRasterMinMax算出的最大最小值，是排除了无心义值后计算出来的最大最小值。
Toggle line numbers
1 >>> band.GetRasterColorInterpretation()
2 2
3 >>> gdalconst.GCI_PaletteIndex
4 2
5 >>> colormap = band.GetRasterColorTable()
6 >>> dir(colormap)
7 ['Clone', 'GetColorEntry', 'GetColorEntryAsRGB', 'GetCount', 'GetPaletteInterpre
8 tation', 'SetColorEntry', '__del__', '__doc__', '__init__', '__module__', '__str
9 __', '_o', 'own_o', 'serialize']
10 >>> colormap.GetCount()
11 256
12 >>> colormap.GetPaletteInterpretation()
13 1
14 >>> gdalconst.GPI_RGB
15 1
16 >>> for i in range(colormap.GetCount()):
17 ... print colormap.GetColorEntry(i),
18 ...
19 (0, 0, 0, 255) (0, 0, 28, 255) (0, 0, 56, 255) (0, 0, 85, 255) (0, 0, 113, 255)
20 (0, 0, 142, 255) (0, 0, 170, 255) (0, 0, 199, 255) (0, 0, 227, 255) (0, 0, 255,
21 255) (0, 28, 0, 255) (0, 28, 28, 255) (0, 28, 56, 255) (0, 28, 85, 255) (0, 28,
22 113, 255) (0, 28, 142, 255) (0, 28, 170, 255) (0, 28, 199, 255) (0, 28, 227, 255
23 ) (0, 28, 255, 255) (0, 56, 0, 255) (0, 56, 28, 255) (0, 56, 56, 255) (0, 56, 85
24 , 255) (0, 56, 113, 255) (0, 56, 142, 255) (0, 56, 170, 255) (0, 56, 199, 255) (
25 0, 56, 227, 255) (0, 56, 255, 255) (0, 85, 0, 255) (0, 85, 28, 255) (0, 85, 56,
26 255) (0, 85, 85, 255) (0, 85, 113, 255) (0, 85, 142, 255) (0, 85, 170, 255) (0,
27 85, 199, 255) (0, 85, 227, 255) (0, 85, 255, 255) (0, 113, 0, 255) (0, 113, 28,
28 255) (0, 113, 56, 255) (0, 113, 85, 255) (0, 113, 113, 255) (0, 113, 142, 255) (
29 0, 113, 170, 255) (0, 113, 199, 255) (0, 113, 227, 255) (0, 113, 255, 255) (0, 1
30 42, 0, 255) (0, 142, 28, 255) (0, 142, 56, 255) (0, 142, 85, 255) (0, 142, 113,
31 255) (0, 142, 142, 255) (0, 142, 170, 255) (0, 142, 199, 255) (0, 142, 227, 255)
32 (0, 142, 255, 255) (0, 170, 0, 255) (0, 170, 28, 255) (0, 170, 56, 255) (0, 170
33 , 85, 255) (0, 170, 113, 255) (0, 170, 142, 255) (0, 170, 170, 255) (0, 170, 199
34 , 255) (0, 170, 227, 255) (0, 170, 255, 255) (0, 199, 0, 255) (0, 199, 28, 255)
35 (0, 199, 56, 255) (0, 199, 85, 255) (0, 199, 113, 255) (0, 199, 142, 255) (0, 19
36 9, 170, 255) (0, 199, 199, 255) (0, 199, 227, 255) (0, 199, 255, 255) (0, 227, 0
37 , 255) (0, 227, 28, 255) (0, 227, 56, 255) (0, 227, 85, 255) (0, 227, 113, 255)
38 (0, 227, 142, 255) (0, 227, 170, 255) (0, 227, 199, 255) (0, 227, 227, 255) (0,
39 227, 255, 255) (0, 255, 0, 255) (0, 255, 28, 255) (0, 255, 56, 255) (0, 255, 85,
40 255) (0, 255, 113, 255) (0, 255, 142, 255) (0, 255, 170, 255) (0, 255, 199, 255
41 ) (0, 255, 227, 255) (0, 255, 255, 255) (28, 0, 0, 255) (28, 0, 28, 255) (28, 0,
42 56, 255) (28, 0, 85, 255) (28, 0, 113, 255) (28, 0, 142, 255) (28, 0, 170, 255)
43 (28, 0, 199, 255) (28, 0, 227, 255) (28, 0, 255, 255) (28, 28, 0, 255) (28, 28,
44 28, 255) (28, 28, 56, 255) (28, 28, 85, 255) (28, 28, 113, 255) (28, 28, 142, 2
45 55) (28, 28, 170, 255) (28, 28, 199, 255) (28, 28, 227, 255) (28, 28, 255, 255)
46 (28, 56, 0, 255) (28, 56, 28, 255) (28, 56, 56, 255) (28, 56, 85, 255) (28, 56,
47 113, 255) (28, 56, 142, 255) (28, 56, 170, 255) (28, 56, 199, 255) (28, 56, 227,
48 255) (28, 56, 255, 255) (28, 85, 0, 255) (28, 85, 28, 255) (28, 85, 56, 255) (2
49 8, 85, 85, 255) (28, 85, 113, 255) (28, 85, 142, 255) (28, 85, 170, 255) (28, 85
50>>>
经过GetRasterColorInterpretation，咱们知道咱们的图像是一个颜色表索引的图像而不是纯粹的黑白灰度图像(PaletteIndex，其余的颜色模型，能够察看gdalconst模块中GCI打头的枚举值)。这意味着咱们读出的数据有可能不是真实的数据。这些数据只是一个个实际数据的索引。真实数据存储在另外一个表中。咱们经过ReadRaster读出的数据值只是对应到这个表的一个索引而已。咱们须要经过读出这些数据，并在真实数据表中找出真实数据，从新组织成一个RGB表才能用来绘制。若是咱们不通过对应，咱们绘制出来的东西可能什么东西都不是。
用GetRasterColorTable得到了颜色表，经过GetPaletteInterpretation咱们知道咱们得到的颜色表是一个RGB颜色表。GDAL支持多种颜色表，具体能够参考gdalconst模块中GPI打头的枚举值。而后咱们能够经过GetCount得到颜色的数量。经过GetColorEntry得到颜色表中的值。这里的颜色值都是一个4值的元组。里面有意义的只有前三个(若是颜色模型是GPI_RGB, GPI_HLS，都使用前3个，若是采用GPI_CMYK，则4个值都有意义了)。
Toggle line numbers
1 >>> help(band.ReadAsArray)
2 Help on method ReadAsArray in module gdal:
3 ReadAsArray(self, xoff=0, yoff=0, win_xsize=None, win_ysize=None, buf_xsize=None
4 , buf_ysize=None, buf_obj=None) method of gdal.Band instance
5 >>> help(band.ReadRaster)
6 Help on method ReadRaster in module gdal:
7 ReadRaster(self, xoff, yoff, xsize, ysize, buf_xsize=None, buf_ysize=None, buf_t
8 ype=None) method of gdal.Band instance
9 >>>
显然，band里的ReadAsArray参数显然比dataset里面的要好用，而ReadRaster则差很少。可是ReadAsArray读出的是数组，能够用Numeric模块进行矩阵魔法。ReadRaster读出的是二进制，虽然能够直接绘制，可是对于一些绘图API来讲，对[[RRR...][GGG...][BBB...]]表的处理明显不如[[RGB][RGB]...]，有的甚至不支持。虽然能够用struct.unpack来拆封，可效率上就差不少(并且拆封出来仍是数组)。数组在矩阵魔法的控制之下则会显得十分方便快捷，最后用tostring直接转化称为二进制绘制，速度也至关快。
好了，快速开始已经使咱们能够初步看清楚了gdal中图像的组织。下面用一句话总结一下：波段组成图像，波段指挥颜色。