2_空间数据文件访问.pdf

Python 空间数据处理（Geoprocessing With Python），国科大 （宋现锋） 第二章 空间数据文件访问（GDAL/OGR 包） 1.图像文件的读写 栅格数据模型 （像元，如：正方形）实质上是对连续空间分布现象的离散化，采用了一 个二维矩阵、使用行号与列号来确定每个象元点 cell）的空间位置。其数据模型主要包括 两部分 –数据矩阵和大地坐标。 栅格矩阵：像元值用于表示地物某种属性的数值，像元大小表示栅格数据的空间分辨率， 波段则表示地物属性有多类，每类属性用一个波段表示。 大地坐标：空间数据参考表示地图投影信息，仿 变换表示行列坐标到平面坐标的映 。 #读GeoTiff 文件示例1 from osgeo import gdal import numpy as np dataset=gdal.Open(fdem.tif) #打开文件 im_width = dataset.RasterXSize #栅格矩阵的列数 im_height = dataset.RasterYSize #栅格矩阵的行数 im_bands = dataset.RasterCount #波段数 im_geotrans = dataset.GetGeoTransform() #仿 矩阵，左上角像素的大地坐标和像素分辨率 im_proj = dataset.GetProj ection() #地图投影信息，字符串表示 im_data = dataset.ReadAsArray(0,0,im_width,im_height) type(im_data),im_data.shape del dataset #显示栅格数据 #numpy 的用法，详见第四章科学计算 print im_data[0, 10:15,20:25] #查看波段0,10~ 14 行和20~25 列的数据 print im_geotrans #(409294.88697, 27.37648, 0.0, 4423871.08338, 0.0, -27.37648) Python 空间数据处理（Geoprocessing With Python），国科大 （宋现锋） #写GeoTiff 文件示例2 #gdal 数据类型 #gdal.GDT_Byte, #gdal .GDT_UInt 16, gdal.GDT_Int 16, gdal.GDT_UInt32, gdal.GDT_Int32, # gdal.GDT_Float32, gdal.GDT_Float64 #判断栅格数据的数据类型 if int8 in im_data.dtype.name: datatype = gdal.GDT_Byte elif int 16 in im_data.dtype.name: datatype = gdal.GDT_UInt 16 else: datatype = gdal.GDT_Float32 #判读数组维数 if len(data.shape) == 3 : im_bands, im_height, im_width = im_data.shape else: im_bands, (im_height, im_width) = 1, im_data.shape #创建文件 driver = gdal.GetDriverByName(GTiff) dataset = driver.Create(fdem_new.tif, im_width, im_height, im_bands, datatype) dataset.SetGeoTransform(im_geotrans) #写入仿 变换参数 dataset.SetProj ection(im_proj) #写入投影 if im_bands == 1 : dataset.GetRasterBand(1).WriteArray(im_data) #写入数组数据 else: for i in range(im_bands): dataset.GetRasterBand(i+ 1).WriteArray(im_data[i]) del dataset Python 空间数据处理（Geoprocessing