王金鑫，趙光成，張廣周，于百順，羅蔚然，張成才，李 穎

(1.鄭州大學(xué)水利與環(huán)境學(xué)院，鄭州 450001；2. 中國氣象局·河南省農(nóng)業(yè)氣象保障與應(yīng)用技術(shù)重點(diǎn)實(shí)驗(yàn)室，鄭州 450003；3. 河南省氣象科學(xué)研究所，鄭州 450003)

0 引 言

遙感墑情監(jiān)測是目前遙感技術(shù)應(yīng)用的前沿領(lǐng)域[1,2]，也是公認(rèn)的世界性研究難題之一[1-4]。究其原因，影響墑情的因素包括土壤類型、降水、太陽高度、大氣、植物長勢、地表狀況以及田間管理等等，使土壤墑情遙感反演異常困難[1,5,6]?；谥脖恢笖?shù)的墑情監(jiān)測或反演是遙感墑情監(jiān)測的常用方法。在一定區(qū)域、一定時(shí)間內(nèi)，氣溫、土壤、地形等因素相對穩(wěn)定，只降雨量變化會對作物生長造成較顯著影響[7]，因而，植被指數(shù)的變化與土壤墑情之間有較顯著相關(guān)性。但植被受缺水脅迫在短期內(nèi)不會在生長狀態(tài)上明顯表現(xiàn)出來，具有一定的時(shí)滯性[8]，無法及時(shí)地反映出土壤水分狀況(后面有詳細(xì)討論)，所以植被指數(shù)是一個靈敏度較差的水分脅迫指標(biāo)[9]，但在農(nóng)作物的整個生育期內(nèi)，墑情的波動與植被指數(shù)的變化頻率與趨勢應(yīng)該是一致的。除早期的基于歸一化植被指數(shù)(NDVI)的距平植被指數(shù)(AVI)、條件植被指數(shù)(VCI)[10]外，考慮地表溫度影響(其是十分復(fù)雜的[11])，人們又提出了植被供水指數(shù)(VSWI)[12,13]、溫度植被旱情指數(shù)(TVDI)[14]、條件植被溫度指數(shù)(VTCI)[15,16]以及特征空間法[14]等新的概念和方法用于土壤墑情監(jiān)測，相應(yīng)于研究案例都取得了較好的效果。經(jīng)過多年的研究積累，墑情遙感監(jiān)測也達(dá)成一些基本共識：①從土壤深度角度，以10～20 cm左右效果最好[1-3,17]，表層和深層土壤水分與遙感資料的相關(guān)關(guān)系往往不佳[3]；②從地表覆蓋角度，裸土和低植被以熱慣量模型和微波遙感效果最好[3,4,18]；完全植被覆蓋以作物缺水指數(shù)(CWSI)、距平植被指數(shù)、條件植被指數(shù)、特征空間模型、能量指數(shù)法和植被供水指數(shù)法等效果較好[1-3,19]；部分植被覆蓋采用農(nóng)田蒸散雙層模型(但模型復(fù)雜)、改進(jìn)求解的熱慣量模型和溫度調(diào)節(jié)指數(shù)效果較好[3,18]。概而言之，裸土?xí)r，土壤的熱特性起決定作用，基于溫度的方法較好；植被茂密時(shí)，植被的生長狀態(tài)是重要的指示因子，指數(shù)法占優(yōu)；半植被覆蓋時(shí)，要同時(shí)考慮土壤和植被熱特性差異的雙重效應(yīng)，綜合指數(shù)或方法最好。充分證明了地表溫度和植被指數(shù)是農(nóng)業(yè)墑情光學(xué)遙感監(jiān)測中的主要指示因子[4]。

國內(nèi)外已有研究表明，任何單一的指數(shù)或方法都有自己的特點(diǎn)和適用條件，地表生態(tài)的復(fù)雜性和遙感技術(shù)的特點(diǎn)也決定了不存在一種普適的墑情監(jiān)測指數(shù)和方法。大范圍農(nóng)田墑情遙感監(jiān)測更是存在遙感數(shù)據(jù)的適用性選擇以及光熱條件的區(qū)域分異問題。本文擬提出一種基于高時(shí)間分辨率、中空間分辨率的遙感時(shí)序數(shù)據(jù)，考慮地域和物候期差異，利用距平植被指數(shù)，實(shí)現(xiàn)大區(qū)域冬小麥全生育期土壤相對濕度的監(jiān)測方法。以河南省中東部黃淮海平原冬小麥主產(chǎn)區(qū)為研究對象，以MODIS NDVI為基礎(chǔ)墑情指數(shù)，針對不同地區(qū)、不同生長階段引進(jìn)相應(yīng)因子進(jìn)行自適應(yīng)改進(jìn)，實(shí)現(xiàn)對農(nóng)田墑情及其時(shí)空分布規(guī)律的監(jiān)測。

