不同數(shù)據(jù)處理方法在高密度電法水庫滲漏成像中的應(yīng)用

李 進(jìn)，王小明，譚 磊

(1.南京市水利規(guī)劃設(shè)計(jì)院股份有限公司，江蘇 南京 210022；2.蘇交科集團(tuán)股份有限公司，江蘇 南京 210000；3.浙江省水利河口研究院(浙江省海洋規(guī)劃設(shè)計(jì)研究院)，浙江 杭州 310020；4.浙江省水利防災(zāi)減災(zāi)重點(diǎn)實(shí)驗(yàn)室，浙江 杭州 310020)

1 引 言

地球物理勘探是通過采集到的多物理場數(shù)據(jù)體來認(rèn)識(shí)地下空間的地層構(gòu)造、礦產(chǎn)賦存、災(zāi)害分布以及工程地質(zhì)環(huán)境等[1-4]，具有高效、無損、透視、大數(shù)據(jù)的特點(diǎn)，在宏觀調(diào)查和精細(xì)勘探中都有著突出的優(yōu)勢，與鉆探手段相比，二維、三維的探測成果能連續(xù)的反映出勘探區(qū)的地質(zhì)情況，極大地提高了勘探密度[5,6]。但是，受限于野外工作環(huán)境、地形條件以及時(shí)間成本等因素的制約，在權(quán)衡精度與效率利弊的情況下，提交的物探圖像往往是適度特征點(diǎn)數(shù)據(jù)經(jīng)網(wǎng)格插值處理后的成果。目前，不規(guī)則離散數(shù)據(jù)網(wǎng)格化運(yùn)算的插值技術(shù)較多，如美國Golden Software公司的商業(yè)Surfer系列軟件中就有12種不同類型的插值方式[7]，豐富的插值功能為數(shù)據(jù)顯示提供了便利，但具體到每種插值算法在數(shù)據(jù)可視化表達(dá)方面也存在一定的差異性。近年來，一些物探學(xué)者也針對(duì)不用插值算法的應(yīng)用及特點(diǎn)方面作出許多有益的工作。吳衛(wèi)國[8]利用Surfer軟件的克里金法對(duì)水系沉積物數(shù)據(jù)進(jìn)行數(shù)據(jù)網(wǎng)格化，并結(jié)合白化功能減少離散數(shù)據(jù)的外推；劉曉霞等[9]推導(dǎo)了Kriging球面應(yīng)變計(jì)算公式，交叉檢驗(yàn)及殘差分析結(jié)果表明，速度場濾波和網(wǎng)格插值是可效的；姚文等[10]以某磁鐵礦區(qū)高精度磁測數(shù)據(jù)為例，在比較了5種網(wǎng)格化平面圖的基礎(chǔ)上，認(rèn)識(shí)到在保證網(wǎng)格化的精度和等值線圓滑的前提下，克里金插值法和最小曲率法網(wǎng)格化效果最好；王咸彬等[11]在總結(jié)了常見散亂數(shù)據(jù)插值方法的基礎(chǔ)上，采用克里金插值方法把實(shí)測測井?dāng)?shù)據(jù)插值成反演的三維速度場初始模型；郭崇光[12]采用Surfer軟件中的克里格方法對(duì)電測深數(shù)據(jù)進(jìn)行網(wǎng)格化，并指出具體插值中的注意點(diǎn)。

