不同地表覆蓋對InSAR技術的相干性影響研究

周鑫城,左小清,李勇發,楊 栩,鄧云龍,周家厚

(昆明理工大學 國土資源工程學院,云南 昆明 650093)

為了解地表的沉降情況,眾多學者采用持續性的形變監測技術對地表進行觀測,從而揭示地表的活動狀態,預測與防范地質災害的發生[1]。當前,常規的地表形變監測技術主要有:水準測量[2]、全球衛星導航定位技術(global positioning system,GPS)[3]等。近些年來,隨著地表監測技術水平和科學技術的發展,諸如全球導航衛星系統(global navigation satellite system,GNSS)[4]、三維激光掃描[5]、地理信息系統(geographic information system,GIS)[6]等高新技術的引入,使得地表形變監測的工作效率有顯著的提高。這些技術的應用,為區域的沉降情況提供了可靠的數據基礎,但存在一定的缺陷,如空間分辨率低、觀測周期長和監測成本高等問題。隨著雷達衛星(synthetic aperture radar,SAR)逐步發育成熟,具有高時空分辨率、不受云雨條件限制、覆蓋范圍廣、數據獲取周期較短、安全性高等優勢的干涉合成孔徑雷達(interferometric synthetic aperture radar,InSAR)成為地表形變監測的主要手段之一[7-9]。相干性[10-11]是干涉測量的基礎,不僅直接影響形變監測數據的處理難度,而且失相干也是地表形變誤差的主要來源。InSAR失相干[12]包括時間失相干和空間失相干,而失相干主要是由時間基線或空間基線過長引起的。目前,多數學者從定性的角度進行研究,發現InSAR相干性與時空基線的長度呈負相關[13-14]。歐洲航天局的哨兵雙星系統的出現極大地克服了大部分空間失相干所帶來的誤差,而引起時間失相干的主要原因是研究區分辨單元內目標的物理和化學特性隨時間的變化而變化。

當前,研究不同地表對相干性的影響較少,余祥偉等[15]分析植被及坡度對InSAR相干性的影響,發現相較于坡度,植被是導致高覆蓋山區InSAR影像時間失相干的主要因素,相干性和植被覆蓋度之間存在一定的相關性。佀文娜等[16]研究典型地物的時間失相干影響,發現風速對相干性沒有影響,而降雨量是導致地物失相干的重要因素。王天祥等[17]使用D-InSAR技術對地表進行形變監測,在獲取大量數據后進行仿真試驗,發現空間基線和時間基線與InSAR相干性之間存在指數關系,而且當空間基線的長度不超過臨界基線的三分之一時,可以保證得到較好的InSAR相干性。

現有研究主要對植被、坡度、降雨量等多方面進行綜合分析,尚未針對不同地表覆蓋對InSAR的相干性影響進行深入研究。因此,本文選取昆明市滇池及周邊環境作為研究區,采用小基線集(small baselines subset InSAR,SBAS-InSAR)技術對歐洲航天局哥白尼計劃(global monitoring for environment and security,GMES)中的地球觀測衛星哨兵1號(Sentinel-1)衛星數據進行處理。從定量角度分析不同地表覆蓋對相干性的影響,深入分析在不同月份、不同季節和不同降雨量下,地表覆蓋類型對相干性的影響,為提高InSAR形變精度提供有力的理論依據和數據支撐。

1 研究方法

1.1 SBAS-InSAR技術

SBAS-InSAR技術是在InSAR技術基礎上發展起來的一種多主影像的時序InSAR方法。該技術的測量原理是利用兩幅具有相同觀測幾何的影像構成干涉像對,結合衛星和觀測目標之間的幾何關系,對干涉像對所產生的相位差進行相位濾波、相位解纏、地理配準等步驟,獲取觀測目標的三維信息和形變信息[18]。

