考虑土壤水分层的SW模型对果园蒸散发的过程模拟

黄可静，王纪超，肖 俊，高学睿，赵西宁2，

（1.西北农林科技大学水利与建筑工程学院，陕西 杨凌 712100；2.西北农林科技大学旱区农业水土工程教育部重点实验室，陕西 杨凌 712100；3.中国水利水电科学研究院，北京 100038；4.西北农林科技大学 水土保持研究所，陕西杨凌 712100）

0 引 言

黄土高原处于半干旱、半湿润地区，常年降雨量稀少，蒸发大，且是世界上最大的黄土堆积区，面临水资源短缺和水土流失的生态问题[1]。中国政府在20世纪90年代大力推动退耕还林（草）项目以缓解黄土高原的生态环境问题[2]。同时，为了保证当地农民的经济收益，大力扩大黄土高原地区红枣，苹果等经济作物的种植面积。在过去十几年中，陕西省苹果的产量和种植面积占到了全国的1∕4[3]。苹果产业由于其较高的经济效益，已成为黄土高原地区应对大规模减产和促进区域经济有效发展的主要产业[4]。由此可见，苹果产业的可持续发展对于农民经济收入和区域生态可持续发展至关重要。黄土高原的部分苹果树种植区地形复杂、灌溉条件差，雨水是果树生产最主要的水分来源，缺水状态下的果树产量不稳定，品质不高。缺水已成为限制果园苹果产量和品质的重要因素[5]。因此，及时、准确的定量估算果树的蒸腾量和土壤的蒸发量，提高农业用水效率，对于果园的用水管理具有重要的生产意义。

近些年，基于遥感数据的蒸散发估算已经越来越成熟，在利用遥感数据来重建、估算作物蒸散发时会使用到各种方法以及模型，例如：水量平衡法[6]、趋势分析法[7]、以及利用机器学习方法进行估算等[8,9]。SW模型由于其参数较少，模拟结果较好已被广泛的应用[10]。值得一提的是，研究者们在选择SW模型时多是由于它可以单独模拟蒸发和蒸腾过程，并且对于不同区域、不同种类的树木等方面都得到了验证与实际运用[11,12]。同时，SW模型也在被不断地进行改进，Hu等人[13]在SW模型中引入Ball-Berry气孔导度模型[14]以及GPP模型[15]形成了新的SWH模型，以期更加准确的估算植物冠层阻力及蒸散发，新的SWH模型在实际运用中也得到了良好的模拟结果[16]。然而，SW模型在实际的应用过程中是将土壤水视为一个输入变量，土壤水数据较难获取的特点会使得模型的应用也受到限制；同时，由于SW模型并未对土壤水进行分层，对于水文过程的刻画上不够细致，没有考虑不同根深处的土壤水分变化，这会在一定程度上影响模型的准确性和可靠性。因此，本研究改进了SW模型，在SW模型中加入了土壤水分层模块，以降雨作为输入数据，通过模拟水文过程将土壤水进行分层计算，以更准确、及时的估算田块尺度蒸散发；并利用增加了土壤水分层模块的SW模型对黄土高原台塬区田块尺度的果树蒸腾量（T）和土壤蒸发量（E）模拟，探索其对苹果树蒸散发模拟的适用性，为果园的用水管理提供理论支撑。

1 材料与方法

1.1 研究区概况

研究区域位于陕西省咸阳市淳化县马家镇，该位置居于陕西省中部偏西，咸阳市北部，渭北黄土高原沟壑区南缘。该研究区属温带大陆性季风气候，2019年平均气温11.6 ℃，年降水量545.8 mm，降雨日数为97 d，降雨集中在每年7-10月。该区域光照时间长，降雨适中，昼夜温差大，是黄土高原中南部地区苹果树栽培的理想位置。作为黄土高原地区苹果园的典型代表，对该区域果树蒸散发及土壤含水量的研究可以为黄土高原地区苹果的种植提供有效的参考。

1.2 数据来源

气象数据：本研究在研究区域内安装了气象站，自动记录的气象数据包括温度、湿度、风速、降雨、辐射强度，数据时间分辨率为1 h。

遥感数据：本研究主要使用了下载于欧空局数据中心（https:∕∕scihub.copernicus.eu∕）的Sentinel-2A遥 感 影像数据。选取了研究区域从2019年1-12月的24景遥感影像，对哨兵2A L1C级数据进行辐射校正、几何校正和大气校正，重采样为10m分辨率，并基于ENVI对数据进行波段融合和剪裁。

