李 旭,王 志 慧,李 春 意,劉 楊,馬 力
(1.河南理工大學(xué)測繪與國土信息工程學(xué)院,河南 焦作 454003;2.黃河水利委員會(huì)黃河水利科學(xué)研究院,河南 鄭州 450003;3.水利部黃土高原水土保持重點(diǎn)實(shí)驗(yàn)室,河南 鄭州 450003)
地表反射率表征不同物體對太陽輻射的反射和吸收能力[1],受傳感器、大氣、地形和光照等多種因素影響。山區(qū)地形起伏變化改變了太陽、傳感器和地表三者之間的位置關(guān)系,從而產(chǎn)生“同物異譜,異物同譜”現(xiàn)象[2],因此,地形校正是山區(qū)遙感影像預(yù)處理中不可缺少的一步。地形校正指通過一系列的變換,將斜坡上的像元輻射值轉(zhuǎn)化為理想表面(通常取水平面)上的像元輻射值,以此減弱或消除因地形效應(yīng)引起的輻射值變化,從而獲取目標(biāo)像元的真實(shí)光譜特征[3]。現(xiàn)有地形校正模型可分為三類[4]:①經(jīng)驗(yàn)?zāi)P?應(yīng)用簡單、適用性強(qiáng),且不需要額外的物理參數(shù),但模型輸出的結(jié)果不具備實(shí)際的物理意義,如波段比值法[5]和經(jīng)驗(yàn)統(tǒng)計(jì)模型[6]等;②物理模型,通過使用大量物理參數(shù)獲得精確的結(jié)果,但模型復(fù)雜且需要較多參數(shù),如Proy校正模型[7]、Sandmeier校正模型[8]和ACTOR校正模型[9];③半經(jīng)驗(yàn)?zāi)P?結(jié)合了經(jīng)驗(yàn)?zāi)P秃臀锢砟P?利用DEM數(shù)據(jù)計(jì)算坡度、坡向和太陽入射角等地形參數(shù)并應(yīng)用于校正模型中,計(jì)算簡便且輸出結(jié)果較準(zhǔn)確,如C校正(C-correction)模型[6]、SCS+C校正(Sun-Canopy-Sensor+C-correction)模型[10]和Minnaert(MIN)校正模型[11]。
地形校正結(jié)果受DEM的分辨率與精度[12]、校正方法[13]、大氣條件[14]等多因素的綜合影響,國內(nèi)外學(xué)者對不同校正模型的校正效果展開了對比分析。例如:Hantson等[15]對余弦校正、經(jīng)驗(yàn)統(tǒng)計(jì)校正(Statistical-Empirical correction,SE)[16]、C模型、Minnaert模型和改進(jìn)Minnaert模型進(jìn)行對比分析,結(jié)果表明SE模型校正效果最好;Richter等[17]對C模型、改進(jìn)Minnaert模型和Gamma模型進(jìn)行對比實(shí)驗(yàn),結(jié)果表明改進(jìn)Minnaert模型在多數(shù)情況下校正效果優(yōu)秀;王颯等[18]使用余弦模型、C模型、SCS+C模型和Minnaert模型進(jìn)行地形校正,指出C模型的應(yīng)用范圍最廣泛。隨著AW3D30、SRTM、GDEM、MERIT[19]和ALOS等開源DEM的免費(fèi)開放和廣泛應(yīng)用[20],基于開源DEM的半經(jīng)驗(yàn)地形校正模型成為地形校正中常用的方法,然而針對不同開源DEM最合適的校正模型以及最佳的組合方式仍缺乏研究[21]。為此,本文以黃土丘陵區(qū)為例,基于Landsat8 OLI影像,利用3種開源DEM(AW3D、SRTM、GDEM)和4種地形校正模型(C模型、SCS+C模型、SE模型、MIN模型)開展地形校正,并通過構(gòu)建地形校正效果評價(jià)指標(biāo),對由局地入射角引起的光譜異常、坡向引起的光譜異常和同地類光譜變異性的校正效果進(jìn)行定量評估,以期為黃土丘陵區(qū)地形校正模型選擇提供參考,對山區(qū)地表定量遙感應(yīng)用具有重要意義。
研究區(qū)為陜西省榆林市綏德縣和子洲縣(圖1),平均海拔為1 028 m,坡度在0.02°~50.7°之間,平均坡度為20.02°,草地、農(nóng)田、裸地面積占比分別為72.3%、22.7%和4.95%。

