六分量地震觀測在工程勘查中的應(yīng)用試驗

秦林鵬, 王赟*, 張東明, 郭高源, 陳暢, 陳鋒,錢榮毅, 馬振寧, 陳永陽, 汪超

1 中國地質(zhì)大學(xué)(北京)地球物理與信息技術(shù)學(xué)院, 北京 100083 2 北京自動化控制設(shè)備研究所, 北京 100074 3 中國科學(xué)院地球化學(xué)研究所, 貴陽 550081

0 引言

三個平動自由度和三個旋轉(zhuǎn)自由度的組合完整刻畫了空間中任意一點(diǎn)的運(yùn)動狀態(tài)(胡德綏,1989),一般稱為六分量地震觀測.旋轉(zhuǎn)運(yùn)動的直接測量提供了對地震波場的額外獨(dú)立信息,補(bǔ)充了地震波場的空間梯度測量,對近地表面波層析成像等有重要意義(Lee et al., 2009; Schmelzbach et al., 2018; 王赟等,2021).

面波層析成像的關(guān)鍵點(diǎn)之一是頻散曲線的提取.20世紀(jì)90年代以來,MASW方法(Multichannel analysis of surface wave,面波多通道分析方法)(Park et al., 1999)得到廣泛運(yùn)用.對主動源激發(fā)、線性陣列布設(shè)多道臺站記錄的地震數(shù)據(jù),或由背景噪聲互相關(guān)方法(NCFs)得到的虛擬炮集(Bensen et al., 2007),使用一定的數(shù)學(xué)變換方法就可以得到頻散曲線;被動源的空間自相關(guān)法(SPAC)(Aki, 1957)也可以得到頻散曲線,但上述方法均需要多臺站同時記錄地面振動.

Igel等(2005)的研究揭示單臺站記錄的六分量地震數(shù)據(jù)中,長周期(1～150s)Love波信號的水平切向加速度和垂向旋轉(zhuǎn)速率擁有高波形相似性,幅度比能準(zhǔn)確反映Love波相速度,從而可以利用單物理點(diǎn)上觀測的平動加速度與旋轉(zhuǎn)速率實(shí)現(xiàn)橫波速度的一維反演.這一方法的有效性得到多個天然地震實(shí)例的證實(shí)(Cochard et al., 2006; Igel et al., 2007; Kurrle et al., 2010; Yuan et al., 2020; 操玉文等,2022);在微震定位和利用微震實(shí)現(xiàn)面波相速度頻散曲線提取方面也獲得了應(yīng)用(Hadziioannou et al., 2012; Tanimoto et al., 2015; Taylor et al., 2021; Yuan et al., 2021).

Wassermann等(2016)首次給出了利用頻率域單臺站六分量地震數(shù)據(jù)提取面波頻散曲線的流程與誤差估計及質(zhì)控方法.該方法對面波后方位角和相速度進(jìn)行掃描,使用正交回歸的算法進(jìn)行誤差評估;同時引入經(jīng)驗系數(shù)對采集質(zhì)量較差的數(shù)據(jù)進(jìn)行質(zhì)控.但是在主動源六分量地震勘探中,該方法是否可行值得探究.因此,本文針對人工源采集六分量地震數(shù)據(jù)具有高信噪比特點(diǎn),我們避免使用不確定的經(jīng)驗系數(shù)質(zhì)控,嘗試通過簡單的窄帶濾波,在時間域?qū)崿F(xiàn)基階面波頻散曲線的提取.

1 試驗區(qū)概況與地震觀測

如圖1所示,試驗區(qū)位于河北省唐山市開平區(qū)東湖生態(tài)區(qū),城區(qū)環(huán)境復(fù)雜,背景噪聲較強(qiáng),如圖1a所示.試驗區(qū)淺表為古灤河中-晚期沖積平原地貌(胡云壯等, 2014),生態(tài)恢復(fù)工程使得現(xiàn)階段全場地內(nèi)均為人工回填,地形起伏較小,如圖1b所示.由于歷史上試驗區(qū)地下煤礦開采時間長,存在嚴(yán)重的采空塌陷隱患(劉曉燕等, 2013),需要進(jìn)行地球物理探測.在該試驗區(qū)的淺層地質(zhì)勘查中,我們對主動源條件下的單臺站六分量地震與主被動源條件下多道提取面波頻散曲線的方法進(jìn)行了試驗對比.

