基于遥感反演的南方丘陵灌区农田灌溉水量估算
2023-10-17和玉璞胡梦阳麦紫君史一平
付 静,和玉璞,胡梦阳,孙 浩,麦紫君,史一平
(1.南京水利科学研究院水文水资源与水利工程科学国家重点实验室,江苏南京 210029;2.河海大学农业科学与工程学院,江苏南京 210029;3.江苏省农村水利科技发展中心,江苏南京 210029;4.高淳区水务局,江苏南京 211301)
农业用水占我国用水总量的62.1%,其中灌溉水量占农业用水90%以上[1]。灌区作为我国粮食生产及经济作物种植的重要区域,区域内水资源状况的精准评估及科学调度对提高我国农业水资源利用效率具有重要意义[2-4]。南方低山丘陵灌区内分布有大量小型水库、塘坝等,调蓄水源用于拦蓄区域降雨及地表径流,再向周边区域农田灌溉供水。由于灌区调蓄水源数量众多,现有的水位-流量监测点位仅能覆盖小部分规模较大的水库及塘坝,无法全面反映灌区实际灌溉面积及用水情况,限制了灌区水资源利用与配置的进一步优化[5]。
目前农田灌溉用水量估算方法分为3 类:一是采用典型调查的灌溉水量或现行灌溉定额与实灌面积数据进行估算,灌溉面积与灌水定额通过抽样调查与层层上报的方法获取,因其时效性、宏观性和准确性等方面的限制,统计的精度需要进一步提高[6-11]。二是基于数值模型计算灌区农田灌溉水量,众多学者以水文模型为工具开展了灌溉水量及变化规律的研究[12-13],基于验证后的数值模型得到的灌溉用水量,虽然弥补了数据精度的短板,但模型存在着诸多假设和不确定性,只适合特定区域,其空间异质性很难推广应用。三是基于遥感数据计算灌溉水量,随着各类高时间、空间、光谱分辨率卫星的出现,遥感技术以其观测面积大、周期短和数据具有较强的综合性等特点广泛应用于农业[14]。基于遥感影像光谱特征和植被指数差异的作物种植结构提取、灌溉面积识别[15-16]为灌区用水量估算提供了新路径。目前我国灌溉面积识别主要针对北方旱作物,如通过垂直干旱指数识别灌溉面积,旱作物灌水频次低、灌溉水量少,而南方灌区稻田灌溉频繁且灌水量大,灌溉前后土壤水分处于饱和的状态,基于垂直干旱指数等识别灌溉面积不适用于南方地区。
基于此,本文以淳东灌区为研究区,融合Landsat-8、MODIS 遥感影像,利用机器学习算法提取水稻种植面积,结合实测数据进行农田灌溉水量估算,分析灌溉前后水稻植被水分指数及其差值的变化特征,确定灌溉面积的植被水分指数临界阈值,实现稻田实际灌溉面积的识别,预期研究成果可为高效、精准获取灌区灌溉用水数据提供技术支撑。
1 灌区概况及监测点布置
1.1 区域概况
淳东灌区位于南京市高淳区的东部,设计灌溉面积2.07 万hm2,有效灌溉面积1.92 万hm2,主要灌溉水源为胥河。主要种植水稻、小麦、油菜和蔬菜等作物。灌区气候四季分明,寒暑显著,降水充沛,无霜期共有247.5 d,多年年平均气温16.4 ℃,平均蒸发量1 406.6 mm,平均降水量1 218.1 mm。灌区属于南方低山丘陵区域典型的“长藤结瓜”灌区。由于灌区旱作物生育期内的自然降水满足作物需水要求,故本研究以水稻为典型作物,开展基于遥感影像的灌区农田灌溉用水量估算。
1.2 监测点位布置
2021年,在淳东灌区共布设2个监测区域,监测区域1位于高淳区义保村,面积5.07 hm2,用于分析、获取区域内植被水分指数阈值;监测区域2 位于高淳区马家,面积8.46 hm2,用于对比分析不同种植结构提取估算的灌溉水量与实测水量差值。
将监测区域1 内的渠道作为典型渠道,布置Odyssey 电容式水位计,水位记录间隔为10 min;采用LB50-1C 型旋杯式流速仪测定典型断面不同时间、不同水位时的流速,采用Origin 2018 拟合得到渠道内水位-流量函数关系,通过典型断面的测量水位推算流量,从而得到灌溉水量。
2 研究方法
2.1 遥感数据及预处理
本研究将重现周期为1 d的MODIS影像与空间分辨率30 m的Landsat-8影像融合,采用ESTARFM模型实现融合。为提高后续的提取精度,对选用的遥感影像进行预处理,主要包括辐射定标和大气校正,大气校正采用FLAASH模块。
2.2 农田灌溉用水量遥感监测方法
利用灌区范围内的遥感影像数据及机器学习算法,实现灌区种植结构的提取,进而,结合观测数据建立农田灌溉取用水量动态分布数据集,融合遥感监测数据与原位实测数据进行灌区农田灌溉水量计算。将监测区域1的实测数据作为典型调查水量,监测区域2 基于水稻种植结构的提取结果与实测数据结合估算灌溉水量,并与区域2 实测的灌溉水量进行对比分析,以验证农田灌溉用水量遥感监测精度。
2.3 基于植被水分指数的灌溉面积识别方法
本研究利用植被水分指数变化反演实际稻田灌溉面积,采用IMSI和IMSII两种植被水分指数进行计算,两者的计算公式为
式中,rSwir1为短波红外1 的反射率;rSwir2为短波红外2 的反射率;rNir为近红外的反射率。
计算灌水前后的植被水分指数与差值,结合监测区域的实际灌水记录,对植被水分指数的变化值进行分析,得到适合的灌溉临界值。当指数差值大于灌溉临界值时,表明此时作物含水量较之前变低,该区域属于未灌溉状态;当小于临界值时,表明此时作物含水量较高,该处已进行了灌溉。将验证后的临界值用于其他灌溉过程的面积识别,对比反演灌溉面积与实际灌溉面积,验证该方法的精度。
3 结果与分析
3.1 农田灌溉用水量遥感监测精度
通过布设的水位计和流速仪计算7 月31 日监测区域1 灌溉总量为861.90 m3,结合监测区域1 的灌溉面积,可得此次灌溉水层厚度为1.70 cm。同时在监测区域2 测量多点的灌溉水深,得到平均值为1.62 cm,测得该区域灌溉水量为1 370.44 m3。
利用基于目视解译结合现状调研的方法获取监测区域2 水稻的种植面积,将监测区域1 的灌溉水量作为典型调查水量,应用于监测区域2,通过水稻种植面积与灌溉水深估算监测区域2 灌溉用水量为1 400.80 m3。监测区域2 的实测灌水量1 370.44 m3与基于遥感提取水稻种植面积估算灌溉用水量1 400.80 m3相差30.36 m3,精度较高,表明将监测区域1实测灌水量作为典型区域估算灌溉用水量较为合理。
3.2 灌溉面积识别精度
3.2.1 植被水分指数变化特征
选择7月21日、8月1日、8月6日及8月30作为典型稻田灌溉时段,监测区域1 在上述4 个时段均进行了全域灌溉,灌水深度分别为16、5、12、50 mm。通过融合后的遥感影像计算监测区域1 水稻MSI、MSII 植被水分指数,4 个典型灌溉时段监测区域1 MSI、MSII指数的差值见图1~2。
图1 MSI指数差值分布
图2 MSII指数差值分布
由图1~2 可以看出,MSI 指数差值在-0.22 至0.11 范围内变化,MSII 指数差值在-0.09 至0.2 之间波动,除极个别像元外,灌水后的作物植被水分指数MSI、MSII数值均小于灌水前,与实际规律相符。
3.2.2 实际灌溉面积反演与精度分析
水稻生长阶段共选择了4次灌溉过程,利用8月1日、8月30日灌溉过程确定临界值阈值,利用剩余2次灌溉过程验证反演面积的精度。目前临界值确定的方法有两种:一是选取灌溉点和非灌溉点,通过分析比较不同临界值产生的结果确定阈值[17];二是分析不同临界值时灌溉面积的变化情况,当临界值达到某一值时,灌溉面积出现明显变化,通过试算确定临界值[4]。由于农户高频次灌水,无法获取较好的非灌水时段,故采用第二种方法确定临界值。
根据选取的典型灌溉过程,试算不同的临界值反演得到的淳东灌区灌溉面积,经统计MSI、MSII差值在-0.01至0处分布较多,由于灌水后水稻水分胁迫指数下降,故以-0.001 为起算点向下试算,当反演的灌溉面积出现较大变动时,则认为此时的临界值为最优临界值。MSI、MSII 差值分布统计见图3,不同MSI、MSII 差值反演灌区水稻灌溉面积见表1~2,不同MSI、MSII 临界值反演水稻灌溉面积变化曲线见图4。
表1 不同MSI差值反演灌区水稻灌溉面积
表2 不同MSII差值反演灌区水稻灌溉面积
图3 MSI、MSII差值分布统计
图4 不同MSI、MSII临界值反演水稻灌溉面积变化曲线
从图表中可看出,MSI差值取值大于-0.002时,反演的灌溉面积较为稳定,并没有明显变化,当MSI差值小于-0.002时,各临界值反演的灌溉面积变化较大,因此,确定灌溉临界值阈值-0.002。将MSI临界值-0.002 代入剩余2 次灌溉过程,反演灌溉面积并验证精度,如图5 所示。当MSII 差值取值大于-0.006 时,反演的灌溉面积较为稳定,并没有明显变化,当MSII 差值小于-0.006 时,各临界值反演的灌溉面积出现变化,因此,确定灌溉临界值MSII临界阈值为-0.006。将临界值-0.006代入剩余的2次灌溉过程,反演灌溉面积并验证精度,如图6所示。
图5 MSI反演水稻灌溉面积
图6 MSII反演水稻灌溉面积
基于MSI水分指数计算的7月22日与7月21日稻田灌溉面积为4.45万m2,8月7日与8月6日稻田灌溉面积为4.73 万m2。监测区域1 在这两次灌水中全域灌溉,反演的精度分别为87.82%和93.42%。基于MSII水分指数7月22日与7月21日反演的灌溉面积为2.37万m2,8月7日与8月6日反演灌溉面积为4.64 万m2,反演的精度分别为46.71%和91.66%。已有研究表明,使用遥感手段反演旱作物实际灌溉面积可行,精度基本达到80%以上[18],本研究在南方丘陵灌区应用MSI植被水分指数反演典型时段水稻灌溉面积平均精度为90.62%,与现有研究对旱地灌溉面积识别精度的报道较为接近。
在本研究中,MSI 水分指数反演的灌溉面积精度较好,反演灌溉面积具有一定的可行性,MSII 水分指数7月22日与7月21日灌溉面积精度为46.71%,具有不稳定性。总体MSI水分指数反演精度较MSII高,这主要由于在短波红外波段中,作物冠层水分对SWIR1(1.57~1.65 μm)波段反射率较为敏感,与遥感监测植被冠层水分的波段在1.48~1.759 μm 最佳这一结论保持一致[19-21]。
3.3 淳东灌区应用结果
3.3.1 农田灌溉用水量估算
基于遥感反演的水稻种植面积与灌溉面积均可以较好地估算灌区用水量,精度较高。而利用植被水分指数变化反演实际灌溉面积,由于在研究中无法获取质量较好的逐日遥感影像,只能选取影像质量较高的时间段作为典型灌溉时段,此种方法无法落实到整个稻季。故本研究利用实测数据结合水稻种植面积估算淳东灌区在水稻生育期内的总灌溉用水量。
为提取水稻灌溉面积,首先采用2021 年3 月26 日、5月29日、9月18日和12月7日4个时间段的Landsat-8遥感影像对水稻的种植结构进行提取,分别采用CART 决策树、支持向量机和随机森林分类器进行提取和精度对比分析,经验证CART 决策树提取精度最高,为96.48%。提取的淳东灌区2021年水稻种植面积为6 793.3 hm2。
通过义保村监测区域1 的灌溉过程,结合淳东灌区水稻种植面积计算淳东灌区在水稻生育期内的总用水量。经过调查,淳东灌区在水稻泡田期时的泡田定额约为1 500 m3/hm2,本研究估算的灌溉用水量不包括泡田期。本研究水稻于6 月7 日插秧,并于当日安装水位计,至9 月11 日水稻进入黄熟期,此后水稻未进行灌溉。从6 月7 日至9 与11 日,估算淳东灌区总灌溉用水量5 108.84 万m3。
7—8 月是淳东灌区水稻生育期用水最多的时段,这与实际情况相符。6 月下旬灌溉总量为774.48 万m3,占总灌溉用水量的15.16%。7 月灌溉总量为2 017.72 万m3,占水稻生育期灌溉总量的39.49%,也是占比最多的一个月份。8 月灌溉总量为1 562.55 万m3,仅次于7 月,占水稻生育期灌溉总量的30.59%,其中8 月上旬共灌水658.99 万m3,中旬由于发生连续降雨,故灌溉用水量共计190.22 万m3,8月下旬共灌溉了713.34万m3水量,其中8月30日1 d的灌溉水量为339.68万m3,是水稻生育期灌溉水量最大的1次。9月共灌溉754.10万m3,占灌溉总量的14.76%。
3.3.2 典型灌溉时段实际灌溉面积识别
选取典型时段识别灌区实际灌溉面积,考虑到降雨影响、灌溉次数以及遥感数据质量,最终选择7 月30 日至8 月8 日作为典型灌溉时段,其中8 月3 日至8月5日淳东灌区发生降雨,故不考虑灌溉水量。利用前述方法计算水稻种植区域的MSI水分指数以及灌溉前后指数差异,结合MSI 临界值阈值识别水稻区域的实际灌溉面积。淳东灌区水稻种植区域典型时段灌溉面积见图7。
图7 淳东灌区典型灌溉时段实际灌溉面积
4 结 语
针对南方丘陵灌区灌溉用水特点,构建了基于遥感反演的农田灌溉水量估算方法,在典型区域验证的精度为97.78%。结合卫星遥感与机器学习算法提取灌区种植结构,与实测数据融合计算南方丘陵灌区灌溉用水量。选取淳东灌区典型区域,基于遥感反演的农田灌溉水量计算值为1 400.80 m3,与实测灌水量1 370.44 m3相差仅30.36 m3,精度较高。
基于MSI、MSII 植被水分指数反演稻田实际灌溉面积的精度均值分别为90.62%、69.19%,基于MSI 临界阈值提取的水稻灌溉面积,总体精度优于MSII。结合监测区域的实际灌溉过程,通过试算确定植被水分指数的灌溉临界值分别为-0.002及-0.006。
利用基于遥感反演的农田灌溉水量估算方法明确了淳东灌区稻田灌溉水量。淳东灌区在水稻生育期内总灌溉用水量5 108.84 万m3,其中7 月和8 月的灌溉用水量占比较大,分别占总灌溉用水总量的39.49%和30.59%,与实际情况相符。