楊熙鐳,劉懷山,2
(1.中國海洋大學 海底科學與探測技術教育部重點實驗室,山東 青島 266100;2.海洋國家實驗室 海洋礦產資源評價與探測技術功能實驗室,山東 青島 266071)
地震資料去噪一直是地震資料處理的一個重要環(huán)節(jié),去噪結果的好壞直接影響到地震資料處理整個流程。受限于地震數據采集過程中野外條件以及設備成本的限制,實際采集到的地震資料中經常存在大量的隨機噪聲。低信噪比的地震資料會給后期處理環(huán)節(jié)的速度分析、偏移成像等帶來困難,影響地震資料成像的可靠性。因此地震資料去噪是提高成像分辨率的重要手段,是后續(xù)地震資料處理流程高質量進行的基礎。
在實際地震資料中往往無法直接獲得隨機噪聲的先驗參數,而本文旨在利用曲波變換得到地震資料的噪聲估計作為先驗信息,自適應選擇合適的隨機噪聲標準差得到K-SVD字典學習重構的誤差參數,從而實現在無人為設定參數的情況下,通過K-SVD字典學習方法實現對隨機噪聲的壓制,提高地震資料的信噪比。
K-SVD字典學習的方法是通過從已知的信號中獲得訓練樣本集,在稀疏約束的條件下不斷訓練學習,對信號特征值進行提取,實現信號的最大稀疏化,并獲得經過訓練后的超完備字典和稀疏系數[17]。K-SVD算法的主要數學思想是利用一個包含K個原子dk的字典矩陣D∈Rm×K,來最大化稀疏線性表示原始信號Y∈Rm×n,y∈RM×1對應二維信號中的列向量,稱為訓練原子,則有Y近似于Y1=DX,其中m、n分別對應二維信號的行列數,X∈RK×n為稀疏系數矩陣,xi(i=1,2,…,K)為該矩陣中的行向量,這樣就轉化為了如下的優(yōu)化問題
(1)
(2)
其中,
(3)
此時(1)式的優(yōu)化問題可以轉化為
(4)
(5)
(6)
在坡度小于15°區(qū)域內,農村居民點斑塊面積與數量占比高達91.9%與86.8%。在小于5°區(qū)域農村居民點密度、平均面積與標準差高出其他區(qū)域較多,而5°~15°區(qū)域內斑塊數量與景觀形狀指數大于0°~5°區(qū)域,體現出0°~5°區(qū)域內農村居民點具有連片聚集大規(guī)模的特征,而5°~15°區(qū)域內農村居民點斑塊更加碎小、并具有離散分布特性。隨著坡度增大農村居民點總面積、面積比例、數量等各項指數均減小,且在坡度大于25°時降至極低水平,表明坡度對農村居民點分布與規(guī)模均具有極大的影響。農村居民點形狀指數在0°~5°區(qū)域內最高,表明區(qū)域內農村居民點形狀復雜,而其他區(qū)域差異不明顯。
K-SVD字典學習去噪是通過選取一部分含噪的地震資料作為樣本進行學習,獲得一個與該地震資料相適應的超完備字典,其中地震波有效信號在該字典內能夠稀疏表示,而地震資料中含有的隨機噪聲則是在整個稀疏域均勻分布,通過重構算法便能夠實現有效信號和隨機噪聲的分離。為了在實現較好的隨機噪聲壓制效果的同時,能夠最大程度地保留有效信號不被切除,需要合理設置重構算法的誤差參數。目前大部分主流的重構算法,例如正交匹配追蹤(OMP,Orthogonal Matching Pursuit)、分段正交匹配追蹤(StOMP,Stagewise Orthogonal Matching Pursuit)等,都需要獲取隨機噪聲的先驗信息來控制誤差參數,才能實現較好的去噪效果,但在實際地震資料處理中,無法直接獲取數據中隨機噪聲的先驗參數,只能通過不斷地調試輸入的誤差參數來控制去噪效果。本文將曲波變換中所利用到的曲波噪聲估計[22]單獨提取出來作為一個模塊來獲取地震資料的噪聲先驗信息。由于曲波變換的性質,不同方向上的曲波系數能夠表征圖像在該方向上的特征信息,從而能夠幫助人們在曲波域內獲得噪聲方向,實現在曲波域的噪聲估計。
曲波噪聲估計是對地震資料進行曲波變換,即
(7)
其中C表示曲波系數,j為尺度參數,l為方向參數,k為位置參數,f為含噪地震資料,φj,l,k是在曲波變換的基函數。根據已獲取的含噪地震資料的曲波系數,選取其中某位置k的的曲波系數矩陣,選擇其中尺度參數j最大的一系列曲波系數。設Egj,l為含噪圖像在j個尺度第l個方向上投影得到的平均能量,Efj,l為原始圖像在第j個尺度第l個方向上投影得到的平均能量,Eηj,l為圖像噪聲在第j個尺度第l個方向上投影得到的平均能量,其中
Efj,lmin≤Efj,l
(10)
(11)
說明在lmin方向隨機噪聲能量最大,可以在此方向上獲得更好的噪聲估計。所以地震資料經曲波變換后,選取某點曲波系數中尺度參數最大且符合lmin條件的,則該曲波系數對應該點的最佳噪聲估計,并通過曲波系數得到相應的隨機噪聲標準差,再依此獲得地震資料中每一點的隨機噪聲標準差,并通過求取所有這些點隨機噪聲標準差的平均值,獲得整個地震資料的隨機噪聲標準差。
因此基于曲波噪聲估計的K-SVD字典學習去噪的基本流程步驟是:
1)對含有隨機噪聲的地震數據進行曲波變換,在曲波域利用曲波分解系數進行噪聲估計,即通過獲得符合最優(yōu)噪聲估計的曲波分解系數確定隨機噪聲標準差;
2)選取部分地震數據作為字典學習的訓練樣本集,利用OMP算法得到這部分地震資料的稀疏表示,獲得稀疏矩陣;
3)利用K-SVD字典學習方法對稀疏矩陣和字典進行迭代更新,獲得適應該地震資料的超完備字典及其對應的稀疏矩陣;
4)利用步驟3)中獲得字典和稀疏矩陣,通過OMP算法重構地震資料,在重構過程中,通過步驟1)中獲得的隨機噪聲標準差控制重構的迭代參數,最終重構出去噪后的地震資料。
為檢驗本方法對地震資料中隨機噪聲的去除效果,模擬單炮地震記錄,假設簡單的四層水平均勻介質,以主頻為45 Hz的Rick子波作為地震子波進行激發(fā),實現中間激發(fā)兩端接收,得到512道地震數據,單道采樣點數為512個,采樣間隔為0.001 s,得到模擬地震記錄如圖1(a)所示,加入標準差為4的高斯隨機白噪聲,得到含噪模擬地震記錄如圖1(b)所示。本文圖件縱橫端點所標注的數值,均為端點處刻度所對應的數值。
圖1 模擬地震記錄加噪前后對比Fig.1 Comparison of simulated seismic records before and after noise addition
分別對含噪模擬地震記錄(高斯隨機噪聲標準差為4)使用中值濾波、維納濾波、小波變換、本文方法進行去噪處理,并獲得去噪后殘差剖面。其中本文方法通過曲波噪聲估計得到含噪模擬地震記錄的噪聲標準差估計為3.9。設定基本字典參數,包括字典分塊大小為8,余度系數為4。
圖2展示了應用這四種去噪方法得到的地震記錄以及對應的殘差剖面。圖2縱坐標為時間,橫坐標為道號(CDP號)。從圖中可以看出,中值濾波去噪方法雖然有效地壓制了隨機噪聲,恢復了反射波同相軸,但是通過殘差剖面可以清晰地看到部分有效信息也被切除;維納濾波去噪方法相較于中值濾波更好地壓制了隨機噪聲,但是通過殘差剖面圖像可以看出部分有效信息也被切除,而且應用中值濾波去噪方法和維納濾波去噪方法后所獲得地震記錄中反射波同相軸并不光滑,依然存在噪點;應用小波變換去噪方法后,有效地壓制了隨機噪聲,較好地保留了有效信號,反射波同相軸更為光滑,但是可以看到其對隨機噪聲的壓制效果不及維納濾波,去噪后地震記錄的部分區(qū)域仍可見噪點;采用本文方法進行去噪后,同相軸光滑,有效信息得到保留,隨機噪聲壓制效果較好。
圖2 去噪后的地震記錄及其對應的殘差剖面Fig.2 Seismic records and their corresponding residual profiles after denoising
為了更清晰地對比本文方法與小波變換的去噪效果,提取第250道的單道地震記錄進行對比,如圖3所示,可以看到本文方法對比小波變換去噪方法,壓制隨機噪聲的同時,能更好地還原有效信號的幅值,避免了小波變換后幅值出現增益的情況。
圖3 第250道單道地震記錄去噪前后對比Fig.3 Comparison before and after denoising of the 250 th single channel seismic record
為了定量分析應用各個方法后的去噪效果,本文采用信號的信噪比(SNR,Signal to Noise Ratio)和峰值信噪比(PSNR,Peak Signal to Noise Ratio)作為衡量標準,其中二維信號的SNR和PSNR分別定義為
(12)
(13)
(14)
表1為對含噪模擬地震記錄應用不同去噪方法效果的比較。從表1中可以知道對于同樣的含噪地震資料,采用各種去噪方法壓制隨機噪聲后,信噪比均有所提升,但是通過對比可知,利用本文去噪方法進行處理后的SNR和PSNR最大,說明本文方法的去噪效果優(yōu)于其他幾種方法。
表1 幾種不同去噪方法的效果比較
為了驗證本文方法在實際地震資料去噪中的效果,對某海域含噪拖纜地震資料的疊加剖面進行去噪,實際地震資料如圖4所示。分別利用中值濾波、維納濾波、小波變換以及本文方法這四種方式對實際地震資料進行去噪處理,得到圖5所示的實際地震資料去噪效果對比圖。圖4和圖5的縱坐標為時間,橫坐標為道號(CDP號)。其中本文方法采用同模擬數據測試部分相同的參數設置。
圖4 實際含噪地震資料Fig.4 Actual seismic data with noise
從圖5中可以看出,中值濾波去噪處理后,同相軸基本已不再連續(xù),無法識別出實際地震資料中原本存在的斷層構造;維納濾波去噪后,雖然在淺部的強反射界面得到了較為連續(xù)的同相軸,但是其余深部的同相軸均在去噪過程中被切除,導致有效信號大量損失;小波變換去噪后,得到了大部分有效信號的同相軸信息,且可以明顯依據去噪后的信息推斷出實際資料中存在的斷層構造,但是其去噪后所得反射波同相軸仍然存在部分不連續(xù)的位置;本文方法去噪處理后,同樣恢復出大部分反射波同相軸,并且淺部強反射界面的同相軸連續(xù)性和光滑程度都要比小波變換去噪處理后更好。
結合圖4和圖5可以看出,本文運用的所有去噪方法,均導致650 ms下方兩個弱反射同相軸被壓制,說明運用包括本文方法在內的各種方式進行去噪處理后,都會不同程度地對反射波同相軸產生影響,但是相較于傳統(tǒng)方法,本文方法最大限度地降低了對有效信號的干擾。同時也說明在面對具有較低信噪比的地震數據時,本文方法也存在一定的局限性,雖然能夠準確地估計隨機噪聲,但是當噪聲幅值與弱反射同相軸的幅值過于接近時,進行去噪處理后仍無法有效保證弱反射同相軸的連續(xù)性。
圖5 實際地震記錄去噪前后對比Fig.5 Comparison of actual seismic records before and after denoising
本文提出基于曲波噪聲估計的K-SVD字典學習去噪方法,能夠直接對需要去噪的地震資料進行噪聲估計,使得K-SVD字典學習去噪方法不僅能夠自適應地基于地震資料獲取超完備字典用于重構,同時能夠自適應地獲得該地震資料的隨機噪聲標準差,增強了K-SVD字典學習去噪的泛化能力。對模擬資料和實際數據的測試均表明該方法具有較好的去噪效果,能夠在壓制噪聲的同時,最大限度地保留真實的有效信號。