高可靠度的原始數(shù)據(jù)是地球物理探測數(shù)據(jù)網(wǎng)格化成像的基礎(chǔ)，然而在觀測過程中顯然不能排除外來噪聲的干擾[13]，因此有必要采用一定的手段對(duì)網(wǎng)格數(shù)據(jù)體進(jìn)行濾波降噪處理。李文杰[14]運(yùn)用二維自動(dòng)調(diào)平處理了頻率域航空電磁數(shù)據(jù)的零漂誤差，提高視電阻率和平面剖面圖質(zhì)量；余金煌等[15]采用小子域?yàn)V波技術(shù)對(duì)高密度電法探測的拋石層圖像進(jìn)行增強(qiáng)處理，提高了圖像精度和可讀性；黃基文[16]根據(jù)基階模態(tài)與高階模態(tài)面波之間視速度存在差異，利用二維數(shù)字濾波法進(jìn)行分離；王瑞雪等[17]把均值濾波、中值濾波以及小波變換對(duì)測井電導(dǎo)率數(shù)據(jù)進(jìn)行濾波處理。近年來，高密度電法在水庫大壩滲漏隱患探測中得到廣泛應(yīng)用[18]，但較少研究不同插值方法對(duì)電阻率圖像上隱患的優(yōu)化識(shí)別，并且涉及到具體濾波方法的差異性以及相應(yīng)參數(shù)設(shè)置等內(nèi)容還需進(jìn)一步探討。為此，本文以水庫大壩蟻穴滲漏高密度電法探測數(shù)據(jù)為試驗(yàn)對(duì)象，比較了多種網(wǎng)格插值方法圖像的差異性，并對(duì)自然鄰點(diǎn)法插值圖像中的濾波方法及插值結(jié)果進(jìn)行分析，從而得到適合于土石壩電阻率可視化表達(dá)的成像技術(shù)，以此提高隱患的識(shí)別能力及效果。

2 試驗(yàn)工程概況

浙江某水庫始建于1957年10月，并于1958年7月建成蓄水，壩址以上集雨面積0.301 km2，總庫容為16.395萬m3，正常庫容為12.149萬m3，灌溉面積720畝，是一座以灌溉為主的小(2)型水庫。大壩壩型為類均質(zhì)壩，壩頂高程為55.70 m，壩頂長為103.5 m，最大壩高7.58 m，壩頂寬3.0 m，壩頂鋪廣場磚。

該水庫已服役運(yùn)行60多年，大壩防滲結(jié)構(gòu)存在不同程度的老化，盡管多次對(duì)大壩進(jìn)行加固處理，但在極端天氣作用下，于2019年在大壩左壩段出現(xiàn)滲漏、塌陷的突發(fā)險(xiǎn)情，為盡早排查出大壩滲漏的原因以期保障水庫的安全運(yùn)行，采用高密度電法對(duì)水庫大壩滲漏原因進(jìn)行探測，并對(duì)異常區(qū)進(jìn)行防滲處理，在灌漿過程中大壩迎水坡、背水坡均出現(xiàn)局部的冒漿現(xiàn)象，與探測得到的滲漏通道走向基本吻合。以此探測數(shù)據(jù)為基礎(chǔ)，從不同插值方法、濾波分析等角度，提升對(duì)探測成果的可識(shí)別性。

圖2是地下記錄視電阻率值的散點(diǎn)分布圖，記錄點(diǎn)數(shù)據(jù)共有946個(gè)，最小電阻率值為11.98 Ω·m，最大電阻率值339.15 Ω·m，均值95.93 Ω·m，中值86.39 Ω·m。

圖1 試驗(yàn)現(xiàn)場示意Fig.1 Schematic diagram of the test site

圖2 原始數(shù)據(jù)記錄點(diǎn)分布Fig.2 Distribution of original data recording points

3 電阻率數(shù)據(jù)插值及濾波處理方法

3.1 常用數(shù)據(jù)插值方法的簡介

數(shù)據(jù)的網(wǎng)格化處理工作是空間離散測點(diǎn)原始數(shù)據(jù)可視化成像的基礎(chǔ)，在高密度電法勘探過程中，水平方向上各個(gè)電極測點(diǎn)的最小間距是保持固定的，垂直方向上相鄰深度點(diǎn)的距離又與水平長度有關(guān)，并且隨著深度不斷的加大，水平層間的測點(diǎn)呈等差數(shù)列不斷減少，總體表現(xiàn)為矩陣排列的形式。因此，為把縱、橫向上數(shù)據(jù)點(diǎn)之間的變化趨勢反映出來，對(duì)原始數(shù)據(jù)采用網(wǎng)格化插值是重要的環(huán)節(jié)。在地球物理數(shù)據(jù)網(wǎng)格插值方面應(yīng)用較多的方法[19]，通常有反距離加權(quán)插值法(Inverse Distance to a Power)、克里金插值法(Kriging)、最小曲率插值法(Minimum Curvature)、改進(jìn)謝別德插值法(Modified Shepard’s Method)、自然鄰點(diǎn)插值法(Natural Neighbor)、最近鄰點(diǎn)插值法(Nearest Neighbor)、多項(xiàng)式回歸插值法(Polynomial Regression)、徑向基函數(shù)插值法(Radial Basis Function)、線性插值的三角剖分插值法(Triangulation with Linear Interpolation)、移動(dòng)平均插值法(Moving Average)、數(shù)據(jù)度量插值法(Data Metrics)、局部多項(xiàng)式插值法(Local Polynomial)12種方法。同一組數(shù)據(jù)采用不同的插值方法，圖像中有用信息的顯示程度也存在差異，因此有必要通過比較選擇最優(yōu)的插值方法[20]，從而指導(dǎo)工程應(yīng)用中的數(shù)據(jù)處理。

