GIS使用gdal的gdal_grid实现空间插值

如何通过点数据生成面数据是地理信息领域常见问题之一，GDAL是一个好用的开源GIS项目，基于GDAL可快速实现空间插值。

1.下载gdal_grid.exe文件

在gdal的windows包里面有，如果gdal_grid.exe在cmd里面直接执行没反应很可能是缺少依赖包，双击gdal_grid.exe就知道缺少哪些包了，再按顺序找到相应的包放到同级目录即可。

我是从MapServer中拷出来的包，因此缺少：

expat.dll、freexl.dll、geos.dll、geos_c.dll、iconv.dll、libcrypto-1_1.dll、libcurl.dll、libmysql.dll、libpq.dll、libssl-1_1.dll、libxml2.dll、openjp2.dll、proj.dll、spatialite.dll、sqlite3.dll、xerces-c_3_2.dll、zlib1.dll

这几个dll文件，不同情况可能有所不同，注意版本要对应上。

2.准备点数据datatest.txt

内容如下，要按结构化的格式准备好，后续需要修改为csv格式：
id,x,y,z

0,-4033.176588,-8013.121063,229.4500001,-4100.832155,-8015.092321,229.2600002,-3340.176588,-8013.121064,234.2200003,-3638.176588,-8022.121064,234.0100004,-3638.176588,-7987.121064,233.5000005,-4451.933582,-8013.124110,228.3200006,-4458.433112,-7977.121063,227.8800007,-4707.932199,-7989.121063,226.250000

3.根据数据情况生成dem.vrt文件

dem.vrt内容如下：

<OGRVRTDataSource> <OGRVRTLayer name="datatest"> <SrcDataSource>datatest.csv</SrcDataSource> <GeometryType>wkbPoint</GeometryType> <GeometryField encoding="PointFromColumns" x="x" y="y" z="z"/> </OGRVRTLayer></OGRVRTDataSource>

4.准备执行命令：

./gdal/gdal_grid.exe -a invdist:power=2.0:smoothing=1.0 -of GTiff -ot Float64 -l datatest dem.vrt dem.tiff

按官方的命令行执行后报错：
ERROR 1: Failed to open datasource `datatest.txt'.No point geometry found on layer datatest, skipping.
此错误的主要原因是官方示例采用的csv作为点数据源，我采用的txt作为数据源，txt内肯定是没有layer这个概念的，但是csv中有sheet的概念，对应就是layer。
将datatest.txt文件后缀改为csv，名为datatest.csv，然后修改dem.vrt文件中SrcDataSource的值为新文件名就可以了。
再执行提示：
0...10...20...30...40...50...60...70...80...90...100 - done.
查看文件夹内，生成了dem.tiff文件

5.效果如下：