吳展開(kāi)，王星東, ，王 成

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

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

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

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

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

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

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

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

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

表1 FY-3 MWRI傳感器若干參數(shù)

2 海冰分布探測(cè)方法

2.1 極化差

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

圖1 海冰和海水發(fā)射率與頻率的關(guān)系

2.2 Otsu算法

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

首先假設(shè)閾值具有L級(jí)數(shù)據(jù)的圖像劃分為兩類：C0∈[0,T]及C1∈[T+1,L+1]圖像各數(shù)據(jù)級(jí)對(duì)應(yīng)的概率為

(1)

(2)

(3)

C0和C1類的均值分別為

(4)

(5)

(6)

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

(7)

2.3 極化差結(jié)合Otsu算法的海冰分布探測(cè)

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

(1) 數(shù)據(jù)預(yù)處理：輻射定標(biāo)、掩膜、異常數(shù)據(jù)點(diǎn)處理(當(dāng)某個(gè)數(shù)據(jù)點(diǎn)與以其為中心的窗口內(nèi)的數(shù)據(jù)點(diǎn)的中值相差太大時(shí)，該數(shù)據(jù)點(diǎn)被定義為異常數(shù)據(jù)點(diǎn)，用窗口內(nèi)的數(shù)據(jù)點(diǎn)的中值代替)。

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

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

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

(5) 驗(yàn)證：所得結(jié)果與美國(guó)冰雪數(shù)據(jù)中心結(jié)果進(jìn)行對(duì)比驗(yàn)證。

3 結(jié)果與驗(yàn)證

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

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

圖2 19 GHz數(shù)據(jù)的海冰分布探測(cè)流程

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

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

圖5 NSIDC 2016年1月海冰范圍變化趨勢(shì)

日期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數(shù)據(jù)利用Otsu算法結(jié)合19 GHz頻段極化差得到的2016年1月海冰范圍，與NSIDC提供的結(jié)果在日變化趨勢(shì)上基本一致，其中MWRI得到的結(jié)果，以平均每天0.038 4×106km2的速度增加，NSIDC結(jié)果以平均每天0.041 9×106km2的速度增加，二者的增長(zhǎng)速度僅相差0.003 5×106km2；NSIDC 1月海冰范圍平均為13.493 3×106km2，MWRI數(shù)據(jù)結(jié)果1月海冰范圍平均為12.905 6×106km2，二者同比相差4.355%；NSIDC 1月海冰范圍極差為1.302 5×106km2，MWRI數(shù)據(jù)結(jié)果1月海冰范圍極差1.252 5×106km2，二者同比相差3.838%；NSIDC 1月海冰范圍單日極大值為13.958 7×106km2，出現(xiàn)日期為1月29日，MWRI 1月海冰范圍單日極大值為13.439 6×106km2，出現(xiàn)日期為1月30日，造成二者的差異可能是由于2個(gè)衛(wèi)星過(guò)境時(shí)間不同，成像不同步引起的，但綜合來(lái)看二者的極大值比較接近，同比相差僅為3.718%；NSIDC 1月海冰范圍單日極小值為12.656 2×106km2，出現(xiàn)日期為1月1日，MWRI 1月海冰范圍單日極小值為12.187 1×106km2，出現(xiàn)日期為1月1日，二者較為同步，且同比相差僅為3.706%。因此，結(jié)合以上幾個(gè)指標(biāo)來(lái)看，基于本文算法得到的北極海冰分布范圍與美國(guó)冰雪中心的海冰分布范圍及日變化趨勢(shì)基本一致，這在一定程度上說(shuō)明本文提出的方法可行。

4 結(jié) 語(yǔ)

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