機(jī)載LiDAR和高光譜數(shù)據(jù)融合提取冰川雪線

李光輝，王 成，習(xí)曉環(huán)，鄭照軍，駱社周，岳彩榮

(1.西南林業(yè)大學(xué)林學(xué)院，昆明 650224;2.中國科學(xué)院對地觀測與數(shù)字地球科學(xué)中心，北京 100094;3.國家衛(wèi)星氣象中心，北京 100081)

0 引言

冰川是自然界重要的淡水資源，地球陸地表面10%被冰川覆蓋，80%的淡水儲存在冰川中[1]。同時冰川作為冰凍圈系統(tǒng)的主要組成部分，也是氣候變化的天然指示器[2]。世界各國對冰川的研究非常重視，了解和研究冰川雪線對研究冰川及其動態(tài)變化具有重要的現(xiàn)實(shí)意義。在冰川學(xué)中，雪線是指連接積雪能在夏季全部融化的、平坦的、非陰處地面的最高點(diǎn)所構(gòu)成的面，在這個面上雪的積累與消融量達(dá)到平衡，即冰川上雪的積累量與消融量相等的地方，因此雪線又稱冰川物質(zhì)平衡線。一個地區(qū)的雪線高度(即平衡線高度)能綜合反映該區(qū)氣候、地勢和地形等因素對冰川發(fā)育的影響，對冰川的物質(zhì)平衡以及冰川區(qū)環(huán)境具有較好的指示意義[3]。

最初獲取冰川雪線多采用人工手段，通過設(shè)立野外觀測地面站，由工作人員利用GPS等儀器進(jìn)行測量。這種方法周期長、成本高、獲取樣本數(shù)量有限，而且由于山地冰川區(qū)地形復(fù)雜，人員難以到達(dá)，因此該方法在高原冰川區(qū)應(yīng)用有限，不能滿足宏觀、綜合和快速的變化監(jiān)測要求。隨后遙感方法逐漸取代了人工方法，遙感技術(shù)快速、大范圍、低成本及周期短等特點(diǎn)，使其成為獲取冰川雪線的有效技術(shù)手段[4]。張世強(qiáng)等[5]利用冰川在TM圖像的可見光波段具有強(qiáng)反射和在中紅外波段具有強(qiáng)吸收特性，采用譜間關(guān)系比值法提取了青藏高原喀喇昆侖山區(qū)現(xiàn)代冰川邊界。紀(jì)鵬等[6]依據(jù)西昆侖地區(qū) 1991—2009年Landsat TM/ETM+遙感圖像，獲取近20 a西昆侖地區(qū)冰川及其冰川融水形成的郭扎錯湖和阿克賽欽湖的動態(tài)變化信息，并結(jié)合氣候資料，分析了冰川、湖泊變化與氣候變化的線性關(guān)系;金鑫等[7]利用1989年TM和2000年ETM+遙感圖像提取冰雪信息，結(jié)果表明，運(yùn)用線性混合像元分解模型能很好地監(jiān)測實(shí)驗(yàn)區(qū)的冰雪覆蓋變化;張瑞江等[8-9]通過對青藏高原現(xiàn)代冰川分布區(qū)ETM遙感圖像的分析，提出了利用遙感圖像確定現(xiàn)代雪線高度的方法，基本查明了青藏高原各山系的現(xiàn)代雪線高度，以及近30 a來現(xiàn)代雪線高度的變化狀況;而且又以1960—1970年間1∶10萬比例尺地形圖、1975年MSS和2000年ETM數(shù)據(jù)為信息源，基本查明了青藏高原各山系內(nèi)現(xiàn)有冰川的面積以及近30 a來冰川面積減少的數(shù)量;李邦良等[10]利用TM圖像與瓊西北地區(qū)內(nèi)的地質(zhì)地球物理資料進(jìn)行綜合分析，發(fā)現(xiàn)冰蝕和冰磧地貌之間內(nèi)在的成生規(guī)律;李震等[11]以1973—1994年的RBV，MSS和TM數(shù)據(jù)為信息源，用非監(jiān)督分類法和目視解譯的方法分析了青藏高原布卡塔格峰冰川20 a的變化信息，并提取了冰川的邊界及其變化信息。綜上所述，前人研究多基于星載多光譜數(shù)據(jù)或星載數(shù)據(jù)與其他數(shù)據(jù)的協(xié)同研究，融合機(jī)載激光雷達(dá)和高光譜數(shù)據(jù)提取雪線研究的報道較少。機(jī)載激光雷達(dá)點(diǎn)云數(shù)據(jù)精度高、密度大，不僅具有傳統(tǒng)遙感手段難以獲取的高精度高程信息，而且還具有多次回波和反射強(qiáng)度信息;高光譜數(shù)據(jù)的紋理特征連續(xù)、光譜信息豐富;綜合利用機(jī)載激光雷達(dá)和高光譜數(shù)據(jù)各自的優(yōu)勢不僅可以獲取研究對象的平面信息，而且可以獲取高精度的高程信息。本研究以TerraSolid，ENVI和ArcGIS等軟件為數(shù)據(jù)處理平臺，對冰川區(qū)域的機(jī)載激光雷達(dá)數(shù)據(jù)進(jìn)行處理，提取了冰川的數(shù)字高程模型(digital elevation model，DEM);并利用DEM對機(jī)載高光譜數(shù)據(jù)的分類結(jié)果進(jìn)行正射糾正，獲取冰川的數(shù)字正射影像圖(digital orthophoto map，DOM);最后融合研究區(qū)DOM和激光點(diǎn)云數(shù)據(jù)提取冰川雪線。