3.2 常用網(wǎng)格濾波方法

數(shù)據(jù)網(wǎng)格化處理有效預(yù)測了測點(diǎn)之間的數(shù)據(jù)體，但不能消除原始信號(hào)中的噪聲信息。高密度電法數(shù)據(jù)噪聲來源復(fù)雜，包括天然電場、游散電流、接地不良、地形起伏以及儀器雜音等[21]，外來噪聲往往會(huì)降低數(shù)據(jù)信號(hào)的信噪比，從而降低圖像信息的真實(shí)度以及圖像自動(dòng)識(shí)別的能力。因此，開展網(wǎng)格化數(shù)據(jù)體進(jìn)行濾波工作對(duì)水庫大壩隱患的識(shí)別十分重要，但如何選擇適合于電阻率成像的濾波方法成為圖像可視化的關(guān)鍵。常用的濾波器分為線性濾波器和非線性濾波器，并且線性濾波器中有滑動(dòng)平均濾波、距離加權(quán)平均濾波、反距離濾波、高斯低通法等。

4 結(jié)果分析

4.1 不同插值方法的網(wǎng)格化結(jié)果及分析

圖3是12種插值方法網(wǎng)格化后的電阻率圖像，所有圖像采用相同的色階歸一化。從圖像上可以直觀看出，不同的插值方法得到的電阻率值圖像存在較大的差異，主要表現(xiàn)在經(jīng)插值后的電阻率范圍、分布以及圖像的平滑度等。其中，多項(xiàng)式回歸插值法、移動(dòng)平均插值法、數(shù)據(jù)度量插值法以及局部多項(xiàng)式插值法平滑度較大，不能有效揭示出大壩的電阻率值分布，不適合運(yùn)用在水庫滲漏探測的高密度電法數(shù)據(jù)處理中；反距離插值圖像上出現(xiàn)較多的“牛角眼特征”，不利于對(duì)隱患的識(shí)別；改進(jìn)的謝別德插值對(duì)圖像的邊界處理較差，局部出現(xiàn)負(fù)值，不符合電阻率的物理意義；最近鄰點(diǎn)法插值在圖像上形成連續(xù)的網(wǎng)格，網(wǎng)格之間的連續(xù)性較差；最小曲率插值的視電阻率范圍為8.29～518.94 Ω·m(表1)，比實(shí)際測量最大值339.15 Ω·m高約34.6%；線性插值的三角剖分插值中的電阻率值也出現(xiàn)負(fù)值?？死锝鸱?、自然鄰點(diǎn)法以及徑向基函數(shù)插值法圖像電阻率等值線平滑度較好，插值后的電阻率值與原始值范圍基本吻合，在測線上18～24 m段、核心深度3 m的高阻閉合區(qū)域與隱患位置基本一致，總體上能反映出大壩地下結(jié)構(gòu)的分布特征，也表明高密度電法在探測土石壩隱患中結(jié)果可靠度高。但是，自然鄰點(diǎn)法插值外延性較差，克里金法、徑向基函數(shù)插值能向無值區(qū)外推出電阻率的變化趨勢，但計(jì)算時(shí)間較長。因此，水庫滲漏高密度電法原始數(shù)據(jù)的插值優(yōu)選采用自然鄰點(diǎn)法，其次是克里金插值法、徑向基函數(shù)插值法，在工程應(yīng)用中采用合理的插值算法能更加真實(shí)地反映出大壩的電阻率變化特征。

圖3 不同插值方法的電阻率Fig.3 Resistivity plots of different interpolation methods

