基于深度學習的位場邊界識別方法

張志厚， 姚禹， 石澤玉， 王虎， 喬中坤， 王生仁，覃禮貌， 杜世回， 羅鋒， 劉慰心

1 西南交通大學地球科學與環境工程學院， 成都 611756 2 西南交通大學， 高速鐵路線路工程教育部重點實驗室， 成都 610031 3 吉林大學地球探測科學與信息技術， 長春 130026 4 中鐵一院勘察設計院集團有限公司， 西安 710043

0 引言

高精度的位場(重力、磁力)勘探被廣泛應用于礦產資源勘察(Kamm et al., 2015; Krahenbuhl and Li, 2017; Liu et al., 2017)、巖石結構探測(Prutkin and Saleh, 2009; Astort et al., 2019)、水文地質勘察(Williams and Stevens, 2002)、環境地質調查(Zunino et al., 2009; Gavazzi et al., 2016; Yannah et al., 2019)、地質體定位(Beiki et al., 2012; Tang et al., 2018; Hu et al., 2019)、油氣資源遠景評價(Gadirov et al., 2018)等.位場的邊界識別可為地質體解譯提供更加豐富的信息.地質體的邊界通常是指構造界線，不同巖性不同地質體的接觸線，其物理性質表現為具有一定密度或磁性差異的邊緣位置(王萬銀, 2010).位場邊界識別方法主要有數值計算和數理統計兩大類方法，兩類方法針對不同深度不同姿態目標體的邊界識別(Cordell，1979；王萬銀，2010; Wang et al., 2010; Sertcelik and Kafadar, 2012; 馬國慶等, 2012; Ma, 2013; Ma et al., 2015; Zhou et al., 2013; Beiki et al., 2014; Yuan et al., 2014; Hidalgo-Gato and Barbosa, 2015; Yuan and Yu, 2015; 顏廷杰等, 2016; 周帥等, 2016; 于平等, 2019; 湯井田等, 2019; 郭華等, 2019)都取得了一些較好的應用效果.以上方法的主要思路是對位場數據進行各類型基本運算，然后研究其最大值、零值與地質體邊界之間的對應關系.這些計算方法都屬于無監督式的機器運算，其計算結果精度與目標體空間分布有一定的關聯性，因此在實際應用中，選擇適合的邊界識別方法，才能獲得可靠的結果.

近年來，大數據已滲透到諸多行業里，已成為社會發展的重要生產要素和戰略資源.人工智能(Artificial Intelligence, AI)、機器學習(Machine Learning，ML)和深度學習(Deep Learning,DL)等已成為許多領域科學家和工程師關注的熱點，地球物理領域也不例外.諸多ML方法在地球物理的數據分析與處理中(Guo et al., 2011; Sun and Li, 2015; Jiang et al., 2016; Nawaz and Curtis, 2019; Xie et al., 2019；劉昊楠等；2020；王昊等，2020)都表現出較好的優勢.ML屬于AI的領域，它是通過統計計算從大數據中提取相關特征.其優勢是計算速度快，精度高；其挑戰是需要獲取足夠的樣本數據和優化相關參數(Kim and Nakata, 2018).DL(LeCun et al., 2015; Schmidhuber, 2015)作為ML的一個重要新分支，在語音和圖像的識別與分類中取得了超越性進展(Russakovsky et al., 2015; Greenspan et al., 2016; Davoudi et al., 2019; Pan et al., 2019; Schoppe et al., 2020)，其標志著ML的里程碑.尤其是自2012年ImageNet大賽上，AlexNet架構的DL網絡結構以極大的優勢獲得圖像分類的冠軍(Krizhevsky et al., 2012)，卷積神經網絡(Convolutional neural networks, CNN)模型開始受到計算機視覺領域的廣泛關注.CNN能夠從大量的監督樣本中自動學習到由底層位置特征到高層屬性特征的層次化特征表達，避免了人工特征設計的弊端，極大地增強了模型的泛化能力，已經在地球物理資料自動解譯方面取得了巨大的突破(Xiong et al., 2018; Li et al., 2019; Shi et al., 2019; Wu et al., 2019; Zhang et al., 2020a).

