結合SAR與光學遙感的疑似稀土礦非法開采情況識別方法

樂穎，夏元平，劉媛媛，錢文龍

(1.東華理工大學 測繪工程學院，南昌 330013；2.福州市勘測院，福州 350108)

0 引言

作為我國獨特的一種稀土礦，離子吸附型稀土礦具有儲量多、分布區域廣、種類豐富、用途廣等特點[1]。根據省第三環境保護督察組的反饋可知，贛州市證外稀土礦山的總數大約是證內稀土礦總數的2.5倍，分別占據901座與355座，需要治理的礦山數量達到了522座，同時廢棄礦山的數量占據了很大比例。因此，針對稀土礦山的違法開采現象，如何實時準確監測開采是目前的研究熱點。

合成孔徑雷達差分干涉測量(differential interferometric synthetic aperture radar，DInSAR)是地表形變監測、災害監測、變化監測等諸多應用領域重要信息的獲取手段[2-4]。國內外部分學者對雷達變化監測進行了相關研究[5-6]，Kamila等[7]運用DInSAR和時序短基線集差分干涉測量(smallbaseline subset-insar，SBAS-InSAR)技術反演研究區域完整的形變趨勢，并用水準數據對結果進行了評價；Xia等[8]提出了一種基于概率原理的方法，利用積分法和合成孔徑雷達干涉測量(interferometric synthetic aperture radar，InSAR)技術，成功反演了地下采空區的位置；馬威等[9]對比分析了基于InSAR相干系數圖與DInSAR對地表信息的提取效果，獲得連續變化開采區39處；何曙光等[10]提出了結合InSAR相干系數的差值方法對城市進行變化監測，并證明了該方法提取實驗區變化信息具有較高的準確性。同時，國內很多學者通過影像分類或建立明確的解譯標志，從而研究分析礦產資源開采變化情況[11-14]。邵艷坡等[15]結合紋理與光譜信息進行分類，結果表明,該分類方法用于定南縣稀土礦的總體精度為90.33%。岳建偉等[16]疊加分析遙感與矢量數據，能夠自動提取礦區違法開采信息。然而，高分辨率影像一般為光學影像，既受氣候環境的限制，又受觀測時間的影響，在夜間無法監測。雷達數據具有全天時全天候的優勢，不僅提供了反映地物特性的幅度信息和相位信息，還可以獲取稀土礦開采引起的地表形變信息。高分辨率影像和雷達數據具有較高的互補性，將二者結合起來進行變化監測具有廣泛的前景。

因此，本文結合雷達影像的強度信息和光學影像的光譜信息進行變化監測，并借助已有礦權信息在光學影像上對疑似非法開采區域進行圈定，最后根據地表形變信息和相干性系數圖動態監測疑似非法開采區域在2015—2020年間的開采變化情況。該研究揭示了在稀土礦山開采監測中將雷達影像和光學影像相結合進行非法開采識別的能力，同時能夠給礦山管理部門治理礦山周邊環境提供決策依據。

1 研究方法

1.1 基于相干性系數的 DInSAR 地表形變信息提取方法

利用干涉相干性分析可以提取地表覆蓋的變化信息，借助二次差分方法從干涉相位圖中去除地形及其他因素的影響，從而達到提取形變信息的目的，該方法被稱為DInSAR技術[17-18]。

在DInSAR處理過程中考慮相干系數可以提高結果的準確性[19]，相干系數γ能夠直接反映出影像之間的相干性指標[20]。

1.2 基于單波段比值法的遙感變化監測

圖像比值法是遙感圖像處理過程中最常用的方法之一，就是將不同波段或波段組合所對應的像元灰度值進行比值運算，是對兩景影像做除法的變化監測方法。其基本思想是：通過比值運算計算待檢測影像中每個像元在不同時相亮度值的比值，然后根據結果繪制出比值圖像。檢測原理如式(1)所示[21]。

(1)

式中：Ratio表示所需的影像變化監測結果；pwr1和pwr2分別是兩個時相影像。若待影像中地類在兩個時相未發生變化，其灰度值往往大致相同；若地物類別發生了改變，其灰度值也隨之發生變化。因此，通過分析比值影像上地類變化區域的灰度值與背景值之間的差異，然后對比值結果進行一定的處理，最終可獲得研究區域內地物發生變化的信息。

