唐山馬蘭莊鐵礦露天開采邊坡變形監測的GB-InSAR技術

劉作利 劉景玉 申修強 馬 輝 田衛明 左榮虎

(1.唐山首鋼馬蘭莊鐵礦有限責任公司；2.北京理工大學信息與電子學院；3.北京理工雷科電子信息技術有限公司)

我國有露天礦山約7萬座，尾礦庫1. 2萬座，排土場10萬余座，該類邊坡一旦失穩，不僅會影響企業生產，更會威脅作業人員的生命安全[1]。因此，通過對邊坡變形進行實時監測，并對監測數據進行一系列處理分析，對于準確掌握邊坡災害發生、發展和變化機理，進而實現邊坡災害預警具有重要意義[2-3]。全站儀、測量機器人和水準儀等傳統測量技術在變形監測中已經得到了廣泛應用，但精度易受氣候、時間、通視等條件的限制[4]。

合成孔徑雷達(Synthetic aperture radar, SAR)是一種基于微波傳感器的雷達，具有全天候和一定穿透性等優點。差分干涉合成孔徑雷達技術(Differential interferometric synthetic aperture radar，D-InSAR)是SAR的一個重要應用，在近十余年內得到了快速發展和廣泛應用。星載和地基干涉合成孔徑雷達技術是干涉合成孔徑雷達技術(Interferometric synthetic aperture radar，InSAR)的2種重要形式。星載干涉合成孔徑雷達具有獲取變形范圍大的優勢，可應用于大區域地表沉降監測[5-9]，但由于星載干涉合成孔徑雷達具有較長重返周期、固定成像姿態和低空間分辨率等不足，導致其在邊坡變形監測中效果不理想。地基干涉合成孔徑雷達測量技術(GB-InSAR)是星載干涉合成孔徑雷達干涉測量技術的有效補充，不僅具有最優觀測姿態和連續觀測性能，而且具有靈活多變、分辨率高、平臺穩定、觀測周期短、造價相對低廉等優點，屬于非接觸測量方法范疇，適用于對危險邊坡進行實時變形監測[10-12]。本研究以唐山馬蘭莊鐵礦為例，針對該礦邊坡特點，采用GB-InSAR技術進行邊坡變形監測，結合監測數據處理分析結果，進一步評價利用該技術進行邊坡變形監測的可行性。

1 GB-InSAR技術原理

1.1 SAR技術原理

SAR是利用一個小的真實天線的運動來等效構成一個長天線，從而可以在運動方向上獲得一個等效的大孔徑(合成孔徑)天線，利用該方式，雷達對目標的分辨率將會有所提高。根據合成孔徑成像技術原理，設雷達系統發射的信號帶寬為B，電磁波在真空中的傳播速度為c，則雷達斜距向分辨率為c/2B[13]。當信號帶寬為300 MHz時，系統能達到的最大距離向分辨率為0.5 m。若雷達發射的電磁波波長為λ，雷達滑動的最大距離為L，則雷達方位角分辨率為λ/2L。

本研究采用的GB-InSAR雷達系統信號處于Ku波段(12.5～18 GHz)，可以獲得的最高方位向分辨率約為4 mrad。通過將SAR與調頻連續波技術(Stepped frequency continuous wave,FMCW)相結合，監測區域可以被分割為較多二維小單元(圖1)。

圖1 GB-InSAR監測空間分辨單元

GB-InSAR雷達系統在獲取目標區域數據過程中，由于天線運動軌道固定，對某一區域目標重復監測成像時天線都處于相同的空間位置，因而空間基線為0。雷達信號往返路徑之差即為雷達視線方向目標距離的變化量，具體可以表現為相位差[13]。

1.2 地基雷達差分干涉變形測量原理