以上針對地球物理資料自動解譯的DL模型框架是通過修改CNN網絡結構，實現了“端到端”訓練和預測，其實質上屬于DL語義分割框架的方法.Long等(2015)將CNN中的全連接層替換為卷積層，提出了一種全卷積深度神經網絡(Fully Convolutional Networks, FCN)模型，并成功地應用于醫學圖像語義分割(Russakovsky et al., 2015)，實現了高精度的類別解釋.但是同時該網絡模型也引入了兩個問題，一是池化等下采樣步驟提取深層屬性特征時容易顯著降低特征數據的分辨率，從而導致最終上采樣得到的預測結果嚴重損失局部細節信息，因此不能對異常體邊界等細節做出精細化的表達；二是卷積通過一定大小的感受野只能利用局部信息，缺少全局信息的輔助，容易造成局部目標體預測混淆，最終結果不連續且包含一定的椒鹽噪聲(喬文凡等, 2018).基于以上問題，許多學者提出了金字塔解析網絡結構(Zhao et al., 2017)、SegNet網絡結構(Badrinarayanan et al., 2017)、U-Net網絡結構(Ronneberger et al., 2015)和DeepLab系列網絡結構(Chen et al., 2018)等進行了改善，這些網絡結構及其改善策略也極大地促進了FCN方法在地球物理反演問題求解中的研究與應用，如井地電磁二維反演(Puzyrev, 2019)、直流電阻率反演(Liu et al., 2020)、速度模型重建(Li et al., 2020a)、重力異常和重力張量異常聯合反演(張志厚等, 2021a)、以及磁異常和磁梯度異常聯合反演(張志厚等, 2021b)等.雖然基于經典FCN網絡結構在地球物理模型語義分割中都取得了令人滿意的結果，但仍然存在對細小異常體邊界分割不夠精細的情況，如密集細小的次生斷裂，巖性差異較小的界限等.為了進一步提高FCN的預測精度，在計算機視覺領域，相關研究者分別采用密集連接網絡(Zhou et al., 2020)、多尺度漸進融合網絡(Xian et al., 2020)、自適應閾值多模型融合網絡(Jiang et al., 2020)和注意機制密集連接網絡(Xia et al., 2021)等進行了圖像高層特征的提取，實現了不同場景下圖像的精細化表達.這些研究為本文重磁位場邊界識別帶來了更為廣泛的應用前景.

基于此，本文提出了位場邊界識別的DL方法，DL網絡結構為一種融合多尺度特征的注意機制與密集跳躍連接網絡，稱為PFD-Net(Potential Field Detection-Net).PFD-Net首先是在基于改進的U-Net網絡結構中嵌套標準卷積模塊(Standard Convolution Modules, SCM)和密集跳躍連接模塊(Dense Skip Connections Modules，DSCM)；然后通過注意模塊(Attention Module, AM)融合不同大小的低層、高層特征，形成以U-Net網絡為骨干的DSCM和AM網絡結構；此外，相較于網絡結構的精心設計而言，網絡結構充分學習所依賴的訓練樣本數據同樣值得關注，其多樣性與足量性的快速構建決定了方法的實用性，本文采用張志厚等(2021a, b)提出的基于網格點幾何格架的位場空間域快速正演算法進行高效構建.最后，對本文所提方法進行理論數據和實測數據的檢驗，以此來證實的文中方法的可行性和有效性.

1 邊界識別方法

1.1 CNN網絡結構與語義分割

地下半空間存在的異常體(或邊界)會在觀測面上引起位場響應，不同異常體(或邊界)引起的位場響應也表現出不同特征.即具有一定的空間相關性和局部存在性.總之，位場數據的輸入和輸出有兩個特點：(1)空間的對應性；(2)局部的存在性(張志厚等, 2021a,b).而CNN方法正是通過卷積算子重點學習輸入圖像與輸出標簽的局部性和空間性(Shi et al., 2019)，如提取圖像的邊界信息(Deng et al., 2012)等.此外，小子域濾波器法(Jiang and Gao, 2012)(該濾波器由較多卷積算子的串聯組成)的成功應用也表明了通過卷積算子識別位場邊界的可行性.

CNN網絡結構是由一系列不同尺度的卷積算子和相關濾波器組成，卷積算子的主要功能提取輸入數據或者中間結果的特征，每一層卷積計算公式為(Li et al., 2019)：

(1)

典型的CNN網絡結構包含了卷積運算、非線性激活、層間批處理歸一化、池化操作以及節點丟棄等.每一層的輸出構成了下一層的輸入，重復這些步驟可以堆疊一個較深的網絡結構，因此稱為深度神經網絡.網絡節點超參數(如卷積核，偏置項等)是通過監督學習和反向傳播進行更新.優化算法通常采用隨機梯度下降法、以及適合大數據量且具備計算效率高收斂速度快的RMSProp算法和ADAM算法等.隨著訓練輪數的增加，超參數的更新，使得網絡輸出與標簽的誤差逐漸減小，最終獲得了訓練數據到標簽數據的復雜非線性映射關系.

CNN方法典型應用于分類(預測類別標簽)任務，即輸出是一個單一的熱點向量，表明該圖像屬于某一類的可能性.然而，一些任務中，如圖像分割(Shi et al., 2019; Davoudi et al., 2019; Pan et al., 2019; Schoppe et al., 2020)，不僅需要對圖像進行分類，而且還要進行超像素分類(圖1a).本文也將位場的邊界識別視為語義分割任務，即將位場的數據作為輸入端，輸出端分為邊界位置與非邊界位置兩類進行監督學習(圖1b).FCN(Shi et al., 2019)方法作為CNN典型分支可實現密集像素預測，可生成任意大小的分割圖像，并且計算速度快，成為當下最為先進和流行的語義分割方法.隨后，較為典型的FCN網絡結構，如U-Net(Ronneberger et al., 2015)、SegNet(Badrinarayanan et al., 2017)和DeconvNet(Noh et al., 2015)都出色地完成了相關語義分割任務.并且適用于地球物理數據解譯分割和地球物理數據端到端反演的網絡結構大都是在以上網絡結構基礎上的改進(Li et al., 2019, 2020a; Wu et al., 2019; Xiong et al., 2018; Shi et al., 2019; Puzyrev, 2019; Liu et al., 2020; 張志厚等, 2021a,b)，此類改進網絡都屬于編碼器-解碼器的CNN網絡結構，圖1c為典型的編解碼器CNN網絡結構示意圖.

