1.7. GDAL python教程(6)——更多栅格数据处理函数¶
栅格数据的投影
首先要知道输入投影和输出投影的WKT(Well Known Text),可以通过GetProjection()读到,可以用SpatialReference对象创建。
用下面的语句新建栅格数据集并重新投影,投影结果输出到新数据集之中:
gdal.CreateAndReprojectImage( <source_dataset>, <output_filename>, src_wkt=<source_wkt>, dst_wkt=<output_wkt>, dst_driver=<Driver>, eResampleAlg = <GDALResampleAlg>)
下面关注一下性能问题吧
对比两种方法的执行效率:逐个像素处理,使用内建函数处理
内建的函数效率最高,比逐个像素循环要快很多很多。
1.7.1. 滑动窗口¶
按坐标逐个循环¶
data = inBand.ReadAsArray(0, 0, cols, rows).astype(Numeric.Int)
outData = Numeric.zeros((rows, cols), Numeric.Float)
for i in range(1, rows-1): # skipping first & last
for j in range(1, cols-1):
outData[i,j] = (data[i-1,j-1] + data[i-1,j] + data[i-1,j+1] + data[i,j-1] + data[i,j] + data[i,j+1] + data[i+1,j-1] + data[i+1,j] + data[i+1,j+1]) / 9.0
数组切片¶
data = inBand.ReadAsArray(0, 0, cols, rows).astype(Numeric.Int))
outData = Numeric.zeros((rows, cols), Numeric.Int)
outData[1:rows-1,1:cols-1] = (data[0:rows-2,0:cols-2] + data[0:rows-2,1:cols-1] + data[0:rows-2,2:cols] + data[1:rows-1, 0:cols-2] + data[1:rows-1,1:cols-1] + data[1:rows-1,2:cols] + data[2:rows,0:cols-2] + data[2:rows,1:cols-1] + data[2:rows,2:cols]) / 9
用数组切片的方法相对而言效率较高。
1.7.2. 过滤器filter¶
或者叫滤波器吧。低通滤波Low-pass filter,用于平滑数据。高通滤波High-pass filter,用于边缘增强。上面滑动窗口的两段程序就是低通滤波。