潘建平,趙瑞淇,蔡卓言,袁雨馨,李鵬霞

(重慶交通大學(xué)智慧城市學(xué)院,重慶 400074)

Sentinel-1A數(shù)據(jù)從2014年10月對(duì)全球用戶開(kāi)放以來(lái)[1],已經(jīng)被廣泛應(yīng)用于地震與地質(zhì)構(gòu)造監(jiān)測(cè)[2-3]、火山與冰川監(jiān)測(cè)[4-5],以及其他形變監(jiān)測(cè)領(lǐng)域[6-8],現(xiàn)已成為合成孔徑雷達(dá)干涉測(cè)量(interferometric synthetic aperture radar,InSAR)領(lǐng)域內(nèi)重要的C波段數(shù)據(jù)源(波長(zhǎng)為5.6 cm)。但是其短波長(zhǎng)的特性給其應(yīng)用帶來(lái)了諸多限制,特別是在植被區(qū),C波段的雷達(dá)信號(hào)通常難以穿透植被的表層或頂層(如樹(shù)枝和樹(shù)葉),衛(wèi)星重訪周期內(nèi)植被的理化特性變化,甚至是枝葉的微小運(yùn)動(dòng),都會(huì)導(dǎo)致干涉失相干現(xiàn)象的發(fā)生,并將誤差引入InSAR地表形變監(jiān)測(cè)。

以往的研究表明,植被區(qū)的失相干程度可以在兩景單視復(fù)數(shù)影像(single look complex,SLC)發(fā)生干涉后由相干系數(shù)(coherence)進(jìn)行評(píng)估[9],相干系數(shù)反映兩幅SLC影像的相似程度,可以衡量干涉條紋的質(zhì)量與SLC影像中像素間的相位和振幅信息的變化[10]。但是干涉后評(píng)估必然會(huì)導(dǎo)致SAR數(shù)據(jù)存在選擇的盲區(qū),使得研究者不能選擇最有效的干涉對(duì),從而降低InSAR對(duì)地表形變監(jiān)測(cè)的效率,增加了時(shí)間與經(jīng)濟(jì)成本。對(duì)于短波長(zhǎng)的Sentinel-1A數(shù)據(jù),如果失相干能在干涉前得到精確評(píng)估,將十分有助于克服上述干涉后評(píng)估的缺點(diǎn),從而提高監(jiān)測(cè)效率。

作為應(yīng)用最廣泛的植被指數(shù),歸一化差分植被指數(shù)(normalized difference vegetation index,NDVI)與植被覆蓋度密切相關(guān),可用于監(jiān)測(cè)植被的動(dòng)態(tài)變化,對(duì)于地表覆被的可靠表征為建立其與InSAR失相干之間的定量模型提供了可能。

以往針對(duì)NDVI與Sentinel-1A相干系數(shù)的研究主要集中于定性說(shuō)明兩者的關(guān)系[11-12],其定量關(guān)系尚不明確。基于此,本文選擇位于四川的一個(gè)植被覆蓋區(qū)作為研究區(qū),在同極化(VV)與交叉極化(VH)模式下利用免費(fèi)光學(xué)影像Landsat 8計(jì)算得到的NDVI,分別建立其與Sentinle-1A相干系數(shù)之間的定量模型,并在鄰近區(qū)域?qū)δＰ瓦M(jìn)行驗(yàn)證。

1 研究區(qū)概況與數(shù)據(jù)預(yù)處理

1.1 研究區(qū)概況

本文的研究區(qū)位于四川省,其西部為川西高原,東部緊鄰成都平原,面積約為340 km2,如圖1紅色實(shí)線框所示。研究區(qū)內(nèi)植被覆蓋茂盛,主要分布有林地和耕地(主要作物為水稻),并伴隨有少量村落分布,屬于InSAR地表形變監(jiān)測(cè)中的低相干地區(qū),在該區(qū)域建立模型有助于提高模型的普適性與應(yīng)用價(jià)值;同時(shí)選擇位于研究區(qū)西南的鄰近區(qū)域(圖1黑色實(shí)線框所示)作為驗(yàn)證區(qū),以驗(yàn)證模型的可靠性,以及區(qū)內(nèi)地形和植被覆蓋狀況與研究區(qū)類似。