1 工作平臺及軟件

利用芬蘭TerraSolid公司開發(fā)的TerraSolid軟件的TerraScan和TerraModel模塊對機(jī)載LiDAR數(shù)據(jù)進(jìn)行處理。其中TerraScan的主要功能是根據(jù)激光點(diǎn)的坐標(biāo)、強(qiáng)度和反射強(qiáng)度等信息對點(diǎn)云數(shù)據(jù)進(jìn)行濾波和分類;TerraModeler主要用于生成各種數(shù)字模型(如DEM)及完成體積、面積、等高線、輪廓線和洪水淹沒區(qū)域等計算[12]。機(jī)載高光譜數(shù)據(jù)處理采用ENVI軟件，后續(xù)處理過程還運(yùn)用了ArcGIS軟件中的ArcMap和ArcScene工具以及地圖繪制軟件Global Mapper等。

2 試驗(yàn)區(qū)及數(shù)據(jù)源

2.1 試驗(yàn)區(qū)

選擇西藏那曲縣境內(nèi)的“中習(xí)一號”冰川為研究區(qū)。該區(qū)位于 E 91°24′48.006″～91°25′27.232″，N 30°50′11.356″～30°51′7.686″之間(圖 1)，屬高原亞寒帶季風(fēng)半濕潤氣候區(qū)。冬季長達(dá)半年，年溫差小，晝夜溫差大，年平均氣溫為-2.2℃，1月份平均氣溫-14.5℃;每年10月至次年5月為風(fēng)雪期和土壤凍結(jié)期，6—8月為牧草生長期;年降水量為406.2 mm，多集中在6—9 月;年蒸發(fā)量為1 811.2 mm，相對濕度為51%，氣候特征是干燥和多風(fēng)。

圖1 研究區(qū)地理位置Fig.1 Location of the study area

2.2 數(shù)據(jù)源

2.2.1 機(jī)載LiDAR數(shù)據(jù)

機(jī)載LiDAR數(shù)據(jù)于2011年8月10日獲取，相對航高3 000 m，點(diǎn)云密度約為1點(diǎn)/m2。如圖2所示。

圖2 機(jī)載點(diǎn)云數(shù)據(jù)Fig.2 Airborne point cloud data

2.2.2 機(jī)載高光譜數(shù)據(jù)

機(jī)載高光譜數(shù)據(jù)于2011年8月11日獲取，空間分辨率約2 m，共100個波段(波長范圍400～1 000 nm)。采用的傳感器為中國科學(xué)院上海技術(shù)物理研究所研制的PHI-1309。圖3(a)是一個掃描航帶的部分高光譜數(shù)據(jù)的R(41)，G(26)，B(11)合成圖像，3個波段的中心波長分別為642 nm，549 nm和458 nm;圖3(b)為研究區(qū)高光譜數(shù)據(jù)。

圖3 高光譜數(shù)據(jù)Fig.3 Hyperspectral data

3 數(shù)據(jù)處理與分析

3.1 機(jī)載LiDAR數(shù)據(jù)處理

