何勤業(yè)，盧涵宇*，袁詠儀，盧天健

(1.貴州大學 大數(shù)據(jù)與信息工程學院, 貴州 貴陽 550025；2.貴州六盤水三力達科技有限公司, 貴州 六盤水 532001)

0 引言

土壤水分是氣候系統(tǒng)中的一個關(guān)鍵狀態(tài)變量，對理解陸地-大氣能量、水分和碳循環(huán)通量具有重要意義[1]。土壤水分在洪水和干旱預報、氣候變化、天氣預報、水資源管理和農(nóng)業(yè)生產(chǎn)等方面有著廣泛的應(yīng)用，對人類社會和環(huán)境有著重要的影響[2]，因此，精確度和高空間分辨率的土壤水分全球測繪對于滿足這些應(yīng)用的需求至關(guān)重要。

目前最廣泛應(yīng)用的測量地表土壤水分的方法是被動微波遙感，它能全天候全天時觀測數(shù)據(jù)。由于微波遙感的1.41 GHz(L波段)具有更好的穿透大氣和植被覆蓋能力及受射頻干擾(radio grequency interference，RFI)影響較小特點，因此在大多數(shù)與土壤水分相關(guān)的研究以該頻段為主要研究對象。第一個搭載L波段微波輻射計的衛(wèi)星是在2009年11月2日由歐空局發(fā)射的土壤水分和海洋鹽度(soil moisture and ocean salinity，SMOS)衛(wèi)星。該衛(wèi)星可以提供地表的雙極化和大范圍入射角(5°～60°)的亮溫數(shù)據(jù)，同時這也是第一顆監(jiān)測全球土壤水分和海水鹽度的衛(wèi)星[3]；但是L波段輻射計的空間分辨率只有40 km，并不足以滿足水文氣象、生態(tài)學、水資源管理等應(yīng)用超過空間分辨率為10 km的使用需求，所以美國國家航空航天局(national aeronautics and space administration，NASA)于2015年1月31日在加利福尼亞州發(fā)射土壤水分主動-被動(soil moisture active passive，SMAP)任務(wù)[4]。該衛(wèi)星搭載著空間分辨率為3 km的L波段雷達與空間分辨率為40 km的L波段輻射計，可以通過幾種主被動結(jié)合的降尺度方法，提供空間分辨率為9 km的全球土壤水分[4]，但是，在2015年7月7日SMAP雷達故障，無法提供數(shù)據(jù)。而為了能繼續(xù)提供高分辨率數(shù)據(jù)，SMAP使用了Backus-Gilbert(BG)optimal interpolation對原始的SMAP Level 1B亮溫產(chǎn)品(L1B_TB)進行插值處理從而得到空間分辨率為9 km的亮溫數(shù)據(jù)[5]。

