基于PS-InSAR的曹妃甸填海區(qū)域形變監(jiān)測(cè)

李 星，陳星彤

(華北理工大學(xué)，河北 唐山 063000)

從2008年起，曹妃甸填海區(qū)為了響應(yīng)‘一港雙城’的政策，開始了大規(guī)模的加速發(fā)展。尤其是華北理工大學(xué)等一系列高校的搬入組成了曹妃甸大學(xué)城，京唐首鋼等各種鋼廠、鐵廠的搬入所引起的大規(guī)模的人流擾動(dòng)、基礎(chǔ)設(shè)施的建設(shè)，勢(shì)必會(huì)造成地表的擾動(dòng)。為了更好地建設(shè)曹妃甸填海區(qū)以及保護(hù)人民財(cái)產(chǎn)安全，對(duì)曹妃甸填海區(qū)進(jìn)行地表形變監(jiān)測(cè)是極其有必要的。

以往的水準(zhǔn)測(cè)量等監(jiān)測(cè)地表形變受環(huán)境天氣影響，而合成孔徑雷達(dá)(SAR)獲取數(shù)據(jù)有全天時(shí)、全天候的特點(diǎn)。因此學(xué)者們提出了用合成孔徑雷達(dá)干涉測(cè)量(InSAR)技術(shù)、合成孔徑雷達(dá)差分干涉測(cè)量技術(shù)(D-InSAR)監(jiān)測(cè)地表形變。但該技術(shù)受時(shí)空失相干和大氣延遲問題的影響較為嚴(yán)重，學(xué)者們又提出了永久散射體差分干涉技術(shù)(PS-InSAR)。Kiseleva對(duì)Caucasus的黑海海岸進(jìn)行地表沉降監(jiān)測(cè)，得到了很好的效果[1]。Gokhan對(duì)Istanbul地區(qū)采用PS-InSAR技術(shù)進(jìn)行地面運(yùn)動(dòng)監(jiān)測(cè)，并對(duì)沉降原因進(jìn)行分析[2]。裴媛媛等、劉鵬等、Peng H等采用永久散射體時(shí)序分析(PS-InSAR)技術(shù)分別對(duì)上海市臨港新城主城區(qū)、黃河三角洲東營地區(qū)、杭州灣南部沿海地區(qū)這三個(gè)地區(qū)進(jìn)行地表形變監(jiān)測(cè)，反演了各個(gè)地表形變過程，分析其沉降原因[3-5]。

本文詳細(xì)介紹了 PS-InSAR技術(shù)用于形變監(jiān)測(cè)的技術(shù)原理，在此基礎(chǔ)上對(duì)20幅覆蓋唐山市曹妃甸區(qū)的Sentinel_1A影像進(jìn)行了沉降監(jiān)測(cè)實(shí)驗(yàn)，通過取得的沉降規(guī)律，加以分析其沉降原因。

1 PS-InSAR的原理與關(guān)鍵技術(shù)

1.1 PS-InSAR的原理及技術(shù)流程

PS-InSAR技術(shù)是從地區(qū)一致時(shí)間不同的多幅SAR影像中提取出散射特性穩(wěn)定的永久散射體(persistent scatterers，PS)點(diǎn)，如人工建筑體、裸露的巖石等硬物體，用提取出的PS點(diǎn)構(gòu)建網(wǎng)絡(luò)，并提取其相位信息，從而建立關(guān)于相位的函數(shù)模型，采用合理的算法，精確地反演地表形變。

在K幅SAR影像中，首先對(duì)主影像進(jìn)行選取，將其與剩余影像配準(zhǔn)及干涉處理，生成K-1幅干涉條紋圖；去除用外部參考DEM數(shù)據(jù)模擬的高程相位，生成去平后的干涉圖；在差分干涉處理后得到的平均強(qiáng)度圖像上提取PS點(diǎn)，對(duì)PS點(diǎn)進(jìn)行線性形變相位和DEM誤差的估計(jì)，從差分干涉相位中將其去除；通過加入在時(shí)間域和空間域的高通、低通濾波去除獲得的殘余相位信息，包括非線性形變相位、大氣相位、噪聲相位；最后解纏形變相位，對(duì)解纏結(jié)果時(shí)序分析，得到可靠的地表形變信息。