茎流计数据：本研究在马家研究区内选取了3棵直径在20 cm的苹果树，采用热扩散式茎流计（TDP）测量果树的蒸腾，频率为5 min 1次。

土壤水分数据：本研究在选取的样本树下不同深处分别埋设TDR315土壤水分探头，监测土壤水分变化。

1.3 遥感数据处理

归一化植被指数NDVI（Normalized Difference Vegetation Index）是一种遥感中常用的植被指数，它可以反应植物的生长情况和营养参数。NDVI定义为近红外波段（NIR）与可见光红波段（Red）反射率之差与反射率之和的比值[17,18]，本研究基于Sentinel-2A波段融合后的数据计算，其表达式如下：

式中：NIR代表近红外波段的反射率；Red代表红波段的反射率。

叶面积指数LAI决定植被冠层对净辐射的截获量，也是SW模型重要的输入参数，本研究基于NDVI和LAI的关系反演LAI[12,17]，反演公式为：

式中：Kpar为消光系数，取值为0.3；Fcanopy=NDVI+0.05。

1.4 热扩散式茎流计

本研究采用热扩散式茎流计（TDP）测定苹果树的蒸腾量，该研究区域的苹果树树龄在12年及以上，本实验选取3棵生长状况良好的苹果树，并测得树干面积为207.03 cm2，296.97 cm2，246.42 cm2。该热扩散式茎流计（TDP）记录了树木的树干边材每5 min的平均液流速率，通过下式求出果树蒸腾量。

式中：T为蒸腾量，mm∕d；A为树干面积，m2；V为液流速率，mm∕d；n为研究区域内苹果树棵数；S为研究区果园面积，m2。

本研究区域面积44 470 m2，果树间距3.5 m×4 m，除去田埂、道路，大约估算3 000棵树，由此估算研究区域单位面积的果树蒸腾量（T）。

1.5 模型构建

1.5.1 模型原理

本研究基于改进的SW模型，将研究区域划分为10m×10m的网格，利用气象数据和遥感数据模拟研究区苹果树逐日的植被蒸腾量（T）和土壤蒸发量（E）。该模型所需要输入的气象数据包括降水（Pre）、风速（WS）、水汽压差（VPD）、温度（Ta）和光合有效辐射（PAR），遥感数据主要是叶面积指数（LAI）和归一化植被指数（NDVI），模型模拟的流程图如图1所示。对SW模型的改进部分主要体现在土壤水分层模块，公式（13）～（16）对其原理进行了说明。

图1 模型模拟的流程图Fig.1 Flow chart of model simulation

SW模型基本水平衡原理如下：

式中：i为天数；ΔSWi表示土壤水分变化量，mm；TRi表示降雨量，mm；SFi表示地表径流，mm；Ti表示植被蒸腾量，mm；Ei表示土壤蒸发量，mm；SEPi表示深层渗漏量，mm。

其中，TRi为模型的输入降雨数据，SFi基于SCS曲线法计算，主要公式如下：

式中：S为潜在最大保留量，mm；NRi为第i天的降雨量，mm；CN为曲线数，该曲线数通过蒙特卡洛模拟优化。

T和E分别为模型的植被蒸腾量和土壤蒸发量，主要基于SW双源蒸散发原理计算，具体公式如下：

式中：λ表示为水的汽化潜热，J∕kg；PMc和PMs分别表示植被蒸腾和土壤蒸发，mm；Cc和Cs分别表示植被阻力系数和土壤表层阻力系数；R和Rs分别表示冠层以上和土壤表面以上的可用能量输入，w∕m2；raa、rac和ras分别表示冠层与参考高度之间的空气动力阻力、冠层植被要素的整体边界层阻力和冠层与地表之间的空气动力阻力，s∕m；Δ表示饱和水汽压随温度变化的斜率，kPa∕K；ρ表示空气密度，kg∕m3；cp表示定压比热，J∕(kg·K)；D表示参考高度上的水汽压亏缺，kPa；γ表示干湿球常数，kPa∕K。

土壤表层阻力rss和冠层气孔阻力rsc是该模型中的关键参数，通过引入土壤表面阻抗方程和Ball-Berry气孔导度模型实现对这两个参数的估算[14]，具体方法如下：

土壤表面阻抗rss估算方法：

式中，SW和SWs分别为土壤表层含水量和土壤表层的饱和含水量，m3∕m3；b1、b2和b3为经验参数，s∕m。

