羅耀耀 葛良全 熊 超 徐立鵬 花永濤
(成都理工大學 成都 610059)
航空γ能譜測量是一種快速、經濟、有效的核地球物理方法[1,2],航空γ射線能譜儀,記錄0–3 MeV的多道(256 道或更多道)全譜數據。但這些單次測量的原始能譜的統計漲落明顯,對這些能譜數據進行降噪處理,能降低數據中統計漲落帶來的影響,是提高放射性核素含量精度的關鍵。
小波分析是近20多年發展起來的交叉學科,在信噪分離、提取弱信號等方面應用較廣。小波去噪有硬閾值去噪和軟閾值去噪[3–6],但這兩個函數存在一定的缺陷,如硬閾值函數的不連續性,軟閾值函數估計的小波系數與被處理信號的小波系數之間存在恒定偏差等,限制了該方法的應用效果[7,8]。
本文分析了最小二乘法及小波的軟硬閾值去噪方法,并采用雙小波基與Wiener濾波器相結合的方法,實現了航空γ能譜的降噪新方法。
小波變換能同時在時域和頻域中對信號進行分析,具有自動變焦功能,故能有效區分信號中的突變部分和噪聲,從而實現信號的降噪。
航空γ能譜數據f(n)可用式(1)表示

式中,f(i)為能譜信號,s(i)為純信號,n(i)為高斯白噪聲信號,n(i)=σv(i)為加性隨機噪聲,其強度是σ。噪聲信號n(i)與采樣信號s(i)不相關。由于小波變換是線性的,所以經過小波變換后的小波系數模型可寫成:

其中,y=wf,θ=ws,z=wn分別為信號、純凈信號以及噪聲的小波系數,w為小波變換。
噪聲和信號在小波變換下具有如下表現形式:從能量觀點看,噪聲能量分布在所有小波系數上,而信號能量僅分布在一小部分系數上;噪聲經小波變換后幅值較小,數目較多,而信號經小波變換后幅值較大,數目較少。通過分析,可以構造閾值去噪方法,基于小波的消除噪聲的方法可分為三個步驟[8]:
(1) 選擇小波和小波分解的層次,計算信號s第1層到第N層的小波分解系數。
(2) 高頻系數的閾值選取,對于從第 1層到第N層的每一層,選取1個閾值(也可用同1個閾值),對高頻系數用軟、硬閾值進行處理(用硬閾值處理信號較粗糙,一般采用軟閾值)。
(3) 根據第N層的低頻系數和從第1層到第N層經修改后的高頻系數,重構信號。
Wiener濾波器在最小均方誤差意義上是最優的,其前提是已知信號和噪聲的統計特性。因此,實際應用中需要先對信號和噪聲作相應估計,為此,本文根據Wiener濾波器和小波變換的特點,采用小波域Wiener濾波器對航空γ能譜譜線數據進行處理。

由于小波變換對信號有去相關的特性,對小波系數可采用Wiener濾波器,小波的Wiener濾波器設計值為:其中,σ2是噪聲z(i)的方差,因為θ(i)是未知的,所以,使用其估計值(i) ,于是得到經驗的Wiener濾波器的值為:

在設計經驗Wiener濾波器時,須從信號中估計出θ2(i)和σ2,其中噪聲σ2可由小波系數估計,而θ2(i)則較難估計。故采用兩個小波基W1和W2,W1對信號進行降噪,從而得到原始信號較為準確的估計(i),再用經驗Wiener濾波器在W2小波閾中進行處理,由此獲得較好的降噪效果。采用兩個小波基可保證在 W1閾內由于閾值函數被置為零的系數在W2小波閾中得以恢復,達到降噪效果。根據要求設計的小波域經驗Wiener濾波器如圖1所示。

