利用FY-3衛星MWRI數據探測海冰分布

吳展開，王星東, ，王 成

(1. 河南工業大學信息科學與工程學院，河南 鄭州 450001；2. 中國科學院遙感與數字地球研究所，北京 100094)

北極是全球氣候系統的重要組成部分，是全球氣候變化的指示器。作為全球氣候的冷源，其大氣、海洋、海冰的運動直接或間接影響著全球水汽環流的強度、全球熱平衡和氣候變化，其中海冰是影響全球氣候系統的關鍵因素[1-3]。近幾十年由于全球變暖的影響，北極海冰呈逐年退縮、海冰厚度變薄、海冰范圍的最小記錄不斷被刷新、夏季出現北極航道通航等現象，使得連接亞洲、歐洲和北美洲的北極航線有望全線開通。北極航線作為北美洲、歐洲和東亞地區之間潛在的最快捷的海上通道，具有重要的航運價值和戰略價值。相比傳統航線，北極航線可規避海峽擁堵及海盜等多種問題，且大大縮短我國至歐洲地區的航程，因此對我國具有重要的經濟、軍事和戰略價值，也將有利于國家“一路一帶”倡議的推進[4]。

北極航道的開通受多種因素的制約及影響，其中受海冰的影響極為嚴重[4-5]。另外，海冰動態與全球氣候變化間的相互影響機制十分復雜，海-冰-氣相互作用過程的反饋機理還不十分清楚。但是，通過對歷史海冰觀測資料和同期海冰參數的比較，進行逐旬、月、季和年間距分析可以獲取海冰的動態變化信息；結合氣象觀測數據，可對海冰變化趨勢做出預測。

由于微波輻射計能全天候工作，且受天氣影響較弱，因此常用于海冰的研究。Andersen等[6]總結了針對 SSM/I 數據的7種海冰密集度的反演算法。其中，NASA-Team 算法[7]和Bootstrap算法[8]主要是基于低頻數據對海冰密集度進行反演，這類方法受天氣影響較弱，但空間分辨率也相對較低；NASA-Team2 算法[9-10]、SEA LION 算法[11]和 ARTSIST Sea Ice (ASI)[12]等算法均使用高頻數據，可得到空間分辨率相對較高的結果。ASI算法具有不需要額外數據源的優點，而且與利用其他通道的海冰密集度算法有相似的結果[13]。同樣，近年來國內學者，對海冰分布及海冰密集度的探測也做出了諸多研究成果，目前大多是基于AMSR-E數據及SSM/I進行[5,14]，而利用FY-3 MWRI數據對海冰的研究較少。

利用AMSR-E 數據基于ASI[15-16]、NASA-Team2 和Boot-strap Algorithm (ABA)[17]算法3種算法分別對北極海冰密集度進行反演，并將所得結果與2003年3—4月和2004年7—8月“北極星”號考察船船舶觀測數據進行對比分析，結果表明3種算法得到的結果與觀測數據的相關系數分別為0.80、0.79和 0.81。過去30年來，北極海冰范圍減少明顯，約以每年4.51×104km2的速度消退，且季節平均及月平均海冰范圍均表現出下降趨勢，4個季節中，夏季減少最為明顯[18]。1979—2006年3月平均海冰范圍以每10年2.8%的速度減少，而9月以每10年8.6%的速度減少[19]。通過對1978—2000年的微波數據進行分析研究，結果表明：夏季出現北極海冰范圍最小時的海冰主要由多年冰構成，且前后兩個10年間的平均海冰范圍及平均海冰面積分別減少了6.3%和11.0%[20]。1979—2008年北極平均無冰季節長度表現出了較為明顯的增加趨勢，其中1979—2006年間，以無冰季節長度每年1.1 d的速度增長，2007和2008年無冰季節長度分別達到了168和158 d，且2007年的無冰季節長度為1979—2008年中的最高紀錄[21]。

在分析海冰(一年冰、多年冰)與海水的輻射率變化特征的基礎上，提出基于微波成像儀(MWRI) 19 GHz水平極化和垂直極化的極化差(極化差最大)結合Otsu算法(確定海冰和海水的分類閾值)來探測北極海冰分布的新方法，利用2016年1月共31 d的數據得到北極海冰月平均分布圖。

1 研究區域與數據源

北極地區一般是指北緯66°34′以北的地區，由北冰洋、歐亞大陸的邊緣陸地及眾多的島嶼島礁組成，面積約2.1×107km2。北冰洋面積約1.4×107km2，約占世界大洋總面積的3%，約占北極地區面積的67%，其大部分區域常年被海冰覆蓋。白令海峽是北冰洋與太平洋連通的路徑，在格陵蘭島西側的加拿大北極群島有許多與巴芬灣相連接的通道，北極西北航道就在此片島嶼中穿行。