1.2 關(guān)鍵技術(shù)

(1)公共主影像的選取

影響干涉相位的精度以及準(zhǔn)確性的因素有很多，時(shí)間基線、空間基線和Doppler質(zhì)心頻率基線是最直接的，而這三者會(huì)隨著主影像的不同而變化。因此，主影像的選取是PS-InSAR技術(shù)中的一個(gè)關(guān)鍵問題。多個(gè)影響因素綜合考慮，讓其成為最佳組合，即綜合相關(guān)系數(shù)γm最大：

(1)

(2)PS點(diǎn)的選取

選取永久散射體(PS)點(diǎn)也是PS-InSAR技術(shù)中的關(guān)鍵步驟。該技術(shù)通過分析處理這些高相干點(diǎn)目標(biāo)上的干涉相位，從而得到時(shí)間序列上的某地區(qū)形變信息。在處理過程中，建模的條件是影像上要有足夠數(shù)量和密度的PS點(diǎn)。選取PS點(diǎn)的方法主要有相位離差閾值法、相干系數(shù)閾值法和振幅離差閾值法的單種閾值方法以及兩者結(jié)合或三者結(jié)合的復(fù)合閾值方法。PS點(diǎn)選取的方法中最為經(jīng)典的是振幅離差閾值法。該方法的基本原理是，信噪比高的點(diǎn)像元的振幅離差指數(shù)可以和其相位的噪聲水平等價(jià)衡量。對(duì)于標(biāo)準(zhǔn)差為δ的高斯噪聲的SAR影像而言，其振幅信息和 Rice 分布相似，若點(diǎn)具有高斯信噪比的特征，則該點(diǎn)的分布與高斯分布相近，即振幅離差指數(shù)DA：

(2)

式中，δA為平均強(qiáng)度圖像的振幅標(biāo)準(zhǔn)差，mA為平均強(qiáng)度圖像的振幅平均值[7]。

(3)干涉相位模型及形變相位反演

PS-InSAR技術(shù)中最主要、最關(guān)鍵的一步就是相位模型的建立以及形變相位的反演。每個(gè)提取的PS點(diǎn)的解纏干涉相位φunw：

φunw=φtopo+φdef+φatm+φnoise

(3)

式中，φtopo為DEM高程誤差引起的地形相位；φdef為沿LOS方向的地表形變相位；φatm為大氣影響相位；φnoise為噪聲相位。

在時(shí)間域與空間域上DEM高程誤差引起的地形相位φtopo、大氣影響相位φatm、噪聲相位φnoise有不同的表現(xiàn)形式，可通過采用高通濾波和低通濾波的方式來去除，最后，形變相位φdef：

(4)

式中，φlin是線性的形變相位；φnon-lin是非線性的形變相位，V是線性的形變速率；T是時(shí)間間隔；Dnon-lin是非線性形變。

由此LOS方向線性形變速率得以獲取[8]。

2 曹妃甸填海區(qū)形變監(jiān)測(cè)與分析

2.1 研究區(qū)概況及數(shù)據(jù)選取

唐山市曹妃甸區(qū)位于唐山南部、環(huán)渤海中心地帶，地理坐標(biāo)為北緯39°07′～39°27′，東經(jīng)118°12′～118°43′。曹妃甸區(qū)的地質(zhì)構(gòu)造為填海造陸，土質(zhì)松軟，近幾年曹妃甸區(qū)進(jìn)行了大規(guī)模建設(shè)，引進(jìn)工業(yè)、商業(yè)和學(xué)校。而地表形變會(huì)對(duì)各種基礎(chǔ)設(shè)施的安全造成威脅，這是曹妃甸填海區(qū)地表變形監(jiān)測(cè)的原因。

本文選取20景2018年8月～2019年4月，覆蓋曹妃甸地區(qū)的Sentinel_1A影像，其分辨率為5 m×20 m。對(duì)曹妃甸填海區(qū)進(jìn)行裁剪，裁剪區(qū)域如圖1黑色線框所示，外部的DEM數(shù)據(jù)采用的是分辨率為30 m的STRM數(shù)據(jù)。

圖1 研究區(qū)域

