李昱昊,安士凱,周大偉,詹少奇,高銀貴
(1.中國礦業大學江蘇省資源環境信息工程重點實驗室,江蘇 徐州 221116;2.平安煤炭開采工程技術研究院有限責任公司煤礦生態環境保護國家工程實驗室,安徽 淮南 232001;3.鄂爾多斯市華興能源有限責任公司,內蒙古鄂爾多斯 010399)
我國煤炭資源賦存現狀是“西部多、中部富、東部區域枯竭”。以晉陜蒙為代表的西部地區,煤炭產量已占到全國總產量的70%以上。隨著大量煤炭的產出,所引起的地表沉陷及災害問題也日益突出,如地表塌陷和塌陷區積水[1-2]。因此,對開采沉陷進行快速、準確、全面的監測,合理地建立下沉盆地模型和反演參數,可以定量研究受開采影響的巖層的時空變化規律。現有的GNSS、水準等大地測量方法野外工作量大,且不足以獲取整個沉陷區的變形特征,因此難以支撐大規模西部礦區煤炭開采沉陷監測。InSAR[3-8]等遙感檢測方法在變形劇烈的區域容易出現是失相干的現象,許多高分辨率影像成本過高。無人機攝影測量技術[9-12]以其自動化、智能化、快速化的特點,通過對地面周期性快速航拍獲取高分辨影像,構建密集點云,從而在礦區沉陷監測中發揮巨大優勢;但在產品生產過程中,點云常常出現噪聲,野點的情況,同時因其精度僅達到分米級,對礦區盆地參數反演產生極大的影響[13]。在沉陷參數預計方面,目前有很多求取參數的研究方法[14],但存在或多或少的問題;正交參數[15]試驗當參數較多時,工作量較大,試驗次數較多;模矢法[16]容易把局部最優解當成全局最優解,對參數初值有一定的要求;遺傳算法[17-18]需要確定種群數量等參數,否則會影響收斂速度與參數精度,且存在局部搜索能力差的缺陷。模擬退火算法[19-21](Simulated Annealing,SA)是由Kirkpatrick Z 提出的一種求解組合優化問題的算法,其優點是:根據Metropolis 接受準則,采用一定的概率接受比當前解差的解,從而跳出局部最優解來搜索全局最優解。
針對西部礦區地形起伏且植被較少的地理環境,本研究利用索尼a7r2 數碼相機搭載垂直起降固定翼無人機開展礦區沉陷監測建模研究。以內蒙古唐家會煤礦工作面地表為實驗區,利用無人機獲取2 期影像數據,制作DEM 并獲取下沉盆地,針對存在噪點問題,對全盆地點云去噪方法進行對比擇優。同時,為了使求參結果更準確,融合概率積分與模擬退火算法,選取全盆地下沉值進行盆地參數反演,最后對沉陷參數進行抗差分析[22-23]。
無人機攝影測量技術流程圖如圖1。
在無人機攝影測量技術應用過程中,為了確保礦區監測工作的順利進行,首要任務是相機標定,設置相控點布設方案以及無人機飛行路線的規劃處理。在具體過程中,需要明確監測區域邊界,設定航向重疊率、旁向重疊率、飛行航高,以此為基礎開展礦區監測工作。
在獲取區域影像數據后,結合無人機飛控系統導出的POS 數據和校正過后得到的相機參數數據,進行影像自動內定向,然后依次進行輸出坐標系、快速檢測、導入像控點、刺相控點、高精度處理,生成高密度點云,最終生成測區數字高程模型。通過多期對該區域的重復監測,獲取不同時期工作面推進到不同位置的DEM,將DEM 相減得到監測區域的地表下沉盆地[24]。
根據概率積分法原理[25-26],若煤層是水平的,單元B(s,t)的開采所引起的地表任意一點A 的下沉We(x,y)為:為煤層走向長為,L 為煤層傾向水平長。
模擬退火算法反演概率積分參數步驟如下:
1)給出概率積分的初始參數,確定初始溫度T,最低溫度T′,降溫系數δ 等退火參數,計算初始解。
2)判斷T<T′,若是,則輸出結果X,若不是,執行步驟3)~步驟5)。
3)若達到迭代長度,則降溫T=Tδ,δ 為介于0~1的實數。
4)生成相鄰方案。假設在某一狀態下參數值為X,則它的相鄰狀態X′為X 在一定鄰域隨機變化產生的值,并且該值必須滿足上下限約束。
式中:W0為煤層的頂板下沉量,W0=mqcosα;m為實際開采厚度;q 為下沉系數;α 為煤層傾角;D
唐家會區域影像圖圖2。
圖2 唐家會區域影像圖Fig.2 Image map of Tangjiahui area
唐家會煤礦位于內蒙古自治區準格爾煤田東孔兌普查區的西南部,其地理坐標為:東經:111°10′27″~111°14′34″,北緯:39°52′45″~39°57′22″。選擇61202工作面作為觀測地點,該工作面尺寸為240 m×1 960 m,平均采深538 m,平均采厚19 m,平均傾角2°,計劃開采時間2020 年7 月,計劃結束時間2021 年7月,在測區東南部有矸石山堆積。
用垂直起降固定翼無人機于2020 年8 月采集第1 期無人機影像數據,2021 年3 月采集第2 期影像數據,監測區域寬約為1 000 m,長約為2 900 m,面積約為3 km2。其中航向重疊率、旁向重疊率均設為80%,航高設置為140 m。在工作面四角、航飛邊緣四角布設通視較好的像控點,其他像控點均勻布設在旁向重疊范圍內,共計29 個。第1 次數據采集共獲得2 079 張有效影像,第2 次數據采集共獲得2 531 張有效影像,通過內業處理獲取DEM,無人機影像處理生成DEM 如圖3。用檢查點驗證DEM 精度,經計算,第1 期DEM 中誤差為0.066 m,第2 期DEM 中誤差為0.059 m,用2 期DEM 相減獲得唐家會地區地表下沉盆地,地面下沉盆地如圖4。
圖3 無人機影像處理生成DEMFig.3 The DEM generated by UAV image processing
圖4 地面下沉盆地Fig.4 Surface subsidence basin
經過濾波后的DEM 仍帶有噪點,將DEM 復原為下沉點云,下沉點云圖如圖5,其中噪點區域用圓圈圈出。
圖5 下沉點云圖Fig.5 Subsidence point cloud
產生噪點的主要原因有:
1)點云濾波算法的局限性。植被點以及低矮建筑物不能完全去除。在研究區域里,由于2 期DEM監測時間在夏冬兩季,植被區域覆蓋不同,疊加后的DEM 偏差較大。
2)插值引起的誤差。無人機穿透能力較低。在研究區域西部方位有落差較大的溝壑,其覆蓋植被較多,無人機攝影測量難以穿透此區域。不含點云數據的區域依賴于DEM 插值。在點云密度較低的區域,其差值DEM 存在較多偏差。
3)無人機精度原因。無人機精度為分米級,在保證整體精度較高的情況下,有少數點高程誤差較高,且水平偏移引起的誤差較為明顯。
這些噪點對后續處理或者預計參數效果帶來了很大的不穩定性,因此需要進行去噪處理。BP 神經網絡是一種按誤差反向傳播(簡稱誤差反傳)訓練的多層前饋網絡,其算法稱為BP 算法,它的基本思想是梯度下降法,利用梯度搜索技術,以期使網絡的實際輸出值和期望輸出值的誤差均方差為最小。針對沉陷點云去噪,選用MATLAB 中BP 神經網絡工具箱對原始下沉點云進行去噪處理,設置訓練次數為1 000,學習速率為0.05,設置顯示頻率為500 次顯示1 次,訓練目標最小誤差為10-5。同時選擇中值去噪、最小二乘去噪、最近鄰去噪與BP 神經網絡去噪對比,其中中值濾波窗口設置為20,最小二乘濾波窗口設置為3,最近鄰點數設置為10,離群點閾值設置為1。去噪對比圖如圖6。
圖6 去噪對比圖Fig.6 Comparison of denoise
經過4 種去噪方法的測試對比發現:最小二乘去噪效果最差,去噪窗口增大后容易出現點云紊亂的現象,其次為中值去噪及最近鄰去噪算法,BP 神經網絡去噪效果最好,同時保留了下沉盆地數據的完整性。
61202 工作面走向達到充分采動,傾向方向為非充分開采,因此計算時走向需要乘上傾向采動系數。為使結果更可靠,均勻選取全盆地觀測值作擬合初值,避免選點帶來的參數不穩定性,其中由于東南出矸石山在工作面開采范圍內,因此選點時避開矸石山區域,共計取到368 個點。根據概率積分法,在Matlab 中多次使用模擬退火算法預計,選取目標函數值最小的參數作為最終參數,最終求得:下沉系數q 為0.6,走向主要影響角正切tanβ 為0.6,傾向下山方向主要影響角正切tanβ1為1.337,傾向上山方向主要影響角正切tanβ2為3.43,開采影響傳播角為89°。根據預計參數畫出下沉等值線,并模擬下沉盆地,模擬下沉盆地圖如圖7。
圖7 模擬下沉盆地圖Fig.7 Simulated subsidence basin graph
用配對T 檢驗判斷2 組數據相似程度,在執行T 檢驗之前需要驗證實測值與預計值的誤差是否符合或近似符合正態分布,因此隨機選取40 組數據執行Shapiro-Wilk(W 檢驗),利用本次實驗樣本計算得出檢驗統計量Wα=0.928。取顯著性水平0.05,根據樣本數量查表得到界限值為0.980 131,可知Wα=0.928<0.981 031,從而在該顯著性水平下,接受原假設,即樣本服從正態分布,實測數據與預計數據間的誤差為偶然誤差,對實測值與預計值做配對T 檢驗,得到結果為0.704>0.05,配對數據無顯著性差異,擬合結果較好。
實測數據由于存在誤差,與真實模擬盆地存在一定的差距,而誤差也會對求參結果帶來很大的影響,因此需要對參數進行抗差分析。各參數對誤差的敏感度不同,有些參數敏感度很高,而有些參數在誤差很大時也不會產生很大的影響。選取唐家會模擬盆地,分別加入不同中誤差來對比參數的影響程度,對每一次加入的誤差多次求參,優選出擬合程度最好的參數來作為在該誤差時的參數,全盆地加入隨機誤差后的求參結果見表1。
由表1 可得,當對實測下沉值加入約(1%~10%)W0,即100~700 mm 的測量中誤差時,求取的下沉系數q、主要影響角正切tanβ、均與參數真值吻合度較高,比較穩定。當對實測下沉值加入超過10%W0,即800 mm 中誤差時,下沉系數與主要影響角正切誤差較大,故認為此時的求參結果不可靠,而原始擬合下沉中誤差為589 mm,占最大下沉值8.1%,因此求參結果可靠。
表1 全盆地加入隨機誤差后的求參結果Table 1 Parameters acquired after adding random error to the whole basin
以內蒙古唐家會煤礦為例,用BP 神經網絡算法對下沉點云做去噪處理,并與中值去噪、最小二乘去噪、最近鄰去噪方法做對比,提高了下沉模型精度。針對傳統依靠走向、傾向主斷面求參精度不高的情況,融合模擬退火算法與概率積分模型,選用全盆地特征點進行參數反演,并對參數結果作抗差分析。結果表明,擬合模型測量中誤差占最大下沉值的8.1%,預計結果滿足工程要求,克服了無人機精度不高對求參帶來的影響,為后續采煤工作、災害治理提供技術支持。