目前土壤水分的降尺度方法可以大致分為基于經(jīng)驗?zāi)Ｐ偷姆椒?、基于?shù)理模型的方法和基于微波遙感數(shù)據(jù)融合的方法。其中基于經(jīng)驗?zāi)Ｐ偷姆椒ㄒ全@取的高分辨率的地表參數(shù)，再利用獲取的地表參數(shù)建立與微波土壤水分的模型，并將模型參數(shù)用于微波土壤水分的降尺度處理；但降尺度結(jié)果的準確程度主要取決于此類經(jīng)驗?zāi)Ｐ偷木_性，而且為了能更好的構(gòu)建經(jīng)驗?zāi)Ｐ瓦€需要獲取大量的實地數(shù)據(jù)，因此限制了經(jīng)驗?zāi)Ｐ驮谳^大區(qū)域內(nèi)的應(yīng)用[6]。而基于數(shù)理模型的方法主要是通過數(shù)理統(tǒng)計模型(如譜分析、多重分形和地學統(tǒng)計等空間插值方法)對微波數(shù)據(jù)進行降尺度處理，其缺陷在于這些數(shù)理模型往往只是對單一的微波數(shù)據(jù)進行空間上的處理，缺少額外的輔助數(shù)據(jù)作為支撐，并不能確保降尺度的精確程度[7]。本文使用的是基于微波遙感數(shù)據(jù)融合的方法，這種方法是通過將2種不同分辨率的微波數(shù)據(jù)進行融合以達到降尺度的目的。其中文獻[8]利用傅里葉變換的方法將主被動微波亮溫數(shù)據(jù)進行結(jié)合，但是SMAP的主動微波雷達已經(jīng)損壞，而且傅里葉變換在分析微波信號這類非平穩(wěn)信號時存在局限性，會導致降尺度結(jié)果出現(xiàn)較大的誤差。文獻[9]使用時間序列觀測回歸分析降尺度法研究被動微波中L波段與S波段亮溫的線性關(guān)系，將S波段作為輔助數(shù)據(jù)對L波段的亮溫進行降尺度處理，但是目前還沒有搭載S波段的衛(wèi)星，而且S波段的分辨率較低遠遠不及10 km，無法滿足土壤水分的應(yīng)用需求。本文在上述2種方法的基礎(chǔ)上進行了改進，首先是在選用微波遙感數(shù)據(jù)時選用了被動微波的X波段作為L波段亮溫降尺度的輔助波段，這是因為X波段的分辨率為10 km，可以滿足土壤水分的應(yīng)用需求，而且X波段的亮溫數(shù)據(jù)可以直接通過AMSR2得到；其次在對X波段與L波段進行融合時使用了二維離散小波變換(2-D DWT)的方法，這種方法在分析非平穩(wěn)信號時沒有傅里葉分析的局限，可以更加精準的對微波信號進行分析與變換。具體步驟是通過2-D DWT分辨率為9 km的X波段亮溫與分辨率為36 km的L波段亮溫進行分解與融合，從而得到降尺度后的L波段亮溫數(shù)據(jù)。實驗結(jié)果表明降尺度后的L波段亮溫可以滿足反演土壤水分的需求，且反演的土壤水分的精度與SMAP 9 km亮溫反演的土壤水分近似，具有極高的適用性。

本文主要通過結(jié)合L波段與X波段亮溫的特性來實現(xiàn)亮溫降尺度的研究，并且將小波變換技術(shù)引入到了亮溫降尺度和土壤水分反演這2個領(lǐng)域，從而拓展了亮溫降尺度和土壤水分反演的研究方向。

1 實驗區(qū)域與數(shù)據(jù)