利用D-InSAR技術可以對SAR圖像進行干涉處理并去除干涉圖中的非形變相位，實現對地表進行大范圍、高精度、全天候變形監測。由于地基雷達主機部分始終在線性導軌上滑動，滑軌的起點和終點位置始終未發生變動，SAR圖像數據的軌道是嚴格重合的，即地基雷達并不存在空間基線。如此，地基雷達和星載雷達的干涉測量方式有較大區別，地基雷達差分干涉處理時，無需去平地效應和去地形效應，因此從理論上分析，差分干涉SAR圖像的相位信息僅包含變形信息。GB-InSAR可通過測量視線方向上回波信號的相位變化實現高精度變形測量(圖2)，視線方向的形變量Δd與干涉相位φdef關系可以表示為

(1)

式中，λ為雷達信號的波長。

圖2 差分干涉雷達變形監測原理

由干涉相位解纏理論[9]可知，干涉圖中相鄰2個永久散射體(Permanent scatterer,PS)點位之間的相位差不大于1/4λ，即為邊坡雷達變形監測的最大形變分辨能力。因此，地基差分干涉SAR圖像的變形測量分辨能力與雷達波長密切相關。如本研究采用的Ku波段雷達波長約為1.8 cm，則相鄰2個形變點位的最大形變分辨能力為4.5 mm。根據差分干涉雷達的基本理論，目標點的干涉相位可以表示為

φInt=φdef+Δφatmosphere+Δφtrack+Δφnoise，

(2)

式中，φdef為地表變形相位；Δφatmosphere為大氣相位；Δφtrack為由定軌誤差引入的相位；Δφnoise為由時間去相關等因素引入的隨機相位。

當式(2)中的后3項誤差較大時將無法獲得高精度的地表變形信息。為此，在本研究邊坡雷達系統中采用了永久性散射體技術(PS)。該方法首先在SAR圖像中根據幅度閾值、幅度離差閾值、相關性閾值等指標優選出高相關性的目標點(PS點)；然后對Δφatmosphere、Δφtrack、Δφnoise進行建模并估計補償大氣誤差、軌道誤差和噪聲誤差相位；最后獲取PS點的形變信息。數據解算過程的重要環節為相位解纏[9]，特別是在GB-InSAR礦區邊坡變形監測數據處理中，直接影響了邊坡變形監測精度。進行相位解纏時首先選擇1組永久散射體候選子集，通過構建三角網進行插值處理，而后進行空間維和時間維解纏，最終通過處理分析永久散射體間的位移信息來消除大氣、軌道等噪聲誤差，進而獲取高精度的變形監測數據。

2 GB-InSAR系統組成與工作流程

2.1 系統組成

GB-InSAR系統的主體部分由高分辨調頻連續波(FMCW)雷達、高精度電控位移臺組成(圖3)。整個系統包括主控計算機設備、遠程數據傳輸系統設備等。調頻連續波雷達通過發送線性調頻的電磁波信號可以實現快速成像，最快觀測一次耗時小于2 min。高精度電控位移臺可保障雷達在運行過程中平穩可靠，為雷達合成孔徑計算過程的必備外界條件。實時監測軟件包括二維和三維2個模式，可以有效地從雷達觀測結果中識別出相應的地形位置信息。雷達系統的主控計算機除了安裝實時監測軟件以外，還可以構建數據庫及發布網站，可將雷達觀測數據以網站的方式進行實時發布。當雷達與監控室之間無法進行有線網絡通訊時，數據遠程傳輸設備可以通過無線方式實時發送雷達觀測數據。

圖3 邊坡監測雷達系統主體部分外觀

2.2 系統工作流程

設備工作時，邊坡雷達實時處理軟件可以實時顯示雷達工作狀態，可以實現數據實時處理、處理結果實時存儲與顯示、歷史數據實時回放等功能。雷達觀測數據通過網絡傳輸至服務器，再利用服務器中構建的數據庫對大量的雷達變形監測數據進行統一管理，最終觀測結果可以通過多個客戶端進行實時顯示。本研究雷達系統的變形監測數據可與三維地形信息進行疊加顯示，對變形數據可以進行打點曲線分析和實時預警，此外，用戶還可以在線查詢特定區域在各時期的變形信息。GB-InSAR系統的工作流程如圖4所示。