2.2 形變信息提取與分析

通過公共主影像的選取方法，最終確定2018年11月30日的圖像為公共主影像，并與其他輔影像進(jìn)行配準(zhǔn)、干涉，而干涉數(shù)據(jù)的時(shí)間基線與空間基線的連接圖如圖2所示。從圖中可以看出，空間基線的范圍是在-76.246 1 ～ 96.068 9 m，時(shí)間基線范圍是在-120 ～ 131 d。最后干涉處理完成后得到了19對(duì)干涉條紋圖，如圖3所示。圖中大部分干涉條紋圖的條紋不是很明顯，對(duì)后面PS點(diǎn)的提取有極大的影響。去除干涉圖中利用外部DEM數(shù)據(jù)模擬的高程相位后，得到的去平后的差分干涉圖效果更好，如圖4所示。

圖2 數(shù)據(jù)的連接圖

在生成差分干涉圖之后，開始選取PS點(diǎn)，由于其地形條件，將振幅離差指數(shù)閾值DA設(shè)置為0.8，PS點(diǎn)的密度為25個(gè)/km2，總共選取了1 107 106個(gè)PS點(diǎn)，結(jié)果如圖5所示。PS點(diǎn)大都落在人工建筑，如房屋、道路等硬面地物上，而不在草地、養(yǎng)殖池等散射特性較弱的地方。

圖3 干涉條紋圖

圖4 去平后的差分干涉圖

圖5 永久散射體(PS)分布

PS點(diǎn)選取后，對(duì)其進(jìn)行線性形變相位和DEM誤差的估計(jì)，從差分干涉相位中去除得到殘余相位。由于大氣引起的相位延遲和噪聲屬于高頻信號(hào)，而非線性形變相位屬于低頻信號(hào)，則可以分別采用時(shí)間域的高通濾波和空間域的低通濾波去除。2019年2月22日的影像能更直觀的表達(dá)濾波前后的效果，如圖6所示，可以看出濾波后的影像相較于濾波前降低了相干噪聲，有效地去除了一些由大氣引起的噪聲，效果明顯。

圖6 差分干涉圖(a為濾波前，b為濾波后)

由此，對(duì)濾波后的影像進(jìn)行線性形變相位解纏，得到真實(shí)的相位信息以及各PS點(diǎn)在LOS方向的形變速率信息，如圖7所示。

圖7 曹妃甸填海區(qū)的年平均沉降速率圖

由圖7可以看出，整體的平均沉降速率約為-26.56～22.45 mm/a，可見從2018年8月起該地區(qū)仍然存在著較為嚴(yán)重的沉降。紅圈范圍內(nèi)的A區(qū)和B區(qū)是沉降嚴(yán)重的區(qū)域，C區(qū)為抬升區(qū)。A區(qū)位于曹妃甸區(qū)填海區(qū)，緊靠海堤，該區(qū)域的沉降速率達(dá)到了-23 mm/a以上，沉降尤為顯著。B區(qū)位于曹妃甸區(qū)的柳贊鎮(zhèn)附近，沉降速率最大為-26.532 5 mm/a，是整個(gè)研究區(qū)域沉降最為嚴(yán)重的地區(qū)。C區(qū)為抬升區(qū)最明顯的地區(qū)，最大抬升速率為22.440 8 mm/a。