圖1 (a) 試驗區(qū)位置(紅色方框)及附近功能區(qū)劃分,及測線位置(白色放大圖).白色矩形是試驗區(qū)的放大圖,其中綠線代表測線; (b) 試驗區(qū)地表情況

如圖1a所示,試驗區(qū)內(nèi)容布置了一條近南北向的地震測線.如圖2所示,測線上布設(shè)了24臺SmartsoloIGU-16節(jié)點(diǎn)地震儀,間距2 m,自然頻率為5 Hz;測線開始處放置了1臺光纖旋轉(zhuǎn)地震儀,儀器頻帶為30s～50 Hz.可控震源車在測線不同位置進(jìn)行縱波激發(fā),掃描頻帶為10～100 Hz,掃頻時間為7 s.光纖旋轉(zhuǎn)地震儀及其自身攜帶的加速度計有效記錄了六分量地震數(shù)據(jù),用于后續(xù)單點(diǎn)提取面波頻散曲線.

圖2 地震儀器布局

觀測期間Smartsolo節(jié)點(diǎn)式地震儀設(shè)置采樣頻率為1000 Hz,儀器有效觀測時間350 min;北京自動化控制設(shè)備研究所研制的FOSN-Ⅱ光纖旋轉(zhuǎn)地震儀設(shè)置采樣頻率為1000 Hz,有效觀測時間235 min.

2 原理及方法

2.1 利用六分量地震數(shù)據(jù)提取頻散曲線的原理

假設(shè)已知空間某點(diǎn)r=[x,y,z]T處的平動位移波場u,自由地表旋轉(zhuǎn)位移ωFS可以近似表示為(Cochard et al., 2006; Schmelzbach et al., 2018; 孫麗霞等,2021):

(1)

在二維工程地震勘查中,可認(rèn)為Rayleigh波主要記錄在平動徑向分量R、垂向分量Z、旋轉(zhuǎn)切向分量T上 (圖3),其六分量的波場可表示為(Sollberger et al., 2020):

圖3 Rayleigh波在RTZ坐標(biāo)系下的傳播示意圖

(2a)

(2b)

Rayleigh波的相速度表示為(Ferreira and Igel, 2009; Kurrle et al., 2010; Lin et al., 2011):

(3)

同理Love波的相速度表示為

(4)

時域中上述關(guān)系仍然成立(Igel et al., 2007):

(5)

公式(5)的有效性在天然地震引起的長周期面波信號(Cochard et al., 2006; Igel et al., 2007; Kurrle et al., 2010; Yuan et al., 2020)及較高頻段的微震事件(Hadziioannou et al., 2012; Tanimoto et al., 2015; Taylor et al., 2021; Wassermann et al., 2016; Yuan et al., 2021)中得到證實(shí).

2.2 頻散曲線提取方法

MASW方法自提出以來,因其對體波干擾不敏感、能識別低速層等優(yōu)點(diǎn),被廣泛使用在提取線性陣列主、被動源地震的面波頻散曲線中(Cheng et al., 2016; Park et al., 1999, 2007).常見用于提取多通道數(shù)據(jù)頻散曲線的數(shù)學(xué)方法包括頻率-波數(shù)譜法(f-k法)(Capon, 1969; Lacoss et al., 1969)、相移法(Park et al., 1998),τ-p法(McMechan and Yedlin, 1981)或拉東變換(Luo et al., 2008).

理論研究表明,同時記錄的不同臺站背景噪聲數(shù)據(jù)的互相關(guān)函數(shù)與經(jīng)驗格林函數(shù)只有幅值上的區(qū)別(Snieder, 2004; Weaver and Lobkis, 2001),因而可以獲取面波的頻散特征.Bensen等(2007)給出了時間域的背景噪聲數(shù)據(jù)處理流程,主要步驟包括單臺數(shù)據(jù)預(yù)處理、臺站對的互相關(guān)數(shù)據(jù)計算及疊加、頻散曲線提取.