本研究使用的L波段的亮溫數(shù)據(jù)采用了SMAP-3級輻射計(SPL3SMP版本7)中的亮溫數(shù)據(jù)，SMAP數(shù)據(jù)于“美國國家冰雪數(shù)據(jù)中心”(https://nsidc.org/data)進行下載。X波段的亮溫數(shù)據(jù)則使用了AMSR2中10.65 GHz的L3亮溫數(shù)據(jù),其數(shù)據(jù)來源于“日本宇宙航空研究開發(fā)機構(gòu)”(https://gportal.jaxa.jp/gpr/)，并且選擇了位于帕里站點獲取的2016年1月1日到12月31日的土壤水分數(shù)據(jù)作為驗證數(shù)據(jù)，該數(shù)據(jù)來源于“國家青藏高原科學數(shù)據(jù)中心”(http://data.tpdc.ac.cn)。

1.1 實驗區(qū)域

本次研究區(qū)域選擇為帕里區(qū)域。該區(qū)域的土地覆蓋以草地、裸露地面為主，地形平坦，平均年降水量小于400 mm，且降雨大多集中在5—10月。在帕里區(qū)域內(nèi)一共部署了25個觀測站點，在本研究中只用到了其中的8個觀測站點，并分為了A、B共2個9 km網(wǎng)格，其中黑色三角點表示土壤水分的觀測站點[10]。帕里研究區(qū)土壤水分采樣點分布如圖1所示。

圖1 帕里研究區(qū)土壤水分采樣點分布圖

1.2 SMAP數(shù)據(jù)

土壤水分主動被動(SMAP)衛(wèi)星于2015年1月31日發(fā)射，運行周期為2～3 d，上升軌道時間為18:00，下降軌道時間為上午06：00[4]。SMAP衛(wèi)星分別搭載了L波段(1.41 GHz)輻射計和L波段(1.26 GHz)雷達，入射角度為40°[4]，可以提供空間分辨率為9、36 km的亮溫數(shù)據(jù)。而空間分辨率為9 km的亮溫數(shù)據(jù)在L波段雷達故障后就無法獲取，因此基于對BG 最佳插值的應(yīng)用，將SMAP的原始觀測數(shù)據(jù)進行了修正從而獲得了空間分辨率為9 km的亮溫數(shù)據(jù)。

在本研究中使用了SMAP衛(wèi)星在帕里實驗區(qū)域的分辨率36、9 km的亮溫數(shù)據(jù)，以及用亮溫反演土壤水分需要的輔助數(shù)據(jù)，為了方便小波變換的應(yīng)用選取的范圍為16×16的9 km網(wǎng)格，具體的經(jīng)緯度為27.259 5°～28.534 2°、88.879 7°～90.373 4°。

1.3 AMSR2數(shù)據(jù)

GCOM-W由Japan Aerospace Exploration Agency(JAXA)在2012年5月18日發(fā)射，軌道位于地球上方約700 km處，上升軌道時間為13：30,下降軌道時間為01：30[11]。GCOM-W衛(wèi)星攜帶的無源微波輻射計Advanced Microwave Scanning Radiometer 2(AMSR2)是由AMSR-E改進而來，AMSR2的天線旋轉(zhuǎn)周期為1.5 s[11]，可以獲取1 450 km范圍內(nèi)的數(shù)據(jù)，使得AMSR2能以2 d為一周期獲取全球99%以上的數(shù)據(jù)[11]。AMSR2的入射角為55°，有6個不同頻率的微波波段(6.9、10.7、18.7、23.8、36.5、89 GHz)，可以提供空間分辨率為10、25 km的雙極化亮度溫度。本研究中使用了10.65 GHz(X波段)采集的空間分辨率為10 km亮度溫度作為降尺度的輔助數(shù)據(jù)。

1.4 現(xiàn)場測量

土壤水分是評價亮溫降尺度結(jié)果的一個重要標準，因此本研究使用了帕里區(qū)域中按照布設(shè)的8個觀測站點(如圖1黑色三角點所示)獲取的土壤水分。同時為了方便進行土壤水分的對比驗證還對這8個觀測站點根據(jù)其所在的網(wǎng)格進行了平均計算，從而得到這2個9 km網(wǎng)格的土壤水分。

2 研究方法

雖然L波段對于土壤水分含量的探測更加精準，但是SMAP的L波段的空間分辨率比AMSR2的X波段低，因此為了提高L波段的空間分辨率，使得L波段更加泛用，本文通過二維離散小波變換(2-D DWT)將L波段與X波段的亮溫數(shù)據(jù)融合，從而得到的同時擁有L波段亮溫數(shù)據(jù)特征和X波段空間分辨率的亮溫數(shù)據(jù)。同時為了對這種降尺度方法的結(jié)果進行對比還加入了SMAP 9 km亮溫數(shù)據(jù)(通過Backus-Gilbert(BG)optimal interpolation獲取)。為了方便將降尺度后的亮溫數(shù)據(jù)與SMAP的分辨率為9 km的亮溫數(shù)據(jù)進行對比，還需要對AMSR2的X波段亮溫數(shù)據(jù)進行插值，得到X波段分辨率為9 km的亮溫數(shù)據(jù)。

2.1 Backus-Gilbert(BG)最佳插值

Backus-Gilbert(BG)最佳插值基于Backus-Gilbert理論，是一種對成像微波輻射計數(shù)據(jù)進行插值的技術(shù)，因為保持與微波輻射計相關(guān)的天線增益函數(shù)的空間分辨率，所以其插值過程稱為最優(yōu)的。SMAP 9 km的亮溫數(shù)據(jù)是通過對SMAP 1B級亮溫產(chǎn)品(L1B_TB)中的天線溫度(TA)的測量值進行Backus-Gilbert(BG)最佳插值得到的，其插值生成的TA數(shù)據(jù)在經(jīng)過校準后，就能在9 km網(wǎng)格上生成SMAP 1C級亮溫產(chǎn)品(L1C_TB_E)。

2.2 二維離散小波變換(2-D DWT)

小波變換作為一種分析工具，可以對不同空間分辨率下的亮溫數(shù)據(jù)、土壤水分數(shù)據(jù)等地球物理分布變量數(shù)據(jù)進行特性變化的分析。

連續(xù)小波變換(CWT)可描述為：假設(shè)存在一個能量有限的信號空間為L2(R)，且f(t)∈L2(R)，可以用小波基將f(t)展開[12]，得到f(t)的連續(xù)小波變換Wf(a,τ)，

(1)

式中：a是伸縮因子；τ是平移因子；函數(shù)ψa,τ(t)為小波基函數(shù)，其逆變換為

(2)

離散小波變換(DWT)是通過對CWT進行二進離散化實現(xiàn)的，可以將信號分解低頻與高頻2個部分到離散小波基的方法一般是將平移因子a和伸縮因子τ變成冪級數(shù)結(jié)構(gòu)[13]，即

(3)

式中：a0≠1；τ0是常數(shù)；j,k∈Z對應(yīng)的離散小波基為

(4)

最終f(t)的離散小波變換為

(5)

2-D DWT就是將二維數(shù)據(jù)分別變成行數(shù)據(jù)與列數(shù)據(jù)后再進行一維離散小波變換。二維離散小波變換主要通過Mallat前向金字塔算法的擴展實現(xiàn)，二維小波變換原理如圖2所示。每進行一次一維小波變換，原始數(shù)據(jù)就會分成高頻與低頻2個分量，因此二維數(shù)據(jù)在進行一次二維小波變換后，原始數(shù)據(jù)將會被分成4個分量，分別是原始二維數(shù)據(jù)的低頻分量LL和高頻分量HL、LH、HH。

圖2 二維小波變換原理圖

其中HL包含表示原始數(shù)據(jù)在水平方向上的高頻分量，LH表示原始數(shù)據(jù)在垂直方向上的高頻分量，HH表示原始數(shù)據(jù)在對角處的高頻分量，這些高頻分量表示的是原始信號的邊緣信息，而LL表示原始數(shù)據(jù)的低頻分量，是原始數(shù)據(jù)的中心部分，代表了原始數(shù)據(jù)的主要信息。而為了能夠獲得更加精確的小波系數(shù)，可以將第一次二維小波變換后獲得的近似值(低頻分量)LL作為下一次二維小波變換的輸入數(shù)據(jù)，最多可進行二維小波變換p=log2N次，N為原始數(shù)據(jù)的長度。在實際應(yīng)用中，分解的級數(shù)也與小波基的選取有關(guān)。

本研究使用了db5小波作為小波基，并進行了2次2-D DWT，最后又使用了小波逆變換重構(gòu)出了新的亮溫數(shù)據(jù)。之所以使用db5小波作為小波基是因為db5小波具有較好的正則性，且消失矩較大，比較光滑，具有較好的頻域局化能力，同時還不會使得時域的緊支撐性過低從而增大計算量使得實時性變差。而進行2次二維小波變換的原因是為了在不失真的情況下獲得更準確的小波系數(shù)。其中二維小波的重構(gòu)過程為：先對二維小波變換后獲得的高低頻分量的每一列進行一維小波逆變換，再對逆變換后的結(jié)果的每一行進行一維小波逆變換，即可得到重構(gòu)的數(shù)據(jù)。

本研究的技術(shù)路線如圖3所示。

圖3 技術(shù)路線

2.3 土壤水分反演

由于沒有實地觀測的分辨率為9 km亮溫數(shù)據(jù)，無法對降尺度后的L波段的亮溫數(shù)據(jù)進行直接有效的驗證。因此用單通道檢索算法(single-channel algorithm，SCA)和地表輔助數(shù)據(jù)[14]。該算法主要是利用對土壤水分變換敏感的亮溫和輔助數(shù)據(jù)(地表溫度、地表土壤質(zhì)地、地表覆蓋情況、植被含水量)來實現(xiàn)對土壤水分的反演。

SCA算法主要是基于零階近似模型ω-τ構(gòu)建的，其模型為

(6)

其實現(xiàn)過程具體分為以下5個步驟：

(7)

② 通過地表覆蓋分類圖確定對應(yīng)的植被類型和參數(shù)b。

③ 根據(jù)地表覆蓋分類圖確定的植被類型獲取植被單次反照率ω和植被光學厚度τ，通過公式(8)、(9)計算出植被透過率γP和植被含水量VWC，用來校正植被影響

(8)

τ=b·VWC。

(9)

④ 在將植被的影響去除后，使用Choudhury模型[15]對土壤粗糙度進行校正，獲取光滑表面的發(fā)射率。

⑤ 最后用菲涅爾方程計算土壤介電常數(shù)后，使用Dobson模型[16]估算出土壤水分含量。

該算法因為使用了大量的輔助數(shù)據(jù)來去除土壤溫度、大氣溫度、植被和地表粗糙度的影響，所以反演得到的土壤水分數(shù)據(jù)精度很高，同時也是目前SMAP衛(wèi)星的官方的土壤水分反演算法。在本研究中的土壤水分反演用到的地表輔助數(shù)據(jù)均由SMAP數(shù)據(jù)集得到，并根據(jù)實驗區(qū)域內(nèi)的植被分布選取單次散射反照率ω和植被參數(shù)b，且b=0.13，ω=0.05。在將亮溫數(shù)據(jù)反演得到土壤水分數(shù)據(jù)后，就可以與實地探測土壤水分數(shù)據(jù)進行對比，根據(jù)土壤水分數(shù)據(jù)的對比的結(jié)果就可以確定降尺度后亮溫數(shù)據(jù)的精度。

3 實驗結(jié)果與分析

本次實驗為了在獲取分辨率為9 km的L波段亮溫數(shù)據(jù)的同時證明加入輔助數(shù)據(jù)對L波段亮溫數(shù)據(jù)降尺度處理的重要性，以及驗證使用2-D DWT對亮溫降尺度的有效性，因此計劃借助AMSR2的X波段亮溫數(shù)據(jù)作為輔助數(shù)據(jù)，將SMAP 36 km的亮溫數(shù)據(jù)降尺度到9 km，并與SMAP 9 km的亮溫數(shù)據(jù)進行比較。由于AMSR2的X波段只有10 km的亮溫數(shù)據(jù)，因此需要使用線性插值的方法將其改變?yōu)? km的亮溫數(shù)據(jù)，以此滿足小波變換降尺度方法的要求。