估算冠层阻抗是结合气孔阻抗模型Ball-Berry模型求得：

式中：g0为最小气孔导度，m∕s；a1为经验参数；pn为光合速率，µmol∕(m2·s)；hs为冠层表面空气相对湿度；cs为气孔胞间的CO2浓度，本研究取值500×10-6；此处使用光能利用率模型估算的GPP代替pn。

式中：PAR为光合有效辐射，µmol∕(m2·s)；FPAR为植物冠层吸收的PAR占总射入PAR的比例；ε为根据气温、土壤表层含水量、饱和水气压差标准化后的光能利用率。

土壤分层计算模块原理如下：

如果NRi-SFi-FC(1)+SW(1)i-1≤0

如果NRi-SFi-FC(1)+SW(1)i-1＞0

如果SEP(i1)-FC(2)+＞0

式中：(1)、(2)分别为第一层和第二层土层；SFi为第i天的地表径流，mm；SEP(1)i、SEP(2)i为当前土层到下一土层的渗水量，mm；FC(1)，FC(2)为该土层土壤的田间持水量，mm；，为第i-1天该土层实际土壤含水量，g∕kg；K(1)、K(2)为该土层的土壤导水率，mm∕s；DP(1)、DP(2)为该土层的深度，mm。

本研究使用蒙特卡洛算法[19]估算表1中列出的6个经验参数。在确定6个参数的可能范围时，采用以下3个步骤进行参数校准，以最小化建模的不确定性：①在前人研究成果的基础上，为每个参数指定一个大致范围，并在给定范围内随机抽样参数集进行500万次蒙特卡洛模拟[20]； ②利用线性回归函数对模拟ET和实测数据进行比较，得到决定系数R2和斜率k。③选取R2值最高的20个模拟值，k值在0.95～1.05之间为最佳模拟。取20个参数集的平均值，用最优参数驱动模型，得到最终的仿真数据，确定最优参数集。

表1 6个关键参数的最终取值Tab.1 Final values of the six key parameters

1.5.2 模型验证

本研究通过计算模拟值和观测值的相关系数（R2）和均方根误差（RMSE）对模型结果的准确性进行验证。计算公式如下所示：

式中：Qsimi表示模拟值；Qobsi表示实测值；上划线表示平均值；n是总的实测数量。

R²表示模型模拟值与实测值的相关程度；一般地，R²=0.5被认为是可接受的阈值。RMSE表示模拟值与实测值的偏差程度。

2 结果与分析

2.1 模型验证

本研究测量了研究区域2019年5月17日到2019年8月4日的80 d苹果树液流速率，并基于此计算了果园单位面积的果树蒸腾量T作为实测值，选取前40 d数据利用蒙特卡洛的方法来校准模型，后40 d数据验证模型模拟结果。模型模拟的果树蒸腾量与果园实测值的验证结果如图2所示。

图2 果树蒸腾量模拟结果与实测值的比较Fig.2 Comparison between simulated and measured transpiration of fruit trees

本研究基于5月17到6月25日的实测果树蒸腾量数据对模型进行校准，在校准期，模拟值与实测值的R2为0.83，RMSE为0.39 mm。同时，通过6月26日到8月4日的数据对模型进行验证，在验证期间，相关系数R2=0.76（大于0.5），RMSE=0.50 mm。

由于苹果树根系主要分布在土壤表面下80cm处，因此本研究分别对根深20 cm处（表层）和80 cm（深层）的土壤水进行了模拟，并通过实测数据对模拟的土壤含水量进行验证，验证结果如图3所示。

图3 土壤含水量模拟结果与实测值的比较Fig.3 Comparison between simulated and measured soil water content

本研究基于5月17到8月4日的实测土壤含水量数据对模型进行验证，结果显示20 cm和80 cm处模拟值与实测值的R2与RMSE分别为0.820、0.016和0.840、0.005。该模型模拟结果较好，模拟精度较高，可以准确的量化田块上日尺度的果树蒸腾量及土壤含水量。

2.2 LAI和NDVI的时间变化特征分析

本研究基于研究区域24景的Sentinel-2A数据计算了果园日尺度的NDVI和LAI，结合该研究区域的果树物候信息对遥感输入数据进行评价，物候表见表2。

表2 研究区域果树的物候表Tab.2 Phenology of fruit trees in the study area