1 研究區(qū)域與數(shù)據(jù)源

河南省黃淮海平原位于北緯32°08′～36°21′和東經(jīng)112°51′～116°35′之間，北起豫冀交界，南至淮河一線，西起海拔100 m等高線和豫西北丘陵邊緣，東至豫魯、豫皖分界線(為方便數(shù)據(jù)處理，本文以該區(qū)域的行政界線為準(zhǔn))，面積約8.7 萬km2，占黃淮海平原總面積的25.6%，河南省總面積的52.4%，屬于南暖溫帶半濕潤大陸性季風(fēng)氣候區(qū)[20]。該地區(qū)是河南省的主要農(nóng)業(yè)區(qū)，冬小麥主產(chǎn)區(qū)，年農(nóng)業(yè)總產(chǎn)值占全省農(nóng)業(yè)總產(chǎn)值的70%以上。

本文采用的遙感數(shù)據(jù)包括：研究區(qū)域2011-2015小麥生產(chǎn)年度的b1、b2和b7波段的Terra MODIS數(shù)據(jù)、MODIS的每日陸面溫度產(chǎn)品等，從共享網(wǎng)站(https:∥ladsweb.nascom.nasa.gov/)下載。 b1、b2和b7波段數(shù)據(jù)的數(shù)據(jù)等級是L2級的陸地表面反射率數(shù)據(jù)，已經(jīng)過大氣糾正和幾何糾正。b1和b2波段的空間分辨率為250 m，共下載1367幅；b7波段的分辨率為500 m，根據(jù)需要下載了75幅；MODIS每日陸面溫度數(shù)據(jù)(簡稱遙感地表溫度數(shù)據(jù)，包括每天白天和夜晚兩幅數(shù)據(jù)，這里取其平均值)是從共享網(wǎng)站下載的成品數(shù)據(jù)，空間分辨率是1 000 m，數(shù)量和使用情況同b1、b2波段。主要問題是有一些日期的數(shù)據(jù)“漏洞”很多(沒有溫度反演數(shù)據(jù))。MODIS源數(shù)據(jù)為Sinusoidal projection，利用MRT軟件統(tǒng)一轉(zhuǎn)換為WGS84基準(zhǔn)的Albers等積投影，并將地溫?cái)?shù)據(jù)和b7波段數(shù)據(jù)重采樣為250 m，以與b1、b2波段匹配，進(jìn)行像元級計(jì)算。

本研究所涉及的墑情及其他數(shù)據(jù)：研究區(qū)域2014- 2015年度10和20 cm土壤水分含量實(shí)測數(shù)據(jù)、2011-2014年2月上旬和11月下旬的10和20 cm土壤水分含量實(shí)測數(shù)據(jù)、豫南、中、北地區(qū)的多個監(jiān)測站1981-2010年(30年)逐月平均氣溫?cái)?shù)據(jù)、2014-2015年河南省冬小麥種植區(qū)域數(shù)據(jù)以及河南省行政界線數(shù)據(jù)等，從河南省氣象局調(diào)研獲取。墑情實(shí)測數(shù)據(jù)是研究區(qū)域均勻分布的100多個站點(diǎn)的10和20 cm每日土壤體積含水量的平均值。以各觀測站30年的逐月平均氣溫進(jìn)行內(nèi)插(本文的內(nèi)插方法均采用基于球面函數(shù)的普通克里金插值)，分別得到各區(qū)域每個像元的常年逐月平均氣溫。將2014- 2015年冬小麥種植區(qū)域數(shù)據(jù)、河南省界線數(shù)據(jù)也轉(zhuǎn)換為與上述遙感數(shù)據(jù)統(tǒng)一的基準(zhǔn)與投影，然后對影像數(shù)據(jù)進(jìn)行剪裁和植被指數(shù)計(jì)算。由于冬小麥?zhǔn)呛幽鲜〉闹饕Z食作物，種植范圍年際變化不大，所以本文假設(shè)2014-2015年的小麥種植范圍同前三年的種植范圍相同，并用之對MODIS數(shù)據(jù)做掩膜，得到最終的研究區(qū)域。

