一种基于热红外成像技术估算土壤蒸发强度的方法

1.本发明属于农田水文与水环境领域，涉及土壤蒸发强度的估算，具体来说是一种简便、快速且精度较高的基于热红外成像技术估算土壤蒸发强度的方法。


背景技术：

2.土壤蒸发是指水分从土壤表面以气态的形式向大气逸散的过程。土壤蒸发过程受气象条件、土壤水分状态、地下水位、土壤理化特性等多种因素的影响，涉及陆地和大气之间的质量和热量交换，并直接影响地下水资源和地上水资源管理，是全球水文循环的重要过程之一。
3.在我国西北干旱灌区，由于大量引黄灌溉，导致地下水埋深普遍较浅，而强烈的蒸发过程导致地下水大量流失，造成了水资源的大量浪费。盐随水去，大量的盐分逐渐积累在地表，使土壤盐碱化更加严重，导致灌区大面积减产，并对当地生态环境造成了不可估量的影响。因此，土壤蒸发的精确监测及估算对灌区水资源高效利用和土壤盐渍化控制有着十分重要的作用。
4.当前的土壤蒸发测量方法主要有直接测量法、遥感反演法和水量平衡法三种。直接测量法多使用微型蒸发器进行测量，该方法具有简便、灵活等特点，但频繁更换土体难以确保测量的代表性，且存在以点代面的局限性，不适用于空间异质性较高的地区。遥感反演法是一种较为成熟的区域尺度的蒸散发估算方法，但不同模型之间反演的土壤蒸发量存在较大差异，在土壤蒸发的连续性监测和精确量化分析方面还存在一定的差距。水量平衡法则由于各平衡要素的精确获取存在一定的难度、所需参数多、计算过程相对较为复杂，因此也难以实现大范围的推广应用。
5.因此，有必要提出一种快速准确估算农田土壤蒸发强度的方法，以便指导农田灌溉，提高水分利用效率及降低土壤盐渍化风险。


技术实现要素：