本文使用2016年1月的FY-3 衛星的MWRI數據對北極地區的海冰分布進行分析研究。FY-3氣象衛星能夠全天候連續獲取全球范圍內的大氣、地表及海洋表面多種特征參數，并廣泛應用于三維大氣探測、地表遙感特征提取等各方面研究。其中，MWRI是裝載于FY-3衛星上的重要遙感儀器。MWRI可用于地表熱輻射的監測，通過使用天線繞軸旋轉形成圓錐形跨軌的掃描方式，在10.6～89 GHz間的多個固定頻段上對來自大氣及地表的雙極化電磁輻射進行接收，并通過對數據的加工及分析，用于降雨反演、水汽監測、地表及海表上各類目標物的研究，且可為天氣預報及環境監測提供參考數據。MWRI傳感器具體參數見表1。

表1 FY-3 MWRI傳感器若干參數

2 海冰分布探測方法

2.1 極化差

地物波譜特征是遙感反演的基礎，不同地物的波譜特征曲線反映出了不同地物對電磁波反射的差異。其中，一年冰、多年冰和開闊水的光譜特征曲線如圖1所示[22-23]。由亮溫的定義可知，當物體的物理溫度一定時，一定頻率下亮溫大小只與物體的輻射率有關。因為垂直極化和水平極化的電磁波在海冰表面或海水表面同時發出時，對同一物體的物理溫度相同，因此極化差只受輻射率大小的影響。通過圖1中的仿真結果可知，在微波成像儀(MWRI)的這幾個頻段中，19 GHz的一年冰和多年冰中的極化差值相差不多且與開闊水面相比差別最大，即可利用19 GHz的極化差P=19V-19H識別冰和水進而得到海冰的分布。即如果P>T，那么此像素為開闊水，否則為冰。

圖1 海冰和海水發射率與頻率的關系

2.2 Otsu算法

Otsu算法，又被稱為大津法。該方法是根據圖像數據值之間的方差選取閾值。其原理如下：

首先假設閾值具有L級數據的圖像劃分為兩類：C0∈[0,T]及C1∈[T+1,L+1]圖像各數據級對應的概率為

(1)

(2)

(3)

C0和C1類的均值分別為

(4)

(5)

(6)

海冰和海水的最佳分類閾值T*應使類間方差最大，即

(7)

2.3 極化差結合Otsu算法的海冰分布探測

利用Otsu算法對19 GHz極化差P進行處理，得到海冰和海水的分類閾值，進而得到北極海冰的分布信息。圖2為海冰分布探測方法的流程圖，基本步驟如下：

(1) 數據預處理：輻射定標、掩膜、異常數據點處理(當某個數據點與以其為中心的窗口內的數據點的中值相差太大時，該數據點被定義為異常數據點，用窗口內的數據點的中值代替)。

(2) 計算極化差：利用19 GHz的垂直極化和水平極化計算極化差P=19V-19H。

(3) 確定海冰和海水分類的最優閾值：選擇樣本點，基于Otsu算法對極化差P進行處理得到海冰和海水分類的最優閾值T。

(4) 得到海冰分布圖：根據閾值T對極化差P圖像進行分類，得到北極海冰分布圖。

(5) 驗證：所得結果與美國冰雪數據中心結果進行對比驗證。

3 結果與驗證

以2016年1月共31 d的FY-3 MWRI 19 GHz水平極化數據和垂直極化數據進行海冰分布探測，對19 GHz水平極化數據和垂直極化數據進行上述處理，得到19 GHz下的極化差異P，再選取樣本區域，應用Otsu算法得到閾值為48.79 K。根據此閾值對圖像進行分類得到2016年1月平均北極海冰分布圖(如圖3(a)所示)。獲取同期1月NSIDC發布的北極海冰分布圖(如圖3(b)所示)。對比圖3 可以看出，2個結果海冰外緣線分布狀況基本一致；其中，MWRI數據反演得到的結果在高緯度地區的陸地邊緣的海冰數量要多于NSIDC結果，由于陸地邊緣在一定程度上會對海冰的判別造成誤差，初步判斷可能是美國冰雪數據中心為了消除陸地邊緣造成的誤差，對所有的陸地邊緣進行了處理，但本文認為處于冬季的高緯度陸地邊緣存在一些數量的海冰也是合理的，因此沒有對其做進一步的處理；在鄂霍次克海區域(圖幅的左上偏右區域)NSIDC結果海冰數量多于MWRI數據得到的結果，造成這種差異的原因可能如下：①由于微波數據是單軌成像，進行數據拼接，2種衛星過境時間不同，處于鄂霍茨克海區域軌道成像時間可能差異較大，因此造成這種誤差；②該區域位于西北太平洋，與開闊水域鏈接暢通，可能受季風及洋流的影響較大，因此可能引起海冰移動，進而造成兩者的誤差；③2種微波數據的本身特性存在某種差異。當然，為了進一步地分析驗證FY-3 MWRI數據及本文提出的算法應用到海冰識別方面的精確度，對于以上幾個問題的深入研究，將會是筆者接下來努力的方向。

