王立偉 金濤勇 張勝軍 李大煒
1 武漢大學測繪學院,武漢市珞瑜路129號,430079
CryoSat-2衛(wèi)星是繼ERS-1、ERS-2、Envisat及ICESat后首個針對極地進行觀測的衛(wèi)星,在高緯區(qū)域具有30d的子周期[1]。衛(wèi)星搭載Ku波段雷達高度計SIRAL,有3種測量模式:LRM(低分辨率測量模式)、SAR(合成孔徑雷達測量模式)、InSAR(干涉合成孔徑雷達測量模式)[2]??紤]到數(shù)據(jù)源及數(shù)據(jù)處理方式的一致性,本文選擇基于SAR 模式的觀測數(shù)據(jù)。為顧及數(shù)據(jù)的現(xiàn)勢性,使計算出的海冰干舷高能更好地反映某一時刻某一區(qū)域的特征,計算過程中采用沿軌內插方法,而沒有采用常用的區(qū)域格網(wǎng)化方法。將兩種方法的計算結果進行比較發(fā)現(xiàn),二者在多年冰區(qū)域存在略微差異,但整體趨勢較為一致。
在近岸區(qū)域或海冰區(qū)域,衛(wèi)星測高波形通常會受到較為嚴重的污染,在數(shù)據(jù)處理之前,首先需要進行波形重跟蹤,計算出星載預設門與重跟蹤門之間的距離,進而對測量得到的衛(wèi)星到星下點之間的距離進行改正。常見的波形重跟蹤方法包括重心偏移法(OCOG 法)、閾值法、改進閾值法、神經(jīng)網(wǎng)絡算法、Beta-5參數(shù)算法、Beta-9參數(shù)算法以及海洋算法等[3-9]。
歐空局所發(fā)布的SAR 數(shù)據(jù)分為3 個級別:L1B數(shù)據(jù),即波形數(shù)據(jù);L2數(shù)據(jù),直接給出相對于WGS84橢球的表面高;L2I數(shù)據(jù),在L2數(shù)據(jù)基礎上增加海冰區(qū)域的表面識別信息,以便計算海冰干舷高(海冰露出水面的高度)。由L1B 數(shù)據(jù)得到的L2 數(shù)據(jù)和L2I數(shù)據(jù)的處理方式統(tǒng)一采用OCOG 波形重跟蹤算法,而OCOG 算法最初是針對海洋波形提出的,并不適合海冰的返回波形,而閾值法更適合于海冰區(qū)域的波形重跟蹤[1,8-11]。
閾值法的思路是:首先計算或者直接利用OCOG 方法獲得波形的振幅,然后對振幅乘以一個閾值因子,獲得重跟蹤點的幅值,最后線性內插獲得重跟蹤點的采樣位置(圖1)。計算公式如下:
圖1 閾值法重跟蹤示意圖Fig.1 Sketch map of threshold retracking method
式中,p(i)為第i個采樣點所對應的回波功率;Amax為回波功率的最大幅值,即波形振幅;DC 為利用前c個采樣點所計算的噪聲水平;TL為重跟蹤點所對應的幅值;t為幅值大于TL值點所對應的最近采樣門;LEG 為通過線性插值所獲得的波形前緣中點,即波形重跟蹤點所對應的采樣門[12]。需要說明的是,很多學者在計算波形的噪聲水平時通常采用波形的前幾個采樣值,一般默認從第一個采樣值開始。有的測高衛(wèi)星由于第一個采樣值或前幾個采樣值存在較大偏差,實際計算時是利用后面的一段采樣值對噪聲水平進行估計[8]。本文在實際計算時利用每個波形的前5個采樣值來估計噪聲水平,即a=0,c=5。
本文采用閾值法對CryoSat-2衛(wèi)星在海冰區(qū)域的SAR 波形進行重跟蹤。由于雪密度和雪濕度等的影響,實際中雷達波形并不能完全穿透雪深。一些研究者將閾值系數(shù)直接選為50%,是假設雷達波形能夠完全穿透雪深,使得計算出的冰層厚度存在較大誤差。本文將閾值系數(shù)選為40%,即假設雷達波形只能穿透一定深度的雪深[1]。
計算海冰干舷高的關鍵是要獲得高精度海冰表面高和可靠的瞬時海面高,通常采用格網(wǎng)化的方法,計算格網(wǎng)內的平均海冰高度和瞬時海面高度[8]。高精度海冰表面高可通過重跟蹤方法獲得,平均瞬時海面高一般取為格網(wǎng)內浮冰之間的海面高度。因此,對于海冰干舷高反演來說,最重要的就是進行浮冰之間的海水面(Lead)識別,準確測量出浮冰之間的局部瞬時海面高,然后內插計算海冰處的局部瞬時海面高。
格網(wǎng)法的缺點是,一個格網(wǎng)內的數(shù)據(jù)可能會含有多條軌跡,幾條軌跡間的時間跨度可能很長,通過取中值或平均值所得結果并不一定代表某一時刻的瞬時海面高,從而影響海冰干舷高的計算精度。
本文采用單條軌跡的沿軌內插方法計算海冰位置處的瞬時海面高。先將DTU10平均海平面高內插到衛(wèi)星的每一個星下點,然后通過波形特征參數(shù)識別出海冰的位置。以每一個海冰位置為節(jié)點,向軌跡的兩端同時搜索,識別出與該海冰節(jié)點距離最近的Lead。如果兩端均搜索到Lead,則同時計算兩個Lead處海平面異常,并采用距離加權方法將其內插到海冰位置處,再加上海冰處的平均海平面高得到海冰處瞬時海平面高;如果僅一端搜索到Lead,則直接將該Lead處的海平面異常視為海冰位置處的海平面異常。事實上,僅一端能搜索到Lead的情況非常少見,可以忽略。
為了進行Lead和海冰的識別,將OSI SAF(海洋與海冰監(jiān)測)高緯處理中心發(fā)布的海冰濃度格網(wǎng)數(shù)據(jù)內插到CryoSat-2數(shù)據(jù)中,刪除那些海冰濃度小于70%的數(shù)據(jù)(主要是海洋數(shù)據(jù)),然后再利用脈沖峰值、波束標準差、波峰陡峭度這3個波形特征參數(shù)對Lead和海冰進行區(qū)分。各參數(shù)的門限值設置情況如表1。
表1 波形識別參數(shù)設置Tab.1 Set the waveform parameter identification
脈沖峰值是波形回波功率的函數(shù)[1]:
pp值越大,說明波形越尖銳。波束標準差用來描述波形后向反射角的變化[2],而波峰陡峭度用來描述波形頂端的陡峭程度[13]。
德國AWI利用40%閾值法分析2011-01~2013-03共27個月的數(shù)據(jù),利用格網(wǎng)化方法計算得到海冰干舷高。選取其最后3個月的數(shù)據(jù)與本文結果進行比較,如圖2。可以看出,二者得到的北冰洋海冰干舷高整體分布趨勢較為一致,AWI在多年冰區(qū)域的海冰干舷高要更高一些。圖3給出2013-01~03 二者的差值直方圖??梢钥闯觯@3個月二者的差值直方圖非常接近于正態(tài)分布。最后利用OSI SAF發(fā)布的每天多年冰和一年冰格網(wǎng)數(shù)據(jù),分別計算AWI和本文數(shù)據(jù)的多年冰區(qū)域和一年冰區(qū)域的平均海冰干舷高(圖4)。
從圖4看出,本文計算的2011~2013年北冰洋多年冰和一年冰區(qū)域的平均海冰干舷高結果與AWI較為一致,海冰干舷高呈現(xiàn)明顯的周期性變化。由于夏季海冰融化,使得波形識別受到影響,利用返回波形特征識別Lead和海冰的方法不再可用,所以沒有對6~8月的海冰進行處理。AWI利用波形識別方法得到的6~8月北冰洋海冰干舷高很不理想。對于夏季海冰的識別,還需采用其他更為有效的方法。
圖5為2011~2013年北冰洋各年同期數(shù)據(jù)對比圖,無論是本文結果還是AWI結果,各年同期海冰干舷高結果均較為接近。但同時也要注意到,在2~5月,本文得到的2013年數(shù)據(jù)無論是在多年冰區(qū)域還是在一年冰區(qū)域都較同期其他數(shù)據(jù)偏低,而在9~12月略微偏高,尤其是2013-09的數(shù)據(jù) 較2011-09和2012-09呈現(xiàn)明顯的增長趨勢。由于AWI數(shù)據(jù)只公布到2013-03,所以為了驗證本文結果,還需要采用實地測量數(shù)據(jù)或其他衛(wèi)星測量數(shù)據(jù)進行更好的解釋和說明。
圖2 本文海冰干舷高計算結果與AWI計算結果的比較Fig.2 Comparison between the results of sea ice freeboard calculation in this paper and the results of AWI
圖3 2013-01~03AWI與本文海冰干舷高差值直方圖Fig.3 Histogram of the sea ice freeboard height difference(2013-01~03)
本文采用沿軌內插方法計算北冰洋區(qū)域的平均海冰干舷高,與常用的格網(wǎng)化方法計算結果整體趨勢較為一致,說明沿軌內插方法完全可以用來進行海冰干舷高的反演。在計算中發(fā)現(xiàn),2013-09以后北冰洋區(qū)域無論是多年冰區(qū)域還是一年冰區(qū)域,平均海冰干舷高均呈現(xiàn)明顯上升趨勢,這還需要通過其他數(shù)據(jù)進行驗證。為了對本文結果和AWI結果進行更好的比較,下一步將采用實測數(shù)據(jù)和機載數(shù)據(jù),對二者的計算結果進行更詳細的比較,進而評價格網(wǎng)化方法和沿軌內插方法在計算海冰干舷高中的優(yōu)劣,最終將兩種方法合理組合,更好地反演海冰干舷高。
圖4 2011~2013年北冰洋平均海冰干舷高變化趨勢Fig.4 Variation trend of the average sea ice freeboard height from 2011to 2013year
圖5 2011~2013年北冰洋平均海冰干舷高同期數(shù)據(jù)比較Fig.5 Arctic sea ice freeboard period average data comparison from 2011to 2013
[1]Wingham D J,F(xiàn)rancis C R,Baker S,et al.CryoSat:A Mission to Determine the Fluctuations in Earth’s Land and Marine Ice Fields[J].Adv Space Res,2006,37(4):841-871
[2]Ricker R,Hendricks S,Helm V,et al.Sensitivity of CryoSat-2Arctic Sea-Ice Volume Trends on Radar-Waveform,Interpretation[J].The Cryosphere Discuss,2014,8(2):1 831-1 871
[3]高永剛,郭金運,黃金維,等.湖泊范圍內TOPEX/Poseidon衛(wèi)星測高波形分析[J].大地測量與地球動力學,2008,28(1):45-49(Gao Yonggang,Guo Jinyun,Huang Jinwei,et al.Waveform Analysis of Topex/Poseidon Satellite Altimeter over Lake[J].Journal of Geodesy and Geodynamics,2008,28(1):45-49)
[4]褚永海,李建成,張燕,等.ENVISAT 測高數(shù)據(jù)波形重跟蹤分析研究[J].大地測量與地球動力學,2005,25(1):76-80(Chu Yonghai,Li Jiancheng,Zhang Yan,et al.Analysis and Investigation of Waveform Retracking Data of Envisat[J].Journal of Geodesy and Geodynamics,2005,25(1):76-80)
[5]褚永海.衛(wèi)星測高波形處理理論研究及應用[D].武漢:武漢大學,2007(Chu Yonghai.The Research of Satellite Altimetry Waveform Theory and Application[D].Wuhan:Wuhan University,2007)
[6]常曉濤,李建成,郭金運,等.一種多前緣多閾值的波形重構算法[J].地球物理學報,2006,49(6):1 629-1 634(Chang Xiaotao,Li Jiancheng,Guo Jinyun,et al.A Multileading Edge and Multi-threshold Waveform Retrackef[J].Chinese J Geophys,2006,49(6):1 629-1 634)
[7]European Space Research Institute(ESRI N).European Space Agency and Mullard Space Science Laboratory-University College London[EB/OL].http://earth.esa.int/web/guest/missions/cryosat/ipf-baseline,2013
[8]Rose S K.Measurements of Sea Ice by Satellite and Airborne Altimetry[D].Denmark:DTU Space,National Space Institute,2013
[9]Davis C H.A Robust Threshold Retracking Algorithm for Measuring Ice-Sheet Surface Elevation Change from Satellite Radar Altimeters[J].IEEE Trans Geosci Remote Sens,1997,35(4):974-979
[10]Helm V,Humbert A,Miller H.Elevation Change of Greenland and Antarctic Derived from CryoSat-2[J].The Cryosphere Discuss,2014,(8):1 673-1 721
[11]Stenseng L,Andersen O B.Preliminary Gravity Recovery from CryoSat-2Data in the Baffin Bay[J].Adv Space Res,2012,50(8):1 158-1 163
[12]阮海林.衛(wèi)星測高波形重跟蹤算法在臺灣周邊海域的應用研究[D].廈門:國家海洋局第三海洋研究所,2010(Ruan Hailin.Research and Application of Algorithm in the Waters Surrounding Taiwan Satellite Altimetry Waveform Retracking[D].Xiamen:The Third Oceanography Research Institute,2010)
[13]Brenner A C,Zwally H J,Bentley C R,et al.Derivation of Range and Range Distributions from Laser Pulse Waveform Analysis for Surface Elevations,Roughness,Slope,and Vegetation Heights[EB/OL].http//glas.wff.nasa.govl,2000