王 勃，曾林峰，張 衍，劉盛東，章 俊，陳泓云

(1.中國礦業(yè)大學(xué) 深部巖土力學(xué)與地下工程國家重點實驗室，江蘇 徐州 221116；2.中國礦業(yè)大學(xué) 資源與地球科學(xué)學(xué)院，江蘇 徐州 221116；3.中國礦業(yè)大學(xué) 力學(xué)與土木工程學(xué)院，江蘇 徐州 221116；4.深地科學(xué)與工程云龍湖實驗室，江蘇 徐州 221116)

“雙碳”目標(biāo)下我國能源資源稟賦特征及當(dāng)前復(fù)雜國際能源形勢，決定了煤炭在能源中的基礎(chǔ)和兜底保障作用。煤炭綠色開采、智能精準(zhǔn)開采等對煤礦安全高效開采地質(zhì)保障系統(tǒng)提出了更高的要求，礦井地質(zhì)透明化是當(dāng)前礦井地質(zhì)保障系統(tǒng)發(fā)展的努力方向。礦井地震勘探具有不受復(fù)雜地面條件影響、距離目標(biāo)體近、地震波能量和高頻成分衰減少、分辨率高、探采對比易驗證等優(yōu)點，進而被廣泛運用于地質(zhì)精細(xì)探測及透明地質(zhì)建模。

不同于地面半空間地震勘探，礦井地震屬于全空間條件下勘探。震源在井下激發(fā)后，地震波向四周傳播，檢波點會接收到全方位的地震信號。單分量信號難以確定地震波的方向，雙分量地震信號僅能算出二維平面上地震波傳播方向，三維空間條件下準(zhǔn)確判定地震波傳播方向則依賴三分量信號的矢量特性，并通過極化分析實現(xiàn)。極化分析在地震信號處理中主要用于波場分離，根據(jù)各類地震波的極化屬性差異設(shè)置濾波器來分離提取有效波。極化濾波方法在天然地震、石油地震勘探領(lǐng)域研究較多，DIALLO 等對多分量地震數(shù)據(jù)進行小波變換，將自適應(yīng)瞬時極化濾波引入時頻域，進行面波壓制及轉(zhuǎn)換波分離。KULESH 等提出了基于自適應(yīng)協(xié)方差矩陣的小波域時頻極化分析方法，在時頻域?qū)崿F(xiàn)了波場分離。PINNEGAR采用S變換在時頻域計算三分量地震信號極化參數(shù)，設(shè)置極化率濾波器壓制極化橢圓性地震波。程冰潔等在小波域進行能量分類約束從而實現(xiàn)極化濾波。

在煤炭領(lǐng)域，張平松等在時間域構(gòu)建協(xié)方差矩陣，采用極化濾波方法實現(xiàn)了巷道超前探測多波有效分離。王勃利用三分量信號的矢量特征，提出了一種全空間條件下集波場分離、偏移成像于一體的極化偏移方法。胡澤安等在時間域?qū)ΦV井槽波地震數(shù)據(jù)進行了極化分析，實現(xiàn)了波場分離與噪聲壓制。金丹等利用槽波信號與干擾波在偏振度上的差異，改進頻率域的極化濾波權(quán)函數(shù)，將頻率域極化濾波用于槽波記錄的噪聲壓制。馮磊等對二分量槽波數(shù)據(jù)進行S變換，采用時頻域自適應(yīng)協(xié)方差矩陣極化濾波方法實現(xiàn)了兩類槽波分離。劉盛東等利用三分量地震記錄，通過時窗自適應(yīng)的極化分析方法獲取極化率，提取了縱橫波及勒夫型槽波等線性極化波。