1.3 結合光學和雷達影像的稀土礦非法開采識別方法

根據以上理論，本文提出一種結合SAR與光學遙感的疑似稀土礦非法開采情況識別方法，具體流程如圖1所示。根據雷達影像的強度信息和光學影像的光譜信息，將贛州市定南縣作為研究區域，對離子型吸附型稀土礦的開采情況進行變化監測。本文變化監測結果可通過虛警率(false alarm，FA)、總體錯誤率(overall error，TE)、總體準確率(overall accuracy，OA)和Kappa系數進行精度評價[22]。具體操作步驟如下。

步驟1：分別對SAR和光學遙感影像進行數據預處理。雷達影像經過影像配準、濾波、地理編碼和輻射定標可得到強度圖，再利用基于相干性系數的 DInSAR 地表形變信息提取方法得到相干性系數圖和形變干涉圖。光學影像需經過大氣校正、重采樣、波段合成、影像鑲嵌和裁剪等數據預處理，并對其進行真彩色圖像增強。

步驟2：對同一年份的兩景影像，根據雷達影像的強度信息和光學影像的光譜信息，分別進行基于單波段比值法的遙感變化監測，并對差異圖進行真彩色合成。

步驟3：依據變化標準判斷出可能變化區域。若在兩個變化監測結果中都未變化，該區域沒變化；當在雷達變化監測結果圖上變化而光學變化監測結果影像上不變，則該區域沒變化；只有當雷達變化監測結果圖變化，光學變化監測結果影像也變化才能說明該區域發生了變化。

步驟4：識別出開采區域。通過上一步的變化監測結果可以得到研究區域內的變化區域，結合礦區典型地物(沉淀池、稀土礦采場、礦山建筑和引水管道等)和形變干涉圖判斷是否為稀土礦開采變化區域，借助已有礦權信息圈定疑似非法開采區域。

步驟5：結合形變圖和相干系數圖對比分析圈定區域的稀土礦非法開采情況，探究相干性系數和地物之間的關系。

2 研究區域與數據處理

2.1 研究區域概況

贛州市地處江西省南部，占據整個江西省總面積的23.6%，是中國重要的“稀土王國”[23]。定南縣作為贛州市內離子型稀土主產縣之一，位于114°47′49″E～115°22′48″E，24°33′37″N～25°03′21″N，地處南嶺與武夷成礦帶的中間部位，擁有得天獨厚的成礦條件。定南縣全縣總面積達到1 321 km2，礦產資源非常豐富，其中包含了木子山、甲子背、長坑尾等11個稀土礦區。近幾年以來，定南縣稀土礦存在違法開采現象，累積了大量的礦山環境問題。因此，為了幫助當地實現可持續發展和恢復生態功能，精確、可靠地識別非法開采區域顯得尤為重要。

2.2 數據處理

光學遙感影像數據本文使用了分辨率為10 m的哨兵2號數據，可在歐空局官網(https：//scihub.copernicus.eu/)免費下載。哨兵2號是高分辨率多光譜成像衛星，從可見光和近紅外到短波紅外，具有不同的空間分辨率，用于陸地監測，可提供植被、土壤和水覆蓋、內陸水路及海岸區域等圖像，還可用于緊急救援服務[24]。由于數據質量和資源的限制，選取了2017—2020年10月份的8景影像。首先需要對影像進行預處理，可根據不同地物在光學影像上表現的顏色、大小、形狀等建立解譯標志。

由于雷達數據與光學影像不同，屬于主動遙感，影像包含的是振幅、相位和極化等多種信息。本文使用的雷達影像數據為C波段Sentinel-1A衛星的SLC數據，分辨率為30 m。為了探尋地表形變與稀土礦區開采之間的關系，選擇了時間跨度為2015—2020年的影像數據，具體的雷達影像數據基本參數如表1所示。雷達數據經過預處理后可得到強度圖，1)雷達影像變化監測。一般的變化監測方法只考慮處理兩幅SAR圖像，忽略了這些SAR圖像的相位信息。作為一個相干成像系統，SAR具有更多的潛力，其相位信息應該被挖掘出來[25]。應用比值法將每一年的主影像與副影像做除法運算，可得到2015—2020年間的比值結果。圖2是以疑似非法開采區域C為例的比值結果，圖中的深紅色區域為可能減少的信息，深綠色區域為可能增加的信息。由于地物在兩景影像之間有所變化，從而造成該區域在兩個時相的SAR影像的強度有所差異，同時由于地物的散射單元結構發生了變化，使得其兩個時相的SAR影像相干性較弱，因此，導致了該區域的比值法結果圖和相干性系數圖發生變化。