圖1 CNN端到端示意圖 (a) 圖像語義分割； (b) 位場數據通過FCN預測邊界位置； (c) 典型的編解碼CNN網絡結構.Fig.1 Schematic diagram of CNN end-to-end (a) Image semantic segmentation; (b) Edge detection based on potential field data in FCN; (c) Typical CODEC in CNN.

1.2 PFD-Net網絡結構

為了進一步提高位場邊界細節的識別精度，本文提出了一種名為PFD-Net的網絡結構，其流程圖如圖2所示，該網絡結構采用改進的U-Net作為骨干網，包含了SCM、DSCM、AM等模塊，以及用于信息融合的相關連接.主干網絡通過4個尺度的下采樣來提取特征，而密集跳躍連接則充分利用不同尺度的特征進行組合.注意模塊通過高級特征加權初始低級特征，以對密集邊界進行精準識別.

圖2 位場邊界識別PFD-Net總體結構圖Fig.2 Schematic diagram of potential field edge detection PFD-Net

本文基于改進的U-Net(Badrinarayanan et al., 2017)結構作為PFD-Net結構的骨干網絡，如圖3所示.該骨干網絡是U-Net網絡結構的一個簡化，其包含一個獲取地質特征的收縮路徑和一個用于向上采樣的對稱擴展路徑，以此實現精準定位.并采用基于最大池化和反卷積的收縮-擴展結構.當給定一個固定大小的卷積核(大小為3×3)，網絡結構的有效感受野隨著網絡層數的增加而增大.隨著網絡深度的增加，左側路徑的通道數從16個變為256.通過跳躍層，將右側路徑上的局部、淺層特征映射與左側路徑上的全局、深層特征映射相結合.SCM各層間采用批處理歸一化(Puzyrev, 2019)增加網絡訓練的穩定性，激活函數選用DL類型最常用的線性整流函數ReLU(Xu et al., 2015)來提高訓練的收斂性，編碼階段的各卷積層之間加入Dropout層(Hinton et al., 2012)來增強網絡結構的泛化能力，最后通過1×1的卷積之后再采用Sigmoid(Badrinarayanan et al., 2017)激活函數進行輸出.

骨干網絡結構是對原始U-Net網絡結構進行了三個方面的修改：一是輸入數據不是RGB三色通道的像素值，而是每一個觀測點的位場響應值；二是簡化了網絡結構的深度及卷積通道數，這是由于本文的目標僅僅是兩個類別的語義分割，其復雜度遠低于多類別分割任務，所以不太適用太深的網格數和太多的卷積通道數；三是地質體邊界所對應的標簽不全是封閉的曲線，也不是語義分割區塊，而是與地質體邊界對應的曲線段，其可能不封閉.

圖3 骨干網絡結構—改進的U-Net 由4個下采樣層和4個上采樣層組成，下采樣和左側卷積層對應編碼階段，上采樣和右側卷積層對應解碼階段，圖中黑色數字表示特征圖的通道數，紅色數字表示單通道數據的大小，Conv為卷積，ReLU為激活函數，Sigmoid為輸出激活函數，MaxPool為最大池化，Up-Conv為上采樣和卷積 Dropout為丟棄節點，BN為批處理 歸一化，Copy & Concatenate為復制跳躍鏈接.Fig.3 Backbone network structure——the improved U-Net It consists of four down-sampling layers and four up-sampling layers. The down-sampling layers and the left convolution layers correspond to the coding stage, and the up-sampling layers and the right convolution layers correspond to the decodingstage.The black numbers in the Figure indicate the number of channels in the feature map, and the red numbers indicate the size of single channel data, Conv is the convolution, ReLU is the activation function, Sigmoid is the output activation function, MaxPool is the maximum pooling, Up-Conv is the up-sampling and convolution, Dropout is discarding node, BN is batch normalization, Copy & Concatenate is the copy skip link.

圖4為PFD-Net網絡的拓撲結構.文中將網絡結構中包含1層深度U-Net和1個AM的網絡結構稱為PFD-Net-L1(圖4a)，類似的有PFD-Net-L2(圖4a)和PFD-Net-L4(圖4b)等.為了便于繪圖，用Xij代替圖4中的標準卷積模塊.圖4b中所示的PFD-Net-L4是以U-Net網絡結構為基礎骨干網絡，其中節點X1,1到節點X5,1為編碼結構，節點X5,1到節點X1,5為解碼結構.首先在改進的U-Net基礎上，增加了6個節點(X1,2、X1,3、X1,4、X2,2、X2,3、X3,2)放在不同的尺度上；其次，對所有節點輸出的特征進行密集跳躍連接；第三，在同一層編解碼器上添加注意機制模塊，其輸入為相鄰層的節點Xi,1和Xi+1,5-i；最后，對第一層輸出的各個節點與AM的合成結果進行輸出，從而得到最終的語義分割結果.