由图4遥感数据的时间变化特征来看，NDVI和LAI在4月到8月底持续增长，该研究区域的果树在3月20日进入萌芽期后，果树逐渐经历萌芽期-开花期-幼果期，叶面积指数和植被覆盖度均在逐渐增加，8月下旬果树的NDVI和LAI均达到最大，此时果树已经处于果实的膨大期。10月上旬果实成熟，采摘后果树树叶开始掉落，叶面积指数和植被覆盖度大幅降低，果树也逐渐进入休眠期。总的来说，该研究区域的NDVI和LAI变化与果树的物候信息基本一致，基于Sentinel-2A遥感影像计算的NDVI和LAI日尺度数据能正确的反映果树的生长变化，也能满足模型的输入要求。

2.3 模型模拟结果分析

2.3.1 果树蒸腾量的时间变化特征

本研究基于改进的SW模型模拟了2019年马家研究区域田块尺度的苹果树蒸腾量，结果如图5所示。

模型模拟的果树蒸腾量的时间变化特征和实测果树蒸腾数据如图5所示，实测值与模型模拟结果耦合较好，苹果树的生育期主要集中在4-10月份，模拟的果树蒸腾从3月下旬开始逐渐增加，在7、8月份达到顶峰，主要是由于果树在3月底进入萌芽期，植被耗水逐渐增大，7月份和8月份果树进入果实膨大期，同时，随着降雨量的增加，耗水达到了最大值。果树蒸腾在10月中旬之后逐渐降低，主要是由于苹果采摘后，随着气温降低，果树逐渐进入休眠期，果树的耗水逐渐降低。总的来说，该模型的模拟结果与实测耦合良好，也能正确的反映果树的生理过程。

图5 苹果树蒸腾量的时间变化特征Fig.5 Temporal variation characteristics of apple tree transpiration

2.3.2 研究区域果树蒸散分布特征

SW模型可以分离土壤蒸发和植被蒸腾量，马家研究区域的蒸散组分分布特征如图6所示。结果表明，果树从4-10月份蒸腾量T先增大后减少，7月蒸腾量最大为46.64 mm。7-9月份是果树蒸腾量的高峰期，果树正经历膨大期和成熟期。该模拟结果与李青华[21]实测数据一致。果园总蒸散发ET也是在4-10月份呈现先增大后减少的趋势，ET在7月份最大，为77.22 mm。该区域降雨主要发生在5-10月份，10月份的降雨达到最大150.8 mm，但10月份的果树蒸腾量T仍在减少，主要是由于10月上旬苹果采摘后，果树叶子也逐渐凋落，果树逐渐进入休眠期，果树蒸腾耗水逐渐减少。总的来说，该研究区果树在发育期4-10月，蒸腾量T呈现先增加后减少的趋势。2019年总蒸散发ET为478.43 mm，果树蒸腾量T为265.17 mm，土壤蒸发E为213.26 mm，果树蒸腾T占总蒸发的55.4%。

图6 蒸散发组分分布特征Fig.6 Distribution characteristics of evapotranspiration components

2.4 水分利用效率计算与评价

用水效率（WUE）是指单位用水量所产生的同化物质的量，基本上反映了用水量与干物质生产之间的关系，是评价生长适宜性的综合生理生态指标。在干旱和半干旱地区的极端环境中，植物是否能够适应当地的极端环境条件，最重要的是看它们是否能够很好地协调碳同化和水耗散之间的关系，即植物用水效率（WUE）是其生存的关键因素之一。对于果树，水分利用效率（WUE）是每单位水蒸腾和蒸发产生的干物质质量（GPP）。

本研究利用加入了土壤水分层模块的SW模型模拟了研究区域果树的蒸散发ET和干物质质量GPP，进而计算了WUE，如图7所示，GPP从3月到10月份呈现先增大后减少的趋势，在9月达到最大值32.11 g∕m2，7-9月的GPP值较大，占全年的36.58%。与GPP相比，WUE的年内分布在4月份达到最大，水分利用效率为0.56 g∕(m2·mm)，果树关键生育期为4-9月，除4月用水效率较高外，其他月份用水效率低于0.5 g∕(m2·mm)，用水效率较低。该果树种植地位于黄土高原台塬区，属于干旱半干旱地区，该地区降雨较少，为了保证果实的品质，应在关键生育期保证果树水分的供应，提高用水效率。

图7 苹果树GPP和WUE的时间变化图Fig.7 Time changes of GPP and WUE in apple trees

3 讨 论