基于單臺站六分量地震數(shù)據(jù)提取淺層介質(zhì)的面波頻散曲線已被理論和諸多數(shù)值模擬結(jié)果驗證.陳鋒等(2023)從數(shù)值模擬的角度驗證了不同模型下利用單臺站六分量數(shù)據(jù)估算Rayleigh波相速度的可靠性;Wassermann等提出了在頻率域提取面波頻散曲線的方法,但在實(shí)際數(shù)據(jù)應(yīng)用中發(fā)現(xiàn),高階面波的干擾很難避免(Kurrle et al., 2010; Wassermann et al., 2016).特別是在淺層工程勘查中,高階面波的能量強(qiáng),頻帶范圍寬,頻域的算法會導(dǎo)致一定誤差.本文不再論述數(shù)值模擬測試結(jié)果,而直接針對實(shí)際觀測數(shù)據(jù)測試方法的適用性和精度.針對淺層工程勘查中主動源六分量數(shù)據(jù)能量強(qiáng)、信噪比高,二維探測方向唯一的特點(diǎn),我們提出了在頻域進(jìn)行窄帶濾波,并在時域進(jìn)行相速度的計算及誤差估計,具體技術(shù)流程如下.

(1)數(shù)據(jù)預(yù)處理

去儀器響應(yīng)、去均值、去趨勢.去儀器響應(yīng)的操作必不可少,地震儀記錄的數(shù)據(jù)等于震源項、路徑項與儀器響應(yīng)的卷積結(jié)果.去儀器響應(yīng)操作對應(yīng)著地震儀記錄數(shù)據(jù)與儀器響應(yīng)的反卷積,操作后得到臺站記錄的地面運(yùn)動的絕對振幅值.地震六分量數(shù)據(jù)求解后方位角時需要用到絕對振幅信息.

(2)后方位角校正

根據(jù)實(shí)際后方位角或計算出的后方位角得到圖3坐標(biāo)系下的六分量數(shù)據(jù),計算方法可以采用六分量極化分析方法(Sollberger et al., 2020).

(3)相速度及其誤差估計

受噪聲影響,第二步驟得到的數(shù)據(jù)并不能直接代入公式(5)中.設(shè)定合理長度的時窗,通過平動垂向分量與旋轉(zhuǎn)切向分量(Rayleigh波)或平動切向分量與垂向旋轉(zhuǎn)分量(Love波)的相關(guān)性掃描,需要利用相關(guān)系數(shù)判定數(shù)據(jù)是否適用于相速度計算.掃描前,合適的帶通濾波有助于在時域上初略分離不同階的面波信號.誤差估計使用標(biāo)準(zhǔn)偏差進(jìn)行估計,代表帶通濾波中心頻率相速度值的波動性.

3 數(shù)據(jù)分析與試驗結(jié)果

3.1 數(shù)據(jù)分析

如圖4顯示了一次震源車縱波激發(fā)的結(jié)果, Smartsolo儀器記錄顯示激發(fā)能量充分;由于吸收衰減,遠(yuǎn)炮檢距接收到的面波信號優(yōu)勢頻帶為20～50 Hz,后文數(shù)據(jù)不做特殊說明則均為此次激發(fā)結(jié)果.

圖4 Smartsolo記錄的波形圖(a)及其連續(xù)小波變換(b1,炮檢距0 m,b2炮檢距46 m)

3.2 傳統(tǒng)方法提取頻散曲線

首先,以零偏移距臺站為參考臺站,使用互相關(guān)方法進(jìn)行震源車掃頻結(jié)果的處理,其結(jié)果如圖5a所示,顯然面波發(fā)育.然后利用商業(yè)軟件SurfWave中的MASW方法進(jìn)行分析,相移法得到頻散曲線如圖5b所示,基階面波頻散曲線能量強(qiáng)、連續(xù),信噪比較高.

圖5 (a), (b)為互相關(guān)后的多道地震Z分量記錄及面波頻散曲線;(c), (d)為背景噪聲互相關(guān)方法得到的虛擬炮集及面波頻散曲線

對Smartsolo測線記錄的350分鐘被動源記錄進(jìn)行背景噪聲互相關(guān)得到虛擬炮集及面波頻散曲線,如圖5c和d所示.顯然,該頻散曲線與主動源Z分量記錄中提取的面波頻散曲線(圖5a和5b)具有很好的相似性.