圖1 小波閾Wiener濾波器結構圖Fig.1 Structure of the wavelet domain Wiener filter.
為研究該方法在航空能譜數據處理中的有效性和準確性,將內蒙某測區航空γ能譜原始譜線數據按工作流程分為早晚校正數據集和測區數據集。早晚校正數據集是指航空γ能譜儀在飛行前和飛行后對儀器穩定性、譜線峰位等參數檢查過程中獲得的數據;測區數據是指航空γ能譜儀在測區工作時獲得的原始數據。
取航空γ能譜儀位于內蒙某機場地面靜止測量時的數據,即早晚校正數據進行分析,比較了不同降噪方法后的鉀窗、鈾窗、釷窗內的計數率均值和相應的標準差(表1);與原始譜線相比,采用最小二乘法和小波軟閾值濾波法及小波域 Wiener濾波器方法降噪后,各個窗口的計數率平均值基本不變,但標準差發生了變化,經小波域Wiener濾波器方法處理后誤差小于其他兩種方法。可見,采用小波方法應用于單次測量取得了比較理想的效果。
圖2顯示飛行高度位于142 m航空γ能譜儀獲得的原始能譜數據與經小波域 Wiener濾波器去噪后的譜線對比圖。可以看出,能譜經小波域Wiener濾波器后,降低了噪聲能譜數據統計漲落帶來的影響,特征峰均能清晰顯露,如208Tl的2.62 MeV、214Bi的1.76 MeV和1.12 MeV、40K的1.46 MeV、228Ac的0.911 MeV的特征γ射線峰以及由214Bi的0.609 MeV和208Tl的0.583 MeV構成的重疊峰。
圖3為采用小波域Wiener濾波器和5點3次多項式擬合法去噪后得到的譜線。觀察208Tl的0.609 MeV峰和40K的1.46 MeV以及208Tl的2.62 MeV的峰可以看出,小波域Wiener濾波方法的降噪效果明顯優于傳統的能譜降噪方法(多項式光滑方法),糾正了傳統方法導致的譜線畸變。

表1 采用不同降噪方法處理全能譜數據前后鉀窗、鈾窗、釷窗內的計數率均值和相應的標準差Table 1 Counts of K, U and Th windows before and after denoising by different methods.
本文介紹了采用小波域Wiener濾波器對航空γ能譜數據進行降噪的方法,并應用于某測區的數據中,與常規方法進行了對比分析。通過對航空γ能譜儀在測區工作時早晚校正數據三個窗口數據的平均值和標準差分析,對窗口數據的標準差有所改善;航空γ能譜儀在高空中(142 m)的譜線數據表明,該方法能有效降低航空γ能譜數據的統計漲落,凸顯能譜的特征峰。

圖2 采用小波域Wiener濾波器降噪前(*)后(●)譜線比較Fig.2 Energy spectra before(*) and after (●)denoising.

圖3 采用多項式擬合法(*)和小波域Wiener濾波器(▲)去噪后的譜線比較Fig.3 Denoising result comparison between conventional and wavelet domain Wiener filter.* Smoothed by Polynomial fitting, ▲ Smoothed by wavelet domain Wiener filter
1 Bruce L D. Recent advances in aerial gamma-ray surveying[J]. Journal of Environmental Radioactivity,2004, 76(1–2): 225–236
2 International atomic energy agency[S]. Guidelines for radioelement mapping using gamma-ray spectrometry data, 2003: 52–71
3 于國梁, 張一云, 賀梁國. 小波分析法檢測譜特征峰[J].原子能科學技術, 2008, 42(11): 984–987 YU Guoliang, ZHANG Yiyun, HE Liangguo. Characteristic peak detection by using wavelet analysis in γ-ray spectra[J]. Atomic Energy Science and Technology, 2008,42(11): 984–987
4 鄒永祥, 吳建平. 基于多分辨率分析的 X 熒光能譜的除噪與弱峰檢測[J]. 核技術, 2008, 31(4): 260–264 ZOU Yongxiang,WU Jianping. Multi-resolution analysis for noise reduction and weak peak detection in X-ray fluorescence spectra[J]. Nucl Tech, 2008, 31(4): 260–264
5 肖 剛, 周春林, 李天柁, 等. 一種自適應的射線能譜數據光滑方法[J]. 原子能科學技術, 2001, 35(1): 69–72 XIAO Gang, ZHOU Chunlin, LI Tiantuo,et al. An adaptive method for γ-spectra smoothing[J]. Atomic Energy Science and Technology, 2001, 35(1): 69–72
6 史東生, 弟宇鳴, 周春林. 小波變換與傅里葉變換在能譜降噪處理中的比較研究[J]. 核電子學與探測技術,2006, 26(6): 134–137, 107 SHI Dongsheng, DI Yuming, ZHOU Chunlin. Compara-tive study on Cenergy spectrum denoise by fourier and wavelet transforms[J]. Nuclear Electr onics &Detection Technology, 2006, 26(6): 134–137, 107
7 Mallat S. A wavelet tour of signal processing[M]. London:Academic Press, 1998
8 Donoho D L. De-noising by soft-thresholding[J]. IEEE Trans on Information Theory, 1995, 41(3): 613–627