ICESat-2與多源光學影像的潮灘地形反演方法

王子赫,康彥彥,王敏京

(1. 河海大學海洋學院,江蘇 南京 210098; 2. 自然資源部華東海岸帶野外科學觀測研究站,江蘇 南京 210007)

潮灘是由潮間帶淤泥帶部分與海底淺灘共同形成的沿海特殊地形,是海岸帶的重要組成部分,受潮汐的作用,高潮時被水淹沒,低潮時露出水面或變淺[1]。作為海陸相互作用的前沿地帶,近岸潮灘在海岸災害防護、土地資源利用及保護海洋生態環境等方面發揮著重要作用;但受風浪、潮流、風暴潮等因素的影響,其位置、形狀等不斷變化,這對沿海地形的監測與更新提出了較高要求,因此對高精度潮灘地形反演方法的研究具有重要意義[2]。

考慮到潮灘普遍惡劣的自然環境,大規模的人工跑灘現場測量存在極大危險性和困難性,需要消耗大量人力物力[3]。而遙感測量技術憑借其時效性強、覆蓋范圍廣、經濟耗費少等優勢,受到大量學者的關注。InSAR與LiDAR技術的發展對于提取潮灘敏感地形信息、描述潮灘的形態特征具有很大潛力[4-5]。但受潮灘沉積物和殘留水體的影響,輻射信號受阻減弱,精度受限。同時InSAR技術的研究機理復雜、數據處理難度較大,很少被應用于實際的潮灘反演處理中[6]。雖然機載LiDAR技術的效率和精度都可以保證,但其受自然因素影響會產生數據缺失,同時獲取數據成本過高,對于中等精度高頻次大面積、連續、長期的更新監測并不經濟可行。

在眾多遙感監測方法中,水邊線法以其適應性強、成本低及支持歷史地形構建等優勢,成為目前潮灘地形監測最具操作性的方法[7-8]。文獻[9]將水邊線視作等高線,利用多潮位賦值生成潮灘DEM模型。文獻[10]通過潮位測站基線計算水邊線高度值,構建了江蘇沿海輻射沙脊群DEM。文獻[11]利用二維水動力潮汐數值模型,模擬水邊線瞬時潮位,將高度值賦值于水線對應點,建立了江蘇沿海大面積潮灘歷史地形反演。水邊線法依托于潮位資料,賦值精度難以保證,同時在提取水邊線時存在必然的時間跨度,因此水邊線法反演地形存在局部區域差異。

2018年9月ICESat-2(ice,cloud and land elevation Satellite-2)激光衛星成功發射,為潮灘地形構建帶來新的技術手段。該衛星配備了高分辨率高級地形激光高度計系統(ATLAS),用于測量地球表面高程[12]。文獻[13]利用ICESat-2點云數據提高了航天飛機雷達地形任務(SRTM)高程產品衍生DEM的精度。文獻[14]通過構建ICESat-2高度數據與研究區域的淹沒頻率關系模型,獲得潮灘地形。但同一淹沒頻率下的實際高度同樣存在誤差,且使用遙感影像的時間跨度會更長,僅根據淹沒頻率很難反演真實的潮灘地形情況。

本文通過提取ICESat-2激光點云剖面數據高度值,結合多源光學影像與潮位進行反距離權重插值,提高水邊線賦值精度,得到一種新的潮灘地形反演方法,并與傳統水邊法DEM進行對比分析。

1 研究區與數據

1.1 研究區概況

選擇江蘇岸外輻射沙脊群南側水域的兩處大型沙洲——條子泥與高泥為研究對象。條子泥位于輻射沙脊群的中心地帶,東起東大港,西側與高泥相鄰,寬度約100 km,面積分別約為504.9、350 km2(如圖1所示)。受西大港、東大港、條魚港潮波系統及黃沙洋水域的共同作用,該區域泥沙交換活躍,潮灘地貌復雜多變,最大潮差超過9 m,具有輻射沙脊群水域最寬闊的潮面,灘面寬約14 km[15]。潮灘系統十分脆弱,其形態和動力地貌過程極易受人類活動的影響[16]。條子泥墾區一期建設總面積為6 095.62 hm2,大規模圍填海活動后,條子泥及其毗鄰近岸沉積動力環境的均衡態被打破,新的地貌格局逐漸形成[17]。

1.2 數據及預處理

研究初期收集2021年下半年的42組衛星數據,主要包括Landsat 8、GF-1、GF-6、HJ-1A衛星,經篩選獲取了10景清晰且含云量低于40%的光學圖像,其空間分辨率、幅寬等參數為適合水邊線法的應用,見表1。