本次研究選取的時間段為2016-01-01—12-31日，共計366 d；但是，由于SMAP和AMSR2的運行周期不一樣，因此為了保證在實驗區(qū)域內(nèi)SMAP和AMSR2都能提供亮溫數(shù)據(jù)，在上午時間段選取了55 d作為研究樣本。

3.1 亮溫降尺度結(jié)果

本研究中使用了二維離散小波變換對L波段亮溫進行了降尺度，其降尺度的亮溫圖像如圖4所示，顯示了在2016-06-22—09-26，期間，以32 d為間隔選取了4 d的L波段亮溫在上午時間段V極化和H極化的降尺度結(jié)果。其中第一行是SMAP 36 km原始的L波段亮溫圖像，第二行是AMSR2的9 km X波段亮溫圖像，第三行是二維離散小波變換降尺度的L波段亮溫圖像, 第四行是為了驗證算法加入的SMAP 9 km的L波段亮溫圖像。通過比較可以發(fā)現(xiàn)，第三行的降尺度圖像基本上可以呈現(xiàn)出SMAP 36 km亮溫圖像的空間分布情況與AMSR2 9 km亮溫圖像的空間結(jié)構(gòu)，并與SMAP 9 km的亮溫圖像相比，二維離散小波變換的降尺度圖像的異質(zhì)性更高，與原始觀測的36 km亮溫圖像的關(guān)系更緊密，空間上的細節(jié)更多。