圖1 研究區(qū)位置與影像覆蓋范圍

1.2 數(shù)據(jù)預(yù)處理

本文使用的光學(xué)影像為L(zhǎng)andsat 8 OLI影像,成像時(shí)間為2019年8月11日,云量為0.89%,空間分辨率為30 m,覆蓋范圍如圖1中黑色虛線框所示。對(duì)光學(xué)影像進(jìn)行大氣校正和輻射校正兩步處理,并按照研究區(qū)和驗(yàn)證區(qū)矢量范圍進(jìn)行裁剪,最后計(jì)算NDVI,公式為

(1)

式中,PNIR、PRED分別代表近紅外和紅光波段的反射率。NDVI的值介于[-1,1],隨著地表植被覆蓋度的增加而增大。NDVI為負(fù)值時(shí)代表地表覆蓋有水或雪,NDVI為0表示地表為巖石或裸土,NDVI為正值表示地表有植被覆蓋。研究區(qū)和驗(yàn)證區(qū)的NDVI分布如圖2所示。

圖2 研究區(qū)和驗(yàn)證區(qū)NDVI分布

Sentinel-1A影像分別為成像時(shí)間為2019年7月3日和8月20日的干涉寬幅模式下的SLC影像,時(shí)間基線為48 d,分別作為干涉對(duì)的主影像和輔影像,其覆蓋范圍如圖1中紅色虛線框所示。SLC影像的預(yù)處理包括3部分: ①使用雙過(guò)差分對(duì)主、輔影像進(jìn)行差分干涉,并分別使用精密軌道數(shù)據(jù)和SRTM DEM(30 m)消除軌道誤差和模擬地形相位; ②對(duì)干涉對(duì)進(jìn)行多視處理以消除斑點(diǎn)噪聲并使SAR影像與Landsat 8影像保持相同的分辨率;③使用文獻(xiàn)[13]提出的公式消除空間失相干的影響,并計(jì)算研究區(qū)和驗(yàn)證區(qū)的相干系數(shù),公式為

(2)

式中,γ代表相干系數(shù);E[·]代表數(shù)學(xué)期望;μ*代表SLC影像的共軛復(fù)數(shù)。相干系數(shù)的值介于[0,1],當(dāng)主輔影像完全相干時(shí),相干系數(shù)的值為1;當(dāng)兩者完全失相干時(shí)(如干涉圖完全被噪聲污染),其值為0。

最后使用SRTM DEM(30 m)對(duì)NDVI分布圖和相干系數(shù)圖進(jìn)行配準(zhǔn),配準(zhǔn)后完成數(shù)據(jù)預(yù)處理,數(shù)據(jù)預(yù)處理流程如圖3所示。

圖3 數(shù)據(jù)預(yù)處理流程

2 試驗(yàn)方法

2.1 時(shí)間失相干對(duì)相干系數(shù)的影響

文獻(xiàn)[14]提出了對(duì)相干系數(shù)產(chǎn)生影響的各個(gè)因素,主要有空間失相干、時(shí)間失相干及熱噪聲失相干。其中,熱噪聲失相干可以忽略不計(jì);空間失相干由于研究區(qū)較為平坦的地形(平均坡度小于15°)并不明顯,且已通過(guò)文獻(xiàn)[13]提出的公式消除;通過(guò)建立時(shí)間基線與相干系數(shù)的模型削弱時(shí)間失相干對(duì)相干系數(shù)的影響。

基于多主影像策略,利用覆蓋研究區(qū)的Sentinel-1A影像生成了若干干涉對(duì),干涉對(duì)的時(shí)間基線為12～696 d,時(shí)間基線與相干系數(shù)的關(guān)系如圖4所示。