上述方法主要用于地面半空間或煤礦井下全空間的極化濾波，針對礦井全空間條件下三分量地震波傳播方向研究較少，特別在礦井地震近場勘探多類型波場混疊條件下。筆者提出了一種礦井全空間三分量地震波時頻域極化分析方法，首先通過時間域混疊合成信號驗證地震波優(yōu)勢方位角、傾角準(zhǔn)確度，然后開展全空間條件下中心激發(fā)-全方位接收的三維三分量數(shù)值模擬實驗對比縱波、橫波定向精度，最后通過礦井地震勘探常用的反射波現(xiàn)場試驗驗證極化參數(shù)求解的可靠性。

1 全空間地震波時頻域極化分析方法

礦井地震屬于全空間條件下近場勘探，多類地震波混疊，需從時間、頻率2個維度聯(lián)合計算分析。為準(zhǔn)確計算地震波傳播方向，筆者首先構(gòu)建基于基本小波函數(shù)形態(tài)不固定的廣義S變換的時頻域復(fù)協(xié)方差矩陣，然后求解三分量地震信號的優(yōu)勢能量方向。同時，為進一步驗證方向計算的準(zhǔn)確性，開展已知震源點聚焦定位對比研究，具體理論方法如下。

1.1 基于廣義S變換的時頻域復(fù)協(xié)方差矩陣構(gòu)建

設(shè)時間域信號()的傅里葉變換為

(1)

式中，j為虛數(shù)單位；為頻率；為時間。

對時間序列()乘以一個窗函數(shù)()得

(2)

設(shè)()為歸一化的高斯窗，且利用參數(shù)控制其時窗寬度，利用參數(shù)控制其時窗位置，則有

(3)

將式(3)代入式(2)，可得到對時間序列()在時刻加高斯窗的譜為

(4)

將窗寬控制參數(shù)設(shè)置成與頻率成反比的關(guān)系，以此讓高斯窗的寬度自適應(yīng)于頻率，則有

(5)

和為控制時窗寬度變化的2個參數(shù)，進一步獲得廣義S變換的表達式為

(6)

通過廣義S變換計算時頻譜，再利用Hilbert變換構(gòu)建時頻譜的解析信號。在時頻域內(nèi)時刻，頻率處的復(fù)協(xié)方差矩陣(,)可描述為

(7)

矩陣中各元素定義為

(8)

式中，|SC(,)|和|SC(,)|為2組解析信號的瞬時振幅；和為瞬時頻率；arg為瞬時相位；(,)為均值，其定義為

(,=,,)

(9)

式中，R(·)為復(fù)數(shù)的實部，函數(shù)sin()為辛格函數(shù)，其定義為

(10)

(,)為自適應(yīng)時窗長度，其定義為

(11)

式中，為整數(shù)，是一個經(jīng)驗參數(shù)，用于刻畫不同極化屬性，取較大值時可刻畫三維復(fù)雜極化屬性,一般取1或2。

當(dāng)3個分量時頻譜中對應(yīng)時頻點的瞬時頻率相等時，(,)=(,)=(,)，(,)可簡化為2π(,)，復(fù)協(xié)方差矩陣可以化簡為

(,)=|SC(,)||SC(,)|×

cos[arg SC(,)-arg SC(,)]=

(12)

式中，為復(fù)共軛。

1.2 極化參數(shù)計算

地震波矢量特征可以用極化參數(shù)描述，這些參數(shù)可以通過在時頻域中求取復(fù)協(xié)方差矩陣的最大特征值及其對應(yīng)的歸一化特征向量(,,)獲得。

優(yōu)勢極化方位角計算公式為

()=arctan[R()R()]

(13)

其中-90°≤()≤90°，為極化主軸在面的投影與軸的夾角，當(dāng)極化主軸偏向軸正方向時，()>0°；當(dāng)極化主軸偏向軸負(fù)方向時，()<0°。

優(yōu)勢極化傾角計算公式為

(14)

其中，-90°≤()≤90°，為極化主軸與面的夾角，當(dāng)極化主軸偏向軸正方向時，()>0°；當(dāng)極化主軸偏向軸負(fù)方向時，()<0°。