圖1 研究區(qū)概況
從美國地質(zhì)勘探局(United States Geological Survey,USGS)下載2020年9月16日Landsat8 OLI地表反射率產(chǎn)品(Path=127,Row=34)并進(jìn)行大氣校正(圖1a),空間分辨率30 m,包括可見光、近紅外和短波紅外范圍的7個(gè)光譜波段,太陽高度角和太陽方位角分別為39.22°和149.13°。
研究選用3種常用的開源DEM數(shù)據(jù)集:日本宇宙航空研究開發(fā)機(jī)構(gòu)(JAXA)發(fā)布的AW3D30 DEM數(shù)據(jù)集(http://www.eorc.jaxa.jp/alos/en/aw3d30/index.htm),其高程基準(zhǔn)和參考橢球面分別為EGM96高程和WGS-84橢球;美國航空航天局(NASA)發(fā)布的ASTER GDEMV3和SRTM-1 DEM數(shù)據(jù)集(http://earthexplorer.usgs.gov/),高程基準(zhǔn)和參考橢球均為EGM96高程和WGS84橢球,空間分辨率均為30 m。以Landsat8 OLI影像為基準(zhǔn),利用ENVI軟件對3種DEM數(shù)據(jù)分別進(jìn)行幾何校正處理,幾何配準(zhǔn)誤差控制在0.5個(gè)像元以內(nèi)。另外,地表覆蓋數(shù)據(jù)采用2020年全球地表覆蓋GlobeLand30[22](http://www.globallandcover.com/)。
1.3.1 技術(shù)路線 首先,通過分析遙感影像地表反射率與基于3種開源DEM計(jì)算得到的太陽局地入射角余弦值(cosi)之間的相關(guān)性,篩選出精度最高的開源DEM作為基準(zhǔn)DEM;其次,利用3種開源DEM和4種地形校正模型(C、SCS+C、SE、MIN)構(gòu)建不同組合對研究區(qū)遙感影像進(jìn)行地形校正;最后,通過構(gòu)建地形校正評價(jià)指標(biāo)對不同DEM和校正模型組合的校正效果進(jìn)行定量評價(jià),篩選出適用于黃土丘陵區(qū)遙感影像地表反射率地形校正的最佳組合方式(圖2)。

圖2 技術(shù)路線
1.3.2 地形校正模型 地形校正的前提是計(jì)算太陽局地入射角i的余弦值(cosi)(式(1))[23],據(jù)此運(yùn)行C模型(式(2))[6]、SCS+C模型(式(3))[10]、SE模型(式(4))[16]和MIN模型(式(5))[11]。
cosi=cosφZcosθS-sinφZsinθScos(φA-θA)
(1)
式中:φZ和φA分別為太陽天頂角和太陽方位角,θS和θA分別為坡度和坡向,由DEM經(jīng)空間計(jì)算得到。
(2)
(3)
(4)

(5)

1.3.3 DEM刻畫局地入射角精度評價(jià) 本文通過對原始地表反射率和不同開源DEM的cosi值進(jìn)行皮爾遜相關(guān)性分析,相關(guān)性最高的DEM即為能最精準(zhǔn)刻畫局地入射角、最準(zhǔn)確反映地表反射率隨地形變化規(guī)律的DEM數(shù)據(jù)集,并將該數(shù)據(jù)集作為評價(jià)地形校正效果的基準(zhǔn)DEM。
1.3.4 光譜地形校正效果評價(jià)
1)目視效果對比法。目視效果對比是評價(jià)地形校正效果常用的方法之一,通過觀察影像校正前后的整體色彩差異和陰影區(qū)域信息恢復(fù)判斷基于不同開源DEM和不同校正模型的效果[24]。
2)局地入射角引起光譜異常校正效果評價(jià)。由于植被對近紅外波段的變化更敏感,本文以草地為例,選用近紅外波段進(jìn)行比較。首先,基于地表覆蓋分類圖,在草地類型區(qū)域隨機(jī)選取10 000個(gè)像元,分別計(jì)算地形校正前和校正后的地表反射率與基準(zhǔn)DEM計(jì)算的局地太陽入射角余弦值(cosi)之間的皮爾遜相關(guān)性系數(shù)R,其絕對值越小,說明局地入射角對地表反射率的影響越小[13]。校正前后R的相對變化率RΔ(式(6))越大,表明對局地入射角所引起的光譜異常的校正效果越好。
RΔ=(|R|-|Rcorr|)/|R|×100%
(6)
式中:R、Rcorr分別為校正前、后地表反射率與cosi的相關(guān)系數(shù)。
3)同地類光譜變異性校正效果評價(jià)。基于地表覆蓋分類圖,在草地類型區(qū)域隨機(jī)選取10 000個(gè)像元,分別統(tǒng)計(jì)校正前和校正后地表反射率的四分位數(shù),然后使用上下四分位數(shù)之差(IQR)[4]評價(jià)草地類型的光譜變異性,其值越小,說明光譜變異越小,該方法可消除由影像光譜異常值所引起的地形校正效果偏差。校正前后IQR的相對變化率QΔ(式(7))越大,表明對同地類光譜穩(wěn)定性校正效果越好。
QΔ=(Q-Qcorr)/Q×100%
(7)
式中:Q、Qcorr分別為校正前、后的IQR。
4)坡向引起光譜異常校正效果評價(jià)。基于地表覆蓋分類圖和坡向數(shù)據(jù),在草地類型的陰坡區(qū)域(坡向=太陽方位角+180°±10°)和陽坡區(qū)域(坡向=太陽方位角±10°)分別隨機(jī)選取5 000個(gè)像元,然后計(jì)算陰坡(陽坡)校正前、后地表反射率之差D,差值越接近0,表明坡向?qū)Φ乇矸瓷渎实挠绊懺叫 PU啊⒑驞的相對變化率DΔ(式(8))越大,表明對坡向引起光譜異常的校正效果越好。
DΔ=(|D|-|Dcorr|)/|D|×100%
(8)
式中:D、Dcorr分別為校正前、后地表反射率的差值。
5)綜合評價(jià)方法。為客觀評價(jià)基于不同開源DEM和不同校正模型的地形校正效果,將上述3種校正效果評價(jià)指標(biāo)進(jìn)行算術(shù)平均,值最高的開源DEM和地形校正模型組合則為最佳地形校正方法。
從圖1a選取局部區(qū)域(紅色框線范圍),對原始的假彩色合成影像與AW3D、SRTM和GDEM的cosi進(jìn)行視覺對比(圖3a),可以看出 AW3D的紋理信息更精細(xì),與影像地表反射率明暗起伏更一致,SRTM次之,GDEM最模糊且缺失大量細(xì)節(jié)信息;提取圖3a中黃線上近紅外反射率和不同DEM的cosi,分析二者隨距離的變化規(guī)律,可明顯看出,AW3D的cosi變化趨勢與近紅外波段反射率最吻合(圖3b);對局部區(qū)域可見光—近紅外波段的地表反射率和cosi的相關(guān)系數(shù)進(jìn)行統(tǒng)計(jì)可知,AW3D的相關(guān)性在各波段均最高,SRTM次之,GDEM最低(圖3c)。綜上,AW3D刻畫局地入射角的精度最高,表征地表反射率對地形變化響應(yīng)規(guī)律的能力最強(qiáng),與Liu等[25]的結(jié)論一致,因此以AW3D作為基準(zhǔn)DEM,對后續(xù)不同校正模型的校正效果進(jìn)行評價(jià)。