圖4 時(shí)間基線與相干系數(shù)的關(guān)系

圖4(a)顯示,VV極化下相干系數(shù)在時(shí)間基線≤216 d時(shí)隨著時(shí)間基線的增加呈現(xiàn)指數(shù)衰減,其擬合曲線為y=A1·exp(-x/t1),其中A1=0.743,t1=206 d。在216 d之后沒(méi)有明顯的衰減趨勢(shì),相干系數(shù)在0.370左右浮動(dòng)。圖4(b)顯示,VH極化下相干系數(shù)在時(shí)間基線小于168 d時(shí)呈現(xiàn)指數(shù)衰減,其擬合曲線為y=A2·exp(-x/t2),其中A2=0.560,t2=222 d。在168 d之后沒(méi)有明顯的衰減趨勢(shì),在0.325左右浮動(dòng)。本文所使用的干涉對(duì)時(shí)間基線為48 d,因此有必要考慮時(shí)間失相干對(duì)相干系數(shù)指數(shù)衰減的作用。

2.2 基于皮爾遜相關(guān)系數(shù)絕對(duì)值的建模方法

試驗(yàn)表明,將配準(zhǔn)后的NDVI分布圖與相干系數(shù)圖的所有像素進(jìn)行擬合分析會(huì)引入非常大的全局誤差,從而無(wú)法建立起NDVI和相干系數(shù)之間的可靠模型,因此本文考慮使用文獻(xiàn)[15]中提出的窗口采樣方法。該方法包括3個(gè)步驟。

(1)對(duì)NDVI分布圖設(shè)置一個(gè)采樣窗口,同時(shí)對(duì)相干系數(shù)圖設(shè)置一個(gè)相同大小的采樣窗口,使用采樣窗口對(duì)兩張圖像同時(shí)進(jìn)行采樣,直至窗口遍歷整張圖像。

(2)對(duì)兩個(gè)窗口中的NDVI像素值和相干系數(shù)像素值計(jì)算皮爾遜相關(guān)系數(shù)[16],公式為

(3)

(3)設(shè)置一個(gè)閾值,若兩個(gè)采樣窗口中的像素值皮爾遜相關(guān)系數(shù)大于該閾值,則將窗口中的像素值保留,參與后續(xù)分析;否則,將窗口中的像素值舍棄。

該方法被用于建立NDVI與L波段ALOS數(shù)據(jù)相干系數(shù)之間的關(guān)系,但是ALOS數(shù)據(jù)由于其長(zhǎng)波段特性失相干現(xiàn)象較弱,NDVI與相干系數(shù)之間的負(fù)相關(guān)關(guān)系并不明顯。初步試驗(yàn)表明,Sentinel-1A數(shù)據(jù)的相干系數(shù)與NDVI之間呈現(xiàn)較強(qiáng)的負(fù)線性相關(guān)關(guān)系,直接使用該方法將會(huì)忽略掉大量像素值,導(dǎo)致模型精度大幅降低。因此,本文對(duì)該方法的相似性測(cè)度進(jìn)行了改進(jìn),在步驟(3)中計(jì)算皮爾遜相關(guān)系數(shù)的絕對(duì)值,以顧及采樣窗口中NDVI與Sentinel-1A相干系數(shù)大量呈負(fù)線性相關(guān)的像素值。最佳采樣窗口大小和預(yù)設(shè)閾值大小的設(shè)定需要通過(guò)試驗(yàn)確定,改進(jìn)后的窗口采樣方法如圖5所示。

圖5 改進(jìn)相似性測(cè)度的窗口采樣方法

圖5上下兩張圖像分別為研究區(qū)NDVI分布圖和相干系數(shù)圖,T為皮爾遜相關(guān)系數(shù)的絕對(duì)值,T0為預(yù)設(shè)的閾值。大量試驗(yàn)表明,對(duì)于VV極化,最佳的采樣窗口大小為5×5,閾值大小為0.7,VH極化下的最佳采樣窗口和閾值大小為9×9和0.6。將滿足閾值條件的保留像素值進(jìn)行分析,得到了VV和VH極化下NDVI和相干系數(shù)的關(guān)系,如圖6所示。