為了更好的探討其沉降原因，首先對(duì)大部分沉降文獻(xiàn)進(jìn)行了研究，發(fā)現(xiàn)除了地震等大型自然災(zāi)害的直接因素外，地表沉降的主要原因是地下水的開采和地上人工建筑物建造這兩個(gè)間接因素。其次，查找該區(qū)域的水文地質(zhì)資料，該區(qū)域的含水層形成原因是河流沖積和海積，屬第四系松散巖類，其空隙水一般分為Ⅰ、Ⅱ、Ⅲ組，第Ⅰ、第Ⅱ含水組屬潛水及微承壓水，一般不作為開采層。第Ⅲ含水組為深層承壓水，底界深度一般為160～500 m，水位埋深約為l5 m。20世紀(jì)80年代由于超采深層地下水形成以南堡開發(fā)區(qū)和唐?？h城為中心的兩個(gè)沉降漏斗，受這兩個(gè)沉降漏斗的影響，該區(qū)地下水位呈下降趨勢(shì)。這也是該區(qū)地面整體下降的主要原因。就A區(qū)而言，其位于曹妃甸工業(yè)區(qū)，我們對(duì)該區(qū)域進(jìn)行了實(shí)地考察，發(fā)現(xiàn)其附近有原油碼頭和鋼廠(首鋼京唐等)。2018年首鋼京唐共生產(chǎn)815萬t，每天至少需要近10萬t淡水(首鋼京唐相關(guān)報(bào)道)。由此可見，該區(qū)域地表沉降的原因主要是地下水的開采。B區(qū)位于柳贊鎮(zhèn)—一個(gè)以養(yǎng)殖海產(chǎn)品為主要經(jīng)濟(jì)產(chǎn)業(yè)的地區(qū)，其附近有許多養(yǎng)殖池，如蝦池、蟹池。根據(jù)資料顯示，2019年3月柳贊鎮(zhèn)南部5萬多畝露天海水養(yǎng)殖區(qū)域已全數(shù)建成并進(jìn)入海水養(yǎng)殖準(zhǔn)備階段。由此可見，B區(qū)沉降嚴(yán)重的主要原因?yàn)轲B(yǎng)殖池的建造與開挖。C區(qū)位于唐山灣生態(tài)城，曹妃甸東站附近。生態(tài)城復(fù)合功能用地規(guī)劃用地面積387.90 ha，占城市建設(shè)用地總面積的20.28%左右[9]，曹妃甸東站于2017年開始建設(shè)，2019年運(yùn)營。唐山灣生態(tài)城即將被打造為臨海宜居城市，建設(shè)面積在不斷地增加，可見養(yǎng)殖池的平填與地上建設(shè)用地的建造為C區(qū)地表提升的主要原因。

為了更加直觀的了解PS-InSAR技術(shù)所處理的沉降結(jié)果，分別在沉降尤為嚴(yán)重的A區(qū)、B區(qū)選定了兩個(gè)PS點(diǎn)——Dot1(38°56′10.88′′N、118°31′44.12′′E)、Dot2(39°13′39.32′′N、118°37′25.82′′E)，其位置在圖7中標(biāo)記顯示。這兩個(gè)點(diǎn)的沉降結(jié)果如圖8所示，在時(shí)間區(qū)域內(nèi)，兩者的沉降過程以非線性的方式沉降，且其沉降坡度逐漸增大，由此可見，A、B兩區(qū)域的沉降極為嚴(yán)重，也間接印證了圖7沉降速率的準(zhǔn)確性。Dot1的最大累積沉降量為21.856 mm，Dot2的最大累積沉降量為15.625 mm。從這張折線圖可以看出，在沉降過程中，雖然總體是呈現(xiàn)沉降趨勢(shì)，但是某一時(shí)間區(qū)域內(nèi)會(huì)有些許抬升，這主要是與該區(qū)域的地質(zhì)條件、地下水開采及地上建設(shè)有極大關(guān)系。

圖8 A、B區(qū)域中Dot1、Dot2的累積沉降量

3 結(jié) 語

本文采用PS-InSAR技術(shù)對(duì)曹妃甸填海區(qū)進(jìn)行地表形變監(jiān)測(cè)，得到了很好的監(jiān)測(cè)效果。從結(jié)果看出該地區(qū)的確有著極大的地表變形，最大年平均沉降率達(dá)到了-26.558 76 mm/a，沉降區(qū)域最為嚴(yán)重的地方在鋼廠附近，年平均沉降率達(dá)到了-23 mm/a。該區(qū)域的沉降與本身填海造陸且松軟的地質(zhì)條件有一定的關(guān)聯(lián)，同時(shí)，工業(yè)區(qū)建設(shè)運(yùn)營、長期的地下水開采、道路的開通運(yùn)營會(huì)對(duì)這個(gè)本來就是填海造陸的地區(qū)產(chǎn)生了一定的負(fù)擔(dān)。本文也證實(shí)了采用PS-InSAR技術(shù)進(jìn)行地表沉降監(jiān)測(cè)的有效性，為曹妃甸的開發(fā)提供了參考。