為了更為直觀地對比2種海冰面積結果在1月份的變化趨勢及差異，統計得到海冰范圍統計表(見表2)及變化趨勢圖(如圖4、圖5所示)。

圖2 19 GHz數據的海冰分布探測流程

圖3 2016年1月北極海冰平均分布

圖4 基于FY-3 2016年1月海冰范圍變化趨勢

圖5 NSIDC 2016年1月海冰范圍變化趨勢

日期NSIDC海冰范圍(×106)/km2FY-3海冰范圍(×106)/km2FY-3與NSIDC海冰范圍同比差值/(%)2016-01-0112.656212.18713.70612016-01-0212.746212.33093.25832016-01-0312.724312.37322.75922016-01-0412.820612.47002.73482016-01-0512.836812.39403.44952016-01-0613.011212.41044.61742016-01-0713.143112.43905.35572016-01-0813.015612.59453.23522016-01-0913.135012.57814.23962016-01-1013.348012.67205.06952016-01-1113.505012.65936.26152016-01-1213.392512.76714.66912016-01-1313.513712.95954.10112016-01-1413.455613.05732.95992016-01-1513.460613.07072.89612016-01-1613.521813.01603.74042016-01-1713.577513.02844.04392016-01-1813.663713.04074.55922016-01-1913.650613.08954.11032016-01-2013.625613.12073.70512016-01-2113.715613.19983.76052016-01-2213.704313.17463.86512016-01-2313.646213.20463.23572016-01-2413.748713.27573.44002016-01-2513.742513.27343.41322016-01-2613.866213.22854.59862016-01-2713.863113.21004.71122016-01-2813.861213.22144.61602016-01-2913.958713.25265.05842016-01-3013.928713.43962.79322016-01-3113.933713.33604.2892平均值13.493312.90564.3550極大值13.958713.43963.7185極小值12.656212.18713.7061極差1.30251.25253.8387

由表2及圖4、圖5可知，基于FY-3 MWRI數據利用Otsu算法結合19 GHz頻段極化差得到的2016年1月海冰范圍，與NSIDC提供的結果在日變化趨勢上基本一致，其中MWRI得到的結果，以平均每天0.038 4×106km2的速度增加，NSIDC結果以平均每天0.041 9×106km2的速度增加，二者的增長速度僅相差0.003 5×106km2；NSIDC 1月海冰范圍平均為13.493 3×106km2，MWRI數據結果1月海冰范圍平均為12.905 6×106km2，二者同比相差4.355%；NSIDC 1月海冰范圍極差為1.302 5×106km2，MWRI數據結果1月海冰范圍極差1.252 5×106km2，二者同比相差3.838%；NSIDC 1月海冰范圍單日極大值為13.958 7×106km2，出現日期為1月29日，MWRI 1月海冰范圍單日極大值為13.439 6×106km2，出現日期為1月30日，造成二者的差異可能是由于2個衛星過境時間不同，成像不同步引起的，但綜合來看二者的極大值比較接近，同比相差僅為3.718%；NSIDC 1月海冰范圍單日極小值為12.656 2×106km2，出現日期為1月1日，MWRI 1月海冰范圍單日極小值為12.187 1×106km2，出現日期為1月1日，二者較為同步，且同比相差僅為3.706%。因此，結合以上幾個指標來看，基于本文算法得到的北極海冰分布范圍與美國冰雪中心的海冰分布范圍及日變化趨勢基本一致，這在一定程度上說明本文提出的方法可行。

4 結 語

以2016年1月份共31 d FY-3 MWRI數據為例，基于海水與海冰在不同頻段的極化差異的特性不同，選取19 GHz頻段的數據，通過計算極化差，然后利用Otsu算法得到了海冰和海水的分類閾值，進而獲取了北極的海冰分布信息。并與美國冰雪數據中心海冰分布的業務化產品進行對比分析，結果表明二者海冰外緣線分布基本一致，海冰面積也相差不大，海冰范圍日變化趨勢也較為同步。相對于其他光學及高頻數據海冰分布的反演算法，本文提出的基于19 GHz頻段極化差異特性的算法，更能用于大范圍長時間序列的海冰監測與研究；同樣在算法的實現過程中僅需確定一個分類閾值即可，具有較高的運算效率，對于即時性的大量數據的海冰信息發布與監測研究，具有明顯的優勢。因此表明本文提出的這種算法可行，尤其是在更高效率的運用國產衛星數據的層面具有一定價值與意義。