1.3 密集跳躍連接

經典U-Net是將編碼階段高分辨率特征與解碼階段上采樣特征進行融合，該方式是一種淺層信息融合，未將深層信息充分利用，造成編解碼器之間存在語義鴻溝(邢妍妍等，2020).因此，有必要將深層信息上采樣融合淺層信息或者將淺層信息下采樣融合深層信息，從而可以加強編解碼器之間的語義聯系，提升模型的性能.本文重新設計的PFD-Net是在編碼器的高分辨率特征與解碼器相應的豐富語義特征合并之前捕獲背景對象的細粒度細節.同時，密集跳躍連接路徑重用不同層的特征作為輸入，并連接不同深度的網絡特征層，其目的是縮小語義鴻溝，從而使得編碼特征和相應解碼特征在語義上相似時，分割效果得到了極大的改善.

PFD-Net的密集跳躍連接由多個SCM和上采樣模塊組成，SCM和跳躍路徑是核心之一.通過嵌套的密集跳躍連接，每個SCM的輸入為同一層先前所有卷積輸出和下一層上采樣到同一級聯層輸出的串聯結果，這種方式有助于恢復下采樣引起的信息丟失，并使所有層信息在網絡中充分傳播.圖5所示為密集跳躍連接的結構，可以看出，在節點X1,1和X1,5之間加入了三個SCM，每一個SCM都有一個連接層，它都融合了同層以前的卷積輸出和下一層SCM的上采樣輸出.融合后模塊輸出繼續與下一層相對應模塊上采樣特征再融合，如此迭代直至下一層沒有對應的上采樣模塊.每一模塊的輸出結果為：

圖4 PFD-Net網絡結構 (a) PFD-Net-L1和PFD-Net-L2； (b) PFD-Net-L4.Fig.4 PFD-Net network structure

圖5 密集跳躍連接結構 Conv為卷積，BN為批處理歸一化，ReLU為激活函數，Dropout為 丟棄節點，Up-Conv為上采樣和卷積，Skip Connections為跳躍連接.Fig.5 Dense skip connection structure Conv is convolution, BN is batch normalization, ReLU is activation function, Dropout is discard node, Up-Conv is up-sampling and convolution.

(2)

式中，Xi,j表示當前特征提取模塊的輸出，i(i=1、2…5)表示編碼器下采樣的索引，j為同一層模塊序號，也是跳躍連接路徑的索引序號，C(·)表示特征提取的卷積操作，U(·)表示上采樣操作，[·]表示特征通道連接.

因此，PFD-Net整體特征融合結構呈倒金字塔形，每一層都包含了密集模塊，通過密集跳躍連接(圖5)將編碼階段提取的特征向后傳遞并深度融合，縮短編解碼語義鴻溝，提升模型的學習能力.同時，密集跳躍連接可以充分有效組合網絡模型中各深度層的大量特征，提高預測精度.此外，密集跳躍連接可以有效地減少梯度消失(Xia et al., 2021).

1.4 注意力機制

密集跳躍連接通過編碼器的淺層信息和下一層的高級特征來重建細節信息.然而，這種設計沒有考慮特征的多樣性和全局的信息，可能會導致邊緣細節信息不完全一致(Xia et al., 2021).注意力機制是一種能夠突出關鍵特征點的算法，機器學習的注意模塊可以對數據中的特征信息進行重點提取，區分主次信息，削弱消除背景噪聲等無用信息，使得輸入數據的特征得以有效利用(Zhang et al., 2020b).圖6給出了本文注意模塊網絡結構，該模塊將全局平均池化應用于高層特征層，并利用高層特征為低層特征提供加權指導(Li et al., 2018).因此，注意模塊能夠增強解碼器上采樣過程中高層特征權重系數的影響，并有效地利用金字塔結構的多尺度特征映射，從而提高目標區域的像素級分類精度.

圖6 全局注意力的上采樣模塊結構 Conv為卷積，ReLU為激活函數，Global Pooling 為全局池化，BN為批處理歸一化.Fig.6 Global attention upsample module structure Conv is convolution, ReLU is activation function and BN is batch normalization.

2 模型試驗

2.1 樣本數據集及訓練