圖6 NDVI和相干系數(shù)的關(guān)系

圖6中藍(lán)色星星為真實(shí)像素值,兩條直線分別為VV和VH極化下的擬合曲線。由于時(shí)間失相干對(duì)相干系數(shù)的衰減作用,因此需要在模型中加入時(shí)間失相干指數(shù)衰減因子,在95%的置信區(qū)間下采用最小二乘擬合得到了VV和VH極化模式下NDVI和Sentinel-1A相干系數(shù)的定量模型,分別為

(4)

(5)

3 試驗(yàn)結(jié)果與討論

3.1 試驗(yàn)結(jié)果與誤差分析

基于式(4)和式(5),利用NDVI分布圖模擬得到了VV和VH極化下的研究區(qū)相干系數(shù)圖,與研究區(qū)真實(shí)相干系數(shù)的對(duì)照如圖7所示。

圖7 研究區(qū)相干系數(shù)

圖7(a)和圖7(c)分別為研究區(qū)VV和VH極化下的真實(shí)相干系數(shù),圖7(b)和圖7(d)分別為VV和VH極化下的模擬相干系數(shù)。利用真實(shí)相干系數(shù)減去模擬相干系數(shù)得到了圖7(e)和圖7(f)所示的誤差分布圖。由誤差分布圖可以看出,VV極化下誤差分布除研究區(qū)內(nèi)河流分布區(qū)域外其他區(qū)域不明顯,VH極化下誤差分布除了河流區(qū)域外主要集中在圖7(f)中的村落分布區(qū)域。對(duì)誤差分布進(jìn)行統(tǒng)計(jì),如圖8所示。研究區(qū)VV極化下誤差均值為-0.037,標(biāo)準(zhǔn)差為0.205,大部分誤差集中在-0.3～0.3之間,誤差整體基本服從正態(tài)分布。研究區(qū)VH極化下誤差均值為0.045,標(biāo)準(zhǔn)差為0.200,大部分誤差同樣集中分布在-0.3～0.3之間,雖然誤差整體高于VV極化,但是仍然基本服從正態(tài)分布。

圖8 研究區(qū)誤差統(tǒng)計(jì)

3.2 模型驗(yàn)證

為了驗(yàn)證模型的可靠性,在驗(yàn)證區(qū)進(jìn)行了相同的試驗(yàn),數(shù)據(jù)處理方式與研究區(qū)保持一致,驗(yàn)證區(qū)的相干系數(shù)如圖9所示。圖9(a)和圖9(c)分別為驗(yàn)證區(qū)VV和VH極化下的真實(shí)相干系數(shù),圖9(b)和圖9(d)分別為VV和VH極化下的模擬相干系數(shù),利用真實(shí)相干系數(shù)減去模擬相干系數(shù)得到了圖9(e)和圖9(f)所示的誤差分布。由于驗(yàn)證區(qū)分布有兩條明顯的河流,因此由誤差分布圖可以看出,驗(yàn)證區(qū)的誤差主要分布在河流分布區(qū)域及區(qū)域南部的村落分布區(qū)。對(duì)驗(yàn)證區(qū)誤差分布進(jìn)行統(tǒng)計(jì),如圖10所示。驗(yàn)證區(qū)VV極化下的誤差均值為-0.067,標(biāo)準(zhǔn)差為0.256,VH極化下的誤差均值為-0.065,標(biāo)準(zhǔn)差為0.230。兩種極化下的誤差均集中分布在-0.3～0.3之間,且服從正態(tài)分布。綜合圖7(e)、圖7(f)、圖9(e)和圖9(f)的差值及圖8和圖10的誤差分布情況來(lái)看,利用式(4)和式(5)模擬得到的相干系數(shù)沒(méi)有隨機(jī)誤差分布,模型可靠。

圖9 驗(yàn)證區(qū)相干系數(shù)