首先利用TerraSolid軟件中的TerraScan模塊對研究區(qū)LiDAR點(diǎn)云數(shù)據(jù)進(jìn)行濾波分類，然后利用TerraModel模塊生成研究區(qū) DEM，并在 Global Mapper軟件中顯示，如圖4所示。

圖4 冰川區(qū)DEMFig.4 Glacier area DEM

3.2 機(jī)載高光譜數(shù)據(jù)處理

利用ENVI4.7軟件，對機(jī)載高光譜數(shù)據(jù)進(jìn)行預(yù)處理，包括輻射定標(biāo)、大氣校正和幾何糾正。然后采用主成分分析對高光譜數(shù)據(jù)進(jìn)行降維，并通過最大似然法得到冰川分類圖像。

3.3 三維地形模擬

利用ArcGIS Desktop中的ArcScene模塊將冰川DEM數(shù)據(jù)分別同高光譜的柵格數(shù)據(jù)和矢量數(shù)據(jù)進(jìn)行三維地形模擬，如圖5所示。

圖5 三維地形模擬Fig.5 3D terrain simulation

3.4 影像校正

激光雷達(dá)和高光譜數(shù)字影像的地形輻射校正過程是:利用激光雷達(dá)的測距數(shù)據(jù)建立地面的數(shù)字地形圖和高程模型，研究成像瞬間激光雷達(dá)和高光譜數(shù)字影像與地面地形之間的定量關(guān)系，分別建立地形輻射校正模型，將激光雷達(dá)和高光譜數(shù)字影像上各點(diǎn)的影像值校正為地面水平時激光和太陽光照射角度相同的影像。經(jīng)過地形糾正處理的激光雷達(dá)和高光譜影像可視為各點(diǎn)在地面激光和太陽光照度相同。大氣和地形狀況一致的最佳成像條件下得到的激光雷達(dá)和高光譜數(shù)字圖像，即正射影像。正射影像最大程度地消除了因距離、地形和大氣狀況對激光雷達(dá)和高光譜數(shù)字影像的影響，充分表達(dá)了地物的波譜特性，提高了激光雷達(dá)和高光譜對地球表面物質(zhì)的系統(tǒng)探測能力，改善了激光雷達(dá)和高光譜數(shù)字圖像的質(zhì)量以及數(shù)字仿真能力[13]。

選擇Global Mapper菜單中影像校正功能，打開高光譜最大似然法分類數(shù)據(jù)，進(jìn)入幾何糾正界面，然后以DEM作為地理參考坐標(biāo)對高光譜遙感分類圖像進(jìn)行幾何精糾正和坐標(biāo)歸一化處理。在DEM和高光譜遙感分類圖像上分別選取了10個控制點(diǎn)(GCP)，配準(zhǔn)精度控制在1個像元之內(nèi)(圖6)。

圖6 研究區(qū)正射影像Fig.6 Orthophoto of the study area

3.5 冰川雪線提取

首先，在TerraSolid軟件的TerraPhoto模塊中，同時顯示激光點(diǎn)云與高光譜正射影像(圖7);然后，用TerraPhoto模塊中的Model工具，提取正射影像的冰川邊界點(diǎn)云。

圖7 正射影像和點(diǎn)云數(shù)據(jù)疊加Fig.7 Orthophoto and point cloud data overlay

由于所獲取的高光譜數(shù)據(jù)是條帶數(shù)據(jù)，沒有完整包含“中習(xí)一號”冰川，因此所提取出的雪線不是該冰川全部的雪線。結(jié)果如圖8所示。

圖8 “中習(xí)一號”冰川部分雪線及其高程信息Fig.8 Snow line of“Zhongxi－1”glacier part and its elevation information

4 精度分析

4.1 LiDAR DEM的精度分析

DEM精度評定可以采用2種不同的方式:①平面精度和高程精度分開評定;②2種精度同時評定。但在實(shí)際應(yīng)用中，一般只討論DEM高程精度的評定問題[14]。由于“中習(xí)一號”沒有實(shí)測數(shù)據(jù)，無法用實(shí)測數(shù)據(jù)來檢測LiDAR DEM，因此利用了“中習(xí)一號”的ASTER DEM和LiDAR DEM進(jìn)行比較，選取了36個坐標(biāo)點(diǎn)進(jìn)行高程值的比較(表1)。

表1 LiDAR DEM和ASTER DEM的高程比較Tab.1 Comparing LiDAR DEM with ASTER DEM elevation (m)