6.针对灌区土壤蒸发的测量存在精度低，代表性差以及计算过程复杂等问题，本发明的目的是提供一种基于热红外成像技术估算土壤蒸发强度的方法，该方法选用高分辨率的热红外成像仪作为主要工具，并依据该仪器所拍摄数据进行计算和定量分析。该方法可根据实际情况快速识别并计算出所测区域的土壤蒸发强度，具有测量精度高，范围广等特点，可为灌区水分管理及作物灌溉时间提供一定的技术支持。
7.为了达到上述目的，本发明提供了如下技术方案：
8.一种基于热红外成像技术估算土壤蒸发强度的方法，所述方法包括如下步骤：
9.步骤1、获取地表温度图像；
10.使用便携式热红外成像仪对所测区域进行拍摄，得到热红外图像；
11.步骤2、数据预处理；
12.将热红外图像导入软件fluke中；软件fluke依据实际气温矫正背景温度，并将拍
摄时的气温记为t0；将图像温度数据导出为txt格式，并从中提取温度矩阵p；对温度矩阵p进行预处理，使用小波去噪原理去除温度矩阵p中的高频部分，得到新的温度矩阵p1；
13.步骤3、去除温度矩阵干扰项；
14.土壤蒸发干扰项主要为植物及地面石块，根据干扰项的分布特点，使用两种聚类法依次进行聚类；
15.3a.温度矩阵聚类
16.首先采用基于密度的聚类方法进行聚类，得出若干个聚类，并分别计算各个小类的温度平均值；
17.其次，根据各个温度平均值进行基于距离的二次聚类，将聚类结果分为三类，温度最低的一类为植物温度矩阵，温度最高的一类为地面石块温度矩阵；
18.3b.将植物温度矩阵从温度矩阵p1原始矩阵中删除，得到温度矩阵p2；
19.3c.将地面石块温度矩阵从温度矩阵p2中删除，得到温度矩阵p3；
20.3d.计算温度矩阵p3的霍普金斯常量，若温度矩阵p3的霍普金斯常量为0.75～1.0之间，则判定仍然具有明显的聚类趋势，重复步骤3a到3d，直到温度矩阵p3的霍普金斯常量为0.5～0.75之间；若温度矩阵p3的霍普金斯常量为0.5～0.75之间，则判定分布较为均匀，进行下一步；
21.步骤4、估算蒸发强度；
22.4a.使用软件spss遍历出温度矩阵p3中的最大温度px，将最大温度px与空气温度t0之差记作m；
23.4b.估算土壤蒸发强度，包括如下步骤：
24.4b.1：将温度矩阵p3各项减去空气温度t0，得到温度矩阵p4；
25.4b.2：将所得温度矩阵p4各项均除以m，得到温度矩阵p5；
26.4b.3：根据当地气象数据获取净辐射值rn，将净辐射值rn除以2.45
×
106，记作r；
27.4b.4：将温度矩阵p5各项均除以r，所得矩阵即为蒸发强度矩阵p6；
28.步骤5、绘制土壤蒸发强度等值线图；
29.首先，将蒸发强度矩阵p6导入软件matlab中，绘制土壤蒸发强度热力图；其次，利用工具包contour绘制等值线图，调整图像的大小和显示效果。
30.步骤2中，对温度矩阵p进行预处理，使用小波去噪原理去除温度矩阵p中的高频部分，得到新的温度矩阵p1，包括如下步骤：
31.步骤2b.1：将温度矩阵p导入软件matlab中；
32.步骤2b.2：对温度矩阵p进行小波分解；设置小波分解层数为3层，并获得低频分量wt1和高频分量wt2；
33.步骤2b.3：分尺度去噪；对高频分量wt2进行阈值处理，当小波系数小于阈值时，将其置为零；当小波系数大于阈值时，将其减去阈值；阈值为：其中，σ为温度矩阵p的方差，m为温度矩阵p的行数，n为温度矩阵的列数；
34.步骤2b.4：小波逆变换；对处理后的分量进行小波重构，得到新的温度矩阵p1。
35.步骤3a中，包括如下步骤：
36.步骤3a.1：首先进行基于密度的聚类，即只要临近区域的密度超过某一阈值，就进
行聚类；根据聚类结果，将矩阵分别记为n1、n2、n3
…
，并计算每一分类的温度平均值，分别记为a1、a2、a3
…
；
37.步骤3a.2：将计算所得的a1、a2、a3
…
进行基于距离的k-means聚类，将聚类结果设置为3类，记为t1、t2、t3；
38.具体聚类过程为：
39.1、随机从a1、a2、a3
…
选择其中三个，其值记作s1、s2、s3；
40.2、分别计算a1、a2、a3
…
到s1、s2、s3之间的距离；
41.3、分别判断a1、a2、a3
…
到s1、s2、s3的距离，若距离s1近，则归为t1；若距离s2近，则归为t2；若距离s3近，则归为t3；
42.4、计算t1、t2、t3的平均值ta1、ta2、ta3，计算与s1、s2、s3与ta1、ta2、ta3的距离，若距离小于1，则结束本步骤；若大于1，则将s1、s2、s3的值替换为t1、t2、t3，重复步骤3a.2；
43.步骤3a.3：分别计算t1、t2、t3分类的平均温度值，记作at1、at2、at3。
44.与现有技术相比，本发明的有益效果在于：
45.1、本发明应用广泛，可运用于大部分地区土壤蒸发强度的估算。
46.2、地表温度可采用便携式热红外成像仪获取，对实验器材的限制较小，可在野外实现现场观测，省时省力。
47.3、本发明克服了微型蒸发器测量土壤蒸发的局限性，避免了土壤异质性对测量值的干扰，克服了频繁更换土壤蒸发器所带来的测量不确定性。
48.4、本发明避免了土壤蒸发推算过程中的复杂计算，提高了土壤蒸发强度的估算效率。
49.5、本发明填补了点尺度和区域尺度之间的估算方法的空白，可为农业灌溉提供一定的技术指导，提高灌溉水利用效率。
50.6、本发明考虑了植物和地面石块对土壤蒸发强度估算的影响，并利用两步聚类对两类干扰项温度进行了剔除，提高了土壤蒸发强度的估算精度，并为后续的进一步改进提供了方向。
51.7、本发明的估算过程所需参数较少，核心代码模块化设计、使用简便。
附图说明
52.图1为坡面及拍摄过程示意图；
53.图2为本发明一种基于热红外成像技术估算土壤蒸发的方法的流程图；
54.图3为干扰项温度去除过程的示意图；
55.图4为估算土壤蒸发强度等值线图。
56.其中的附图标记为：
57.1、坡面
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
2、地下水埋深
58.3、水位观测井
ꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀꢀ
4、便携式热红外成像仪
59.5、土壤蒸发测量区域
具体实施方式
60.下面结合附图对本发明的具体实施方式作进一步说明。
61.如图2所示，本发明提供了一种基于热红外成像技术估算土壤蒸发强度的方法，该方法可对地势平坦或有一定坡降的地表土壤蒸发强度进行估算，包括如下步骤：
62.步骤1、获取地表温度图像。
63.使用便携式热红外成像仪对所测区域进行拍摄，得到热红外图像。优选地，所测区域应尽量避免杂草及多石块区域。
64.步骤2、数据预处理。
65.将热红外图像导入软件fluke中；软件fluke依据实际气温矫正背景温度，并将拍摄时的气温记为t0；将图像温度数据导出为txt格式，并从中提取温度矩阵p；对温度矩阵p进行预处理，使用小波去噪原理去除温度矩阵p中的高频部分，得到新的温度矩阵p1。
66.步骤3、去除温度矩阵干扰项。
67.土壤蒸发干扰项主要为植物及地面石块，根据干扰项的分布特点，使用两种聚类法依次进行聚类。
68.3a.温度矩阵聚类
69.首先采用基于密度的聚类方法进行聚类，得出若干个聚类，并分别计算各个小类的温度平均值；
70.其次，根据各个温度平均值进行基于距离的二次聚类，将聚类结果分为三类，温度最低的一类为植物温度矩阵，温度最高的一类为地面石块温度矩阵；
71.3b.将植物温度矩阵从温度矩阵p1原始矩阵中删除，得到温度矩阵p2；
72.3c.将地面石块温度矩阵从温度矩阵p2中删除，得到温度矩阵p3；
73.3d.计算温度矩阵p3的霍普金斯常量，若温度矩阵p3的霍普金斯常量为0.75～1.0之间，则判定仍然具有明显的聚类趋势，重复步骤3a到3d，直到温度矩阵p3的霍普金斯常量为0.5～0.75之间；若温度矩阵p3的霍普金斯常量为0.5～0.75之间，则判定分布较为均匀，进行下一步。
74.步骤4、估算蒸发强度。
75.4a.使用软件spss遍历出温度矩阵p3中的最大温度px，将最大温度px与空气温度t0之差记作m；
76.4b.估算土壤蒸发强度，包括如下步骤：
77.4b.1：将温度矩阵p3各项减去空气温度t0，得到温度矩阵p4；
78.4b.2：将所得温度矩阵p4各项均除以m，得到温度矩阵p5；
79.4b.3：根据当地气象数据获取净辐射值rn，将净辐射值rn除以2.45
×
106，记作r；
80.4b.4：将温度矩阵p5各项均除以r，所得矩阵即为蒸发强度矩阵p6。
81.步骤5、绘制土壤蒸发强度等值线图。
82.首先，将蒸发强度矩阵p6导入软件matlab中，绘制土壤蒸发强度热力图；其次，利用工具包contour绘制等值线图，调整图像的大小和显示效果。
83.如图1所示，坡面1为人工堆砌的斜坡；地下水埋深2为坡面不同位置与地下水位之间的距离；观测井3位于坡面顶部，镶嵌于坡面1内，用于观测地下水位；便携式热红外成像仪4位于坡面1底部，用于拍摄土壤蒸发测量区域5的热红外图像；土壤蒸发测量区域5为位于坡面处的矩形区域。本发明的实施例来自于一人工改造的坡面1，坡面距离地下水位的不同而形成了不同的地下水埋深2，使土壤蒸发测量区域5内形成不同的土壤蒸发条件，使发
明的估算结果更有代表性。改造后的斜坡与地下水位具有较好的水力联系，因而具有广泛的代表性。使用水位观测井3对地下水位进行观测，确保拍摄过程中地下水位的稳定，提高土壤蒸发强度的估算精度。包括如下步骤：
84.步骤1、获取地表温度图像
85.步骤1的目的是获取土壤蒸发测量区域5的地表温度。使用便携式热红外成像仪4获取地表温度，所述便携式热红外成像仪4(例如为fluke tis60)，温度数据拍摄于当日14时，拍摄时应确保无大的云块遮挡阳光。便携式热红外成像仪4与土壤蒸发测量区域5之间的距离为1.0～2.5米之间，便携式热红外成像仪4与地面的水平夹角为45
°
左右，最大限度减少距离对热辐射的衰减作用。根据实际情况将反射率调整为0.95-0.97之间。坡面1平整，无大的杂草、多石块区域。
86.步骤2、数据预处理
87.步骤2a、获取温度矩阵
88.步骤2a的目的是获取热红外图像的温度矩阵。包括如下步骤为：
89.步骤2a.1：将热红外图像导入软件fluke中；
90.步骤2a.2：软件fluke依据实际气温矫正背景温度，并将拍摄时的气温记为t0；
91.步骤2a.3：将温度数据导出为txt格式，得到温度矩阵，将温度矩阵记为p。
92.步骤2b、温度矩阵去噪
93.本步骤的主要目的是将温度矩阵中的高频部分去除。采用小波分析作为主要方法进行去噪，得到新的温度矩阵p1，包括如下步骤：
94.步骤2b.1：将温度矩阵p导入软件matlab中；
95.步骤2b.2：对温度矩阵p进行小波分解。设置小波分解层数为3层，并获得低频分量wt1和高频分量wt2；
96.步骤2b.3：分尺度去噪。对高频分量wt2进行阈值处理，当小波系数小于阈值时，将其置为零；当小波系数大于阈值时，将其减去阈值。阈值为：其中，σ为温度矩阵p的方差，m为温度矩阵p的行数，n为温度矩阵的列数；
97.步骤2b.4：小波逆变换。对处理后的分量进行小波重构，得到新的温度矩阵p1。
98.步骤3、去除温度矩阵干扰项
99.如图3所示，植物及地面石块的表面温度是估算土壤蒸发的主要干扰项。本步骤的目的主要是去除这两项干扰项，提高土壤蒸发估算精度。
100.步骤3a、温度矩阵聚类
101.步骤3a的目的是对温度矩阵进行聚类。基于植物和石块在地面温度分布的特点(植物分布多呈条状或羽状；石块分布多呈现点状)，本发明主要采用两种聚类方法依次分别进行聚类，包括如下步骤：
102.步骤3a.1：首先进行基于密度的聚类，即只要临近区域的密度超过某一阈值，就进行聚类。该阈值的选取不应大于温度矩阵的极大极小值之差的5％。根据聚类结果，将矩阵分别记为n1、n2、n3
…
，并计算每一分类的温度平均值，分别记为a1、a2、a3
…
；
103.步骤3a.2：将计算所得的a1、a2、a3
…
进行基于距离的k-means聚类，将聚类结果设置为3类，记为t1、t2、t3；
104.具体聚类过程为：
105.1、随机从a1、a2、a3
…
选择其中三个，其值记作s1、s2、s3；
106.2、分别计算a1、a2、a3
…
到s1、s2、s3之间的距离；
107.3、分别判断a1、a2、a3
…
到s1、s2、s3的距离，若距离s1近，则归为t1；若距离s2近，则归为t2；若距离s3近，则归为t3；
108.4、计算t1、t2、t3的平均值ta1、ta2、ta3，计算与s1、s2、s3与ta1、ta2、ta3的距离，若距离小于1，则结束本步骤；若大于1，则将s1、s2、s3的值替换为t1、t2、t3，重复步骤3a.2。
109.步骤3a.3：分别计算t1、t2、t3分类的平均温度值，记作at1、at2、at3。
110.步骤3b、去除植物干扰项
111.本步骤的目的是去除温度矩阵中的植物温度。所拍摄时，植物表面温度普遍低于土壤表面温度，比较at1、at2、at3中的最小值所在的分类，并从温度矩阵p1中删除，得到温度矩阵p2。
112.步骤3c、去除地面石块干扰项
113.本步骤的目的是去除温度矩阵中的地面石块温度。所拍摄时，地面石块表面温度普遍高于土壤表面温度，比较at1、at2、at3中的最大值所在的分类，并从温度矩阵p2中删除，得到温度矩阵p3。
114.步骤3d、效果检验
115.本步骤的主要目的是检验干扰项的去除效果。使用软件matlab计算分类t1、t2、t3的霍普金斯常量，并记作h，若h值为0.5～0.75之间时，则判定样本已接近随机分布，已经充分排除了以上两类干扰项的影响，则进行下一步。若h值为0.75～1.0之间，则判定样本仍然具有较为明显的聚类趋势，即温度矩阵p3中仍含有一定的干扰项数据，则重复步骤3a到3d，直到温度矩阵p3的霍普金斯常量为0.5～0.75之间。
116.步骤4、估算土壤蒸发强度
117.步骤4a、计算地表温度最大值
118.本步骤的主要目的是计算温度矩阵p3的最大值。其步骤为：将温度矩阵p3导入软件spss中，使用工具包遍历温度矩阵p3的最大值，并将其记作px，将最大温度值px与大气温度t0之差记为m。
119.步骤4b、估算土壤蒸发强度
120.本步骤的目的是利用温度矩阵3估算所拍摄区域的土壤蒸发强度。包括如下步骤：
121.步骤4b.1：将温度矩阵p3导入到软件matlab中；将温度矩阵p3的各项数值分别减去t0，得到温度矩阵p4；
122.步骤4b.2：将温度矩阵p4的各项数值分别减去m，得到温度矩阵p5；
123.步骤4b.3：根据气象数据获取净辐射rn，将净辐射rn除以2.45
×
106，记作r；
124.步骤4b.4：将温度矩阵p5除以r，所得新矩阵即为蒸发强度矩阵，记作p6。
125.步骤5、绘制土壤蒸发强度等值线图
126.本步骤的主要目的是绘制所拍摄区域的土壤蒸发强度示意图。包括如下步骤：
127.步骤5.1：将蒸发强度矩阵p6导入到软件matlab中，利用工具包heatmap绘制土壤蒸发强度热力图；
128.步骤5.2：利用工具包contour绘制等值线图，调整图像的大小和显示效果，清晰地
表达土壤蒸发测量区域5的土壤蒸发状况，如图4所示。