采空區無人機影像地表裂縫活動特征研究

訾應昆　張文君　耿詩畫　郭峻麟　徐夢輝

摘要：地表裂縫是一種典型地質災害。以龍首礦西二采區為研究區，利用2020年8月和2021年6月兩期沉降區的無人機正射影像，通過圖像增強和多尺度分割處理得到用于分類的多邊形對象，結合地表裂縫的光譜特征和形狀特征構建模糊規則集，實現對地表裂縫的提??；選取覆蓋該區域2020年8月19日-2021年6月27日52景Sentinel-1A升降軌影像，采用小基線集（SBAS）技術，獲取地表二維形變信息，并采用經驗貝葉斯克里金法（EBK）和徑向基函數法（RBF）對形變數據進行空間加密；將獲取的裂縫分布信息和加密后的地表形變信息疊加，結合采空區、塌陷坑等從垂直向和水平向上揭示研究期間礦區地表裂縫的形變特征和分布變化規律。研究表明：（1）不同裂縫群距離采空區越近，垂直向和水平向形變速率越大，裂縫群整體向采空區側偏移；（2）采空區周圍的地表裂縫群主要受垂直向不均勻形變的影響，采空區上方的地表裂縫群主要受水平不均勻形變影響；（3）受周圍塌陷坑的影響，地表裂縫群兩側會產生形變速率差異，且靠近塌陷坑的一側形變較大，導致新裂縫的產生和裂縫群的規模擴大。

關鍵詞：采空區 地表裂縫 地表形變 裂縫活動 無人機影像

中圖分類號：TP751;TD325.2文獻標志碼：A文章編號：1671-8755（2023）04-0063-10

Research on the Characteristics of Surface Crack Activity of UAV Image in Goaf

ZI Yingkun1，2， ZHANG Wenjun1，2，3， GENG Shihua1，2， GUO Junlin1，2， XU Menghui1，2

Abstract：? Surface crack is a typical geological hazard. Taking the West 2 mining area of Longshou mine as the study area， the UAV orthophotos of the subsidence area in two phases in August 2020 and June 2021 were used to obtain polygonal objects for classification by image enhancement and multi-scale segmentation processing， and fuzzy rule sets were constructed by combining the spectral features and shape features of surface cracks to realize the extraction of surface cracks; Sentinel-1A elevation track images covering the area from August 19， 2020 to June 27， 2021 were selected， the small baseline set （SBAS） technique was used to obtain two-dimensional surface deformation information， and the deformation data were spatially encrypted by using the empirical Bayesian kriging （EBK） and radial basis function （RBF） methods; the acquired crack distribution information and the encrypted surface deformation information were superimposed to reveal the deformation characteristics and distribution of surface cracks in the study area from the vertical and horizontal directions by combining the mining areas and collapse pits. The study shows that： （1） the closer the different crack groups are to the mining area， the larger the vertical and horizontal deformation rates are， and the crack groups are shifted to the side of the mining area as a whole; （2） the surface crack groups around the goaf area are mainly affected by the vertical uneven deformation， and the surface crack groups above the goaf area are mainly affected by the horizontal uneven deformation; （3） influenced by the surrounding collapse pits， the two sides of the surface crack group will have different deformation rates， and the deformation on the side close to the collapse pits is larger， leading to the generation of new cracks and the expansion of the fracture groups.

Keywords：? Goaf; Surface cracks; Extraction zones; Surface deformation; Fracture activity

礦區地表裂縫是礦山在開采過程中引發的一種典型地質災害，是采區巖層和地表移動變形的結果［1-2］。礦區地表裂縫的產生會破壞地表結構，裂縫的發育會引發次生地質災害，對當地的基礎設施和生活造成嚴重的安全隱患［3-5］。地表裂縫的研究很大程度上受制于影像的分辨率。傳統的實地野外勘查及量取的裂縫數據精度雖高但耗時耗力，甚至受地表形變的影響不能持續監測。近年來，無人機科技飛速發展，其高效率、高分辨率和不受復雜地形影響的特點，為動態監測地表裂縫創造了條件。國內對裂縫影像數據的提取做了大量研究。魏長婧等［6］利用無人機影像和TM影像，根據裂縫的特征使用ERDAS軟件建立提取模型，較好地提取出地表裂縫；肖瑤等［7］采用面向對象自動提取技術，結合目視解譯提取了研究區地面塌陷和地表裂縫信息，精度較高；楊豪杰等［8］根據地表裂縫的影像特征構建了地表裂縫的提取規則集，精確提取出地表裂縫；湯伏全等［9］針對影像初提取的特征集采用隨機森林算法制作掩膜實現采動裂縫的精確提??；楊奇讓等［10］以神南礦區地表的無人機影像為數據源建立地表裂縫的提取特征集，構建K最近鄰法、支持向量機、樸素貝葉斯、隨機森林4種監督學習分類模型對影像進行分類，探討不同方法對地表裂縫的提取效果。