本文使用的激光點云剖面數據是由ICESat-2衛星搭載的ATLAS測高系統探測所得的。ATLAS是一套非掃描推帚式、高重頻、微脈沖、光子計數式激光測高系統。其沿軌道有3束激光束,相鄰光束之間的交叉軌道距離為3.3 km。每束激光被分為2個子束(兩子束激光距離較近,為90 m),其中一子束激光脈沖能量較強,約是較弱子束的4倍。地面探測波束結構如圖2所示,6束激光沿飛行方向同時向地面進行推掃[18]。

圖2 ICESat-2衛星工作示意

激光點云高度數據采用NSIDC提供的ICESat-2二級產品——全球地理定位光子數據(ATL03),選取2021-07-09與2021-12-18兩組共12條激光軌道數據(如圖1所示)。其中黑色實線為2021-07-09軌道,由左至右分別為gt3r、gt3l、gt2r、gt2l、gt1r、gt1l的軌道;白色實線為2021-12-18軌道,由于ICESat-2在航行一段時間后會顛倒航線方向,因此該航道順序相反。原則上,在對數據進行篩選構建模型時,應選取軌道覆蓋潮灘面更為廣泛、能夠與水邊線產生更多交叉點,且激光脈沖能量較強的數據。因此選取2021-07-09中gt1l與2021-12-18中gt3r的高度值數據用于驗證,其余軌道數據用于構建DEM模型。

圖1中分別標注了東大港(DDG)、新條魚港(TYG)、洋口港(YKG)、大豐港(DFG)、東沙港(DSG)5處圍繞整個江蘇近海范圍的潮汐觀測站,能夠提供2021全年潮位數據?？堪兑粋葹镚PS RTK測量儀探測到的2021年10月的2組斷面實測高度數據,每組數據按照約20 m的漸進步長獲取高度值,用于DEM模型的驗證。

2 研究方法

2.1 ICESat-2斷面數據處理

根據強光波束對不同表面散射回光束的能量參數值的不同,ATL03產品將目標地物分類為陸冰、海冰、海洋、陸地、內陸水。產品根據每一光子事件的置信度值,將光子分為可能的信號光子點或背景光子點,通過信噪比與閾值的關系對光子事件進行歸類,以0、1、2、3、4為編號,分別對應噪聲、占位,以及低、中、高的置信度信號[12]。在進一步去噪前,選取高置信度信號(編號4)中被分類為陸地的信號光子事件,剔除其余光子事件,完成初步降噪,對應圖3中的降噪過程。

圖3 激光點云高度值提取過程

采用低通濾波窗口對其余光子點云剖面進行迭代運算。濾波降噪的主要思想是將信號中特定波段頻率濾除,是抑制和防止干擾的一項重要措施。本文采用的低通濾波能夠通過動態變化的方式改變頻率數值,且能夠適應不斷變化的潮灘地形高度,從而對激光點云數據進行降噪處理。其平滑與降噪的效果主要取決于截止頻率的設定,即

Wc=2·fc/fs

(1)

式中,Wc為3 dB截止頻率;fc為邊界頻率;fs為采樣頻率。去除背景噪聲后,激光點云剖面數據獲得最大限度的平滑,高度剖面線起伏變化位置周邊的噪聲得到削弱。對應圖3中的第4條降噪結果。

獲得最終潮灘高度數據前,需注意兩點:

(1)考慮到選取的光學影像分辨率為30 m,并以此基準提取水邊線數據,為簡化數值計算過程、方便處理數據;同時,為反映真實的地面高度情況,將濾波后的激光點云數據根據衛星航行的距離(即光子事件間隔),選取每隔30 m范圍內的光子事件為一組,提取每組數據中的高度中值,作為該區域的地形高度值。

(2)ATLAS后向散射的綠色激光束對清澈的Ⅰ類水體具有良好的穿透性,但對于輻射沙脊群水域渾濁的Ⅱ類水體,透射率極低,因此激光衛星僅能獲取條子泥地區陸地與海表面的高度信息。而條子泥水域受潮汐作用的影響,海表面信息在不斷變化,本文選取的兩組ICESat-2衛星過境軌道數據獲取時刻的潮位值均略高于該時間范圍的最低水位。因此去除多余的水平面數據,獲得潮灘斷面的高度值。如圖3中最后獲得的潮灘地形高度數據。

2.2 批量獲取光學影像水邊線圖像數據