該系統在形變測量的基礎上，設計了振動觀測環節，通過對目標回波進行計算，從而可以計算出強散射體的振幅和頻率值。

3 工程應用

3.1 邊坡概況

馬蘭莊鐵礦資源儲量為8 912.15萬t，SFe平均品位為29.64%，目前露天采礦規模達到 300萬t/a，選礦處理原礦能力達180萬t/a，精礦產量為80萬t/a。該礦位于河北省遷安市馬蘭莊鎮，西距北京約200 km，西南距唐山約80 km。礦區中心的地理坐標為東經118°36'、北緯40°06'。該礦露天礦區呈橢圓狀，整個采場長半軸約1 100 m，短半軸約900 m，采場邊坡為典型巖質邊坡，無植被覆蓋，邊坡高度落差大，最大開采深度超過200 m；邊坡角為38°～47°，在如此陡峭的邊坡上安裝監測設備難度較大；此外，采場存在霧霾、粉塵，能見度不高。隨著該礦露天開采深度進一步增加，邊坡穩定性成為制約采場安全生產的重要因素。

圖4 邊坡雷達系統工作流程

3.2 邊坡變形監測方案

根據該礦露天采場的實際條件，結合被監測邊坡危險區域的分布特點，將邊坡雷達安置于礦坑東坡和西坡地基較穩定處。為使邊坡變形監測工作不受雨雪等不利天氣條件影響，滿足惡劣氣象條件下不間斷監測要求，在礦坑東側邊坡監測位置建立了活動房(圖5)，將GB-InSAR安置于該活動房中進行邊坡變形觀測，可以有效避免惡劣天氣對儀器的影響(圖6)。將邊坡監測系統服務器設置于礦區機房，監控人員通過客戶端可以完成對雷達工作的實時控制以及對監測結果的分析(圖7)。GB-InSAR采用24 h不間斷監測方式，測量頻率為2～ 10 min/次，采用實時處理方式實現對整個邊坡表面的連續變形進行實時監測。

圖5 雷達活動房

3.3 監測數據分析

礦坑西坡監測時間為2016年9月12日—2016年12月20日。礦坑東邊坡監測時間為2016年12月10日—2016年12月20日，監測期分為2段：第1段監測期為2016年12月10日—2016年12月16日；第2段監測期為2016年12月17日—2016年12月20日。礦坑東邊坡2次監測的區域完全一致，雷達觀測距離為300～900 m，觀測范圍恰好覆蓋整個露天采場的東邊坡。GB-InSAR通過主動發射調頻連續波信號并接收成像，而后對成像結果進行PS點選擇、相位解纏及誤差改正，單軌變形結果計算的最短耗時小于2 min，從而實時獲取邊坡變形信息。礦坑西坡的監測期超過3個月，相對于東坡，盡管監測期較長，但未監測出明顯的變形信息，表明在監測期內該區域的邊坡基本處于穩定狀態。礦坑東坡的監測期僅有11 d，雖然監測期較短，但2個監測時段均監測出礦坑東坡在同一區域內出現形變，該區域為90 m×40 m(長半軸×短半軸)橢圓形，標高由高至低，變形量逐漸增大(圖8)。

圖6 安裝于活動房內的邊坡雷達

圖7 監控室遠程邊坡雷達控制

分析圖8可知：第1個監測段內的最大變形量達到-150 mm(符號為負，表現為接近雷達方向的視線向形變)；第2個監測段內雷達重新啟動，變形結果重新清零，邊坡最大變形量約-120 mm，也表現為接近雷達方向的視線向形變。2次觀測結果均顯示出邊坡變形區下部變形速度快、變形量大的特點，2次監測期積累的雷達視線向最大變形量達到-270 mm。