已有的裂縫研究主要集中在裂縫的識別和提取，對于地下開采引起的地表裂縫和地表裂縫變化的研究鮮有涉及。礦區地表裂縫所處環境更為復雜，裂縫的形成與發育受到開采擾動的影響［11-13］。本文以甘肅省金川銅鎳礦的龍首礦西二采區為研究區，以無人機影像、InSAR和GNSS數據為數據源，將研究時段內提取的裂縫變化信息與采空區形變信息進行疊加，實現對礦山地表裂縫的動態監測。

1 研究區概況及數據采集

1.1 研究區概況

研究區位于甘肅省金昌市金川銅鎳礦的龍首礦區（如圖1所示），該區域是典型的干旱半干旱大陸性氣候，區域地表為厚風積沙覆蓋，缺乏地表水系，植被類別稀少，覆蓋度較低，生態環境比較脆弱。地表受地下長時間開采的影響，使該區域地形有明顯的起伏，存在大量裂縫群以及塌陷坑。

1.2 數據采集

采集實驗數據選用大疆精靈4 RTK旋翼無人機，搭載的是大疆FC6310R相機。將地面分辨率設置為2 cm，無人機的相對航高為73 m。依據規范要求，將航向重疊度和旁向重疊度均設置為80%，為保證影像清晰度，將飛行速度控制在10 m/s 。無人機主要參數如下：最大水平飛行速度13.89 m/s；最大上升速度 6 m/s；工作環境溫度 0～40 ℃；最大承受風速10 m/s；衛星定位系統為單頻高靈敏度 GNSS/多頻多系統高精度 RTK GNSS；障礙物感知范圍 0.2～7 m。 傳感器參數如下：傳感器尺寸13.2 mm×8.8 mm；分辨率5472×3648像素；像元大小2.41 μm；焦距8.8 mm。采集時間分別為2020年8月19日和2021年6月27日兩期，采用三維建模軟件Context Captue進行影像的匹配和拼接，獲取研究區地表高分辨率影像。

2 面向對象的地表裂縫提取

面向對象的信息提取技術將具有相同屬性的多個相鄰像元作為目標進行整體提取，打破了傳統的僅限于像元提取的局限性［14-15］。由于研究區域采空區裂縫比較多，為了高效快速提取地表裂縫信息，對無人機影像進行增強處理，并針對裂縫的光譜特征、空間結構等特征，采用隸屬度函數模糊分類法構建地表裂縫的模糊規則集，提高裂縫的提取精度。

2.1 改進的非銳化掩蔽進行圖像增強

為消除無人機采集數據中的噪聲，需對獲取的無人機影像數據進行圖像增強處理。圖像增強流程如圖2所示，采用改進的非銳化掩蔽進行圖像增強處理。首先對無人機原始影像進行Laplace算子銳化，將結果與GHP濾波器處理Sobel梯度影像結果相乘得到鈍化模板，把鈍化模板與原影像相加得到增強影像。經過非銳化掩蔽處理后的無人機影像減弱了影像的噪聲，增加了地表裂縫的對比度。最后使用Gamma變換增強無人機影像的灰度值。如圖3 所示，利用普通非銳化掩蔽與改進的非銳化掩蔽均可增加無人機影像的對比度，但經過改進的非銳化掩蔽的增強處理裂縫邊緣更加清晰。

2.2 多尺度影像分割