DL學習樣本數據主要來自于真實數據和人工合成數據.相關學者探索了基于精細化標注真實數據訓練的語義分割模型(Bai et al., 2017)，然而真實樣本數據標注耗時、昂貴，且對人為錯誤敏感；人工合成方法可以獲取大量的數樣本數據，足量代表性數據可極大提高訓練模型的泛化性能，如在深度學習的重磁反演(張志厚等, 2021a，b)、電磁反演(Puzyrev, 2019; Li et al., 2020b)、直流電阻率反演(Liu et al., 2020)、地震反演(Li et al., 2020b)中取得了較理想的成果.因此，本文通過正演計算來獲取樣本數據集.模型試驗是將地下半空間劃分為31×31×15個長方體單元，其中長方體單元大小為100 m×100 m×50 m，模擬觀測數據區大小為3.1 km×3.1 km.坐標系采用笛卡爾坐標系，z軸向下，計算水平面位場數據的個數為32×32個，即采樣數據點為地下半空間網格點在水平面的投影.因此，樣本輸入數據的大小為32×32個，樣本標簽數據個數也為32×32個.每一個標簽數據是由多個長方體單元組成模型邊界在水平面上的投影線，標注為“1”，非邊界位置標注為“0”.文中將長方體位置設置于地下半空間的非邊緣區域，與文獻(張志厚等, 2021a, b)中長方體單元組合類似，其具體如下：I型：1×1×1塊體(均勻抽樣4761對)；II型：2×2×2塊體(3388對)；III型：4×4×4塊體(4000對)；IV型：6×6×6塊體(1944對)；V型：8×8×8塊體(512對).各模型塊體剩余密度為1 g·cm-3，或者剩余磁化強度大小為1.0 A·m-1.為了便于計算及分析，磁場正演計算所采用的地磁場與剩磁異常體的磁化傾角和偏角都一致，其中磁化傾角為90°，偏角為0°.通過正演獲得的位場(重力異常或磁異常)樣本數據對個數為14605對.由于重磁位場的正演為線性關系(張志厚等, 2021a,b)，故可利用已正演的樣本數據進行不同物性大小的樣本數據擴充，擴充后的樣本數據為29210對.將樣本數據劃分為訓練數據、測試數據以及驗證數據，其對應比為18∶1∶1，即訓練數據集包含26452個樣本數據，驗證集和測試集各有1460個樣本數據.訓練輪數為100，每次訓練后進行一次驗證，訓練誤差為L2范數，其表達式為：

(3)

圖7 PFD-Net訓練誤差及驗證誤差Fig.7 Training error and validation error of PFD-Net

2.2 組合模型邊界識別

采用了四組長方體模型(圖8a、b、c、d)對本文所提方法進行檢驗.模型一、模型二和模型三在地下空間位置及物性參數如表1所示，模型四為圖8h中白色折線與x軸、y軸所圍成的近似上三角塊體，塊體厚度為0.5 km，塊體的中心深度為0.55 km，塊體密度為1 g·cm-3.通過正演獲得模型一、二、三、四的重力異常或磁異常分別如圖8e、f、g、h所示.

表1 長方體模型參數Table 1 Cuboid model parameters

利用文中已訓練完成的PFD-Net網絡分別對四組模型對應的位場數據進行邊界識別，其結果分別如圖9a、b、c、d所示.同時為了更加貼合實際資料處理，對位場數據分別加入其最大值10%的隨機噪聲(圖8i、j、k、l)，再對四組模型進行邊界識別，結果分別如圖9e、f、g、h所示.此外，對四組模型的正演數據及其含噪聲數據分別采用總水平導數法(Total Horizontal Derivative, THDR)(Cordell, 1979)、加強解析信號傾斜角法(Enhanced Analytic Signal Tilt Angle, EASTA)(顏廷杰等，2016)以及基于三維構造張量的兩種歸一化反正切濾波器方法(Normalization Structural Tensor Edge Detection, NSED; Normalization Vertical Derivative Structural Tensor Edge Detection, NVSED)(周帥等，2016)進行邊界識別，以此來對比本文所提方法的效果.THDR、EASTA、NSED及NVSED的識別結果分別如圖10、11、12、13所示.

圖8 三維模型透視圖及正演結果 (a)、(b)、(c)、(d)分別為四種模型的三維透視圖；(e)、(f)、(g)、(h)分別為(a)、(b)、(c)、(d)圖中模型對應的正演數據； (i)、(j)、(k)、(l)分別為(e)、(f)、(g)、(h)正演的加噪聲數據；圖中白線為模型對應的水平邊界位置.Fig.8 Perspective drawing of three-dimensional model and forward calculation (a),(b),(c) and (d) Correspond to perspective drawing of four kinds of three-dimensional model respectively; (e),(f),(g) and (h) Correspond to forward calculation data of (a),(b),(c) and (d) respectively; (i),(j),(k) and (l) are the noise-contained forward calculation data of (e),(f),(g) and (h) respectively; White line in the Figure is the horizontal boundary position of the model.

圖9 四種模型邊界識別結果 (a)、(b)、(c)、(d)分別為四種模型正演數據采用文中方法的邊界識別結果；(e)、(f)、(g)、(h)分別為相對應加噪聲數據 采用文中方法的邊界識別結果.Fig.9 Forward data and edge detection result of four kinds of models (a),(b),(c) and (d) Correspond to edge detection result of forward calculation data of model one, two, three and four with PFD-Net; (e),(f),(g) and (h) Correspond to edge detection result of noise-contained forward calculation data of model one, two, three and four with PFD-Net.

圖10 THDR識別結果 (a)、(b)、(c)、(d)分別為模型一、二、三、四正演數據的THDR邊界識別結果； (e)、(f)、(g)、(h)分別為四種模型正演加噪數據THDR的 邊界識別結果；圖中白線為模型對應的水平邊界位置.Fig.10 Edge detection result of THDR (a),(b),(c) and (d) Correspond to edge detection result of forward calculation data of model one, two, three and four with THDR; (e),(f),(g) and (h) Correspond to edge detection result of noise-contained forward calculation data of model one, two, three and four with THDR; White line in the Figure is the horizontal boundary position of the model.