1.3 聚焦定位方法

與微震定位不同，地震勘探震源點位置已知，為了驗證上述極化參數(shù)準(zhǔn)確性，利用檢波點位置及其極化傾角、方位角反向求解震源點位置，其計算方法如下：設(shè)檢波點的空間坐標(biāo)和主極化軸歸一化矢量分別為=(,,)和=(,,)，設(shè)空間任意一點=(,,)。向量與向量的向量積代表向量和向量共起點的情況下所構(gòu)成平行四邊形的面積。對該向量積除以空間矢量的?？色@得點到空間矢量的最短距離：

(15)

(16)

求取的最小值所對應(yīng)的(,,)為震源點在三維空間內(nèi)的位置。

2 三維數(shù)值模擬及算法驗證

2.1 合成信號及極化特征分析

為了驗證極化分析對三分量信號在時頻域上區(qū)分多類型波和計算極化參數(shù)的準(zhǔn)確性，開展正弦波合成信號分析實驗。合成信號由A,B,C,D,E,F共6個信號組成，合成信號具體參數(shù)見表1。其中B，C信號在時間域混疊，E，F(xiàn)信號在時間域混疊，在時間域可分為4段信號，如圖1所示。其中每組信號均為501個采樣點，采樣頻率為2.5 kHz，對合成信號進行時頻域極化分析，結(jié)果如圖2所示。時頻譜上可清晰區(qū)分A，B，C，D，E，F(xiàn)信號，圖2(a)顯示的6個信號的方位角與理論方位角一致，正、負(fù)方位角均無誤差；圖2(b)顯示的6個信號的傾角與理論傾角一致，正、負(fù)傾角均無誤差。

表1 合成信號基本參數(shù)

圖1 合成三分量信號

圖2 合成三分量信號方位角及傾角分布

2.2 正演模擬信號極化特征分析及聚焦定位驗證

為研究全空間條件下地震波矢量特征，建立三維數(shù)值模型。模型中心為震源點，在震源點的左側(cè)、右側(cè)、前方、后方、頂部及底部布置檢波器，形成中心激發(fā)全方位接收的三維三分量地震觀測系統(tǒng)，如圖3所示。三維模型在,,方向的大小分別為400 m×400 m×400 m，震源點位于(0，0，0)原點處，三分量檢波器在,,面上以原點為中心直徑為200 m的圓形測線上以15°間隔布置，累計66個檢波器。模型在,,方向上進行網(wǎng)格化，網(wǎng)格間距均為0.3 m。模擬采用主頻為125 Hz的零相位雷克子波，采樣頻率為2.5 kHz。模型添加PML吸收邊界，采用三維時空域高階有限差分法進行數(shù)值模擬。

圖3 三維觀測系統(tǒng)

在數(shù)值模擬地震數(shù)據(jù)中選擇一道(第64道)數(shù)據(jù)進行時頻域極化分析，通過對比理論、計算的極化參數(shù)驗證上述方法效果。利用直達波的理論到達時間、震源主頻確定直達波的時頻范圍，在圖4中可分析直達縱波、橫波的極化參數(shù)，從圖5可見直達縱波、橫波質(zhì)點振動軌跡。圖4為第64道數(shù)據(jù)求得的極化方位角、極化傾角。由圖4(a)可看出，在時間0.03 s 和0.06 s附近各有一團能量，分別為縱波、橫波能量。直達縱波、橫波的方位角均為3.3°，與理論值0°存在部分偏差。由圖4(b)可看出，時間0.03 s附近的直達縱波傾角為45°，與理論值45°無偏差。時間0.06 s附近處的直達橫波傾角為-45°，與理論值-45°無偏差。計算其他65道地震波時頻域極化參數(shù)，并利用直達縱波的極化信息對震源點進行反向聚焦定位，結(jié)果如圖6所示，在200 m直徑范圍內(nèi)，求解震源點與已知震源點的直線距離誤差僅為 0.548 6 m。