多尺度分割算法［16］基于每個對象異質性的值，通過遍歷對象將小于設定同質性尺度閾值的對象進行合并 ，當異質性值大于設定的同質性尺度閾值時則停止合并，影像的最終結果受分割尺度、波段權重和均值因子的影響［17］。在ECognition軟件中，默認各波段權重為1，通過多次調試參數，依次將分割尺度、形狀因子、緊致度因子設置為40，0.2和0.6，得到的圖像分割效果如圖 4所示。

2.3 無人機影像地表裂縫的特征分析

（1）光譜特征。根據不同地物在遙感影像上呈現出光譜信息的差異，結合野外實地考察和目視解譯，將研究區域地表地物劃分為裸地、沙土、植被和地表裂縫。典型地物類型解譯標志如表1所示。

（2）形狀特征。在部分區域，其他地物與地表裂縫光譜特征相近，難以通過單一的光譜特征將裂縫準確提取出來。通過對無人機影像裂縫進行統計分析，區內裂縫在形狀上表現為直線狀、折線狀和分叉狀，且面積較小、寬度較窄，具有良好的線性特征。沙土、植被和塌陷坑多呈現為不規則的面狀，面積相對較大，具有良好的面狀特征。因此，本文在ECognition軟件中采用不對稱性、密度、圓度、橢圓擬合和形狀指數5個參數對區內裂縫進行統計分析，得到地表裂縫最佳形狀特征閾值。

2.4 構建裂縫提取的模糊規則集

通過對影像的增強處理和多尺度分割，得到用于分類的多邊形對象，并對地物對象進行特征選取，包含光譜特征和形狀特征。考慮到無人機影像中存在混合像元，會存在異物同譜或同物異譜現象，在ECognition軟件中采用隸屬度函數分類器，構建研究區內裂縫提取的模糊規則集（表2）。

2.5 提取精度分析

本研究對礦區兩期無人機影像進行處理，提取結果如圖 5所示。從圖5可以看出，裂縫的提取效果相對較好，但也存在部分顏色較深的植被和土路邊界錯分為裂縫。

為了驗證面向對象自動提取的精度，結合長度統計和位置統計信息等（表3），以提取條數、條數提取率、自動提取總長度、長度提取率和位置重疊率作為評價指標，對裂縫的提取精度進行評價。

從表3可以看出，自動提取條數和自動提取總長度略低于實際值，主要原因是研究區內信息較為復雜，不同地物間存在相似特征，一些細小的裂縫在影像分割時被包含在其他地物中。

3 獲取采空區二維形變場

3.1 獲取采空區二維形變信息

基于無人機影像時間跨度，選取研究區2020年8月19日-2021年6月27日時段內Sentinel-1A共52景升降軌影像，利用SBAS-InSAR技術進行形變分析。采用基于ENVI遙感影像處理平臺的SARscape雷達圖像處理軟件進行數據處理。在數據預處理階段，將Sentinel-1A數據POD精密定軌星歷數據導入SARscape中進行TOPSAR模式處理，并將結果根據DEM地理坐標系裁剪文件。

研究過程中使用25景升軌影像和27景降軌影像，通過設置臨界基線閾值為45%和時間基線閾值為120 d，得到升軌187個干涉對和降軌215個干涉對，生成3D解纏連接圖。在干涉處理部分，多視處理采用距離向視數和方位向視數的比值為4∶1，干涉圖去平采用Goldstrin濾波方法，相位解纏采用最小費用流法，閾值設置為 0.2。為消除還存在的相位誤差，進行軌道精煉和重去平，以干涉處理較好的升降軌影像的干涉圖為地圖，根據地面控制點的選點要求，最終選取參與升降軌處理的地面控制點有44個和48個。在SBAS反演處理過程中，第一次反演采用線性模型獲得較好效果，第二次反演利用Goldstein濾波方法，在上一步基礎上進行大氣濾波得到研究區的形變速率。最后，利用DEM數據將兩次反演得到的結果轉換到常用的WGS-84地理坐標系下便于后續研究分析。根據二維形變的分解原理，將向東偏移定義為正方向，獲取研究區內地表二維形變信息。

由于InSAR獲取的地面形變數據空間分辨率較低，需對解算的二維形變信息進行空間插值以分析地表裂縫的形變特征。采用經驗貝葉斯克里金法（EBK）和徑向基函數法（RBF）對地面形變數據進行空間插值，然后對插值結果進行對比分析。