圖11 EASTA識別結果 (a)、(b)、(c)、(d)分別為模型一、二、三、四正演數據的EASTA邊界識別結果； (e)、(f)、(g)、(h)分別為四種模型正演加噪數據EASTA的 邊界識別結果；圖中白線為模型對應的水平邊界位置.Fig.11 Edge detection result of EASTA (a),(b),(c) and (d) Correspond to edge detection result of forward calculation data of model one, two, three and four with EASTA; (e),(f),(g) and (h) Correspond to edge detection result of noise-contained forward calculation data of model one, two, three and four with EASTA; White line in the Figure is the horizontal boundary position of the model.

圖12 NSED識別結果 (a)、(b)、(c)、(d)分別為模型一、二、三、四正演數據的NSED邊界識別結果； (e)、(f)、(g)、(h)分別為四種模型正演加噪數據NSED的 邊界識別結果；圖中白線為模型對應的水平邊界位置.Fig.12 Edge detection result of NSED (a),(b),(c) and (d) Correspond to edge detection result of forward calculation data of model one, two, three and four with NSED; (e),(f),(g) and (h) Correspond to edge detection result of noise-contained forward calculation data of model one, two, three and four with NSED; White line in the Figure is the horizontal boundary position of the model.

圖13 NVSED識別結果 (a)、(b)、(c)、(d)分別為模型一、二、三、四正演數據的NVSED邊界識別結果； (e)、(f)、(g)、(h)分別為四種模型正演加噪數據NVSED的 邊界識別結果；圖中白線為模型對應的水平邊界位置.Fig.13 Edge detection result of NVSED (a),(b),(c) and (d) Correspond to edge detection result of forward calculation data of model one, two, three and four with NVSED; (e),(f),(g) and (h) Correspond to edge detection result of noise-contained forward calculation data of model one, two, three and four with NVSED; White line in the Figure is the horizontal boundary position of the model.

(1)模型一

模型一由四個大小相同、埋深相同、剩余密度相同的異常體組成(圖8a)，其在水平面投影的邊界如圖8e中白色線段所示.從圖9a、e可以看出無論是無噪聲重力異常還是加噪數據的PFD-Net識別結果都與地質體邊界相吻合，相比THDR(圖10a、e)、EASTA(圖11a、e)、NSED(圖12a、e)及NVSED(圖13a、e)的識別結果，本文方法精度更高，其計算結果都“聚焦”于真實邊界處，尤其是在含噪聲情況下，EASTA(圖11e)和NVSED(圖13e)受噪聲影響較大，本文方法的魯棒性更強.表明本文PFD-Net網絡不僅可以有效的對地質體邊界位置進行準確識別，而且對單一塊體模型所構建的樣本數據集進行訓練，可對類似塊體的組合模型進行邊界位置預測，相應數值檢驗結果也驗證了理論上CNN卷積核具有局部性和權值共享的特點和優勢，也表明了本文PFD-Net網絡所具備的泛化性能.

(2)模型二

模型二由兩個大小相同、埋深不同、剩余磁化強度不同的磁性異常體組成(圖8b)，兩個異常體在水平方向上有一定距離，其在水平面投影的邊界如圖8f中白色線段所示.該模型無噪聲數據的PFD-Net識別結果(圖9b)與地質體的邊界一致，淺部異常體識別結果相比深部異常體識別結果較好，兩異常體臨邊位置識別效果稍差.但相比其他四種方法(THDR、EASTA、NSED及NVSED)對不同深度組合模型的識別結果(圖10b、11b、12b、13b)，PFD-Net方法識別精度最優，尤其是深部異常體水平邊界的刻畫也更加精確.

含噪聲數據的識別結果(圖9f)相比無噪聲數據的識別結果較差，尤其是深部異常體的結果不夠清晰，具體表現在某些邊界位置變寬，但并不影響識別效果.相比其他四種方法(THDR、EASTA、NSED及NVSED的識別結果(圖10f、11f、12f、13f))，PFD-Net方法仍然具備識別的準確性.表明本文PFD-Net網絡的抑制噪聲能力強.

此外，本文PFD-Net對淺部異常體的識別精度相比深部異常體識別精度高的原因可能如下：大小相同、物性相同的模型在不同埋深處的位場異常存在較大的差異，而模型在水平面投影的邊界位置相同，即樣本數據存在位場響應相差較大，但其對應的標簽卻一樣.因此，訓練的網絡在具備泛化性能的同時也必須具備一定的容錯性，該容錯性是導致深部異常體識別精度稍遜的原因之一.

(3)模型三

模型三為兩個大小不同、埋深不同、剩余密度不同，但水平中心位置相同的異常體組成(圖8c)，兩個異常體在垂直方向上相連，其在水平面投影的邊界如圖8g中白色線段所示.從圖9c可以看出，在無噪聲情況下，除去淺部異常體附近的干擾，其余部分識別效果較好；含噪聲情況下，深部異常體的三個邊部雖有影響(圖9g)，但并不影響整體的識別效果，其在異常體邊界部分的定位是極為準確的，并且淺部異常體的識別精度未降低.說明了本文PFD-Net對深大異常體及其上部的淺小異常體都具有一定的識別能力.

