0 引言
水分可以帮助作物更好地吸收利用养分，作物生物量与土壤含水率息息相关。降雨和灌溉需要转化为土壤水分，才能被作物根系直接吸收和利用。降雨量不足、土壤含水率过低以及缺乏科学合理的灌溉均会导致农作物减产。因此，如何使用遥感技术实时准确监测反演大面积区域土壤的湿度情况具有重要意义[1]。研究表明，改进的表观热惯量计算模型监测土壤含水率是可行的[2]；利用增强型植被指数（EVI）构建的遥感植被供水指数（VSWI）总体上更适用于监测河南省的春季干旱[3]；利用叶面积指数（LAI）代替归一化植被指数（NDVI）反演土壤水分具有较高的监测精度[4]；以归一化植被指数（NDVI）和比值植被指数（RVI）计算得到的温度植被干旱指数（TVDI）与表层土壤水分相关性较好，且NDVI更适用于中、低植被覆盖时的土壤水分反演[5]。
基于地表温度（Ts）和NDVI之间的特征空间，综合利用二者之间的关系是干旱监测领域中常见的方法。NDVI和Ts的相对变化与土壤含水率的关系在大多数气候条件和地表覆盖条件下是相对稳定的[6]。温度植被干旱指数通过比较监测期内NDVI所对应的地表温度与干边地表温度（最干旱）和湿边地表温度（最湿润）之间的关系确定区域旱情。当前对于区域土壤湿度的监测反演研究中，大多数研究采用的是MODIS影像进行反演，MODIS影像具有范围广、综合性强、易于操作的优点，但是MODIS影像是由多期影像合成且分辨率较低，时间和空间精度略低。Landsat 8 OLI影像为单日影像，分辨率达到了30 m，可有效提升土壤湿度监测反演中时间和空间精度。为此，使用Landsat 8 OLI产品数据，采用温度植被干旱指数（TVDI）反演太谷县任村乡5、7、8、10月的土壤湿度，并利用野外实地采样的数据验证遥感影像反演结果，分析太谷县任村乡土壤湿度的时空分布变化情况。
1 材料与方法
1.1 研究区概况
任村乡位于山西省晋中市太谷县（ 112°41′8.25″E—112°48′6.15″E，37°28′55.62″N—37°33′57.87″N），属于温带季风性气候，雨热同期。研究区总面积46.67 km2，耕地面积3000 hm2，年平均气温9.7℃，10℃以上积温3517～3796℃，年降雨量460 mm。
1.2 数据来源与修正
数据来源于美国地质调查局（USGS）提供的Landsat 8 OLI数据。影像覆盖研究全区域，符合研究要求，影像日期为2017年5月5日、7月8日、8月9日、10月12日。
数据下载后由于传感器本身和大气的一些影响，会对辐射的计算造成一定误差，所以需要对影像进行辐射校正以便后期的相关参数计算。由于影像在下载前已经进行过基本的几何校正，故直接使用ENVI软件对影像进行辐射定标和大气校正，并利用任村乡行政区域矢量图对影像进行裁剪，从而修正及消除影像畸变可能对研究参数产生的影响。
1.3 研究原理与方法
有学者提出了通过温度植被干旱指数（TVDI）估算土壤水分的方法[7]，该方法利用全遥感数据，不需要借助其他任何地面实测资料，通过归一化植被指数（NDVI）和地表温度（Ts）计算温度植被干旱指数（TVDI）。
1）归一化植被指数（NDVI）。其计算式为:
式中:ρNIR为近红外波段反射率；ρR为红光波段反射率。
2）地表温度（Ts）。采用单通道反演法计算Ts，结合Jiménez-Muñoz等[7]在原有的算法基础上增加的针对Landsat 8 OLI的大气参数，具体算法为:
式中:br为常数，Landsat 8的TIRS热红外传感器具有2个波段，分别编号为TIRS 10和TIRS 11，TIRS 10对应的br为1324 K，TIRS 11对应的br为1199 K；γ和δ是根据普朗克函数求得的常数，可以通过亮温T、大气顶部辐射L和br计算得到[10]，即γ≈T2/(brL),δ≈T-T2/br；ε是地表比辐射率；ψ1、ψ2和ψ3是大气参数，ψ1=1/τ,ψ2=-L↓-L↑/τ,ψ3=L↓，其中τ为大气透过率，L↑和L↓是大气上行和下行辐射强度，τ、L↑和L↓均可以从美国航空航天局(NASA)的网站（https://atmcorr.gsfc.nasa.gov/）获得。分别计算出TIRS 10和TIRS 11对应的地表温度T10和T11后，使用二者的均值Td以缩小地表温度的反演误差。
地面高程影响地表温度反演精度，故需要结合高程数据对温度反演结果进行修正，计算式为:
式中:Ts为修正后的地表温度；Td为修正前的地表温度；H为地面高程；a为修正系数，表示温度随地面高度的增加而降低，根据相关参考文献[8]，a一般取-0.6℃/100 m。
3）温度植被干旱指数（TVDI）。其计算式为:
式中:Ts是影像中任意像元对应的地表温度，文中为结合高程数据修正后的地表温度，与式（3）得到的Ts一致；Tsmin是影像中同一NDVI值对应的最低地表温度，Tsmin=a1+b1·NDVI，其中a1和b1为TVDI湿边方程的系数；Tsmax是影像中同一NDVI值对应的最大地表温度，Tsmax=a2+b2·NDVI，其中a2和b2为TVDI干边方程的系数。TVDI的取值范围为0～1，TVDI越小，代表土壤含水率越大，反之，代表土壤含水率越小[9]。
2 结果与分析
2.1 不同时期的Ts-NDVI特征空间及干湿边方程拟合
利用Ts-NDVI特征空间中的归一化植被指数、地表最高温度和最低温度可以分别拟合出干、湿边方程。研究[13]发现，当一个地区的植被覆盖率在15%以下，NDVI不能表明该地区的植物生物量，当该区域的植被覆盖度在15%～80%之间时，NDVI能够随着植被覆盖度的增加而增加，而当该区域的植被覆盖度大于80%时，由于植被覆盖密集，NDVI对植被检测的灵敏度将降低。因此，NDVI适用于作物生长发育的中期阶段或者中等植被覆盖度地区的植被监测。
利用反演计算得到的NDVI和Ts，使用ENVI软件统计影像中同一NDVI分别对应的最大、最小地表温度，通过像元直方图确定拟合范围，绘制不同时期的Ts-NDVI特征分布图，同时利用归一化植被指数及其对应的最高、最低地表温度，结合实际情况，在删除掉大于0.80和小于0.12的NDVI后，拟合出干、湿边方程，结果如图1所示。由图1可以看出，随着NDVI的增大，最高温和最低温的温差逐渐减小，最终呈现出一个三角特征空间。
图1 研究区不同时期Ts-NDVI特征空间
2.2 TVDI等级图分布
根据TVDI计算公式，结合干湿边方程，计算不同时期影像各像元的TVDI，以TVDI作为土壤湿度分级标准，将土壤湿度分为5个层次:极干旱（0.8＜TVDI≤1）、干旱（0.6＜TVDI≤0.8）、正常（0.4＜TVDI≤0.6）、湿润（0.2＜TVDI≤0.4）和极湿润（0＜TVDI≤0.2），由此得到研究区不同时期土壤湿度分布图，如图2所示。
图2 研究区不同时期土壤干旱等级分布图
从图2可以看出，5月，西北部土壤湿度比较大，干旱区集中在南部和东部；7月，南部和东部土壤湿度增大，西部土壤湿度减小；8月，研究区土壤湿度整体增大，土壤整体较为湿润；10月，土壤湿润度整体减小，整个研究区的土壤状况偏干旱。
2.3 TVDI与实测含水率的相关性验证
2017年5月5日在任村乡进行卫星同步的样点土壤体积含水率的测定，综合考虑地理位置、种植作物、土壤质地等影响因素，选取位置平坦、无明显遮挡的100个均匀分布的样点（图3），使用TDR水分速测仪分别测量各样点0～10、10～20、20～30 cm土层土壤体积含水率，同时使用RTK测量仪器测定记录每个样点的位置，根据仪器测定的样点的位置，找出Landsat 8 OLI影像图上对应的像元，计算像元的TVDI。对测量样点的体积含水率和样点对应像元计算出的TVDI进行线性回归分析，并构建以土壤体积含水率为横轴，TVDI为纵轴的散点图，结果见图4。
图3 土壤含水率实测采样点分布图
图4 TVDI与不同土层实测含水率的相关性
从图4可以看出，不同土层温度植被干旱指数（TVDI）均随着土壤含水率增大而明显减小。对TVDI和土壤体积含水率的线性回归结果进行F检验，发现线性回归方程都达到显著水平（α=0.05或P＜0.05），这表明TVDI能够体现土壤表层的水分状况，可以作为评价土壤表层水分状况和干旱程度的指标。从不同深度土层和TVDI的相关性看，10～20 cm土层土壤含水率与TVDI相关程度最高，20～30 cm土层次之，表层0～10 cm相关性最低，表明TVDI更能稳定反映和指示10～20 cm土层的水分状况。但0～10 cm土层水分数据之间相对离散，这是因为Landsat 8 OLI影像的空间分辨率为30 m，TVDI反映的是单位像元内表层土壤水分的总体状况，而地面实测含水率为点数据，遥感反演的数值和地面实测值之间存在误差。
2.4 任村乡土壤湿润程度分析
通过分析研究试验区域不同时期TVDI反演结果（图2）可知，5月，西北部土壤湿度比较大，干旱区集中在北部、南部和东部；7月，南部、东部和北部土壤湿度增大，西北部土壤湿度减小；8月，研究区土壤湿度整体增大，土壤整体较为湿润；10月，土壤湿润度整体减小，其中东部土壤湿度减小幅度较大，东部则从湿润将为正常。
结合作物分布、生长和气候变化等因素可知，研究区种植作物多为玉米，东北部每年4月中旬播种，10月上旬收获，南部和西北部每年5月下旬播种，9月下旬收获。任村乡北部和南部土壤质地以砂壤土为主，东部和西部土壤质地则以黏壤土为主。5月，降水较少，东北部玉米生长旺盛，需水量较大，土壤偏干旱，南部作物尚未播种，土壤质地为砂壤土，保水能力较弱，所以南部为极干旱；西北部土壤质地为黏壤土，保水能力较强，土壤湿润；7月和8月，降水量增大，灌溉量增加，使得土壤含水率增大，居民区为人工下垫面，对太阳辐射的吸收率大，再加上植被较少，导致居民区成为干旱区，而耕作区土壤整体较为湿润；10月，北部和南部玉米基本上都已收获，裸土面积大大增加，土壤水分减少，所以研究区西部和南部土壤湿度降低，东部土壤则为正常偏湿润。
3 讨论
目前，关于区域土壤湿度的遥感反演研究主要采用MODIS数据和Landsat数据。MODIS数据提供16 d合成的植被指数数据以及8 d合成的地表温度/反射率数据，研究时能够避免使用原始数据繁杂的处理过程，有研究[23]使用MODIS数据对区域土壤水分进行监测，并指出MODIS数据空间分辨率较低，适用于长时间、大面积干旱监测。前期针对大区域的干旱遥感监测，多是采用MODIS数据和Landsat TM/ETM数据，而Landsat OLI/TIRS系列数据的开放，为区域土壤水分监测提供了新的数据源，但以Landsat 8 OLI为数据源进行的干旱监测研究较少[24]；Ts-NDVI特征空间方法是一个接近于实时的干旱监测方法，可以尝试使用更短的观测周期跟踪干旱发生发展的动态演变过程[25]。本研究在前人研究基础上以Landsat 8 OLI为数据源，利用温度植被干旱指数法对区域进行土壤湿度动态监测，结合实测数据对监测结果进行检验，结果表明以温度植被干旱指数法反演监测区域土壤湿度是切实可行的，这与已有的研究结果[13,26]均一致。同时，本研究指出TVDI与10～20 cm土层处的土壤湿度相关性最高，这与已有相关研究结果[1,27]一致，但是这2个研究并没有考虑NDVI对植被覆盖程度高低的指示灵敏度不同。本研究虽然对NDVI的指示灵敏度有所考虑，但是，本研究没有对裸土和高植被覆盖区的土壤水分采用相应的模型进行反演，在今后的研究中可以结合其他方法对全年的土壤水分监测进行深入研究。
4 结论
1）随着土壤湿度的增大，TVDI逐渐减小，TVDI和0～30 cm各层土壤湿度的线性回归结果经过F检验后都能达到显著水平，TVDI可以作为指示地表土壤湿度的指标和干旱评价指标。
2）在0～30 cm土层中，TVDI与10～20 cm土层土壤湿度相关性最高，TVDI更能稳定反映和指示10～20 cm土层土壤的水分状况。
3）遥感反演的研究区不同时期土壤湿度分布变化特征与研究区作物分布生长情况及气候变化规律基本吻合，说明温度植被干旱指数法可以用来动态监测大面积区域的土壤湿度状况。
总体上，运用GIS技术和遥感影像数据，结合合理的模型，使用温度植被干旱指数法进行大面积长时效的土壤湿度反演监测，为农业生产提供科学有效的指导是切实可行的。
[1]王纯枝,毛留喜,何延波,等.温度植被干旱指数法(TVDI)在黄淮海平原土壤湿度反演中的应用研究[J].土壤通报,2009,40(5):998-1005.
[2]吴黎,张有智,解文欢,等.改进的表观热惯量法反演土壤含水量[J].国土资源遥感,2013,25(1):44-49.
[3]李琪,孙晓宇,王连喜,等.基于不同植被指数的VSWI在河南省春季干旱监测中的应用分析[J].作物杂志,2016(1):162-168.
[4]高中灵,郑小坡,孙越君,等.利用地表温度与LAI的新型土壤湿度监测方法[J].光谱学与光谱分析,2015,35(11):3129-3133.
[5]张喆,丁建丽,李鑫,等.TVDI用于干旱区农业旱情监测的适宜性[J].中国沙漠,2015,35(1):220-227.
[6]CARLSON T N,PERRY E M,SCHMUGGE T J.Remote estimation of soil moisture availability and fractional vegetation cover for agricultural fields[J].Agricultural&Forest Meteorology,1990,52(1):45-69.
[7]JIMÉNEZ-MUÑOZ J C,SOBRINO J A,SKOKOVIĆ D,et al.Land Surface Temperature Retrieval Methods From Landsat-8 Thermal Infrared Sensor Data[J].IEEE Geoscience&Remote Sensing Letters,2014,11(10):1840-1843.
[8]袁志杰.基于ATI模型和TVDI模型的晋中土壤水分遥感反演研究[D].太谷:山西农业大学,2015.
[9]SANDHOLT I,RASMUSSEN K,ANDERSEN J.A simple interpretation of the surface temperature vegetation index space for assessment of surface moisture status[J].Remote Sensing of Environment,2002,79(2):213-224
[10]JIMENEZMUNOZ J C,CRISTOBAL J,SOBRINO J A,et al.Revision of the single-channel algorithm for land surface temperature retrieval from Land sat thermal-infrared data[J].IEEE Transactions on Geoscience&Remote Sensing,2009,47(1):339-349.
[11]QIN Z,KARNIELI A,BERLINER P.A mono-window algorithm for retrieving land surface temperature from Landsat TM data and its application to the Israel-Egypt border region[J].International Journal of Remote Sensing,2001,22(18):3719-3746
[12]ROZENSTEIN O,QIN Z,DERIMIAN Y,et al.Derivation of land surface temperature for Landsat-8 TIRS using a split window algorithm[J].Sensors,2014,14(4):5768.
[13]赵杰鹏,张显峰,廖春华,等.基于TVDI的大范围干旱区土壤水分遥感反演模型研究[J].遥感技术与应用,2011,26(6):742-750.
[14]徐涵秋,黄绍霖.Landsat 8 TIRS热红外光谱数据定标准确性的分析[J].光谱学与光谱分析,2016,36(6):1941-1948.
[15]姚春生,张增祥,汪潇.使用温度植被干旱指数法(TVDI)反演新疆土壤湿度[J].遥感技术与应用,2004(6):473-478.
[16]徐涵秋.新型Landsat 8卫星影像的反射率和地表温度反演[J].地球物理学报,2015,58(3):741-747.
[17]高懋芳,覃志豪,刘三超.MODIS数据反演地表温度的参数敏感性分析[J].遥感信息,2005(6):3-6.
[18]张成才,吴泽宁,余弘婧.遥感计算土壤含水量方法的比较研究[J].灌溉排水学报,2004,23(2):69-72.
[19]李斌,秦明周,张鹏岩.Landsat 8遥感影像NDVI的对比研究[J].河南大学学报(自然科学版),2017,47(2):155-161.
[20]魏强,张吴平,吴亚楠,等.基于SEBAL模型的小麦水分生产率研究[J].灌溉排水学报,2017,36(7):38-46.
[21]韩丽娟,王鹏新,王锦地,等.植被指数-地表温度构成的特征空间研究[J].中国科学(D辑:地球科学),2005(4):371-377.
[22]沙莎,郭铌,李耀辉,等.温度植被干旱指数(TVDI)在陇东土壤水分监测中的适用性[J].中国沙漠,2017,37(1):132-139.
[23]孙玉龙,郝振纯.TDR技术及其在土壤水分及土壤溶质测定方面的应用[J].灌溉排水,2000,19(1):37-41.
[24]季国华,胡德勇,王兴玲,等.基于Landsat 8数据和温度-植被指数的干旱监测[J].自然灾害学报,2016,25(2):43-52.
[25]于敏,程明虎,刘辉.地表温度-归一化植被指数特征空间干旱监测方法的改进及应用研究[J].气象学报,2011,69(5):922-931.
[26]范辽生,姜纪红,盛晖,等.利用温度植被干旱指数(TVDI)方法反演杭州伏旱期土壤水分[J].中国农业气象,2009,30(2):230-234.
[27]熊世为,李卫国,贾天山,等.基于HJ卫星数据的土壤含水量反演及其旱情预测[J].江苏农业学报,2014,30(5):1044-1050.