苹果树的植被蒸腾量在整个生育期内变化规律明显，呈现出先增加后减少的趋势，从3月份开始逐渐增加，7月份达到最高值，然后开始下降直到休眠期，这与徐巧等人[22]的研究基本吻合。植被蒸腾量与植物叶片密切相关，原因在于其蒸腾主要依靠叶片气孔进行，而其气孔的开闭主要受太阳辐射和气温的影响[23]，本研究使用茎流计估算果树蒸腾量，而影响茎流的主要气象因素是辐射量[24]，植被蒸腾量必然与叶面积指数LAI相关。果树在3月底进入萌芽期，叶面积指数LAI也从这个时间开始增加，植被蒸腾量T也随之增加；从6月份开始进入果实膨大期，耗水量急剧上升，在成熟期之前的7、8月份达到最大值；随后，植被蒸腾量在9、10月份苹果成熟期间快速下降。与植被蒸腾量不同的是，LAI的最大值并非出现在7、8月份，而是在9月份。两者峰值错开的原因可能是叶片在蒸腾量最高峰时便开始逐渐老化以及受辐射条件影响，植被蒸腾能力减弱，所以在LAI略有增长的情况下植被蒸腾量却有所下降。

苹果树的植被蒸腾量T在生育期内的绝大部分时间要大于土壤蒸发量E，10月下旬之后才转变为E大于T，转变原因是10月下旬苹果树叶片开始脱落，极大降低了植被蒸腾量T。生育期内大部分时间T大于E的结论与李青华等人的研究结果一致，但值得注意的是，也有学者的研究出现了不同的结果，Wang等人[25]的研究中出现了苹果树绝大部分生长期内E大于T的情况，其主要原因可能在于树龄，随着树木的年龄增长，茎流也会随之增加[26,27]。本研究中苹果树的树龄在12年以上，前文提到E大于T的情况出现在树龄为7年的苹果树上，所以，苹果树的年龄较小时可能会有E大于T的情况出现。

SW模型在ET建模方面具有良好的性能[28]，本文的的验证结果也表明该模型可以用来量化田块尺度的果树蒸发量。通过蒙特卡罗方法模拟校准模型的关键参数，并将土壤水分纳入气孔导度的计算可以使SW模型的可靠性大大提高[29,30]。本研究中模型模拟的ET结果略低于实际观测值。可能的一个原因是模型没有考虑冠层截获水分的蒸发，这也是模型不确定性的来源之一，可能导致的结果是在降雨较多的时期降低模型对ET模拟的精度。但也有研究在对冠层蒸发单独计算之后，认为冠层截流对模型模拟精度的影响比较有限[20]。此外，由于模型是固定参数运行，这也可能导致模型对ET的模拟结果偏低[31]。多项研究表明，植被和气象条件的季节性变化会使得模型一些参数产生相应的改变[13,23]，因此，可以通过使用动态参数建立模型，以实现参数随植被和气象条件的变化，从而进一步提高模型的准确性。土壤水分层模块的SW模型对于果树的蒸散发及土壤含水量的模拟结果较好，可以较为准确的量化田块上日尺度的果树蒸腾量，同时果树蒸腾随时间的变化特征也正确的反应了果树的生理过程。本研究选取了典型的黄土高原苹果种植区作为研究区域，改进的SW模型在该区域得到了良好的验证，可以认为改进的SW模型在黄土高原苹果种植区具有一定的适用性。

本研究在实验设计阶段主要考虑了果树的蒸腾量，关于土壤蒸发的考虑不够充分并未布设相关的测量仪器，这可能会增加本研究的不确定性，在后续的实验中将会进一步增加试验站点；此外，由于侧重点在于对模型进行改进，并验证其适用性，因此，未来的研究可以通过网格化模拟等手段，将研究区从田块尺度扩大，进一步将SW模型应用于灌区尺度、流域尺度等更大的区域，实现更大面积的果树蒸腾及土壤蒸发的预测，为黄土高原苹果种植区提供更加行之有效的参考。

4 结 论

（1）改进的SW模型对于苹果树蒸散发的模拟结果和实测值的相关系数R2为0.78，RMSE为0.62 mm，模拟结果较好，该模型可用于黄土高原地区田块尺度的果园蒸散发模拟计算。

（2）模拟结果表明，果树蒸腾T和果园蒸散发ET在4-10月份呈现先增大后减小的趋势，7月份T达到最大值58.14 mm。

（3）综合考虑苹果树生长期内水分利用效率WUE和干物质质量GPP的变化，应当在果实膨大期6-8月份保证水分供应，并提高水分的利用效率。