利用Anusplin制作降雨数据插值

一、准备工作

降雨量数据（年、月、日），本文为月平均降雨量数据

研究区dem（鄂尔多斯市及周边地区）

下载Aunsplin软件，将splina.exe、lapgrdexe放到后续实验数据所在文件中。

Spss软件

Arcgis

Notepat++

二、数据转换

1.DEM数据

2.降雨量数据

三、CMD文件编写

1.将Anusplin软件包中的splina.exe、lapgrdexe放到后续实验数据所在文件中

2.splina_rain.cmd

3.Lapgrd_rain.cmd

4.运行两个CMD文件

四、成果展示

1.grd文件转tif

2.剪裁

3.重采样

4.去除异常值（降雨量小于0）

5.渲染

一、准备工作

降雨量数据（年、月、日），本文为月平均降雨量数据
研究区dem（鄂尔多斯市及周边地区）
下载Aunsplin软件，将splina.exe、lapgrdexe放到后续实验数据所在文件中。
Spss软件
Arcgis
Notepat++

二、数据转换

1.DEM数据

①确定经纬度范围，挑选dem压缩包

ASTGTMV003_N36E105

N36 条带号：36 E105 行编号105

鄂尔多斯市范围：N36E105-N41E112 ，扩大了dem范围，选则了N36E105-N42E113范围dem数据。

②合并鄂尔多斯市及周边地区dem

③统一投影 WGS-1984（XY范围单位为度，当单位为米时，需定义投影坐标系）

④重采样 1000m

⑤导出.dat

地图处理-ArcToolBox-转换工具-由栅格转出-栅格转ASCII

2.降雨量数据

①在arcgis中添加excel数据，添加数据→添加x、y数据，文件定义为WGS1984坐标系

②经纬度转换打开属性表，编辑，添加X、Y字段，计算几何

以度为单位（地理坐标系），可不计算；以米为单位（投影坐标系），需计算。

③提取高程（由上文dem文件提取）

Spatial Analyst分析→提取分析→多值提取至点

④导出excel

地图处理-ArcToolBox→转换工具→表转Excel箭头得到降雨量数据excel

⑤利用spss软件将降雨量数据转为.dat格式

添加变量：站点、X、Y、DEM（高程）、RAIN（降雨量数据）。需注意，站点类型为字符串，宽度一般为5，其他设置如下图：

添加站点、X、Y、高程、降雨量数据（从上文获得的降雨量数据excel中粘贴至SPSS软件）

文件→另存为：保存类型固定ASCII，编码选择本地编码，将文件保存至数据文件中。

三、CMD文件编写

1.将Anusplin软件包中的splina.exe、lapgrdexe放到后续实验数据所在文件中

2.splina_rain.cmd

①用notepad++软件创建一个cmd文件，命名（本文命名为splina_rain.cmd）

②CMD文件编写：

splina_rain    //文件名
1      //单位：米
2     //自变量个数
1     //协变量（高程）
0     //表面协变量个数
0     //表面的样条个数
104.99986111   113.99914111   0 5     //x最小值、最大值 无转换 单位：度
35.99869905   43.00013889   0 5     //y最小值、最大值 无转换 单位：度
-800 3000 1 1     //高程最小值、最大值 A转换 单位：米
0.00899928     //分辨率（本例dem分辨率为1000米）
0     //独立变量转换参数：不转换
2     //样条次数
1     //输出面个数
0     //表面权重：不设置
1     //优化方式：对所有值进行普通平滑
1     //平滑方式：GCV；2-MSE；3-固定值；4-GML
ordos_rain_202201.dat     //进行插值的降雨数据文件
500     //数据最大数量，需多余市级个数
5     //站点位数（多为5位）
(A5,F10.4,F12.4,F12.4,F8.2)     //降雨dat数据格式，在spss中设置、查看
rain11.res     //输出文件
rain11.opt
rain11.sur
rain11.lis
rain11.cov     //后续需空六行

③X、Y最小值最大值查看方法：

经纬度最小值计算方法，从上文dem.txt获得：

ncols 1000

nrows 778

xllcorner 104.999861111

yllcorner 35.998699049

cellsize 0.00899928

NODATA_value -9999

其中：ncols:行数 nrows:列数 xllcorner:X最小值 yllcorner:Y最小值 cellsize: 像元大小

X最大值计算方法：X最小值+行数×像元大小；

Y最大值计算方法：Y最小值+列数×像元大小。

104.99986111 113.99914111 0 5 //x最小值、最大值无转换单位：度

35.99869905 43.00013889 0 5 //y最小值、最大值无转换单位：度

范围可以适当选择大一点

④降雨dat数据格式，在spss中设置、查看Format。

(A5,F10.4,F12.4,F12.4,F8.2)

3.Lapgrd_rain.cmd

①用notepad++软件创建一个cmd文件，命名（本文命名为Lapgrd_rain.cmd）

②CMD文件编写：

rain11.sur     //加载文件，由上文运行splina_rain.cmd获得
1     //输出面个数
1     //对表面所有值进行计算
rain11.cov     //加载文件，由上文运行splina_rain.cmd获得
2     //误差类型，1-模型标准误差；2-预测标准误差//标准差最大值，可不填，需空行
1     //网格采样位置：像元中心
1     //第一栅格变量指标
0.10499986E+03  0.11399914E+03  0.00899928     //X范围
2     //第二栅格变量指标
0.35998699E+02  0.43000139E+02 0.00899928     //Y范围
0     //掩膜文件，本例无
2     //独立协变量数据格式：输入输出一致
ordos_dem_1000m.txt     //协变量文件（高程文件，需转换为1000m分辨率）
2     //输出的栅格格式：ARC/INFO GRID
-9999.0     //空值数据
RAIN2204.grd     //输出栅格文件名
(f10.4)     //输出文件格式
2     //误差表面格式
-9999.0     //空值数据
RAIN5202204.grd     //标准差文件，后面空6行
(f10.4)     //后续需空六行

③X、Y（纬度范围）可从splina产出的.sur文件中获得：

ANUSPLIN 4.4 UNITS: 1 TITLE: splina_rain

1 2 75 2 3 0 0 0 0

0.10499986E+03 0.11399914E+03 0 0.00000000E+00 0.00000000E+00 5

0.35998699E+02 0.43000139E+02 0 0.00000000E+00 0.00000000E+00 5

-0.80000000E+03 0.30000000E+04 1 0.00000000E+00 0.00000000E+00 1

0.89992800E-02

三、四行标蓝地方分别为X、Y范围

4.运行两个CMD文件

①打开“运行”，输入CMD，调用CMD程序

②利用cd命令，将CMD的运行文件夹转至数据文件夹。

C:\Users\PC>cd d:\rain2

③执行splina.cmd

d:\rain2>splina<splina_rain.cmd>rain202201.log

④执行lapgrd.cmd

d:\rain2>lapgrd<lapgrd_rain.cmd>rain202201lp.log

四、成果展示

1.grd文件转tif

将运行lapgrd.cmd得到的rain.grd文件拖至acegis中，另存为tif文件

2.剪裁

根据鄂尔多斯市shp，将降雨量插值.tif文件进行裁剪

3.重采样

将1公里分辨率重采样为10米分辨率

4.去除异常值（降雨量小于0）

在acrgis 的栅格计算器中输入

Con("ordos_rain_202201_se.tif"<0,0,"ordos_rain_202201_se.tif")

5.渲染