2 研究方法與過程

2.1 研究思想闡述

中原地區(qū)的冬小麥生長周期較長，從頭年十月到來年六月，跨深秋、冬季、春季和初夏幾個季節(jié)，長達(dá)八個月，溫度起伏大，生態(tài)因素多變。要進(jìn)行大尺度麥田墑情監(jiān)測，既要考慮地域差異，又要考慮冬小麥的物候特征。根據(jù)相關(guān)研究成果，黃河和駐(馬店)漯(河)交界是河南省重要的綜合地理分區(qū)和小麥生態(tài)分區(qū)界線[21,22]。本文以33°N緯線(駐馬店市與漯河市交界處)和黃河為界，把河南省黃淮海平原冬小麥主產(chǎn)區(qū)劃分為豫南、豫中、豫北三個地區(qū)(圖1)；大致以返青和乳熟期為基點(diǎn)(咨詢河南農(nóng)業(yè)科學(xué)院和國家小麥工程中心專家)把小麥全生育期劃分為：前期(從播種到返青，為裸土到半植被覆蓋期)、中期(從返青到乳熟，全植被覆蓋期)和后期(從乳熟到收割，為麥子變黃時(shí)期)，見表1。

表1 研究區(qū)冬小麥生長期劃分Tab.1 Division of growing period of winter wheat in research

圖1 河南省冬小麥主產(chǎn)區(qū)區(qū)域劃分Fig.1 Regional division of main wheat producing areas in Henan Province

在小麥生長前期，麥苗尚不能完全覆蓋裸土，NDVI不能完全反映墑情狀況，必須考慮土壤溫度因子。根據(jù)土壤的熱特性，在同一(季節(jié))時(shí)間(光熱條件下)，若土壤含水量高，其溫度要降低；若土壤含水量低，其溫度必升高?？紤]到NDVI為一個無量綱的量，這里對原始NDVI進(jìn)行以下修正：

(1)

在小麥生長中期，麥子生長茂盛，基本全部覆蓋地表，NDVI可以指示根部墑情。但正如前面所述，NDVI相對于墑情波動具有滯后性。圖2是豫中龍城站中期NDVI與10和20cm土壤含水量的曲線圖。NDVI與10、20cm土壤含水量的相關(guān)系數(shù)分別為-0.218，-0.213；將植被指數(shù)曲線前移5d，NDVI與10、20cm土壤含水量的相關(guān)系數(shù)分別為0.664，0.467。可見，NDVI滯后于墑情波動約4～6d。此外，眾所周知，在高植被覆蓋期，NDVI具有飽和性。需注意的是，這里的NDVI是利用遙感波段計(jì)算的原始值，由于受到太陽光照角度、觀測角度、云、氣溶膠等因素的影響，其往往包含著很多噪聲，并不都是植被指數(shù)的真實(shí)值(一般比實(shí)際值要小)，尤其是低值部分。

圖2 豫中龍城監(jiān)測站NDVI相對于墑情變化的滯后性Fig.2 The lag of NDVI relative to soil moisture change in Longcheng monitoring station

基于上述考慮，對于中期NDVI修改如下：

(2)

式中：KMNDVI指中期修正NDVI(MediumNDVI)；β為修正加系數(shù)；Tn為某像元某天的實(shí)際地表溫度；Tn-1為相應(yīng)像元該天的前一天的實(shí)際地表溫度；其他符號同前。

這里修改的前半部分(分子第一項(xiàng))主要作用是將NDVI的滯后性提前。其理由是如果當(dāng)天下雨或大面積灌溉，地表溫度肯定低于前一天的溫度，溫差較大時(shí)，修改效果明顯；如果不下雨或沒有灌溉，那么二者的溫度接近，修改效果可以忽略。后一部分(分子第二項(xiàng))的作用與前期類似，這里主要用于改善NDVI的飽和性。