圖10 驗(yàn)證區(qū)誤差統(tǒng)計(jì)

3.3 結(jié)果討論

從研究區(qū)和驗(yàn)證區(qū)的角度來(lái)看,兩種極化模式下驗(yàn)證區(qū)模擬相干系數(shù)的誤差分布更離散且整體大于研究區(qū),主要原因有兩點(diǎn):①盡管驗(yàn)證區(qū)的地形與研究區(qū)相似,但仍存在微小差異,雖已去除空間失相干的影響,但是微小地形差異引起的空間失相干會(huì)將誤差引入模型中;②相較于研究區(qū),驗(yàn)證區(qū)分布有兩條明顯的河流,河流區(qū)域的誤差成為驗(yàn)證區(qū)的主要誤差來(lái)源。

從極化方式來(lái)看,無(wú)論研究區(qū)還是驗(yàn)證區(qū),VV極化下的相干系數(shù)整體均高于VH極化下的相干系數(shù),這是由于研究區(qū)和驗(yàn)證區(qū)的地表覆被主要為耕地(主要作物為水稻),其次是林地,已有研究表明水稻的理化參數(shù)與反射特性與VH極化的關(guān)系比VV極化更強(qiáng)[16],導(dǎo)致其對(duì)VH極化的雷達(dá)波更敏感,相同時(shí)間基線下的失相干相較于VV極化更嚴(yán)重。同時(shí)研究區(qū)和驗(yàn)證區(qū)VV極化下的整體誤差均低于VH極化下的整體誤差,這是由于盡管經(jīng)過(guò)大量試驗(yàn)得到了兩種極化方式下的最佳采樣窗口與閾值,但是VH極化下的像素值相較于VV極化分布更離散(如圖6所示),整體擬合精度低于VV極化,從而導(dǎo)致VH極化下模型精度低于VV極化下模型精度。

綜合上述模型的誤差源與缺點(diǎn),在未來(lái)的研究中還需考慮如何完全去除空間失相干對(duì)模型的影響,另外水系不可避免地會(huì)出現(xiàn)在植被覆蓋區(qū),因此應(yīng)考慮如何消除水系對(duì)模型的影響,如建立光學(xué)遙感水體指數(shù)與相干系數(shù)的定量模型并融入已建立的模型),從而提高模型的普適性與可靠性。其他諸如研究區(qū)的天氣狀況、地表植被類型及影像不同的成像視角(升軌與降軌)等因素也應(yīng)當(dāng)在模型中加以考慮。除了本文建立的NDVI與C波段Sentinel-1A相干系數(shù)之間的定量模型之外,NDVI與同極化L波段的ALOS-1/PALSAR-1相干系數(shù)之間的定量模型已經(jīng)建立[15],因此在未來(lái)的工作中需要建立NDVI與交叉極化的L波段數(shù)據(jù)的相干系數(shù)及與X波段數(shù)據(jù)的相干系數(shù)之間的定量模型。

4 結(jié) 語(yǔ)

本文建立了兩個(gè)基于Landsat 8影像計(jì)算得到的NDVI與Sentinel-1A不同極化方式下干涉相干系數(shù)之間的二階線性模型,從而得到了InSAR失相干與地表植被覆蓋之間的定量關(guān)系。相干性隨著植被覆蓋度的增加而呈線性減小,且線性減小的趨勢(shì)因極化方式的不同而有所區(qū)別,同時(shí)這兩個(gè)模型在驗(yàn)證區(qū)得到了驗(yàn)證,表明了模型的可靠性。基于本文建立的模型,可以利用免費(fèi)光學(xué)影像計(jì)算得到的NDVI在進(jìn)行InSAR處理之前(即SLC影像干涉之前)對(duì)研究區(qū)的InSAR失相干現(xiàn)象做出定量評(píng)估;同時(shí)在干涉前對(duì)SAR數(shù)據(jù)選擇提供指導(dǎo),從而提升研究者的效率并減少時(shí)間與經(jīng)濟(jì)成本。