3.2 精度分析

為驗證時空插值后的地表二維形變場的精度，采用離散的GNSS監測點數據作為特征點，與對應點二維形變信息對比分析。GNSS數據由礦山地測科提供。根據研究區范圍選取的有效監測點34個，位置如圖6所示，其中點號為A-1，A-3，A-5，A-8，A-9，A-13，B-1，B-2，B-3，B-4的10個點按裂縫群走勢分布，其余24個點在采空區上方大致呈5行6列均勻分布。

精度評價基于已知點的形變量和預測的形變量的平均絕對誤差（MAE）和均方根誤差（RMSE）兩個指標評估兩種插值法的適用性。其中，MAE可以反映插值誤差的大小，而RMSE能反映出插值結果的偏差。兩種插值法精度結果如表4所示。

從表4可以看出，經驗貝葉斯克里金法在垂直向和水平向的MAE和RMSE兩個指標值均明顯低于徑向基函數法在垂直向和水平向的MAE和RMSE兩個指標值。因此，利用經驗貝葉斯克里金插值更有利于本實驗數據加密。

4 采空區地表裂縫活動監測結果與分析

通過兩期裂縫提取結果的對比，選取出新增的主裂縫和寬度變化較大的裂縫作為代表性特征裂縫，將提取結果與加密后的形變結果進行疊加分析，得到研究區地表裂縫的分布如圖 7所示。

為分析地表裂縫活動性，根據地表裂縫群的分布與走向，將其分為裂縫群Ga，Gb和Gc 3個區域。

4.1 地表裂縫群Ga活動性分析

對于地表裂縫群Ga ，選取平行于該裂縫群的兩側剖線Ha1（裂縫群外側）和Ha2（裂縫群內側）以及垂直于該裂縫群的兩條剖線Pa1和Pa2 分析裂縫群活動與地表形變規律的關系。裂縫群Ga形變剖線圖如圖8所示。

平行于裂縫群Ga的剖線Ha1和剖線Ha2以靠右側為起始點，在垂直向上，結合圖9（a），研究區所在區域整體呈下沉趨勢，隨著距離的增加，裂縫群垂直形變速率逐漸增大，且越靠近采空區下沉越嚴重。剖線Ha1 和剖線Ha2 在采空區的右側下沉速率明顯不同，兩者最大沉降率相差為27 mm/a，在采空區上方形變速率基本一致。在水平方向上，結合圖9（b），研究區整體上呈向東偏移的趨勢，隨著距離的增加，裂縫群偏移速率逐漸增加。兩剖線在采空區的左側水平形變速率基本一致，但在采空區上方展現出明顯的差異。垂直于裂縫群Ga的剖線 Pa1和Pa2，以靠上側為起始點，在垂直方向上，結合圖9（c），隨著距離的增加，裂縫群的垂直形變速率緩慢增加，越靠近采空區形變越嚴重，與剖線形變圖9（a）結果一致。在水平方向上，結合圖9（d），Pa2 的形變速率高于Pa1 的形變速率，裂縫越靠近采空區形變速率越快。

綜上，研究區裂縫群Ga在越靠近采空區的位置垂直向和水平向形變速率越大，且裂縫群兩側存在明顯差異。結合形變剖線位置圖可知，采空區右側的地表裂縫群Ga主要受垂直向不均勻形變的影響，采空區上方的地表裂縫群Ga主要受水平不均勻形變影響。采空區上方和右方均有新的較寬裂縫產生，且原裂縫也在逐漸變寬，說明垂直向與水平向的不均勻形變均可導致地表裂縫群的規模擴大。

4.2 地表裂縫群Gb活動性分析

裂縫群Gb在采空區的西側，受附近新生成的塌陷坑影響，產生較多較細的裂縫。對于地表裂縫群Gb，平行于該裂縫群的兩側剖線Hb1（裂縫群外側）和Hb2（裂縫群內側）以及垂直于該裂縫群的剖線Pb的空間位置如圖10所示。