(a)V極化

其亮溫降尺度的圖像分析結(jié)果見表1，分別計算了SMAP 9 km的亮溫圖像與二維離散小波變換的降尺度亮溫圖像在不同極化下的信息熵、平均梯度、標準差、空間頻率。結(jié)果表明通過2-D DWT得來的降尺度圖像在不同極化的條件計算得到的標準差和空間頻率均高于SMAP 9 km的亮溫圖像，說明通過2-D DWT降尺度得到的亮溫圖像的有用信息更多，圖像空間的總體活躍度更高。而在信息熵和平均梯度的計算上2-D DWT圖像和SMAP 9 km圖像各有優(yōu)劣，總體持平，說明這2種亮溫圖像在平均信息量和圖像對比度上有著近似的特征。這些情況表明二維離散小波變換的亮溫圖像有著更多的細節(jié)表達能力，圖像信息更多，圖像空間的活躍程度更高。

表1 亮溫降尺度圖像分析結(jié)果

3.2 土壤水分降尺度結(jié)果和驗證

為了驗證降尺度的結(jié)果，需要利用SCA算法與SMAP提供的輔助數(shù)據(jù)，對2-D DWT降尺度后的9 km亮溫數(shù)據(jù)和SMAP提供的9 km亮溫數(shù)據(jù)進行表面土壤水分反演。