在小麥生長后期，麥子將逐漸變黃，NDVI的值要逐漸變小，但并不一定意味著土壤濕度低。因此也要對后期NDVI進(jìn)行修正。杜曉等[23]依據(jù)對水的吸收率曲線和植被與土壤反射率曲線的分析，指出MODIS的b6、b7波段反映了含水土壤和植被的波譜信息，并構(gòu)造了地表含水量指數(shù)(SWCI)模型。在此基礎(chǔ)上，張紅衛(wèi)等[24]等考慮綠色植被葉綠素對光譜的吸收作用，引入b1波段信息，以b1/b7為修正系數(shù)，構(gòu)造了增強(qiáng)型土壤表層水分含水量指數(shù)(ESWCI)模型，并應(yīng)用于植被茂盛期的墑情監(jiān)測。小麥生長成熟期是一個葉綠素含量逐漸減少，作物由茂盛過渡到黃枯的階段，所以，與上述修正系數(shù)相反，本文對該階段的植被指數(shù)修改如下：

(3)

式中：KLNDVI指后期修改NDVI(LateNDVI)；γ為修正乘系數(shù)；b1、b2、b7為MODIS相應(yīng)波段光譜反射率。

將修正NDVI與實(shí)測墑情數(shù)據(jù)進(jìn)行相關(guān)分析，驗(yàn)證其有效性。在通過驗(yàn)證的基礎(chǔ)上，構(gòu)建修改NDVI時(shí)間序列，并與前三年的平均值進(jìn)行距平分析，得到監(jiān)測年與前三年平均相比的墑情分布時(shí)空規(guī)律。

2.2 技術(shù)路線與過程

綜上所述，本文的技術(shù)路線如圖3所示。

圖3 研究技術(shù)路線Fig.3 Research technical route

本文通過旬最大值合成(MVC)和改進(jìn)S-G濾波法[25]，去除噪聲數(shù)據(jù)，構(gòu)建反映冬小麥生長實(shí)際規(guī)律的、較光滑連續(xù)的修正NDVI生長曲線。

3 結(jié)果與分析

3.1 修正前后像元生長曲線對比

圖4表示豫中淮陽站2014-2015年植被指數(shù)修改前后的情況。

圖4 豫中淮陽監(jiān)測站2014-2015年NDVI修正前后的效果Fig.4 The effect comparisons of NDVI before and after modification of Huaiyang observation station at central Henan province in 2014-2015

可以發(fā)現(xiàn)，修正前后，植被指數(shù)曲線的走勢基本是一樣的，但修正后植被指數(shù)的變化幅度一般較小。

3.2 修正NDVI與實(shí)測墑情的相關(guān)分析

為了驗(yàn)證修正NDVI的有效性，分別在3個地區(qū)均勻選擇8個監(jiān)測站，將修正前后的植被指數(shù)與10和20cm實(shí)測墑情進(jìn)行了相關(guān)分析。其結(jié)果如表2所示。

由表2可知，原始NDVI與實(shí)測墑情的相關(guān)情況，后期比前期和中期要好，中期比前期要好。修正后NDVI與墑情相關(guān)的規(guī)律與原始NDVI基本類似，但修正NDVI在各時(shí)期的相關(guān)系數(shù)都有不同程度提高，說明了修正系數(shù)的有效性，修正后的NDVI更適宜作為農(nóng)田的墑情指數(shù)。

表2 NDVI修正前后與實(shí)測墑情的相關(guān)分析對比Tab.2 Comparison of correlation between NDVI before and after modified and measured soil moisture

續(xù)表2 NDVI修正前后與實(shí)測墑情的相關(guān)分析對比

注：*表示在0.05水平(雙側(cè))上顯著相關(guān)，**表示在0.01水平(雙側(cè))上顯著相關(guān)。

3.3 經(jīng)MVC和S-G濾波后像元生長曲線情況

圖5是豫北某像元修正NDVI曲線情況。

圖5 豫北某像元經(jīng)MVG合成和改進(jìn)S-G濾波后的修正NDVI生長曲線Fig.5 The growth curve of a Pixel modified NDVI after MVC synthesis and improved S-G filtering in the north of Henan province

3.4 基于修正NDVI距平分析的冬小麥墑情監(jiān)測

基于生長曲線得到距平值系列(圖6)，并進(jìn)行距平值等級劃分，最終得到監(jiān)測年(2014-2015年)與前三年均值相比，研究區(qū)域冬小麥田相對濕度的時(shí)空分布規(guī)律(圖7)。

圖6 研究區(qū)像元的距平值時(shí)間序列Fig. 6 The time series of the anomalies of the study area

圖7 2014-2015年河南省冬小麥主產(chǎn)區(qū)土壤濕度相對于前三年均值偏離情況的時(shí)空分布Fig.7 Temporal and spatial distribution of soil moisture relative to the mean of the previous three years in main winter wheat production areas in Henan province from 2014 to 2015