圖3 局部影像與3種DEM的cosi對比
選取地形校正前后假彩色合成影像差異明顯區(qū)域(圖4綠色方框區(qū)域),可以看出,原始影像在未校正前凹凸不平,具有強(qiáng)烈的地形效應(yīng),經(jīng)過地形校正后影像變得平坦,立體感降低,陰影區(qū)域的亮度提升,光譜信息得到恢復(fù)。對比4種校正模型可知,C模型、SCS+C模型和SE模型的整體亮度與校正前基本一致,而MIN模型的亮度發(fā)生明顯改變,表明C、SCS+C和SE模型的校正效果較理想,MIN模型的校正效果較差;對比3種DEM發(fā)現(xiàn),使用AW3D校正后的陰影區(qū)域信息恢復(fù)效果最好,圖像的整體平整度最高,陽坡和陰坡的差異顯著降低,而使用GDEM校正后陰影區(qū)域信息恢復(fù)效果最差,圖像的立體感最強(qiáng),坡度起伏明顯,SRTM介于二者之間。
影像校正前后近紅外波段反射率與cosi的擬合結(jié)果如圖5所示,可以看出原始反射率隨局地入射角變化明顯,擬合斜率較大,表明原始影像存在明顯的地形效應(yīng)。經(jīng)過地形校正后, C模型、SCS+C模型和SE模型的校正結(jié)果接近,而MIN模型的校正效果偏差。在3種DEM中,基于AW3D校正后的散點(diǎn)擬合線斜率接近0,SRTM次之,而使用GDEM的校正效果不明顯,擬合線傾斜嚴(yán)重。

圖5 影像校正前后近紅外波段擬合結(jié)果對比
由圖6可知,基于AW3D校正的RΔ最大,其次為SRTM,GDEM的RΔ最小;4種校正模型中,SE模型和不同DEM組合的RΔ均最大,上述結(jié)果表明,AW3D結(jié)合SE模型對局地入射角引起的光譜異常的校正效果最好,減幅達(dá)98.2%。

圖6 研究區(qū)校正前后各波段地表反射率與cosi的相關(guān)系數(shù)相對變化率對比
由圖7可知,校正后的Q都有一定降低,3種DEM的校正效果排名為AW3D>SRTM>GDEM;4種校正模型中,SE模型的QΔ最大,SCS+C和C模型次之, MIN模型的QΔ最小,表明AW3D結(jié)合SE模型對同地類光譜變異性校正效果最好,減幅達(dá)97.4%。