采用最大似然監督分類方法對光學影像進行批量提取。提取前首先將光學遙感影像進行必要的預處理,包括輻射定標、大氣校正、幾何校正等。由于創建樣本時,同一類地物選擇過多將導致分類精度下降,因此應盡量選擇特征明顯的區域,如圖4(a)所示,將影像分為潮灘、水體兩類地物且作為訓練數據。使用最大似然法監督分類選用多個波段對圖像進行分類,再經濾波、聚類分析等分類后處理操作,得到圖4(b)。圖中高濁度懸沙或淺灘表層殘留水體等復雜地物,需結合潮汐、泥沙運動等因素,目視解譯對影像進行修改,提取影像中的水邊線數據,結果如圖4(c)所示。白色邊緣線即為水邊線,線條連續清晰,能夠較好地展現影像拍攝時刻下潮灘的整體形態。采用提取的潮灘地形高度數據與序列水邊線在二維平面上進行疊加,產生多處交叉點,如圖4(d)所示。

圖4 最大似然法監督分類批量水邊線提取

2.3 構建潮灘DEM模型

考慮到水邊線與激光剖面軌道的交點集中于激光衛星軌道一側,其位置和數值均與距離較遠的潮位值存在一定差距。因此選用反距離權重(IDW)的方式構建海面并高度對水邊線進行賦值。該方法對于距預測位置最近的測量值對預測值的影響將會更大。假定每個測量點都存在局部影響,而這一影響會隨著距離的增大而減小[19]。

將圖像采集時刻對應的潮位與高度交點數據進行IDW計算,得到覆蓋整個研究區域范圍的柵格高度插值數據,將結果賦值于對應水邊線中。為進一步提高模型精度,選取地形高度值中較為顯著的極值點與水邊線共同參與插值,這些點未能與水邊線相交,但能夠展現潮灘地形的部分細節,作為補充修正。利用自然鄰域插值法獲得條子泥和高泥 2021年下半年的DEM。

本文使用的ICESat-2數據,捕捉灘面的高度值范圍為-3～2 m。受潮汐水位的制約,選取最低潮位以上部分為掩膜興趣區,提取該部分DEM后進行分析。

3 結果分析與討論

將預留的兩條軌道用于驗證模型的效果,RMSE分別為0.27、0.29 m。由此可知,IDW算法將潮灘地形高度值與潮汐水位值進行插值,獲得較準確的柵格高度值,模型與ICESat-2地形高度數據結合度較好,軌道處的細節得到了很好的展現。

圖5分別為傳統的潮位賦值水邊線法與本文基于ICESat-2高度值優化的DEM模型。對比兩種DEM模型可知,優化后的DEM模型展示了條子泥內部更多的地形變化,深色部分主要集中于靠近江蘇近岸一側和沙洲系統的中心位置,淺色部分集中于沙洲系統邊緣的較低水位區域;經過地形高度數值優化后的模型能夠較好地識別出復雜多變的潮流、潮溝結構。

圖5 2021年下半年兩種方法獲得的DEM

為評估ICESat-2高度值對模型的優化效果,將兩種方式構建的DEM模型結果與兩組實測斷面地形進行對比,如圖6所示。本文方法與傳統潮位賦值水邊線法獲得的DEM模型相比,R2由0.86提高至0.89,RMSE由0.45 m提高至0.34 m。兩種模型結果均可觀測到清晰的地形起伏變化,并與實測數據大致吻合。但在潮位較高的區域,僅依靠潮位的賦值無法準確反映潮灘上實際的地形。ICESat-2數據的輸入提高了賦值水邊線的精度,使得地形變化的細節被捕捉,進而顯現在模型中。模型與實測斷面存在的微小差異可能是由數據不同的采集日期導致的。因此,受研究區域廣闊、復雜且多變的地形,以及潮灘地區本就存在的時空動態變化等因素的影響,構建模型的數據應該盡可能控制在更短的時間跨度內。

圖6 DEM結果誤差對比

4 結 語

本文通過對ICESat-2激光點云剖面數據進行低通濾波的降噪處理,提取了潮灘地形高度數據,并與潮汐水位數據相結合,使用反距離權重的插值,進而提高水邊線賦值的精度,構建了新的DEM模型?；谲壍榔拭娣秶母叨戎祵φw潮灘地形進行優化,實現了從二維到三維的轉換,結果符合潮灘地形基本規律,進一步驗證了ICESat-2數據應用于潮灘地形反演的可行性,為構建人類活動影響下的潮灘演化理論奠定了數據基礎,并為濕地生態環境保護和海岸防災減災提供了科學支撐。然而,目前的研究仍存在一定的局限性,如ICESat-2發射的綠色激光對于輻射沙脊群水域渾濁水體透射率極低,因此激光衛星無法獲得水面以下的水深信息,潮位情況不可避免地阻礙衛星獲取地形數據。但對于清澈水體綠色激光具有良好的穿透性,展現了ICESat-2在潮灘地形反演、水深勘探等領域的巨大潛力。