圖4 計算所得極化方位角與極化傾角

圖5 第64道三分量信號質(zhì)點振動軌跡

圖6 檢波點反向聚焦定位震源點(右上角為原點位置局部放大)

3 現(xiàn)場試驗

實際礦井地震勘探過程中，除利用透射波之外，還常采用反射波，特別是煤層條件下反射槽波勘探是井下通用方法。為此，開展了現(xiàn)場實測并針對反射槽波信號進一步驗證方法有效性。

安徽某礦1034工作面位于三采區(qū)深部，煤層頂板標(biāo)高-603.8～-466.0 m，煤厚3.0～5.2 m，平均3.9 m，煤層傾角5°～18°，平均10°，頂板為灰白色中粒砂巖，底板為粉砂巖。1034回風(fēng)巷外幫發(fā)育有F13正斷層，與回風(fēng)巷相距57～108 m，走向NE，傾向NW，傾角60°～70°，落差100～200 m，其在礦區(qū)內(nèi)延展長度5.20 km，利用地面、井下鉆探進行斷層探查，斷層控制程度可靠。以R1檢波點為原點，每隔10 m布置1個檢波點，震源點設(shè)置在132.5 m處，三分量地震勘探觀測系統(tǒng)如圖7所示，S1炮激發(fā)R1～R27道接收的三分量地震記錄如圖8所示。

圖7 某礦1034工作面地震觀測系統(tǒng)

圖8 三分量地震記錄

觀察三分量地震記錄，圖中紅圈處(R1～R8道記錄)存在明顯反射槽波信號。選取R4道進行分析，S1炮、F13斷層與R4檢波點的反射路徑距離約213.28 m，R4檢波點接收的反射槽波理論方位角約為27°，傾角則等同于煤層傾角，約10°。S1炮激發(fā)、R4三分量檢波器接收的地震信號進行時頻域極化分析，時頻域方位角及時頻域傾角如圖9所示。根據(jù)已知斷層位置計算反射槽波理論傳播路徑，時間0.18～0.25 s、頻率120～180 Hz區(qū)域為反射槽波。圖9(a)中該區(qū)域的方位角計算結(jié)果可見，反射槽波的方位角較穩(wěn)定，約為28°，分析其主要為Rayleigh型反射槽波，與理論方位角偏差1°。圖9(b)顯示反射槽波的傾角為10°～12°，與實際煤層10°傾角吻合。

圖9 R4道數(shù)據(jù)時頻域極化方位角傾角分布

4 結(jié) 論

(1)對三分量合成信號進行時頻域極化分析，得到信號的時頻位置與理論時頻位置相對應(yīng)，求解的方位角、傾角與理論值吻合；在時間和頻率2個維度上，時頻域極化分析可以精準(zhǔn)確定混疊情況下多類型信號極化方向。

(2)針對中心激發(fā)-全方位接收的三維時空域高階三分量模擬信號，直達波的方位角、傾角誤差分別為3.3°，0°；在200 m直徑范圍內(nèi)，根據(jù)檢波點位置及地震波方位角、傾角聚焦定位的震源點與已知震源點的距離誤差僅為0.548 6 m，驗證了全空間條件下三分量地震波的時頻域極化參數(shù)計算的精度。

(3)針對礦井實測數(shù)據(jù)進行廣義S變換的時頻域極化分析處理，反射槽波的時頻域方位角、傾角參數(shù)與理論值一致，證實了礦井全空間三分量地震波時頻域極化分析方法的有效性。

(4)全空間條件下時頻域極化分析方法確定的地震波傳播方向可為礦井地震勘探精細(xì)成像提供基礎(chǔ)性支撐，并為井下地質(zhì)構(gòu)造的透射波層析成像、反射成像等提供基礎(chǔ)約束條件。