無噪聲情況下，其余四種方法也較好地識別出“淺小深大”的組合模型邊界(圖10c、11c、12c、13c)，雖然NVSED方法也表現出較高的識別精度(圖13c)，但并沒有本文所提方法的精度高(圖9c).含噪聲情況下，THDR(圖10g)和NSED(圖12g)的識別精度降低，EASTA(圖11g)和NVSED(圖13g)受噪聲影響較大，識別能力顯著下降，而本文所提方法具有較強的穩定性(圖9g).

(4)模型四

圖8d所示為模型四的三維透視圖，圖8h為模型四的邊界水平投影及正演位場數據.可以看出無噪聲數據的識別結果在中間部分的識別精度較高，與理論模型輪廓保持高度一致，在白色曲折線段的首尾部分對應的邊界識別結果與理論模型輪廓存在一定的差異(圖9d).該差異并不能說明PFD-Net的識別效果在邊部稍差，恰恰相反能夠說明本文所提方法具備較強泛化性與高精度性.這是由于在異常體角點處能夠識別單一塊體位場響應數據一是不完整，二是整個模型并不是標準的臺階模型，而是近似上三角塊體.同理，在x軸和y軸附近也存在一定識別結果.圖9h為含噪聲數據的識別結果，識別精度沒有明顯下降，表明本文PFD-Net識別精度高且具有較好的穩定性.

相比其他方法，無噪聲情況下NVSED方法的識別結果相對較好(圖13d)，其識別效果與本文所提方法相近，但在含噪聲情況下，NVSED的識別效果較差(圖13h)，EASTA同樣也受到噪聲的影響較大(圖11h)，THDR(圖10d、h)、EASTA(圖11d)和NSED(圖12d、h)的識別結果除去邊部和拐角的影響，也能準確的識別出邊界的位置.但總體而言，本文所提方法無論是在識別精度還是穩定性上都優于其他四種方法.

2.3 未化極磁性體模型邊界識別

為了檢驗本文PFD-Net對未化極磁性體的邊界識別效果，需要先構建不同磁化傾角和磁化偏角的樣本數據集，樣本數據集類型如前文所述的4個類型(I型、II型、III型和IV型)，然后訓練PFD-Net再進行邊界預測.檢驗模型如圖14所示，由三個長方體模型組合而成，長方體相關參加見表2，將磁化傾角分別設計為60°和30°兩種情況進行檢驗.文中正演了兩種不同磁化傾角下組合磁性模型體的理論異常(圖15a、e)，以及加入理論異常最大值5%的隨機噪聲(圖15c、g)，圖15a、c、e、g對應的識別結果如圖15b、d、f、h所示.看出未化極理論磁異常數據以及加噪數據的識別結果都與地質體邊界相吻合，尤其是含噪數據的識別結果也具有較好的穩健性.表明了本文PFD-Net網絡同樣適用于未化極磁異常的邊界識別.

圖14 磁性體三維模型透視圖Fig.14 3D perspective view of the magnetic body model

表2 未化極磁性體模型參數Table 2 Model parameters of unpolarized magnetic body

3 實際資料應用

采用本文的方法對藏東南某鐵路隧道西段的航空磁測數據進行邊界識別.藏東南某隧道西段長達近20 km，大致位于念青唐古拉山脈東段與喜馬拉雅東構造結北部前緣交會區域.該部位為典型的高山峽谷地貌：一是海拔高差大，從平均約2000 m跨度到約4500 m(圖16a)；二是雅魯藏布江及其支流的快速下切和廣泛分布的冰川剝蝕作用共同造就了深切峽谷的發育.構造框架方面，該隧道位于岡底斯巖漿弧地層，隧道西段洞身段巖性主要為燕山期片麻巖、片麻狀二長花崗巖、花崗巖等侵入巖體，局部夾有一些變質巖捕虜體；隧道出口地層為白堊系大理巖夾變鈣質砂巖、變基性火成巖；隧道北側出露地層為泥盆系-石炭系大理巖，以及大理巖夾變鈣質砂巖，該套地層與下覆燕山期片麻狀二長花崗巖呈不整合接觸(圖16a)(羅鋒等，2021).另外，隧道還穿過嘉黎—察隅構造混雜巖帶，可見碎裂巖、斷層角礫巖，局部發育有構造蝕變巖.總體而言，該段落由于印度板塊與歐亞板塊在喜馬拉雅東構造地區的陸-陸碰撞及碰撞后的持續向北推移和楔入作用，導致斷層變形歷史和破裂樣式異常復雜及地層單元和巖性特征多變性，再受到高山峽谷型地表過程的疊加效應，難以在該地區直接有效地開展相關地質單元和構造界限的野外判別.