通過表1的比較，LiDAR DEM的高程值均小于ASTER DEM的高程值，統(tǒng)一偏向一側(cè)，說明LiDAR DEM和ASTER DEM具有很強(qiáng)的系統(tǒng)性，且方向都一致，這一特征符合正常情況。

4.2 機(jī)載高光譜數(shù)據(jù)的精度分析

本文還選擇了“中習(xí)一號”地區(qū)環(huán)境減災(zāi)衛(wèi)星(HJ-1)的多光譜數(shù)據(jù)與機(jī)載高光譜數(shù)據(jù)進(jìn)行精度比較。研究區(qū)2011年8月10日的HJ-1多光譜CCD數(shù)據(jù)(軌道號為29/80)共4個波段(波譜范圍為0.43 ～0.9 μm)，空間分辨率為 30 m。首先，對圖像進(jìn)行輻射定標(biāo)、幾何糾正和大氣校正，并選擇R(3)，G(2)，B(1)的方案進(jìn)行波段組合;其次，采用最大似然法進(jìn)行分類;最后，通過對環(huán)境衛(wèi)星多光譜數(shù)據(jù)目視解譯獲取的地表真實(shí)感興趣區(qū)來驗(yàn)證分類結(jié)果，采用混淆矩陣方法對分類結(jié)果進(jìn)行精度分析，并與高光譜數(shù)據(jù)的主成分變換后的最大似然法分類結(jié)果精度進(jìn)行比較(表2)。

表2 機(jī)載高光譜和HJ－1多光譜數(shù)據(jù)分類精度比較Tab.2 Precision comparison of airborne hyperspectral classification data and HJ－1 multi－spectral data

由表2可得，運(yùn)用高光譜數(shù)據(jù)的分類精度高于同地區(qū)HJ-1多光譜數(shù)據(jù)的分類精度。

5 結(jié)論

本文對機(jī)載激光雷達(dá)數(shù)據(jù)和高光譜影像融合提取冰川雪線進(jìn)行了嘗試，其原理清晰，技術(shù)方法和流程簡單。相對于衛(wèi)星遙感監(jiān)測冰川及雪線變化，雖然機(jī)載數(shù)據(jù)獲取成本高，但機(jī)載數(shù)據(jù)的精度和空間分辨率高，而且數(shù)據(jù)獲取方式及時間靈活，特別適于災(zāi)害的應(yīng)急反應(yīng)監(jiān)測。同時如果能獲取時間序列的機(jī)載LiDAR數(shù)據(jù)和高光譜數(shù)據(jù)，即可提取高精度的冰川及雪線的動態(tài)變化信息。

由于獲取研究區(qū)雪線的實(shí)測數(shù)據(jù)有一定困難，本文未能用實(shí)測數(shù)據(jù)對提取結(jié)果進(jìn)行驗(yàn)證，另外研究僅利用了激光點(diǎn)云數(shù)據(jù)的高程信息，如何有效利用點(diǎn)云數(shù)據(jù)的強(qiáng)度信息進(jìn)行冰川雪線的提取還有待進(jìn)一步研究。

[1]劉志輝，王紅娟，裴 歡，等.基于RS和GIS技術(shù)的近40 a新疆昌吉州冰川變化分析[J].新疆大學(xué)學(xué)報:自然科學(xué)版，2005，22(2):127-133.Liu Z H，Wang H J，Pei H，et al.Analysis of the glacier change near 40 a in Changji Prefecture，Xinjiang based on RS and GIS[J].Journal of Xinjiang University:Natural Science Edition，2005，22(2):127-133.

[2]張國梁，王 杰，潘保田，等.冰川變化遙感監(jiān)測的研究進(jìn)展[J].蘭州大學(xué)學(xué)報:自然科學(xué)版，2010，46(6):1-10.Zhang G L，Wang J，Pan B T，et al.Progress in research on glacier changes using remote sensing[J].Journal of Lanzhou University:Natural Sciences，2010，46(6):1-10.

[3]王利平，謝自楚，丁良福，等.基于GIS的羌塘高原冰川系統(tǒng)雪線場的建立及其空間分布特征[J].干旱區(qū)地理，2010，33(5):692-702.Wang L P，Xie Z C，Ding L F，et al.Establishing the field of equilibrium line altitude of the glacier system in the Chiangtang Plateau based on geographic information system and its distribution characteristics[J].Arid Land Geography，2010，33(5):692-702.