在評價干涉質量的好壞時,分析相位相干性常被作為一個重要的標準。相干性被定義為一個干涉像對中主、副影像是否存在相似,而兩者之間的相似程度目前常使用相干系數進行描述。相干系數(γ)數值分布于[0,1]區間內,相干系數為0表示主、副影像之間無相似程度,相干系數為1表示主、副影像之間完全相同。相干系數(γ)計算模型如下:

式中,N和M為最終獲取相干性像元的尺寸大小;n和m為數據內像元的行列號;μ*為共軛復數;μ1(n,m),μ2(n,m)為主、副影像數據內像元再坐標(n,m)處的復數值;|·|為復數的絕對值。

由于雷達兩次成像期間觀測區域內目標的物理和化學特性隨時間變化而變化,使得相干性在目標條件不穩定的區域通常較低,而在目標條件穩定的區域(建筑區域及工礦用地)相干性通常較高。

1.2 歸一化水指數

歸一化水指數是用遙感影像中的特定波段進行歸一化差值處理,以突顯影像中水體的信息[19-20]。該指數的主要優點之一是能夠很好的區分水體區與非水體區,其利用水體吸收的綠光(492～577 nm)和反射近紅外(780～2 526 nm)的原理建立表達式,計算公式如下:

(1)

式中,B3為遙感影像中的綠波段;B5為遙感影像中的近紅外波段。INDWI的取值范圍為[-1,1]。

1.3 垂直不透水面指數

垂直不透水面指數是利用所選樣本的光譜特性建立較高的關聯表達式,以凸顯影像中不透水面(建筑區域及工礦用地)的信息[21-23]。該指數相較于歸一化建筑指數能夠較好的區分不透水面和透水面,計算公式如下:

(2)

式中,a為選取樣本構造出的關聯表達式的斜率;b為選取樣本構造出的關聯表達式的截距;B2為遙感影像中的藍波段。

1.4 歸一化植被指數

歸一化植被指數是反應植被生長趨勢的重要參數之一,被廣泛應用于干旱季的監測、農作物生產情況的監測和預測等[24-25]。該指數的主要優點之一是能夠較好的區分植被區域與非植被區域,其利用植被吸收紅光(625～740 nm)和反射近紅外(780～2 526 nm)的原理建立表達式,計算公式如下:

(3)

式中,B4為遙感影像中的紅波段。INDVI的取值范圍為[-1,1]。

2 研究區與數據集

2.1 研究區

云南省昆明市總體的地貌趨勢為北部高,南部低,呈階梯式降低,高程主要分布在1 500～2 800 m之間;以亞熱帶高原季風氣候為主,因此,植被以常綠闊葉林為主,年平均氣溫為15℃左右,年日照平均時長為2 200 h左右,年平均降水量為1 035 mm。經全國第三次全國國體調查,昆明市耕地面積約為3 865.73 km2;林地面積約為11 317.67 km2;草地面積約為1 244.40 km2;建筑區域及工礦用地面積約為1 480.53 km2。本次研究區包括云南省昆明市滇池及周邊區域,包括安寧區的東部地區,地處中國西南地域、云貴高原中部地區,其覆蓋研究區覆蓋面積約為2 287 km2,位于東經:102°28′05″E—102°57′40″E,北緯:24°36′34″N—25°13′52″N,如圖1所示。

圖1 研究區位置Fig.1 Location of study area

2.2 數據集