這里的關(guān)鍵是距平值(距平植被指數(shù))的等級劃分。實(shí)驗(yàn)數(shù)據(jù)(圖6)表明，距平值的分布區(qū)間約為[-1,1]，集中在[-0.2,0.2]之內(nèi)。如果等級劃分過粗，那么就看不出墑情的空間分布結(jié)構(gòu)；如果劃分過細(xì)，則會引起數(shù)據(jù)失真，不能反映真實(shí)情況。綜合相關(guān)研究文獻(xiàn)[18,26]和實(shí)驗(yàn)結(jié)果，本文把距平值劃分為7個等級：<-0.5，[-0.5,-0.15), [-0.15,-0.07)，[-0.07,0.07)，[0.07,0.15)，[0.15,0.5)，>0.5，分別表示相對濕度從干到濕的分布序列。其中，[-0.07,0.07)為與往年持平。需要說明的是，圖7所示的情況是監(jiān)測年相對于前三年平均的濕度偏離情況，并不是絕對的墑情等級。

從圖7中可以看出，2014-2015年的10月下旬、3月中旬、4月中旬至5月上旬，研究區(qū)冬小麥田的墑情基本與前三年平均持平；3月下旬、4月上旬和5月中旬至下旬稍偏干；11月中旬至12月上旬偏干；12月下旬至3月上旬偏濕。

為了驗(yàn)證上述結(jié)果的正確性，我們選取了明顯偏干的2014年11月下旬和明顯偏濕的2015年2月上旬作為實(shí)證。將上述兩個時(shí)期前三年的10和20 cm實(shí)測墑情(取最大值)的平均值與2015年相同時(shí)期實(shí)測墑情(取最大值)的值做比較(均值為被減數(shù))，其結(jié)果如圖8所示(等級劃分的相對比例與圖7一致)。從圖8(a)、8(c)與8(b)，8(d)、8(f)與8(e)的比較可以看出，監(jiān)測值與實(shí)測值分布規(guī)律基本是一致的，誤差基本都在一個量級內(nèi)；比較而言，監(jiān)測值與10 cm實(shí)測土壤水分符合度更好。需要說明的是，圖8的色標(biāo)適用于除8(b)、8(e)以外的其他圖。

圖8 相對濕度驗(yàn)證實(shí)例Fig.8 Comparison of relative moisture verification

4 結(jié) 論

土壤墑情是作物旱情診斷和灌溉管理的基礎(chǔ)參數(shù)，是農(nóng)作物產(chǎn)量的限制因子，監(jiān)測土壤墑情對提高水分利用效率和效益具有重要意義[27]。本文充分利用MODIS遙感的高時(shí)間分辨率特性，構(gòu)建基于區(qū)域和作物生育期自適應(yīng)修改的植被指數(shù)生長曲線，通過距平植被指數(shù)，實(shí)現(xiàn)大區(qū)域全生育期的冬小麥墑情監(jiān)測。研究表明：

(1)與NDVI值相比，基于地域和物候期分異的修正NDVI值與土壤水分的相關(guān)性更好，更適合作為農(nóng)田墑情的標(biāo)識指數(shù)；

(2)基于自適應(yīng)修正NDVI作物生長曲線的遙感墑情監(jiān)測是一種成本低、精度較好、適宜性較廣的有效方法；

(3)農(nóng)作物生長曲線是一條隱含諸多農(nóng)情信息的特征曲線，對其進(jìn)行農(nóng)情信息挖掘是一個很有價(jià)值的研究方向。

本文提出對中期NDVI修正的初衷之一是將其滯后性提前，但實(shí)驗(yàn)中我們發(fā)現(xiàn)，恰是在降雨的時(shí)候，往往也是云層遮擋較為嚴(yán)重的時(shí)候，MODIS的每日陸面溫度就存在較為嚴(yán)重的漏洞。為保證結(jié)果的科學(xué)性，這時(shí)我們采取的措施是不修改。從圖4可以看出，滯后性修改的效果不明顯，這在一定程度上影響了實(shí)驗(yàn)結(jié)果。此外，距平植被指數(shù)法雖然不能準(zhǔn)確地劃定墑情的等級，但如果時(shí)間序列積累足夠長，則可直接得到監(jiān)測年相對于常年的偏離情況，具有一定的實(shí)用參考價(jià)值。