空间插值,这是通过离散点表征区域特征值的一种方式,空间插值具体的概念及原理请自己百度。实验中对某一时刻的数据进行插值,我们可以使用ArcMap中的工具进行手动处理。但是项目中使用的空间插值,如对每分钟/小时/天的AQI数据进行插值,就不可能手动进行处理了,这时需要考虑开发插值服务解决问题。

        空间插值属于GIS的范围,GIS中常见的解决方案包括开源和商业,其中GDAL和ArcGIS分别是各自的优秀代表,下面就对GDAL和ArcGIS空间插值进行说明。

1、GDAL空间插值(以python为例https://gdal.org/python/

(1)GDAL中进行空间插值的函数是 Grid,下面是官网对其介绍:

Grid(destName, srcDS, **kwargs)

Create raster from the scattered data.
Arguments are :
  destName --- Output dataset name
  srcDS --- a Dataset object or a filename
Keyword arguments are :
  options --- return of gdal.GridOptions(), string or array of strings
  other keywords arguments of gdal.GridOptions()
If options is provided as a gdal.GridOptions() object, other keywords are ignored. 

它的参数包括:

destName:输出数据集(一般为tif文件)的路径(包括文件名);

srcDS:输入数据集对象或文件名称,可以使用OpenEx打开geojson格式的数据,或使用Open打开shapefile或直接输入shapefile路径;

**kwargs:插值参数,可以使用 key = value 关键字参数的形式进行设置,也可以使用 option = GridOptions()的形式进行一次性设置,

(2)下面对GridOptions的内容进行说明:


GridOptions(options=[], format=None, outputType=GDT_Unknown, width=0, height=0, creationOptions=None, outputBounds=None, outputSRS=None, noData=None, algorithm=None, layers=None, SQLStatement=None, where=None, spatFilter=None, zfield=None, z_increase=None, z_multiply=None, callback=None, callback_data=None)
 
Create a GridOptions() object that can be passed to gdal.Grid()
Keyword arguments are :
  options --- can be be an array of strings, a string or let empty and filled from other keywords.
  format --- output format ("GTiff", etc...)
  outputType --- output type (gdal.GDT_Byte, etc...)
  width --- width of the output raster in pixel
  height --- height of the output raster in pixel
  creationOptions --- list of creation options
  outputBounds --- assigned output bounds: [ulx, uly, lrx, lry]
  outputSRS --- assigned output SRS
  noData --- nodata value
  algorithm --- e.g "invdist:power=2.0:smoothing=0.0:radius1=0.0:radius2=0.0:angle=0.0:max_points=0:min_points=0:nodata=0.0"
  layers --- list of layers to convert
  SQLStatement --- SQL statement to apply to the source dataset
  where --- WHERE clause to apply to source layer(s)
  spatFilter --- spatial filter as (minX, minY, maxX, maxY) bounding box
  zfield --- Identifies an attribute field on the features to be used to get a Z value from. This value overrides Z value read from feature geometry record.
  z_increase --- Addition to the attribute field on the features to be used to get a Z value from. The addition should be the same unit as Z value. The result value will be Z value + Z increase value. The default value is 0.
  z_multiply - Multiplication ratio for Z field. This can be used for shift from e.g. foot to meters or from  elevation to deep. The result value will be (Z value + Z increase value) * Z multiply value.  The default value is 1.
  callback --- callback method
  callback_data --- user data for callback

GridOptions中的参数,需要以 key=value 关键字参数的形式进行设置,基本的参数包括:

format(输出文件格式)、outputType(输出数据类型,如字节GDT_Byte、32位浮点型gdal.GDT_Float32等)、width和height(输出栅格宽度和高度)、outputBounds(插值范围)、outputSrs(输出坐标系)、noData(没有值的标记)、algorithm(插值算法,下面会做详细介绍)、SQLStatement(通过sql语句对输入数据进行过滤)、where(通过条件对参与插值的数据进行过滤)、zfield(插值字段。如aqi、co等)、z_increase和z_multiply()对插值字段进行处理的算子),callback(插值完成后的回调函数)、callback_data(回调函数的用户数据)。

(3)插值算法

最常用的是 反距离权重算法,下面对该算法的参数进行说明:

'invdist:power=3.6:smoothing=0.2:radius1=0.0:radius2=0.0:angle=0.0:max_points=0:min_points=0:nodata=0.0'

invdist表示算法名称 - 反距离权重,即距离越近权重越大,对待计算的网格影像越大;

power 距离对权重的影响程度)(数字表示指数),默认值为2, 值越大表示距离越近的点对插值的影响程指数倍增;对于该参数的设置,若点较密集且均匀分布,则值设置在2以下,若点相对较少且分布不均,则为了保障插值效果,可将参数设置在3以上;

smoothing 平滑系数,默认为0,数值越大表示越平滑,同时精度也会越低。

radius1 表示椭圆x轴方向上的半径;

radius2 表示椭圆y轴方向上的半径;

angle 表示椭圆旋转的弧度;

max_points 表示参与插值的最大点数,默认值0表示搜索到多少就是多少;

min_points 表示参与插值的最小点数,默认值0表示搜索到多少就是多少;

nodata 对nodata的填充值。

 

2、ArcGIS空间插值(以ArcPy为例说明)

下次再进行补充...

Logo

开放原子开发者工作坊旨在鼓励更多人参与开源活动,与志同道合的开发者们相互交流开发经验、分享开发心得、获取前沿技术趋势。工作坊有多种形式的开发者活动,如meetup、训练营等,主打技术交流,干货满满,真诚地邀请各位开发者共同参与!

更多推荐