圖7 研究區(qū)校正前后草地地表反射率的IQR相對變化率對比
由圖8可以看出,地形校正前近紅外波段反射率變化隨坡向變化明顯,兩者呈高度相關(guān),陽坡和陰坡處的反射率受地形效應(yīng)影響最明顯,陽坡反射率變大,陰坡反射率變小,其他坡向受地形影響較小,與Holben等[26]的結(jié)論一致。地形校正后,反射率隨坡向變化的程度減弱,陽坡反射率降低,陰坡反射率得到補(bǔ)償,反射率隨坡向分布更均勻。3種DEM中,AW3D校正后的反射率分布最均勻,校正效果最好,SRTM次之,GDEM的校正效果最差,陽坡和陰坡均出現(xiàn)不完全校正現(xiàn)象;4種校正模型中,MIN模型出現(xiàn)過校正現(xiàn)象,其他模型的校正效果基本一致。綜上,AW3D和SE模型組合對坡向引起的光譜異常校正效果最好,減幅達(dá)38.8%。

圖8 草地校正前后近紅外波段反射率隨坡向的變化
由草地校正前后DΔ的對比(圖9)可知,校正效果排名為AW3D>SRTM>GDEM,基于AW3D校正的DΔ最大,說明校正后陰坡和陽坡的差異幾乎完全消除;4種校正模型均減小了陰坡、陽坡反射率的差異,排名為C>SE>SCS+C>MIN。

圖9 研究區(qū)校正前后草地地表反射率之差的相對變化率對比
由表1可知,基于AW3D的SE模型對地形引起的光譜異常減幅最高(78.16%),基于GDEM的MIN模型減幅最低(20.25%)。此外,AW3D和不同地形校正模型結(jié)合的校正效果均較好,這與Wu等[27]的結(jié)論一致,DEM的質(zhì)量直接決定了地形校正后的效果,精度越高的DEM校正效果越好,AW3D的質(zhì)量最高,所以校正效果最好。在4種校正模型中,SE模型和不同DEM結(jié)合的校正效果均較好,與Sola等[4]的結(jié)論一致,SE模型由于引入了經(jīng)驗(yàn)參數(shù)和平均反射率,可以更好地消除地形效應(yīng)。

表1 不同開源DEM和地形校正方法的評價(jià)指標(biāo)平均值
遙感影像的地表反射率受地形起伏的強(qiáng)烈擾動(dòng),消除地形效應(yīng)至關(guān)重要。本文以黃土丘陵區(qū)為研究區(qū),基于3種DEM(AW3D、SRTM、GDEM)和4種校正模型(C、SCS+C、SE、MIN),開展不同組合試驗(yàn)對研究區(qū)進(jìn)行地形校正,并構(gòu)建地形校正效果評價(jià)指標(biāo)對不同組合校正效果進(jìn)行綜合評判。結(jié)果表明,AW3D刻畫局地入射角的精度最高,能真實(shí)反映衛(wèi)星觀測地表反射率對地形變化的響應(yīng)規(guī)律;AW3D和不同地形校正模型組合的校正效果均較好,SRTM次之,GDEM最差;與其他校正模型相比,SE模型和不同開源DEM組合的校正效果均較好;AW3D和SE模型組合為適用于黃土丘陵區(qū)衛(wèi)星地表反射率地形校正的最佳組合。
不同地形校正模型消除地形效應(yīng)的能力一直存在爭議,需要多種評價(jià)標(biāo)準(zhǔn)來綜合評價(jià)校正模型的校正效果[13],但不同評價(jià)指標(biāo)的物理意義不同,量綱也不一致,不能直接對所有指標(biāo)進(jìn)行綜合計(jì)算,無法體現(xiàn)校正方法的綜合校正效果。針對該問題,本文采用相對變化率的形式對3項(xiàng)指標(biāo)統(tǒng)一單位(%),通過求算術(shù)平均值實(shí)現(xiàn)綜合評價(jià)。然而,本文仍存在如下不足:①黃土丘陵溝壑區(qū)具有塬、墚、峁等特殊地貌,本研究僅在該區(qū)域進(jìn)行試驗(yàn),缺乏不同復(fù)雜山區(qū)地形綜合對比分析,未來將考慮加入對不同山區(qū)地貌的研究分析;②遙感影像地表反射率受地表覆蓋類型、植被參數(shù)、土壤背景、地形特征等多種因素綜合影響,本文隨機(jī)選取樣本點(diǎn)代表研究區(qū)所有可能的植被參數(shù)、土壤背景、地形特征,導(dǎo)致實(shí)驗(yàn)過程中產(chǎn)生不確定性,未來需要考慮更多的影響因素,進(jìn)行更全面的驗(yàn)證。