本文使用的數據集包括從歐空局(https://search.asf.alaska.edu/)獲取2020年1月11日—2021年2月22日的35景降軌Sentinel 1A數據,該數據的波段為C波段,波長為5.6 cm,分辨率為5 m×20 m,觀測模式為干涉寬幅 (interferometric wide swath,IW),極化方式為VV(單極化),單視復數(single look comple,SLC)。土地利用分類數據來源于中國科學院資源環境科學與數據中心數據平臺 (http://www.resdc.cn)。從地理空間數據云(http://www.gscloud.cn/)網站下載得到的2020年2月5日Landsat 8影像。2020年全年MODIS系列NDVI產品的原始數據是從Nasa(https://www.nasa.gov/)網站下載得到的2020年全年MODIS數據。降雨數據是從中國地面氣候資料日值數據集中獲取。

2.3 數據預處理

本文通過使用SBAS-InSAR技術對35景Sentinel 1A降軌數據進行處理,將SBAS-InSAR的時間基線閾值設置為120 d,經過數據導入、研究區裁剪、地理配準、反演等步驟,最終得到204幅相干性圖,并將獲取相干性的兩個日期對相干性圖進行命名,每幅影像的分辨率為17 m×17 m,地理坐標為WGS_1984。

本文對地理空間數據云中2020年全年的影像進行篩選,最終使用2020年2月5日的Landsat 8影像(云量小于5%)作為研究影像,利用光譜特性與三種指數(INDWI、IPII和INDVI)的計算公式,計算出研究區內各像元在2020年2月5日的三種指數值,其分辨率為30 m*30 m,投影坐標為WGS_1984_48N。

為研究不同類型的地表覆蓋與InSAR相干性之間的關系,以三種指數的像元分辨率和坐標為準,對相干性數據進行匹配、投影轉換、重采樣、掩膜提取等步驟使4種數據的行列號保持一致。對各類型的地表覆蓋隨機選擇1 200個點,分析所選點的指數數據與相干性之間的關系,以確定各類型的地表覆蓋是否存在研究的意義,并進一步了解不同地表覆蓋對相干性的影響。

MODIS系列的歸一化植被指數產品在獲取時間序列上的間隔為16 d,因此,需對該數據進行預處理,獲取每月的月度歸一化植被指數數據。降雨量是日值數據,同樣也要對此進行月份劃分,獲取月度降雨量數據。相干性數據是由每兩幅不同時間段的影像計算得到的,而本文所使用的時間基線為120 d,因此,會出現跨月時間影像生成的相干性圖,針對該現象,本文將命名中較早的日期作為分類依據,進行歸納各月的相干性圖,以方便計算各月的相干性平均值(如:2020年3月11日的影像和2020年4月4日的影像所計算出來的相干性圖歸納為三月份的相干性圖,2020年3月11日的影像和2020年6月3日的影像所計算出來的相干性圖也歸納為三月份的相干性圖)。最終結合降雨量數據與季節數據,研究分析在降雨量與季節變化的條件下,不同地表覆蓋對InSAR相干性的影響。

3 結果與分析

3.1 三種指數結果與分析

圖2(a)為求取平均值后的相干性圖,從該圖可知,研究區內的相干性取值范圍為0～1,其中,建筑區域的相干性最高,其次是植被覆蓋區域的相干性,水體的相干性最低。本文的土地利用一級類別有五類,分別為:水域、林地、草地、耕地和建筑區域及工礦用地,其二級分類如圖2(b)所示。統計一級類別的像元,可以發現建筑區域及工礦用地的像元最多,其次為林地。

(a)求取平均值后相干性圖 (b)二級分類圖2 相干性與土地利用分類圖Fig.2 Coherence and land use classification map

圖3(a)為研究區的INDWI指數圖,由圖可知,其INDWI取值范圍為-0.61～1。圖3(b)為研究區的IPII指數圖,IPII指數的獲取與INDWI指數不同。在影像中選取暗色不透水面(居民建筑區等)、亮色不透水面(工廠等)、藍色不透水面(廠房等)、透水面(林地、草地、耕地等)5種類別作為樣本,根據樣本生成的特征空間,計算出(2)式中的a和b分別為1.514和2 118.066,因此,(2)式可以表示為

圖3 INDWI、IPII、INDVI指數圖Fig.3 INDWI, IPII and INDVI index diagram

由圖3(b)可知,IPII數值的取值范圍為-7 299.6～8 065.1。圖3(c)為研究區的INDVI影像圖,其取值范圍為-1～0.82。

圖3中,三種指數可以明顯的區分出水體區域、建筑區域和植被覆蓋區域。為定量分析不同地表覆蓋與相干性之間的關系,利用隨機點繪制散點圖,并通過最小二乘法進行擬合,最終得到三種指數與相干性之間的關聯程度圖,如圖4所示。由于每種地表覆蓋類型的光譜特性不同,用不同的指數表達不同地表覆蓋可以更好的體現地表覆蓋的差異性,因此,本文對不同地表覆蓋的隨機點賦予不同的指數。INDWI指數能夠很好的區分水體區域與非水體區域,所以對水體類型的隨機點賦予INDWI數據;IPII指數主要用于區分透水面和不透水面,而透水面主要是植被覆蓋區域,不透水面主要是建筑區域,因此,對兩者區域賦予IPII數值能夠進一步了解兩者與相干性之間的關系,由于建筑區域及工礦用地的建設需符合國家綠化面積要求,其周邊附著一定的植被及綠化帶,而在農村居民點及城鎮周邊,多為植被覆蓋區與建筑區的交匯,影像中存在混淆像元,因此,為減少此方面的誤差,對建筑區域及工礦用地與植被覆蓋區的隨機點賦予INDVI數值。

圖4 INDWI、IPII、INDVI指數與相干性的關聯程度Fig.4 Correlation degree of INDWI, IPII, INDVI index and coherence

從整體上看,圖4中存在一些離散點,出現離散點的主要原因是三種指數由影像通過波段差計算得到的,而相干性的計算方式是由兩幅Sentinel 1A影像通過干涉處理得到的,在影像源的使用和計算方法上存在差異;各類型地表覆蓋的隨機點是在土地利用的分類基礎上進行選擇的,而土地利用在分類上存在誤差,因此,在數值上會出現異常值,在圖像上出現無規律分布的離散點。如圖4(a)所示,結合光學遙感影像和土地利用類型可以發現,該點的利用類型為水體,而在遙感影像中,該點位于昆明南繞城高速上,且從相干性和INDWI的數值大小可以斷定,該點屬于異常點。因此,研究三種指數與相干性之間的關系前,需對離散點進行剔除。

圖4(a)為INDWI指數相關性圖,圖4(b)為IPII指數相關性圖,圖4(c)為INDVI指數相關性圖。由圖4(a)可知,INDWI數值集中分布在0.75～0.85之間,其對應的相干性分布在0.05～0.58之間,兩者的相關系數R2為0.180 6,呈負相關。且從圖中的點分布可知,INDWI數值與相干性之間無關聯程度,故研究水體對相干性的影響沒有意義,后續將不對水體進行研究分析。

由圖4(b)可知,不透水面(建筑區域及工礦用地)的IPII數值主要分布在0～2 000之間,其對應的相干性主要分布在0.4～0.9之間;透水面(植被覆蓋區)的IPII數值主要分布在-2 000～0之間,其對應的相干性數值主要分布在0.4～0.8之間,兩者的相關系數R2為0.180 3,呈正相關。導致IPII指數與InSAR相干性之間的關聯程度低有兩個原因,一是本文所選取的隨機點多為植被覆蓋類型,所以在圖中點分布較為密集,呈現出上部分稀疏下部分密集的分布情況;二是IPII數值取值范圍遠大于InSAR相干性。雖然在圖中呈現出隨著IPII數值的增加,相干性增加,但是兩者之間的關聯程度低,無法顯著的體現出不同地表對相干性的影響。

從圖4(c)可知,植被覆蓋區域的INDVI數值主要分布在0～0.55之間,相干性主要分布在0.4～0.8之間;建筑區域及工礦用地的INDVI數值主要分布在-0.2～0之間,相干性主要分布在0.6～0.98之間,兩者的相關系數R2為0.399 8,呈負相關。結合INDVI指數特性,不同地表覆蓋的INDVI取值范圍不同,即可以通過INDVI指數區分不同地物。故本文將通過使用INDVI作為主要的研究指數,研究不同地表覆蓋對相干性的影響。

3.2 不同地表覆蓋類型對相干性影響分析

為研究不同地表覆蓋對相干性的影響,本文對4種地表覆蓋類型的1 200個隨機點賦予2020年月度平均INDVI數值和月度平均相干性,將隨機點的INDVI數值以0.2的步長進行劃分,繪制各類型的地表覆蓋在不同INDVI數值下相干性的分布圖,如圖5所示。

圖5 不同地表覆蓋類型的相干性分布圖Fig.5 Coherence distribution map of different surface Types

對比不同地表覆蓋的整體相干性可以發現,建筑區域及工礦用地的整體相干性較高,植被覆蓋區的整體相干性較低。由圖5可知,建筑區域及工礦用地的相干性在INDVI數值為0.34～0.57的區間內下降程度最大,而在0.18～0.34區間內,相干性的衰減程度平穩,由此可見,在建筑區域及工礦用地中,相干性受INDVI數值變化影響顯著,特別是在INDVI數值超過0.34之后,隨著INDVI數值增加,相干性以較大的趨勢下降。林地的相干性在INDVI數值為0.35～0.66的區間內快速下降,而在0.27～0.35和0.66～0.81的區間內下降程度比較平穩,可以確定林地的相干性受INDVI數值的變化影響較大,特別是在INDVI數值為0.35～0.66之間,之后,隨著INDVI數值增加,相干性將以較平穩的趨勢下降。耕地的相干性在INDVI數值為0.36～0.45區間下降趨勢較小,在0.27～0.36和0.45～0.81的區間內下降趨勢程度最大,由此可知,耕地的相干性受INDVI數值變化影響明顯;而草地的相干性受INDVI數值的增加而降低,且以較相同的趨勢下降,由此可見,草地的相干性變化趨勢受INDVI數值的影響較小。

進一步對比4種地表覆蓋類型的相干性變化趨勢,可以發現,在INDVI數值為0.18～0.55之間,不同地表覆蓋的相干性變化趨勢不相同,建筑區域及工礦用地的變化趨勢最大,其次是林地,而耕地與草地的變化趨勢相似;在INDVI數值為0.55～0.65區間內,草地的變化趨勢低于林地與耕地;在INDVI數值為0.65～0.75區間內,草地的變化趨勢高于林地與耕地;INDVI數值大于0.75,三種植被類型的相干性變化趨勢一致,由此可知,當INDVI數值大于0.75之后,相干性不再受地表覆蓋類型的影響,而只受植被覆蓋率影響,且隨著INDVI數值增加而降低。

3.3 基于各月不同地表覆蓋對相干性的影響分析

為研究不同月份下,地表覆蓋類型對相干性的影響,本文提取2020年月度平均INDVI數據和平均相干性數據繪制圖6,圖中的數值刻度為相干性。

圖6 不同地表覆蓋的月季變化Fig.6 The change of Rose in different surface covers

從圖6中可以發現,在各月份中,建筑區域及工礦用地的相干性最高,且不低于0.7,其次是耕地和草地,林地的相干性最低,其主要原因是建筑區域及工礦用地的INDVI最低,林地的INDVI最高。從整體上看,該圖的整體趨勢靠近1月與12月,偏離6月與7月。從整體的相干性變化趨勢上可以發現,在11月至2月,建筑區域及工礦用地的相干性變化趨勢較為平緩,其次是林地,草地與耕地的變化趨勢最大;而在5月至8月,建筑區域及工礦用地、林地和草地的變化穩定,耕地的變化較大。結合不同地表覆蓋類型的植被覆蓋度可以發現,在11月至2月建筑區域及工礦用地的整體植被覆蓋度低,林地的整體覆蓋度高,因此,整體變化趨勢較小,而耕地與草地的植被覆蓋度高于建筑區,低于林地,因此,在該階段內,植被覆蓋率增加較快,相干性變化大;在5月至8月,建筑區、林地與草地的植被生長情況穩定,耕地在該階段屬于豐收期,因此,相干性變化較大。由此可知,植被覆蓋度對相干性影響較大,且對于不同地表覆蓋,對相干性的影響不同。

為進一步分析不同地表覆蓋類型的相干性在各季節中情況,按照節氣的劃分規則,將12個月份劃分成春、夏、秋、冬四個季節,提取對應的相干性求取平均值,并加入研究區內的降雨量數據繪制出圖7。

圖7 不同地表覆蓋的四季相干性變化Fig.7 Seasonal coherence change of different surface covers

圖7中(a)圖為耕地,(b)圖為林地,(c)圖為草地,(d)圖為建筑區域及工礦用地。由圖7可知,建筑區域及工礦用地的相干性變化趨勢較小,耕地、林地和草地的相干性變化趨勢較大。從INDVI數值的季節分布情況,可以發現秋季的INDVI數值最高,其主要原因是研究區以亞熱帶高原季風氣候為主,植被以常綠闊葉林為主[26]。對比不同地表覆蓋的四季相干性分布情況,可以發現不同地表覆蓋的相干性在冬季最高,其次是春季,秋季和夏季最低。其主要原因是冬季的植被處于落葉期,植被覆蓋率低,而且降雨量少,因此,冬季相干性高;春季的植被處于生長期,植被覆蓋率和降雨量均大于冬季的情況,因此,冬季的相干性大于春季的相干性。由此可以發現,相干性受降雨量和植被覆蓋率共同影響。對比夏季和秋季的情況可以發現,夏季的植被覆蓋率低于秋季的植被覆蓋率,而夏季的降雨量高于秋季的降雨量,從圖中可以發現夏季的相干性略低于秋季的相干性。為進一步了解植被覆蓋率和降雨量對相干性的影響,本文使用灰色關聯度進行對比分析,如表1所示。

表1 不同地表覆蓋的INDVI與降雨量之間的灰色關聯度分析Tab.1 Grey correlation analysis between INDVI and rainfall under different surface covers

由表1可以發現,不同地表覆蓋的相干性受INDVI影響較大,降雨量較小,對比4種不同地表覆蓋的降雨量關聯度可知,4種不同地表覆蓋的相干性受降雨量的影響一致,數值在0.4～0.412之間;INDVI對相干性的影響相差較大,林地受INDVI影響較大,其余類型的地表覆蓋類型受INDVI的影響大致一致,數值在0.87左右。

4 結論

相干性是InSAR技術廣泛應用于形變監測與三維重建的基礎,相干性的好壞決定著InSAR技術最終結果的精度。本研究以覆蓋昆明市滇池周邊的Sentinel 1A數據為例,基于不同指數研究地表覆蓋與相干性之間的關聯程度,將效果最好的指數為依據,進一步分析不同地表覆蓋對相干性的影響及在不同月份、季節和降雨量下對相干性的影響,得出以下結論:

1)INDVI指數能夠顯著體現出不同地表覆蓋與相干性之間的關系。隨著INDVI增加,不同地表覆蓋的相干性以不同的趨勢下降。

2)在不同月份下,所有地表覆蓋類型的相干性在1月與12月最高,6月、7月和8月的相干性最低。且在一年的時間內,建筑區域及工礦用地的相干性最高,其次是耕地與草地,林地的相干性最低。

3)在四季的變化中,建筑區域及工礦用地的相干性變化趨勢較平緩,植被覆蓋區域的相干性變化趨勢較大,且春、冬兩季的相干性較高,夏、秋兩季的相干性低。

4)通過灰色關聯度的分析可知,INDVI對不同地表覆蓋的相干性影響較大,降雨量對相干性的影響較小。相比于其他地表覆蓋,INDVI對林地的相干性影響最大。