剖線Hb1 和剖線Hb2 以靠上側為起始點，在垂直向上，結合圖 11 （a），裂縫群Gb所在區域整體呈下沉趨勢，隨著距離的增加，裂縫群垂直形變速率逐漸增大，由于Hb2 在采空區邊界和沉降漏斗邊界處，使其沉降速率較Hb1 的沉降速率大，裂縫群兩側最大相差超過23 mm/a。在水平方向上，結合圖 11（b），Hb1 和Hb2 整體上呈向東偏移，且偏移速率呈先增加后減少的趨勢。剖線Hb1 形變速率相對剖線Hb2 來說變化較大，由圖可知地表裂縫Gb發生在偏移盆地的邊界，該裂縫群兩側最大偏移速率相差超過45 mm/a，裂縫群兩側存在明顯差異性偏移。垂直于裂縫群Gb的剖線Pb，以靠左側為起始點，在垂直方向上，結合圖11（c），隨著距離的增加，Pb所在區域呈整體下沉趨勢，剖線Pb的形變速率呈先增加后減少再增加的特點，并在裂縫群Gb范圍內出現陡降型峰值，裂縫群兩側的沉降速率變化顯著。在水平方向上，結合圖11（d），Pb所在區域整體向東偏移，偏移速率越來越小，裂縫群兩側最大偏移量較小。

綜上，新生地表裂縫群Gb處在垂直向形變與水平向形變盆地的邊界處，是地表裂縫易發區，且塌陷坑的形成對地表裂縫的發育有較大影響。

4.3 地表裂縫群Gc活動性分析

研究區裂縫群Gc處于采空區上方，且周圍伴隨多個塌陷坑，研究期間裂縫群的寬度明顯增大，裂縫兩側落差變化較為明顯。對于地表裂縫群Gc，分別沿其走向和垂直方向上作剖線進行分析。平行于地表裂縫群Gc的兩側剖線Hc1（裂縫群外側）和Hc2（裂縫群內側）以及垂直于該裂縫的兩條剖線Pc1和Pc2 空間位置如圖12所示。由圖12可知剖線Hc1和Hc2 與3個塌陷坑在空間分布上也呈平行分布。

平行于裂縫群Gc的剖線Hc1和剖線Hc2以靠左側為起始點，在垂直向上，結合圖13（a），Hc1與Hc2所在區域呈整體下沉趨勢，隨著距離的增加，沉降速率逐漸減小，在采空區上方由于塌陷坑的影響，兩剖線在平行于3個塌陷坑附近相差較大，最大下沉速率相差超過21 mm/a，且在遠離采空區后兩剖線的垂直形變速率迅速減小，通過剖線的形變速率圖可知周圍塌陷坑影響著地表裂縫的活動性。在水平方向上，結合圖13（b），剖線Hc1與Hc2 所在的區域在采空區上方存在向東偏移，隨著距離的增加，在采空區右側存在向西的偏移。地表裂縫群Gc的兩側剖線Hc1和Hc2 的水平偏移速率相差出現先減小后增大再減小再增大的過程，說明裂縫群Gc受到塌陷坑的影響兩側出現不均勻形變。垂直于裂縫群Gc的剖線Pc1和Pc2，以靠下側為起始點，在垂直方向上，結合圖13（c），剖線Pc1 與Pc2 所在的區域呈下沉趨勢，隨著距離的增加，裂縫群的垂直形變速率緩慢減少，越靠近采空區形變速率越大，但剖線Pc1比剖線Pc2 形變速率大。在水平方向上，結合圖13（d）可以看出剖線Pc1 所在的區域整體水平向東偏移，在經過塌陷坑Ⅱ時，水平形變速率顯著增加，存在明顯的波動。剖線Pc2 經過塌陷坑 Ⅲ 時，水平形變速率也有增加。

綜上，塌陷坑的形成影響著裂縫群的活動，且裂縫群Gc在塌陷坑兩側的形變相差相對較大，靠近塌陷坑一側的形變明顯高于另一側，與實地測量一致。

5 結論