反演的土壤水分和站點土壤水分的長時間序列圖與亮溫反演土壤水分的散點圖如圖5所示。其中區(qū)域A和區(qū)域B分別是選取的2個不同的觀測區(qū)域。為了確保土壤水分反演的精度，只反演了地表溫度大于0 ℃的時期。從長時間序列圖中可知這2種不同的亮溫反演的土壤水分與站點觀測的土壤水分有很高的相似度，特別是V極化下反演的土壤水分相似度更高。從散點圖可知，用降尺度后的亮溫反演的土壤水分和用SMAP 9 km的亮溫反演的土壤水分極為相似，相關(guān)性最高為0.984 8，最低為0.970 3，表明兩者之間相似度極高，而且偏差值最高僅為0.005 6，表明兩者差值極低。其均方根誤差最高為0.008 3 cm3/cm3，表明兩者間的離散度非常低。以上情況說明用2-D DWT進行降尺度的亮溫反演的土壤水分與SMAP 9 km亮溫反演的土壤水分有著極高的相似度。

(a)V極化

反演土壤水分與站點土壤水分的對比統(tǒng)計結(jié)果見表2。將2種亮溫在不同區(qū)域、不同極化的情況下反演的土壤水與相對應(yīng)的站點觀測的土壤水進行了對比，計算了反演的土壤水分與站點土壤水分的相關(guān)系數(shù)(R)、偏差值(BIAS)、均方根誤差(RMSE)、無偏均方根誤差(ubRMSE)。比較結(jié)果表明，反演的土壤水分偏差值均在0.2以內(nèi)，且相關(guān)系數(shù)均大于0.68，均方根誤差在0.2 cm3/cm3以內(nèi)，無偏均方根誤差最高僅為0.021 9 cm3/cm3，表明反演的土壤水分與站點觀測的土壤水分的有著相似的性能，且在V極化情況下能達到最優(yōu)。

表2 反演土壤水分與站點土壤水分對比統(tǒng)計結(jié)果

4 結(jié)語

本文通過二維離散小波變換借助X波段的高分辨率亮溫對L波段亮溫進行降尺度處理，并與SMAP 9 km的亮溫數(shù)據(jù)進行了對比。其結(jié)果表明通過二維離散小波變換得來的亮溫圖像在整體上要優(yōu)于SMAP 9 km的亮溫圖像，且與SMAP 36 km的亮溫圖像的聯(lián)系更加緊密，并且在反演的土壤水分的結(jié)果對比中兩者結(jié)果十分近似，且都與站點觀測的土壤水分有很高的相關(guān)度，無偏均方根誤差都在0.02左右。這些數(shù)據(jù)說明在對L波段亮溫進行降尺度時，可以使用二維離散小波變換對L波段亮溫進行降尺度處理，而且精度與SMAP 9 km的亮溫數(shù)據(jù)相比略有提高。本文的研究只限于帕里區(qū)域，且在反演土壤水分時只在2個像元上進行土壤水分反演結(jié)果驗證，因此今后的工作需要將該方法改進，爭取應(yīng)用到全球區(qū)域的亮溫降尺度，為亮溫降尺度方法提供一種新的參考思路。