表1 雷達影像數據的基本參數

其包含了地面分辨元的雷達后向散射強度信息，然后對強度圖進行比值運算，即將后一時相的強度圖與前一時相的強度圖進行比值。最后，對雷達數據進行DInSAR處理獲得相干性系數圖和形變圖。相干性系數圖反映了兩景影像之間的相干性程度，形變干涉圖可反映出地表形變信息，為后續分析提供參考。

3 實驗結果分析

3.1 變化監測結果與分析

開采行為會造成地物類型發生改變，因此這部分區域在SAR影像上的灰度值會產生較為明顯的變化。由于雷達影像的幾何畸變、斑點噪聲同樣會在變化監測結果中顯示出來，前者是由于雷達成像時角度變化造成的，后者是由一個分辨單元內各個散射點相干回波之間的干涉效應引起的，但是這兩類變化并不屬于真實的稀土礦區開采變化，對于變化結果的準確性有很大影響。本文考慮到這兩個因素在光學遙感影像上并不會發生變化，故采用雷達影像與光學遙感影像相結合的方法，對比分析那些雷達成像時幾何畸變和斑點噪聲引起的變化，最終剔除這些因素的影響。

2)光學影像變化監測。對于光學影像的動態變化監測，大部分研究都需要進行地物分類，但是中低分辨率遙感影像的分類精度并不高。由于同一地物的光譜特征相同，不同地物的光譜特征有所差異，因此，可以根據比值法對光學遙感影像進行變化監測，疑似非法開采區域C在2017—2020年間的光學遙感影像如圖3所示。將變化監測結果與兩景影像進行真彩色合成能夠更好地分析變化區域，其中圖3中的藍色區域為可能減少的信息，黃色區域為可能增加的信息。在圖3中可以看到，疑似非法開采區域C內具有沉淀池這一稀土礦山的典型地物，且在每一時期的影像上都發生了不同程度的變化，則可以證明該區域確實進行了開采活動。圖4為疑似非法開采區域C在2017—2018年、2018—2019年、2019—2020年和2017—2020年的變化情況。

3)變化監測結果分析。在不受外界因素干擾的情況下，選取合適的遙感數據類型和時相，可以發現一些非法開采行為。由于雷達影像和光學影像自身的成像特點，對變化監測都會產生一定的影響，導致出現誤檢的情況。因此，可以充分利用兩種數據的互補性進行變化監測，從而彌補光學影像受天氣影響導致的局限性和改善雷達影像在幾何畸變和斑點噪聲區域的檢測結果，開展雷達數據和光學影像在稀土礦山監測的應用研究。當僅在雷達變化監測和光學變化監測其中一個結果中出現變化而在另一個變化監測結果中表現未變化時，則判定該區域沒變化；只有當雷達變化監測結果圖變化，光學變化監測結果影像也變化才能說明該區域發生了變化，變化監測精度評價如表2所示。相較于雷達影像變化監測和光學影像變化監測，本文方法的準確率分別提高了11.6%和4.46%，虛警率分別降低了38.24%和23.53%。結果表明，和單一影像相比，本文提出的結合雷達和光學影像的變化監測方法結果更準確，利用遙感技術在研究區開展礦業活動開采識別動態監測是行之有效的。

表2 變化監測精度評價

要判斷變化區域是否是由稀土礦開采引起的，可以根據稀土礦山開采的典型地物建立解譯標志。開采的典型地物中比較容易判讀的地物為沉淀池和浸礦池，也就是用來沉淀稀土的池子，主要集中在礦區旁邊，形狀為規則的圓形或長方形。通過對雷達影像和光學影像變化監測結果進行分析，并結合稀土礦山解譯標志和礦權范圍進行目視解譯，圈定了9個疑似非法開采區域，疑似違法圖斑圈定結果如圖5所示，并以疑似非法開采區域C為例說明本文方法的有效性。圖6依次為圖5圈定的9個疑似非法開采區域在2017—2020年間變化監測結果真彩色合成圖，變化監測結果中紅色、綠色和藍色分別代表疑似非法開采區域在2017、2018和2019年間的變化情況。