以金川龍首礦區的西二采區為研究區，采用無人機影像提取的地表裂縫變化信息和同時段InSAR獲取的地面形變數據，在垂直和水平向形變上，結合采空區和塌陷坑對不同位置的裂縫群進行綜合分析，得到如下結論：（1）通過地表裂縫的光譜特征和形狀特征構建模糊規則集，實現對兩期無人機影像地表裂縫的提取，兩期影像的裂縫提取率分別為82.2% 和80.9%，提取效果較好。（2）采用經驗貝葉斯克里金法和徑向基函數法對利用InSAR獲取的形變數據進行空間插值，通過與實測點的對比分析，利用經驗貝葉斯克里金插值更有利于本實驗形變數據加密。（3）將特征裂縫的提取結果與地表二維形變的結果進行疊加，通過在水平向和垂直向兩方面對研究區不同分布位置的裂縫群綜合分析表明地表裂縫群受到采空區和塌陷坑的影響，在其兩側產生不均勻形變，導致了新裂縫的產生和原有裂縫變寬。在采空區上方的裂縫群主要受到水平向形變的影響，采空區附近裂縫群主要受垂直向形變影響。本研究結果對分析地表裂縫的活動性具有一定的參考價值。

參考文獻

［1］ 趙坤陽. 煤礦地下開采誘發地表裂縫與導水裂縫分布規律預測［D］. 北京： 中國地質大學（北京）， 2015.

［2］ 陳超， 胡振琪. 我國采動地裂縫形成機理研究進展［J］. 煤炭學報， 2018， 43（3）： 810-823.

［3］ 范立民， 馬雄德， 李永紅， 等. 西部高強度采煤區礦山地質災害現狀與防控技術［J］. 煤炭學報， 2017， 42（2）： 276-285.

［4］ 閆瑞龍， 屠世浩， 王方田， 等. 高頭窯礦區地表采動裂縫演化規律研究［J］. 煤炭工程， 2016， 48（1）： 95-98.

［5］ 王云廣， 郭文兵. 采空塌陷區地表裂縫發育規律分析［J］. 中國地質災害與防治學報， 2017， 28（1）： 89-95.

［6］ 魏長婧， 汪云甲， 王堅， 等. 無人機影像提取礦區地裂縫信息技術研究［J］. 金屬礦山， 2012（10）： 90-92， 96.

［7］ 肖瑤. 高分系列衛星在煤礦區地質災害監測方面的應用［D］. 合肥： 合肥工業大學， 2017.

［8］ 楊豪杰. 地裂縫高分辨率影像特征提取和紋理分析［D］. 太原： 太原理工大學， 2020.

［9］ 湯伏全， 李林寬， 李小濤， 等. 基于無人機影像的采動地表裂縫特征研究［J］. 煤炭科學技術， 2020， 48（10）： 130-136.

［10］楊奇讓， 胡振琪， 韓佳政， 等. 煤礦區無人機影像采動地裂縫提取方法研究［J］. 煤炭科學技術， 2023， 51（6）： 187-196.

［11］杜江. 煤礦工作中地表裂縫的成因與防治［J］. 當代化工研究， 2021（18）： 101-102.

［12］王子威，蔡毅，嚴家平，等.淮北孫疃煤礦采煤地裂縫發育特征研究［J/OL］.礦業安全與環保：1-8［2022-10-09］. https：//kns.cnki.net/kcms/detail/50.1062.td.20220928.1719.002.html.

［13］李林寬. 基于無人機影像的煤礦區采動裂縫提?。跠］. 西安： 西安科技大學， 2021.

［14］蔡祥， 李琦， 羅言， 等. 面向對象結合深度學習方法的礦區地物提?。跩］. 國土資源遙感， 2021， 33（1）： 63-71.

［15］李琦. 基于面向對象和深度學習的礦區無人機影像地物提取［D］. 北京： 北京林業大學， 2020.

［16］譚衢霖， 劉正軍， 沈偉. 一種面向對象的遙感影像多尺度分割方法［J］. 北京交通大學學報（自然科學版）， 2007， 31（4）： 111-114， 119.

［17］莊喜陽， 趙書河， 陳誠， 等. 面向對象的遙感影像最優分割尺度監督評價［J］. 國土資源遙感， 2016， 28（4）： 49-58.

收稿日期：2022-12-30；修回日期：2023-05-18

基金項目：四川省測繪地理信息學會科技開放基金（CCX202205）

作者簡介：第一作者，訾應昆，碩士研究生，E-mail：363068503@qq.com; 通信作者，張文君，博士，教授，研究方向為遙感信息技術，E-mail： 113999066@qq.com