航磁原始數據經過校正、調平、化極以及網格化后的結果如圖16b所示，圖中觀測區的大小為19.8 km×3.1 km，數據點距為100 m×100 m，共199×32個點.從圖16b中可以看出測區中西部異常呈近似北西向條帶狀分布，有一定的連續性.整個測區在局部都存在磁異常封閉圈，可能表明存在花崗巖等侵入巖體，或者為變質巖捕虜體.

圖15 未化極磁性體模型檢驗效果 (a)、(e)分別對應磁化傾角為60°與30°的正演磁異常；(c)、(g)分別為(a)、(e)對應的加噪聲數據；(a)、(c)、(e)、(g)對應的 識別結果為(b)、(d)、(f)、(h)；圖中白線為模型對應的水平邊界位置.Fig.15 Test effect of unpolarized magnetic model (a) and (e) Correspond to forward calculation magnetic anomalies with magnetization inclination angles of 60°and 30° respectively; (c) and (g) are the noise-contained data correspond to (a) and (e); (b),(d),(f) and (h) are the recognition results corresponding to (a),(c),(e) and (g); White line in the Figure is the horizontal boundary of the model.

圖16 藏東南某鐵路隧道西段地質地形圖(a)、航磁異常(b)和本文方法的識別結果(c)Fig.16 Geological topographic map of the western section of a certain railway tunnel in southeastern Tibet (a), aeromagnetic anomalies (b), and the detection results of this method (c)

PFD-Net網絡輸入與輸出大小為32×32，不能直接對任意大小的實際資料直接預測，需要先將實際資料切割為若干個與訓練網絡大小一致的數據，然后分別預測再合成為最終結果.本文為了避免因邊界的特征提取問題而產生拼接痕跡，繼而影響最終的識別效果，采用了滑動切割方法對網格化原始數據進行裁剪，滑動點距為16個點，切割窗口大小為32×32，共切割13個數據塊.文中將13個數據塊分別采用PFD-Net網絡進行識別，重復區域進行融合處理，最終獲得的結果如圖16c所示.從圖16c可以看出，局部明顯存在封閉塊體，識別封閉線大于0.4的塊體多達13處，這與測區存在的花崗巖侵入體或變質巖捕虜體有較強的聯系.根據識別結果可進一步解譯出兩條東西走向的界線，南側界線與廣義的嘉黎斷裂剪切帶具有較好的一致性，北側界線可能為不同巖性的接觸帶.另外，一些近南北向展布的異常界線，也與該部位發育的系列南北向小規模斷層具有較好的匹配性.由此可見，即便是針對具備復雜地質地貌環境的目標區，本文所采用的PFD-Net算法仍能夠較好地用于解譯構造跡線和巖性跡線，并且解譯結果具有較好的可信度.

4 結論與討論

本文提出了一種融合多尺度全局注意力機制和密集跳躍連接的DL網絡結構—PFD-Net，并用于位場的邊界識別.該網絡以U-Net為骨干網絡，加入了密集跳躍連接和對應的標準卷積模塊，從位場數據中提取異常體的邊界特征.此外，全局注意力模塊對提取的低高層特征進行分離和融合，作為多維度的輸入，使得網絡更具魯棒性.相比經典的U-Net網絡結構，本文所提網絡結構大大提高了預測精度.

文中通過四組模型試驗來檢測本文所提方法效果：(1)對單一塊體模型所構建的樣本數據集進行訓練，可對類似塊體的組合模型進行邊界識別，表明了本文PFD-Net網絡所具備的泛化性能；(2)相同標簽不同輸入的樣本數據對，會增加網絡結構訓練的復雜度，并可能會導致深部異常體預測精度下降，但不影響整體的識別效果；(3)對垂直方向分布大小差異較大的異常體仍具有一定的識別效果；(4)含噪聲情況下，PFD-Net網絡的識別精度并沒有明顯下降，表明本文所提方法具有一定的穩定性.另外，將文中的方法應用于藏東南某鐵路隧道西段的航空磁測數據的邊界識別，其結果與已知的地質情況基本吻合，表明本文的反演方法具有較好的可信度.

本文所提的位場邊界識別方法實質上是利用CNN方法對地下異常體的水平位置進行分類及回歸預測，相比基于DL的重磁異常及其梯度異常反演，有很大的不同之處在于樣本數據的標簽.物性反演所需樣本數據的標簽是異常體的位置向量，不僅僅具有方向性，且具有物性大小值，其有效標注數據量大.而邊界識別所需樣本數據的標簽是異常體在水平面上投影邊界所組成的向量，但該向量只具有方向性，且標注的數據個數少.此外，邊界識別的樣本數據存在不同輸入數據但標簽數據相同，由于地球物理存在多解性，位場邊界識別方法需具有的普適性在于如何構建一個足量的、多樣性的、代表性的高質量樣本數據集.然而，這勢必會引起網絡復雜度與計算量增加.因此，在下一步工作中，我們將進一步增加樣本的多樣性和改進樣本數據的標注方法來提升復雜地質模型的識別效果，同時在保證實驗準確性的情況下，進一步更好地簡化和壓縮模型，提高網絡訓練的速度和效率，從而使得基于DL的位場數據處理方法更加普適化與智能化.

致謝非常感謝匿名審稿專家對論文提出的寶貴修改建議！