[4]王祎婷，陳秀萬，柏延臣，等.多源DEM和多時相遙感影像監(jiān)測冰川體積變化——以青藏高原那木納尼峰地區(qū)為例[J].冰川凍土，2010，32(1):126-132.Wang Y T，Chen X W，Bo Y C，et al.Monitoring glacier volume change based on multi-source DEM and multi-temporal remote sensing images:A case study in the Mount Naimona’nyi Region on the Tibetan Plateau[J].Journal of Glaciology and Geocryology，2010，32(1):126-132.

[5]張世強(qiáng)，盧 健，劉時銀.利用TM高光譜圖像提取青藏高原喀喇昆侖山區(qū)現(xiàn)代冰川邊界[J].武漢大學(xué)學(xué)報:信息科學(xué)版，2001，26(5):435-440.Zhang S Q，Lu J，Liu S Y.Deriving glacier border information on Qinghai Tibet by TM high spectrum image[J].Journal of Wuhan University:Information Science Edition，2001，26(5):435-440.

[6]紀(jì) 鵬，郭華東，張 露.近20年西昆侖地區(qū)冰川動態(tài)變化遙感研究[J].國土資源遙感，2013，25(1):93-98.Ji P，Guo H D，Zhang L.Remote sensing study of glacier dynamic change in West Kunlun Mountains in the past 20 years[J].Remote Sensing for Land and Resources，2013，25(1):93-98.

[7]金 鑫，柯長青.基于混合像元分解的天山典型地區(qū)冰雪變化監(jiān)測[J].國土資源遙感，2012，24(4):146-151.Jin X，Ke C Q.Monitoring of snow cover changes in Tianshan Mountains based on mixed pixel decomposition[J].Remote Sensing for Land and Resources，2012，24(4):146-151.

[8]張瑞江，趙福岳，方洪賓，等.青藏高原近30年現(xiàn)代雪線遙感調(diào)查[J].國土資源遙感，2010，22(s1):59-63.Zhang R J，Zhao F Y，F(xiàn)ang H B，et al.Remote sensing survey of existing snowlines in the past 30 years in Qinghai-Tibet Plateau[J].Remote Sensing for Land and Resources，2010，22(s1):59-63.

[9]張瑞江，方洪賓，趙福岳，等.青藏高原近30年來現(xiàn)代冰川面積的遙感調(diào)查[J].國土資源遙感，2010，22(s1):45-48.Zhang R J，F(xiàn)ang H B，Zhao F Y，et al.Remote sensing survey of existing glaciers in the past 30 years in Qinghai-Tibet Plateau[J].Remote Sensing for Land and Resources，2010，22(s1):45-48.

[10]李邦良.瓊西北地區(qū)冰川地貌陸地衛(wèi)星TM圖像解譯[J].國土資源遙感，1998，10(3):90-94.Li B L.TM image interpretation for glacial geography in the Northwest of Qiong[J].Remote Sensing for Land and Resources，1998，10(3):90-94.

[11]李 震，孫文新，曾群柱.綜合RS與GIS方法提取青藏高原冰川變化信息——以布喀塔格峰為例[J].地理學(xué)報，1999，54(3):263-268.Li Z，Sun W X，Lu Q Z.RS and GIS methods to extract glacier change information to the glaciers of the Qinghai-Tibet Plateau:Cloth Katage Peak as an example[J].Journal of Geographical Science，1999，54(3):263-268.

[12]曹廣強(qiáng)，翟 娜.淺談利用激光點(diǎn)云數(shù)據(jù)制作DEM[J].測繪標(biāo)準(zhǔn)化，2011，27(2):44-46.Cao G Q，Zhai N.On the DEM creating by using radar point cloud data[J].Standardization of Surveying and Mapping，2011，27(2):44-46.

[13]李先華，黃 睿，黃 微，等.激光雷達(dá)與高光譜遙感對地立體探測研究[J].紅外，2009，30(7):5-9.Li X H，Huang R，Huang W，et al.Research on tridimensional ground exploration based on lidar and hyperspectral remote sensing[J].Infrared，2009，30(7):5-9.

[14]肖雁峰.DEM質(zhì)量檢查與精度評定研究[D].成都:西南交通大學(xué)，2008:1-71.Xiao Y F.Quality inspection and accuracy evaluation for DEM[D].Chengdu:Southwest Jiaotong University，2008:1-71.