需要說明的是,在被動源噪聲頻散譜提取中,時域歸一化使用滑動絕對平均法(RAMN),且使用Welch法(Welch, 1967)提高面波的信噪比.Welch法消除了背景噪聲可能對選擇的時間序列起止窗的依賴性,減小了單個高幅度瞬時信號對整體平均值的影響(Lawrence and Prieto, 2011; Prieto et al., 2011).圖6a為不同重疊率及時窗數(shù)量條件下Welch法的改進(jìn)效果.保持時窗數(shù)目不變,重疊率Overlap=50%時,數(shù)據(jù)信噪比更高;保持Overlap不變時,數(shù)據(jù)信噪比隨時窗數(shù)目增加明顯增大;但在時窗數(shù)目從20變?yōu)?0時,信噪比增量不大,因此可以判斷為互相關(guān)結(jié)果逐漸收斂.我們最后使用Overlap=50%,時窗數(shù)目為50作為Welch法的參數(shù),以增強(qiáng)虛擬炮集的面波信噪比.顯然,Welch法使得小炮檢距的數(shù)據(jù)有過大的信噪比,這實(shí)際上是野外試驗期間多種類型的主動源實(shí)驗(錘擊、震源車)以及復(fù)雜微震信號導(dǎo)致的(行車、人員走動等).因此我們使用閾值法將大于噪聲水平的數(shù)據(jù)置零,以去除儀器中大量記錄的主動源及微震信號.如圖6a顯示,去除主動源信號后,近炮檢距的信號信噪比顯著下降.

圖6 (a) 利用Welch法提高互相關(guān)法提取的Rayleigh波信噪比; (b) 閾值法對虛擬炮集數(shù)據(jù)信噪比的影響

3.3 利用單點(diǎn)六分量地震數(shù)據(jù)提取頻散曲線

本次試驗FOSN-Ⅱ光纖旋轉(zhuǎn)地震儀記錄到了高信噪比的六分量地震數(shù)據(jù);采集數(shù)據(jù)為旋轉(zhuǎn)速率,單位rad·s-1,無需進(jìn)行去儀器響應(yīng)處理.由于儀器通頻帶為30 s～50 Hz,因此進(jìn)行50Hz的低通濾波.數(shù)據(jù)預(yù)處理及后方位角矯正后的六分量數(shù)據(jù)波型及其連續(xù)小波變換結(jié)果顯示(圖7):震源車縱波激發(fā)結(jié)果的主要能量集中于平動徑向、垂向和旋轉(zhuǎn)切向的方向上.根據(jù)六分量地震數(shù)據(jù)的對應(yīng)關(guān)系(圖3),顯然該能量主要是Rayleigh波模式.平動切向分量及旋轉(zhuǎn)徑向分量中出現(xiàn)的信號能量是由于儀器放置傾斜,地下介質(zhì)存在非均勻體干擾,或后方位角不完全準(zhǔn)確導(dǎo)致,可以忽略.

圖7 六分量數(shù)據(jù)波形及其連續(xù)小波變換結(jié)果

在天然地震領(lǐng)域,震源位置是未知的,即后方位角是需要求取的參數(shù).基于六分量極化分析可計算后方位角,實(shí)際應(yīng)用已證實(shí)該方法的有效性 (Yuan et al., 2021).假設(shè)本主動源試驗中,震源位置也未知,限定Rayleigh波速在250～600 m·s-1變化,進(jìn)行10～40 Hz頻率掃描,極化分析得到的結(jié)果如圖8所示,計算出的后方位角接近實(shí)際后方位角(偏差為4°),精度可接受.

圖8 六分量極化分析方法估計的后方位角

在時間域通過不同頻段掃描提取頻散曲線需要對時窗長度、不同頻段進(jìn)行試驗.試驗結(jié)果發(fā)現(xiàn)時窗長度選用0.01 s,即所求頻散曲線最短周期的半周期長度進(jìn)行掃描效果最好.頻率窗的大小影響所求頻散點(diǎn)的數(shù)目,窗越短,最終得到的頻散點(diǎn)越多.但當(dāng)頻率窗過小時,一個頻率窗內(nèi)優(yōu)質(zhì)數(shù)據(jù)點(diǎn)的數(shù)量減少,對所求相速度值的誤差估計不利.優(yōu)質(zhì)數(shù)據(jù)點(diǎn)的評判標(biāo)準(zhǔn)為上述時間窗和頻率窗限定下的平動垂向分量和旋轉(zhuǎn)切向分量(Rayleigh波)或平動切向分量和垂向旋轉(zhuǎn)分量(Love波)的相關(guān)系數(shù)是否大于0.75.相速度值使用優(yōu)質(zhì)數(shù)據(jù)點(diǎn)的值代入公式(5)進(jìn)行求解,并用均值代表最終該頻帶范圍中心頻率的相速度值.