表1 不同插值方法的電阻率范圍

4.2 不同濾波方法的等值線結(jié)果及分析

地下空間的地電場信號(hào)頻率較低，因此應(yīng)采用低通數(shù)字濾波器對(duì)電阻率圖像進(jìn)行噪聲處理。根據(jù)多種插值算法網(wǎng)格化成果的比較，采用自然鄰點(diǎn)插值的電阻率值進(jìn)行濾波比較分析。圖4(a)是不進(jìn)行濾波處理的等值線分布圖，在圖像上存在局部極值的小圈閉，并且線條平滑度較差；經(jīng)濾波處理后的圖像，圖中的噪聲信號(hào)明顯得到降低，滑動(dòng)平均、距離加權(quán)法的濾波效果相對(duì)更好。

圖4 不同濾波方法的等值線Fig.4 Contours of different filtering methods

從濾波后的電阻率值范圍(表2)上來看，濾波之前的電阻率范圍為23.01～336.22 Ω·m，經(jīng)濾波之后高、低極值得到削減，其中滑動(dòng)平均濾波電阻率范圍相對(duì)更收斂，有效降低了噪聲對(duì)數(shù)據(jù)的影響。

表2 不同濾波方法的電阻率范圍

4.3 不同過濾器尺寸的濾波結(jié)果及分析

數(shù)據(jù)濾波的效果不僅與濾波器的選擇有關(guān)，而且過濾器尺寸的大小也會(huì)對(duì)濾波效果產(chǎn)生影響。圖5是在垂向-3 m的橫向上以及水平24 m縱向上的電阻率變化圖。由圖5可知，過濾器尺寸越小，濾波效果越差，在曲線上的極值點(diǎn)處變化幅度越大，但隨著尺寸的增大，曲線的平滑度增加，從而導(dǎo)致相鄰異常區(qū)的邊界變得模糊。因此，在對(duì)水庫滲漏的高密度電法數(shù)據(jù)采用滑動(dòng)濾波時(shí)，過濾器的尺寸7×7。

圖5 不同過濾器尺寸的電阻率變化Fig.5 Resistivity variations of different filter sizes

4.4 不同迭代次數(shù)的濾波結(jié)果及分析

濾波的迭代次數(shù)大，則可以提高等值線的平滑程度，但過于采用多次濾波計(jì)算也可能降低對(duì)異常區(qū)的識(shí)別。圖6是在垂向-3 m的橫向上以及水平24 m縱向上的電阻率變化圖。從圖像上可以看出，當(dāng)?shù)螖?shù)為1時(shí)，電阻率曲線在極值點(diǎn)變化強(qiáng)烈，同時(shí)在曲線上有所偏離；當(dāng)?shù)螖?shù)為7時(shí)，電阻率曲線明顯過于圓滑，在拐點(diǎn)位置電阻率的幅度過小，對(duì)異常區(qū)邊界的確定帶來不利影響。綜上所述，在對(duì)水庫滲漏的高密度電法數(shù)據(jù)采用滑動(dòng)濾波時(shí)，迭代次數(shù)采用4次。

圖6 不同迭代次數(shù)的電阻率變化Fig.6 Resistivity variation of different number of iterations

5 結(jié) 論

1)高密度電法數(shù)據(jù)經(jīng)網(wǎng)格化插值處理后能夠把測量點(diǎn)可視化顯示出來，有利于評(píng)價(jià)地質(zhì)體的電阻率分布，建議在網(wǎng)格化插值中優(yōu)先選用自然鄰點(diǎn)法，其次是克里金插值法、徑向基函數(shù)插值法；

2)不同濾波方法、過濾器尺寸以及迭代次數(shù)對(duì)網(wǎng)格化數(shù)據(jù)的去噪效果存在差異，水庫大壩隱患的電阻率濾波優(yōu)先選用滑動(dòng)平均法，過濾器尺寸的行列均為7，迭代次數(shù)為4次；

3)水庫大壩隱患類型復(fù)雜，電阻率等值線形態(tài)及大小也存在差異，在具體使用時(shí)，應(yīng)對(duì)不同電阻率分布特征的可視化方法進(jìn)行綜合比較，從而提高對(duì)不良隱患的識(shí)別能力。