3.2 干涉相干性分析

相干圖表示的是兩幅影像的相關性，是最常用的相位質量圖，也是最直觀的質量評價圖。將光學遙感影像與雷達相干性系數圖進行對比可以發現，二者具有較高的一致性。在這幾個干涉對中，失相干區域多為植被覆蓋區域，相干性高區域多為具有穩定地物的區域。因此，不同地物在相干性系數圖上表現出不同的系數值，且相干系數高低可排序為：建筑物>裸地>植被。在研究區域的相干系數圖中，由于城市區域地物目標和裸地具有較高的穩定性，因此表現出相干性非常高；與此相比，林地、農田等地物覆蓋區域，相干性較低，色調以藍綠色為主；在水體等低散射回波地物區域，則表現為黑色區域。

本文通過對比疑似非法開采區域于2015—2020年間相干系數圖的變化，來探究開采區域與相干性之間的關系。本文研究區域整體相干性系數范圍在0.37～0.71之間，具體如圖7所示。通過該圖可以看出，這幾處疑似非法開采區域相干性系數變化趨勢具有較高的一致性，進一步證實了這些區域進行了開采活動。

3.3 地表形變干涉結果與分析

經過信號采集與數據處理，SAR影像的每一像素既包含地面分辨元的雷達后向散射強度信息，也包含與斜距有關的相位信息，將覆蓋同一地區的兩幅SAR影像對應像素的相位值進行差分，便可得到一個一次差分相位圖，也就是干涉相位圖。通過對12景覆蓋研究區的SAR影像數據進行DInSAR處理，借助相干系數閾值獲取相干性高的像元，并將其作為軌道精煉的控制點引入后續的處理中，最終可得到贛州市定南縣在2015—2020年的地表形變信息。本文通過對比9個疑似非法開采區域在形變圖中的變化，圈定了9個疑似非法開采區域于2015—2020年間的地表形變變化信息(圖8)。通過該圖可以看出，這幾處疑似非法開采區域地表形變趨勢具有較高的一致性，稀土礦的開采能夠引起緩慢的地表形變，其變化范圍為-25～15 mm。

4 結束語

本文在稀土礦非法開采識別過程中，結合SAR和光學遙感影像的特性，給出了一種可更有效識別疑似非法開采區域的方法。該方法綜合利用了Sentinel系列影像的相干性系數、強度信息、形變信息和光譜信息等因素，對雷達和光學影像分別進行變化監測確定變化區域，并依據礦區典型地物、變化標準和已有礦權信息圈定疑似非法開采區域。相較于僅用單一影像進行礦區非法開采識別，結合SAR和光學遙感的疑似非法開采情況識別效果更優。

本文以贛州市定南縣的稀土礦區為例進行實驗，得到以下結論。

1)結合SAR和光學遙感的稀土礦疑似非法開采情況識別方法充分利用兩種數據的互補性，不僅彌補了光學影像受天氣影響的局限性，也改善了雷達影像的幾何畸變和斑點噪聲區域的檢測結果。

2)相較于雷達影像變化監測和光學影像變化監測，本文方法的準確率分別提高了11.6%和4.46%，虛警率分別降低了38.24%和23.53%，可更準確識別出疑似稀土礦開采區域。

3)疑似稀土礦開采區域在相干系數圖和地表形變圖中的變化趨勢具有較高一致性，相干性系數和地表形變范圍分別在0.37～0.71和-25～15 mm之間。

本文揭示了雷達影像和光學影像相結合在稀土礦山開采監測中的疑似非法開采情況識別能力，同時能夠給礦政管理部門治理礦山周邊環境、制定礦產資源分配和確定礦產資源開發秩序等提供決策依據和技術支持。然而，本文采用的影像數據分辨率有限，未來可在本文基礎上使用高分辨率影像進行下一步的研究。