為進一步分析邊坡變形特點，針對第2觀測期內邊坡變形較大的1個點位(圖9)進行了打點曲線分析，如圖9所示。分析圖9可知：該點位變形速率較穩定，觀測期內的最大變形量約-120 mm。為更精確地計算該點位的變形速率，利用邊坡雷達實時監測軟件，在綜合分析邊坡變形特征及監測誤差的基礎上，將礦坑東坡第2監測期內的累計變形量與第2個監測期的總監測時間相除，結果如圖10所示。分析圖10可知：變形速率所顯示的位置與變形區域完全一致，局部區域的最大形變速率超過-1 mm/h，同時可以看到變形區域底部的變形量大于上部。

圖8 開采東坡形變時間序列

圖9 點位變形曲線和變形速率曲線

圖10 觀測期內平均形變速率

依據本研究DB-InSAR邊坡變形監測結果，對邊坡變形區進行了現場考察，發現該變形區頂部的邊坡巖體產生了明顯裂縫，區域上部出現的裂縫顯示局部位置變形量超過180 mm(圖11)，并且下部變形量大于上部。總體上，雷達監測結果與現場實際情況基本一致。礦山工程技術人員根據本研究邊坡變形監測結果，及時對變形區域進行了邊坡處理，有效避免了邊坡災害發生。

圖11 礦坑東坡現場裂縫

4 結 語

針對唐山馬蘭莊鐵礦露天采場邊坡特征，為確保邊坡穩定，采用GB-InSAR技術進行了變形監測，在采坑東部邊坡發現了一橢圓形變形區域，區域最大累計變形量達到-270 mm，并且該區域底部變形量大于上部，與現場實際情況基本吻合。礦山工程技術人員根據本研究監測結果，及時采取了邊坡處理措施，有效避免了邊坡災害發生。研究表明，采用GB-InSAR技術可以對露天采場邊坡進行高精度實時變形監測，在一定程度上可以實現邊坡災害預警。

[1] 佴 磊，徐 燕，代樹林，等.邊坡工程[M].北京：科學出版社，2010.

[2] 鄭穎人，陳祖煜，王恭先，等.邊坡與滑坡工程治理[M].2版.

北京：人民交通出版社，2010.

[3] 王宏偉，王 攀，王永強，等.馬鋼凹山露天采場南幫邊坡監測及穩定性分析[J].金屬礦山，2016(3)：35-37.

[4] 劉 超，高井祥，王 堅，等.GPS/偽衛星技術在露天礦邊坡監測中的應用[J].煤炭學報，2010，35(5)：755-759．

[5] 張 潔，胡光道，羅寧波.InSAR技術在滑坡監測中的應用研究[J].工程地球物理學報，2004(2)：147-152.

[6] BERARDINO P,COSTANTINI M,FRANCESCHETTI G,et al.Use of differential SAR interferometry in monitoring and modelling large slope instability at Maratea (Basilicata, Italy) [J].Engineering Geology, 2003,68(1-2):31-51.

[7] 王桂杰，謝謨文，邱 騁，等.D-InSAR技術在大范圍滑坡監測中的應用[J].巖土力學，2010，31(4)：1337-1344.

[8] 魏海霞，高照忠，葉長斌.融合D-InSAR 技術與ArcGIS 軟件的礦區開采沉陷監測[J].金屬礦山，2016(11)：128-131.

[9] 汪 磊，鄧喀中，薛繼群，等.融合概率積分模型與D-InSAR 的開采沉陷預計[J].金屬礦山，2016(2)：160-163.

[10] 楊紅磊，彭軍還，崔洪曜.GB-InSAR監測大型露天礦邊坡形變[J].地球物理學進展，2012，27(4)：1804-1811.

[11] 李如仁，楊 震，余 博.GB-InSAR集成GIS的露天煤礦邊坡變形監測[J].測繪通報，2017(5)：26-30.

[12] 王志勇，孟慶穎. InSAR相位解纏算法在礦區沉降監測中的對比與分析[J].遙感信息，2014，29(6)：47-51.

[13] 趙東寅，申其鴻，馬海濤，等.國產地基合成孔徑雷達監測預警系統在紫金山金銅礦露天采場邊坡位移監測的應用[J].中國安全生產科學技術，2015(4)：54-58.