如圖9是18 Hz中心頻率下(帶通濾波16～22 Hz)利用掃描法得到的相速度結(jié)果.在Rayleigh波信號范圍時間窗內(nèi)計算垂向加速度和切向旋轉(zhuǎn)速度的相反數(shù)的互相關(guān)系數(shù)大于0.75的相速度,通過平均值得到最終該頻率的相速度值.最終18 Hz中心頻率對應(yīng)的橫波相速度值為282.7 m·s-1,與回填土的速度結(jié)構(gòu)是吻合的(劉培玄等,2019).

圖9 18 Hz中心頻率下的Rayleigh波相速度求取

對不同的時窗進(jìn)行頻域掃描,得到完整的頻散曲線如圖10所示.其中誤差棒代表該點(diǎn)中心頻率對應(yīng)相速度的標(biāo)準(zhǔn)偏差,該值反映了該中心頻率所求相速度值的誤差情況.從圖中可以看出,單點(diǎn)六分量方法所求頻散曲線在30 Hz以下頻段和常規(guī)方法(主動源和被動源的多道分析法)基本一致;但30 Hz以上頻段,求解出的相速度值逐漸偏離.因此,單臺站記錄的地震六分量信號在時間域通過不同時窗分頻濾波的方法,可以提取低頻基階面波頻散曲線,高頻段的影響因素需要進(jìn)一步研究討論.

圖10 不同方法得到的頻散曲線

4 討論

無論主被動源數(shù)據(jù),面波多道分析方法(MASW)都是提取頻散曲線的有力手段.工程勘查領(lǐng)域,被動源數(shù)據(jù)的處理除傳統(tǒng)流程外,Welch法的加入大大提高了被動源數(shù)據(jù)信噪比,為短時觀測的背景噪聲數(shù)據(jù)提取頻散曲線提供了可能.

利用單臺站記錄的六分量地震數(shù)據(jù)提取頻散曲線的原理十分簡潔,但高階面波的干擾消除是重點(diǎn)和難點(diǎn)(Igel et al., 2005; Wassermann et al., 2016; Yuan et al., 2020).在本實(shí)驗中,由單臺站六分量數(shù)據(jù)提取的高頻段的頻散曲線和傳統(tǒng)主被動源方法提取的頻散曲線有明顯差異.這可能是高階面波受體波影響,波場分離不干凈的結(jié)果.由于各階頻散曲線在同一頻率點(diǎn)處的相速度不一,工程地震關(guān)注的地震頻段內(nèi),在時間域進(jìn)行窄帶濾波用以分離提取不同階面波頻散曲線需要進(jìn)一步研究.

在工程地震勘查領(lǐng)域,僅需單臺站儀器便能獲取面波頻散曲線,從而實(shí)現(xiàn)一維的速度結(jié)構(gòu)反演,將會帶來人力、財力、時間等成本的巨大節(jié)約,尤其在城市環(huán)境中應(yīng)用具有明顯的優(yōu)勢.隨著國產(chǎn)光纖旋轉(zhuǎn)地震儀的迅速發(fā)展,其精度及穩(wěn)定性不斷提高(Cao et al., 2021),使用單物理點(diǎn)六分量地震儀提高工程地震勘探效率值得期待.

5 結(jié)論

通過本次的主被動源六分量地震聯(lián)合觀測試驗,可獲得如下結(jié)論.

(1)利用面波多道分析方法 (MASW) 方法,從震源車縱波激發(fā)、Smartsolo單分量節(jié)點(diǎn)地震儀測線記錄中提取主動源面波頻散曲線,和利用背景噪聲虛擬炮集獲得5～50 Hz頻段的Rayleigh波頻散曲線,兩者結(jié)果一致性很好;其中Welch法的使用大大提高了虛擬炮集的面波信噪比.

(2)利用單物理點(diǎn)六分量地震數(shù)據(jù)提取面波頻散曲線,在后方位角缺失的情況下,利用六分量極化分析方法求解后方位角也是可靠的.

(3)唐山城區(qū)觀測試驗與分析說明使用單物理點(diǎn)六分量地震記錄提取頻散曲線同其他兩種傳統(tǒng)方法在低頻5～30 Hz段內(nèi)一致性很好,可用于基階面波提取和反演.

致謝本文部分圖件采用Obspy繪制(https:∥docs.obspy.org/).