γ能譜的數(shù)據(jù)處理

1、能譜的數(shù)據(jù)處理由多道脈沖分析器獲取的譜數(shù)據(jù)需要以一定的數(shù)學(xué)方法進(jìn)行處理才能得到實驗要求的最終結(jié)果。能譜的數(shù)據(jù)處理大致可以分為兩個步驟。首先進(jìn)行峰分析，即由能譜數(shù)據(jù)中找到全部有意義的峰，并計算出扣除本底之后每個峰的凈面積。第二步是放射性核素的活度或樣品中元素濃度的計算，即由峰位所對應(yīng)的能量識別出被測樣品中含有哪些放射性核素或被激發(fā)的元素，并且由峰的凈面積計算出放射性核素的活度或元素在樣品中的濃度。采用不同的物理實驗方法，使用不同的探測器時，能譜的數(shù)據(jù)處理方法也有所不同。在本章中首先討論在各種能譜數(shù)據(jù)處理中經(jīng)常用到的峰分析方法，包括譜數(shù)據(jù)的平滑處理，本底扣除、尋峰、峰凈面積計算和函數(shù)擬合法解譜。然

2、后，以譜分析為例，討論基于計算機(jī)的多道脈沖分析系統(tǒng)中的譜自動分析軟件的工作原理。第一節(jié) 常用的峰分析方法一、譜數(shù)據(jù)的平滑處理由于射線和探測器中固有的統(tǒng)計漲落、電子學(xué)系統(tǒng)的噪聲的影響，譜數(shù)據(jù)有很大的統(tǒng)計漲落。在每道計數(shù)較少時，相對統(tǒng)計漲落更大。譜數(shù)據(jù)的漲落將會使譜數(shù)據(jù)處理產(chǎn)生誤差。其主要表現(xiàn)為在尋峰過程中丟失弱峰或出現(xiàn)假峰、峰凈面積計算的誤差加大等等。譜數(shù)據(jù)的平滑就是以一定的數(shù)學(xué)方法對譜數(shù)據(jù)進(jìn)行處理，減少譜數(shù)據(jù)中的統(tǒng)計漲落，但平滑之后的譜曲線應(yīng)盡可能地保留平滑前譜曲線中有意義的特征，峰的形狀和峰的凈面積不應(yīng)產(chǎn)生很大的變化。對譜數(shù)據(jù)進(jìn)行平滑處理通常使用數(shù)字濾波器。由信號分析理論的觀點(diǎn)出發(fā)，我們可以

3、把原始譜數(shù)據(jù)看成是噪聲（即譜數(shù)據(jù)中的統(tǒng)計漲落）和信號（即峰函數(shù)和本底函數(shù)）的疊加。經(jīng)過數(shù)字濾波器的處理可以提高信號噪聲比。如圖5-1-1所示，令第x道的原始譜數(shù)據(jù)為y(x)，經(jīng)過數(shù)字濾波之后的譜數(shù)據(jù)為()其中，為數(shù)字濾波器的單位沖擊響應(yīng)函數(shù)，并有()y()Y(x)Y(x)圖5-1-1 用數(shù)字濾波器進(jìn)行譜的平滑處理由于譜數(shù)據(jù)是離散量，公式（）、（）的離散量表達(dá)形式為()()只要選擇恰當(dāng)?shù)臄?shù)字濾波器響應(yīng)函數(shù)，就能夠使平滑后的譜既保留了原始譜中的峰和本底的形狀和大小，又得到最佳的信號噪聲比。由頻域的觀點(diǎn)分析，譜中的統(tǒng)計漲落，即噪聲的頻譜分布在整個頻率范圍內(nèi)，而峰函數(shù)和本底函數(shù)的頻譜主要集中在低頻范圍

4、。因此，使用一個低通濾波器進(jìn)行濾波，可以使峰和本底信息都通過濾波器到達(dá)輸出器，而噪聲中的高頻成分被濾波器抑制，從而提高了平滑后譜中的信號噪聲比，減小了譜數(shù)據(jù)的統(tǒng)計漲落。圖5-1-2中畫出了二種常用的數(shù)字濾波器分別在道域和頻域中的響應(yīng)函數(shù)的圖形。1. 最小二乘移動平滑方法1964年A. Savitzky和J.Egolay1提出了一個用于譜數(shù)據(jù)平滑處理的濾波器響應(yīng)函數(shù)。其基本思想是，當(dāng)求平滑之后譜的第m點(diǎn)數(shù)據(jù)時，先在原始譜數(shù)據(jù)第m點(diǎn)的左、右各取K個數(shù)據(jù)點(diǎn)，形成一個共有2K+1個數(shù)據(jù)點(diǎn)的窗口。在這個窗口中用多項式擬合原始譜數(shù)據(jù)，則擬合多項式在m點(diǎn)的值就是平滑后的譜在m點(diǎn)的值。當(dāng)m值沿譜數(shù)據(jù)移動時，就

5、可以得到整個平滑后的譜數(shù)據(jù)。這種方法稱為最小二乘移動平滑法，或最小平方曲線擬合平滑法。原始譜數(shù)據(jù)為，平滑后譜數(shù)據(jù)為，在平滑窗口內(nèi)，用q價多項式逼近原始譜數(shù)據(jù)時，平滑后譜第m點(diǎn)的值為()圖5-1-2 SAVITZKY濾波器的響應(yīng)函數(shù)(a) SAVITZKY濾波器在道域中的響應(yīng)函數(shù)；(b) SAVITZKY濾波器在頻域中的響應(yīng)函數(shù)；1. SAVITZKY5點(diǎn)平滑濾波器；2. SAVITZKY11點(diǎn)平滑濾波器同時還可以把S（x）在m點(diǎn)的各階導(dǎo)數(shù)值作為平滑后的譜在m點(diǎn)的各階導(dǎo)數(shù)值。平滑后的譜在m點(diǎn)的各階導(dǎo)數(shù)值。平滑后的譜在m點(diǎn)的p階導(dǎo)數(shù)值為()根據(jù)上述原理，用最小二乘法函數(shù)擬合可以導(dǎo)出計算平滑后的譜數(shù)

6、據(jù)和其各階導(dǎo)數(shù)值的具體計算公式(5.1.7)規(guī)范化常數(shù)NK和權(quán)因子CK,j的值列在表5.1中。由表中查出NK和CK,j的值就可以寫出平滑譜的計算公式。例如當(dāng)平滑窗口選為5點(diǎn)時（K=2），5點(diǎn)平滑公式為()在文獻(xiàn)【2】中，由數(shù)字濾波器理論可以推導(dǎo)出最小二乘移動平滑公式（）中Ck,j/NK的一般計算公式。當(dāng)平滑窗口為W=2K+1時，()KjK用這個公式計算得出的值，與表5.1中列出的數(shù)值相吻合。表5.1 最小二乘移動平滑計算公式中的NK值和CK,j值2K+1CK,j8-21-7-6-78-67-13-11-518420-36-4278799-21-334122164414

7、-2-2391472169393-3-1421622484546120431672589597171421622484546122391472169393-3334122164414-24278799-21518420-3667-13-117-6-788-21NK32311051434292312135表5.2 最小二乘移動平滑法計算平滑譜的一階導(dǎo)數(shù)公式中的NK值和CK,j值2K+1CK,j8748-7-9812922-6-643-41211133-5-930-14150-660300-4-1002-18334-1578-29486-3-902-17842-1796-

8、532-14222-2-673-13843-1489-503-193-671-1-358-7506-832-296-126-58-80000000013587506832296126588267313843148950319367-13902178421796532142-2241002183341578294-86593014150660-30066434121-1133798-129228748NK-23256334152240245148118825212表5.3 最小二乘移動平滑計算平滑譜的二階導(dǎo)數(shù)公式中的NK值和CK,j值2K+1CK,j840-72591-

9、6125222-51191115-4-8-82628-3-15-29-5-175-2-20-48-10-6-802-1-23-53-13-9-17-3-10-24-56-14-10-20-4-21-23-53-13-9-17-3-12-20-48-10-6-8023-15-29-5-1754-8-8262851191115612522272591840NK397661881001429462427由公式（）也可以計算平滑譜的各階導(dǎo)數(shù)值，只不過權(quán)因子CK,j和規(guī)范化常數(shù)NK的值各不相同。表5.2中列出了采用不同的平滑窗口、用公式（）計算平滑譜的一階導(dǎo)數(shù)時的CK,j與NK的值。表5.3中列出了采用

10、不同的平滑窗口，計算平滑譜的二階導(dǎo)數(shù)時的CK,j與NK的值。根據(jù)表5.2平滑窗口為5點(diǎn)（K=2）時，平滑譜的一階導(dǎo)數(shù)計算公式為()根據(jù)表5.3，平滑窗口為5點(diǎn)時，平滑譜在m點(diǎn)的二階導(dǎo)數(shù)值為：()前面已經(jīng)指出，平滑的本質(zhì)是對譜曲線進(jìn)行低通濾波，去掉高頻成分，保留有用的低頻信息。濾波的效果取決于低通濾波器的頻譜特性。當(dāng)式（）中的權(quán)因子不同時濾波器的頻譜特性不同，濾波的效果也不同。在某些實際應(yīng)用的平滑程序中使用了不同于式（）的權(quán)因子，例如在SPECTRAN-F3程序中使用的平滑公式如下：三點(diǎn)平滑公式為()五點(diǎn)平滑公式為(3)七點(diǎn)平滑公式為(4)這幾個平滑公式的優(yōu)點(diǎn)是權(quán)因子都是正數(shù)，平滑之后的譜數(shù)據(jù)不

11、可能出現(xiàn)負(fù)值，從而提高了平滑之后的譜數(shù)據(jù)的可靠性。這在原始譜數(shù)據(jù)中本底很小、峰很高、而且峰的寬度很窄時是非常重要的。如果平滑之后的譜數(shù)據(jù)出現(xiàn)了負(fù)值（這顯然是不合理的），可能使后續(xù)的計算程序在運(yùn)行時產(chǎn)生錯誤。 2. 采用高斯濾波器的平滑方法如果把譜數(shù)據(jù)中的統(tǒng)計漲落看成是“白噪聲”，當(dāng)使用匹配濾波器進(jìn)行濾波時，可以得到最佳的信號噪聲比。所謂匹配濾波器，就是該濾波器在道域中沖擊響應(yīng)函數(shù)與峰函數(shù)互為鏡象。一般情況下，譜中的峰函數(shù)可以近似為高斯函數(shù)。由于高斯函數(shù)是偶函數(shù)，所以匹配濾波器在道域中的響應(yīng)函數(shù)也應(yīng)該是高斯函數(shù)。實踐表明，僅在平滑窗口比較大的情況下，使用高斯濾波器進(jìn)行平滑處理，才能得到比最小二乘

12、移動平滑方法稍好的平滑效果。因此，目前應(yīng)用最廣泛的還是最小二乘移動平滑方法。3. 譜平滑中的幾個具體問題對譜進(jìn)行平滑處理可以減少譜數(shù)據(jù)的統(tǒng)計漲落，從而減少了尋峰過程中假峰出現(xiàn)的幾率，也可以減小峰凈面積的計算誤差。但是當(dāng)濾波器的參數(shù)選擇不當(dāng)或平滑次數(shù)過多時也會產(chǎn)生某些缺點(diǎn)。例如，在尋峰時可能漏失弱峰，不能分辨距離很近的重峰等等。因此，如何選擇濾波器的參數(shù)和平滑的重復(fù)次數(shù)是很重要的。(1) 平滑窗口的選擇由公式（）可以看出，離散量的卷積運(yùn)算實際上是加權(quán)求和。當(dāng)計算平滑后的譜的第m的數(shù)據(jù)時，需要在原始譜中第m點(diǎn)兩邊各取K個點(diǎn)（共2K+1個點(diǎn)）進(jìn)行運(yùn)算。我們把2K+1叫做平滑窗口。改變平滑窗口的大小對

13、平滑效果有很大的影響。圖5-1-3中 譜平滑處理的效果與平滑窗口大小的關(guān)系曲線圖5-1-3中畫出了經(jīng)過平滑處理之后譜中統(tǒng)計漲落的改善與平滑窗口大小的關(guān)系曲線4。圖中橫坐標(biāo)為最小二乘移動平滑方法的平滑窗口，縱坐標(biāo)是標(biāo)志著平滑效果的統(tǒng)計漲落改善因子。這個改善因子代表了平滑后與平滑前的漲落的比較，即峰的之內(nèi)的面積與本底噪聲的比值或峰高與本底噪聲的比值的變化。由圖5-1-3看出，當(dāng)平滑窗口比較小時，隨著平滑窗口增大，平滑的效果增大。在某一個平滑窗口時，改善因子達(dá)到最大值。但是當(dāng)平滑窗口繼續(xù)增大時，改善因子反而下降。這是因為在平滑窗口較小時，隨著窗口的增大，譜中的統(tǒng)計漲落減小很快，但譜中的峰高和峰的形狀

14、變化不大。當(dāng)平滑窗口超過一定數(shù)值時，經(jīng)過濾波之后峰高將急劇下降，因而平滑效果反而降低。改善因子達(dá)到最大值時的平滑窗口稱為最佳平滑窗口。最佳平滑窗口的大小與譜曲線中峰的寬度有關(guān)。當(dāng)峰的半高寬FWHM比較大時，最佳平滑窗口也較大。因此，為了達(dá)到好的平滑效果，對于寬度不同的峰需要選取不同的平滑窗口。選取平滑窗口大小時應(yīng)該考慮的另一個因素是平滑窗口大小對譜曲線形狀的影響。當(dāng)平滑窗口比峰的FWHM大很多時，平滑之后的譜中的峰將顯著地變寬。這將使譜中原來相互靠得比較近的峰重疊得更加嚴(yán)重，從而使尋峰和峰凈而積的計算更加困難。綜合上面二個因素，在平滑處理時平滑窗口的大小要根據(jù)譜中峰的寬度來選擇。一般的作法是選

15、用的平滑窗口近似等于峰的半高度FWHM（以道為單位）。例如，當(dāng)FWHM7時，取2K+1=5；7FWHM9時，取2K+1=7；9FWHM11時，取2K+1=9，等等。能譜曲線中峰的寬度隨道址的增加而加大。我們可以把整個譜分成若干段，每段采用不同的平滑窗口。(2) 平滑重復(fù)次數(shù)在使用較小的平滑窗口時，對譜數(shù)據(jù)多次重復(fù)地進(jìn)行平滑處理，可以更有效地減小譜數(shù)據(jù)中的統(tǒng)計漲落。一個均值為常數(shù)、服從正態(tài)分布的偽隨機(jī)數(shù)系列，其平滑效果與平滑次數(shù)的關(guān)系曲線如圖5-1-4所示。橫坐標(biāo)MS為平滑的重復(fù)次數(shù)，縱坐標(biāo)是平滑之后該數(shù)列分布的標(biāo)準(zhǔn)偏差的相對值。由圖中可以看出，隨著平滑次數(shù)的增加，譜數(shù)據(jù)的統(tǒng)計漲落逐漸減小。但是

16、在平滑次數(shù)大于3時，曲線下降得很平緩，再增加平滑次數(shù)對統(tǒng)計漲落的改善并不顯著。多次平滑會使譜的形狀產(chǎn)生畸變。這種畸變主要表現(xiàn)為峰高降低，峰寬增大，峰谷被填平。作為一個例子，圖5-1-5和表5.4中給出了平滑次數(shù)對譜中峰形狀的影響。圖5-1-5是用Si（Li）探測器測得的55Fe的X射線譜。表5.4列出了使用五點(diǎn)最小二乘移動平滑方法時，峰位（mp）、峰高（h1）、峰谷（h2）和峰的關(guān)高寬（FWHM）與平滑次數(shù)MS的關(guān)系。由這個例子我們可以看出，當(dāng)平滑次數(shù)增加時，峰位基本上維持不變，但是峰的半高寬增加，峰高降低，峰谷抬高。對于一個單峰來說，雖然峰高下降，但半高寬增加，因而峰的面積變化不大，不會使譜

17、的定量分析產(chǎn)生很大的誤差。在重疊峰的情況下，由于峰的展寬，可以會淹沒位于某一強(qiáng)峰附近的弱峰。在復(fù)雜譜的分析中，這將會造成弱成分的漏失。圖5-1-4 平滑效果與平滑次數(shù)的關(guān)系曲線 圖5-1-5 55Fe的X射線譜表5.4 峰形狀與平滑次數(shù)的關(guān)系MS（次數(shù)）mp（道址）h1（計數(shù)）h2（計數(shù)） FWHM（道數(shù)）01017.31192935024.511017.41182335024.551017.61160435224.9101017.61134735625.4301017.61047737527.5501017.6978541529.41001017.6855151833.4總之，需要對譜數(shù)據(jù)平

18、滑多少次應(yīng)考慮到改善譜的統(tǒng)計漲落、減少譜形畸變兩個因素，根據(jù)譜數(shù)據(jù)的具體情況決定。在譜數(shù)據(jù)中各道計數(shù)較低，統(tǒng)計漲落較大的情況下，平滑次數(shù)可以多些。在譜數(shù)據(jù)中各道計數(shù)較大，或者譜形比較復(fù)雜的情況下，為了減少譜形的畸變，節(jié)省計算時間，平滑次數(shù)應(yīng)當(dāng)少一些，一般不多于3次。下面給出用FORTRAN語言編寫的最小二乘移動平滑法進(jìn)行譜數(shù)據(jù)平滑處理的程序。程序中變量的意義如下：SP：實數(shù)組，存放平滑前后的譜數(shù)據(jù)。A：實數(shù)組，工作單元暫存臨近譜區(qū)左右邊界的譜數(shù)據(jù)。B：實數(shù)組，存放平滑公式中的權(quán)因子。MB：整變量，譜區(qū)的左邊界道址。MF：整變量，譜區(qū)的右邊界道址。IP（1）l：IP數(shù)組的第一個元素，存放平滑次數(shù)

19、。FW：實變量，譜區(qū)內(nèi)的平均FWHM。SUBROUTINE SMTH （MB，MF）DIMENSION A（5），B（6）COMMON/FPDATA/IP（15），RP（10）COMMON/CURRP/FW，SG，CDL，CDR，STEP，NOSP，NPK，NPS，NPF，NN，NMCOMMON/ARRAY/SP（2100），PKL（160），T（128）IF（IP（1）·LE·Ø） RETURNM11DO 1 J1，5J1MBJ61A（J）SP（J1）IF（FW·LT·4）GOTO 3IF（FW·LT·7·

20、16;） GOTO5IF（FW·LT·9·Ø） GOTO7IF（FW·LT·11·Ø） GOTO9B（1）Ø·333333B（2）Ø·27972ØB（3）Ø·13986ØB（4）Ø·Ø2331B（5）Ø·1Ø4895B（6）Ø·Ø41958J5GOTO 103B（1）Ø·5B（2）Ø·25J1GOTO 105

21、B（1）Ø·375B（2）Ø·25B（3）Ø·Ø625J2GOTO 107B（1）Ø·Ø3125B（2）Ø·234375B（3）Ø·Ø9375B（4）Ø·Ø15625J3GOTO 109B（1）Ø·Ø417249B（2）Ø·314685B（3）Ø·Ø 6993 ØB（4）Ø·1282 Ø5B（5）&

22、#216;·Ø34965J410DO 40 IMB，MFSB（1）*SP（I）D0 30 L1，JL1L1L2ILL3IL3 ØS=S（SP（L2）SP（L3）*B（L1）SP（I5）=S40CONTINUEIMF50SP（I）SP（I5）II1IF（I·GE·MB） GOTO 50DO 60 K=1，5J1MBK660SP（J1）=A（K）IF（M1·GE·IP（1）RETURNM1M1GOTO 10END二、尋峰在譜數(shù)據(jù)中精確地計算出各個峰的峰位是能譜分析中的最關(guān)鍵的問題。在譜的定性分析中，只有正確地找到譜中全部峰的位置

23、，才能根據(jù)主峰和各驗證峰的能量來決定在被測樣品中是否存在某種核素。在譜的定量分析中，尤其是用最小二乘法函數(shù)擬合進(jìn)行重峰分析時，一般使用迭代法，峰位作為迭代參數(shù)的初值，如果峰位的誤差很大，或混入了假峰，漏失了真峰，則會造成迭代次數(shù)增多，甚至不收斂，使迭代失敗。由于譜結(jié)構(gòu)的復(fù)雜和統(tǒng)計漲落的影響，從譜中正確地找到全部存在的峰是比較困難的。尤其是找到位于很高本底上的弱峰，分辨出相互靠得很近的重峰更為困難。譜分析對尋峰方法的基本要求如下：(1) 比較高的重峰分辨能力。能確定相互距離很近的峰的峰位。(2) 能識別弱峰，特別是位于高本底上的弱峰。(3) 假峰出現(xiàn)的幾率要小。(4) 不僅能計算出峰位的整數(shù)道址

24、，還能計算出峰位的精確值，某些情況下要求峰位的誤差小于0.2道。很多作者對尋峰方法進(jìn)行了研究，提出了很多有效的尋峰方法。有的方法分辨重峰的能力較強(qiáng)；有的方法適于在很高的本底上尋找弱峰；有的方法計算簡單，尋峰速度快。但是各種不同的尋峰方法都可以概括為譜變換和峰判定二個步驟。譜變換一般采用線性濾波技術(shù)，其目的是減少譜數(shù)據(jù)的統(tǒng)計漲落，消除本底的影響，突出峰位信息。峰判定主要是根據(jù)預(yù)先設(shè)置的條件來識別真峰和剔除假峰。譜變換采用線性濾波技術(shù)，即對譜數(shù)據(jù)進(jìn)行如下的線性變換：（m1，2n）()其中，是譜第m道的數(shù)據(jù)，是變換之后譜的第m道數(shù)據(jù)，Cj是濾波器單位沖擊的響應(yīng)函數(shù)決定的因子，n為譜中的數(shù)據(jù)點(diǎn)數(shù)，2K

25、1為變換窗口。尋峰中使用的數(shù)字濾波器平平滑時使用的濾波器不同。平滑處理時使用低通濾波器，而尋峰中一般使用帶通濾波器。譜中的本底成份隨道址是緩慢變化的，在頻域中屬于低頻成份，不能通過帶通濾波器，因而在變換之后的譜中消除了本底成份的影響。此外，由于帶通濾波器既能抑制高頻噪聲又能抑制低頻噪聲，從而能更有效地減少統(tǒng)計漲落，突出峰信息。圖5-1-6中畫出了變換前后譜的形狀。由曲線（b）可以看出，變換之后的譜數(shù)據(jù)中的統(tǒng)計漲落并不能完全清除。此外，在康普頓邊緣處也可能出現(xiàn)局部極大值。為了避免把較大的統(tǒng)計漲落和康普頓邊緣誤認(rèn)為一個峰，剔除不真實的假峰，需要對變換之后的譜中的峰信息進(jìn)行峰判定。只有峰信息滿足判定

26、條件時才確認(rèn)存在一個真峰。峰判定包括峰高的統(tǒng)計性判定、峰凈面積的統(tǒng)計性判定和峰形判定。幾種不同的判定條件可以用在同一個尋峰方法中，也可以使用其中的一部分。采用不同的譜變換方法和不同的峰判定條件就形成了不同的尋峰方法。下面討論幾種主要尋峰方法的基本原理，并對它們的性能進(jìn)行比較。1. 幾種常用尋峰方法的基本原理(1) 對稱零面積對合法：進(jìn)行譜變換時使用的數(shù)字濾波器的沖擊函數(shù)Cj滿足下列關(guān)系式 （j=K，K）()即沖擊函數(shù)圍繞j0點(diǎn)左右對稱，且對橫軸所包圍的面積為零，因而稱為對稱零面積對合法。由式（），（）可以看出，如果在峰區(qū)范圍內(nèi)本底譜是常數(shù)或按直線分布，則在變換之后的譜中本底的貢獻(xiàn)為零。變換之后

27、譜的形狀完全反映了峰形的變化。下一步使用峰高統(tǒng)計判定條件來確定峰位。沿變換之后的譜進(jìn)行檢索，找出局部極大點(diǎn)，若極大值超過其均方根誤差若干倍時，則認(rèn)為找到了一個峰，該局部極大點(diǎn)所對應(yīng)的道址就是峰位。沿道址在譜中進(jìn)行檢索，找出一系列的局部極大值，并進(jìn)行判定就可以找到一系列的峰。上述的峰高判定條件可以寫為：()其中TRH為尋峰閾值，是預(yù)先給定的常數(shù)。當(dāng)尋峰閾值取得大時，能有效地剔除由統(tǒng)計漲落造成的假峰。但也可能同時漏掉了弱的真峰。當(dāng)設(shè)定較小的尋峰閾值時，能找到全部的強(qiáng)峰和弱峰，但也可能會把由統(tǒng)計漲落造成的假峰當(dāng)成真峰。因此，尋峰閾值的大小應(yīng)當(dāng)根據(jù)譜數(shù)據(jù)的具體情況和物理實驗任務(wù)來確定。在一般情況下，尋

28、峰閾值可選在25之間。利用匹配濾波器可以得到最佳信號噪聲比，導(dǎo)出匹配濾波器尋峰方法。由于譜數(shù)據(jù)中的峰函數(shù)接近于高斯函數(shù)，因而把沖擊函數(shù)為高斯函數(shù)的濾波器作為匹配濾波器。Cj的表達(dá)式如下：()其中，為高斯峰函數(shù)的標(biāo)準(zhǔn)偏差，2K1為變換窗口。實際上，與道址有關(guān)?？梢岳米V儀系統(tǒng)的FWHM刻度公式求出第m道的FWHM值，由FWHM/2.355計算出第i道的值，在變換窗口內(nèi)可以認(rèn)為值是一個常數(shù)。圖5-1-6 變換前后的譜曲線(a) 譜曲線(b) 經(jīng)匹配濾波器變換之后的譜曲線圖5-1-6中曲線（a）畫出了譜曲線中的一個峰，曲線（b）是經(jīng)匹配濾波器變換后，按式（）求出的Rm值與道址m的關(guān)系曲線。在峰位處曲

29、線出現(xiàn)局部極大值，而在峰的兩側(cè)出現(xiàn)局部極小值。兩個極小值之間的距離可以作為峰的寬度信息，大約等于峰的十分之一高全寬度FWTM。下面是用FORTRAN語言編寫的匹配濾波器法尋峰子程序，程序中使用的變量名和子程序名的意義如下：SP：實數(shù)組，存放譜數(shù)據(jù)；PKL：實數(shù)組，存放由尋峰程序找到的峰位的精確道址；C：實數(shù)組，存放匹配濾波器的沖擊響應(yīng)函數(shù)值；MB：整變量，尋峰譜區(qū)左邊界在SP中的下標(biāo)；MF：整變量，尋峰譜區(qū)右邊界在SP中的下標(biāo)；TRH1：實變量，尋峰閾值；FW：峰的關(guān)高寬NP3：在譜區(qū)內(nèi)找到的峰的個數(shù)；NPS：譜區(qū)左邊界道址；GETCRP：子程序名，由峰形刻度參數(shù)計算道址為X時的FW值；FDP

30、S：子程序名，在某一整道址IPK附近計算出精確峰位PS。SUBROUTINE PKFD3（TRH1，TRH2，TRH3，MB，MF，NP3）DIMENSION C（5Ø）COMMON/CURRP/FW，SG，CDL，CDR，STEP，NOSP，NPK，NPS，NPF，NN，NMCOMMON/FPDATA/IP（15），RP（1Ø）COMMON/PSCB/NCSH，CHPK（15），CSIGMA（15），CCDL（15），CCDR（15），CSTEP（15）COMMON/ARRAY/SP（21ØØ），PKL（16Ø），T（128）EQUIVALE

31、NCE（T（1），C（1）UMØ·ØXØ.5*（MBMF）NPSCALL GETCRP（X，IERRi）FW2.355*SGN1=IFIX（1.1*FW）DØ6ØØ5153*FWD2D*DD22*N11DO1ØJ1，N2J1JN11C（J）EXP（J1*J1/D2）SUMSUMC（J）1ØCONTINUES1FLOAT（N2）SUMSUM/S1DO5ØJ1，N25ØC（J）C（J）SUMW2Ø.ØDO 2ØØ IMB，MFDVØ.

32、16;DVVØ.ØM1IN1M2IN1DO 8ØJM1，M2J1JM11C1SP（J）*C（J1）DVDVC18ØDVVDVVC1*C（J1）DVVSQRT（DVV）SDV/DVVIF（W2.LE.S）GOTO 15ØIF（W2.LE.W1）GOTO 15ØIF（W2.LE.TRH1）GOTO 15ØNP3NP31IPKI1CALL FDPS（W1，W2，S，IPK，PS）PKL（NP3）PS15ØW1W2W2S20ØCONTINUERETURNENDSUBROUTINE FDPS（A，B，C，K，PS）

33、W12*BAIF（ABS（W1）.LT.1.ØE1Ø）GOTO 2ØPSKØ.5*（CA）/W1GOTO 3Ø2ØPSK3ØRETURNEND在PKFD3子程序中，首先根據(jù)譜區(qū)內(nèi)平均的FWHM值計算匹配濾波器的沖擊響應(yīng)函數(shù)值C（J），然后逐道計算變換譜，并找出變換譜中出局部極大值的整數(shù)道址。如果在該點(diǎn)由式（5.1.17）計算的Rm值超過了尋峰閾值TRH1，則找到了一個峰，峰位的整數(shù)道址為IPK。在IPK附近，由變換譜用二階差值多項式求出精確峰位PS，并存入PKL數(shù)組中。在匹配濾波器法尋峰的計算過程中，需要重復(fù)地進(jìn)行指數(shù)運(yùn)算，

34、因而運(yùn)算速度較慢。為了簡化運(yùn)算提出了以沖擊函數(shù)為矩形的濾波器代替匹配濾波器，這就是矩形濾波器法5。圖5-1-7 矩形濾波器的沖擊響應(yīng)函數(shù)矩形濾波器的沖擊函數(shù)如圖5-1-7所示，Cj的表達(dá)式為（5.1.19）其中N為正奇數(shù)，等于變換窗口的三分之一。一般N值取近似于峰的FWHM。為了進(jìn)一步提高尋找弱峰的能力，也可以把Cj取為（5.1.20）和匹配濾波器法相同，用矩形濾波器法尋峰時也是按道址順序地進(jìn)行譜交換，在變換后的譜中找出一系列的局部極大值。然后用公式（）進(jìn)行峰高判定。還可以利用變換譜中局部極大值兩側(cè)的局部極小點(diǎn)之間的距離進(jìn)行峰寬判定。當(dāng)滿足判定條件時，確認(rèn)該峰是一個真峰。該局部極大值所對應(yīng)的道

35、址就是峰位。矩形濾波器法尋峰的最大優(yōu)點(diǎn)是計算方法簡單，運(yùn)算速度快，特別適用于運(yùn)算速度較慢的微型計算機(jī)系統(tǒng)。這種方法也能在統(tǒng)計漲落較大的高本底上尋找弱峰。和匹配濾波器法一樣，矩形濾波器尋峰方法的重峰分辨能力也比較差。下面是用FORTRAN語言編寫的矩形濾波器法尋峰子程序PKFD5。程序中使用的變量名和子程序名與前面PKFD3中使用的變量和子程序名相同。SUBROUTINE PKFD5（TRH1，TRH2，TRH3，MB，MF，NP3）COM/MON/CURRP/FW，SG，CDL，CDR，STEP，NOSP，NPK，NPS，NPF，NN，NMCOMMON/ARRAY/SP（2100），PKL（1

36、60），T（128）W20.0DDV0.0XINPSCALL GETCRP（X，IERR1）FW2.355*SGN1IFIX（1.1*FW）MT（N12）/2M1I3*MT2M2I3*MT2DO 2Ø JM1，M2IF（J.GT.IMT）GOTO 105DVDVSP（J）DDVDDVSP（J）GOTO 201ØIF（J.GT.IMT1）GOTO 5DVDV2.0*SP（J）DDVDDV4.0*SP（J）2ØCONTINUEDDVSQRT（DDV）SDV/DVVIF（W2.LE.S）GOTO 50IF（W2.LE.W1） GOTO 50IF（W2.LE.TRH1）

37、GOTO 50NP3NP31IPKI1CALL FDPS（W1，W2，S，IPK，PS）PKL（NP3）PS5ØW1W2W2S1ØØCONTINUERETURNEND(2) 簡單比較法尋峰6；簡單比較法尋峰是最直觀而又快速的一種求峰方法。在譜數(shù)據(jù)中，某一道的數(shù)據(jù)比其鄰近的幾道大很多時，則認(rèn)為該道存在一個峰。在計算中邊疆檢索平滑后的譜數(shù)據(jù)，如在第m道滿足()則第m道附近有一個峰。式中TRH是尋峰閾值。在第m道附近的譜數(shù)據(jù)中用二階差值多項式計算出精確峰位。簡單比較法適于尋找強(qiáng)單峰，在高本底上尋找弱峰和分辨重峰的能力都比較差。但算法簡單，程序運(yùn)行速度較快。下面給出用FO

38、RTRAN語言編寫的簡單比較法尋峰子程序PKFD4。其中調(diào)用了利用二階差值多項式計算精確峰位的子程序FDPS。SUBROUTINE PKFD4（TRH1，TRH2，TRH3，MB，MF，NP3）COMMON/ARRAY/SP（2100），PKL（16Ø），T（128）MØDO 2Ø IMB，MFIF（SP（I）TRH1*SQRT（SP（I）IF（SP（I+2）.GT.W1）GOTO 20IF（SP（I2）.GT.W1） GOTO 20N1IS1SP（I）DO 10 J1，3，2N2N12JIF（SP(N2).LE.S1）GOTO 10N1N12JS1SP（N1）1

39、0CONTINUEIF（N1.EQ.M）GOTO 20MN1NP3NP31CALL FDPS（SP(N11)，SP(N1+1)，N1，PS）PKL（NP3）PS20CONTINUERETURNENDSUBROUTINE FDPS（A，B，C，K，PS）W12*BACIF（ABS(W1)）.LT1.0E10）GOTO 20PSK0.5*(CA)/W1GOTO 3020PSK30RETURNEND(3) 導(dǎo)數(shù)法尋峰：如果我們把譜看成是一個連續(xù)曲線，譜曲線和它的一、二、三階導(dǎo)數(shù)的圖形如圖5-1-8所示。由圖中可見，在峰頂處一階導(dǎo)數(shù)值由正到負(fù)過零點(diǎn)。二階導(dǎo)數(shù)在峰位處出現(xiàn)負(fù)的局部極小值。三階導(dǎo)數(shù)在峰位附

40、近由負(fù)到正過零點(diǎn)。因而利用譜曲線在峰位附近形狀的特征，由譜曲線的斜率，曲率的變化可以準(zhǔn)確地定出峰位。值得指出的是，利用譜曲線的曲率變化可以分辨重峰，確定重疊中各個組分峰的峰位。由圖5-1-9中可以看出，由于兩個組分峰相距很近，疊加而成的重峰只有一個極值。用二階和三階導(dǎo)數(shù)的變化特征仍然可以準(zhǔn)確地定出每個組分峰的峰位。在二個組分峰的峰位處，出現(xiàn)負(fù)的極小值，由負(fù)到正過零點(diǎn)。因此，用二階導(dǎo)數(shù)法和三階導(dǎo)數(shù)法尋峰能夠分辨重疊峰中各個組分峰的峰位。此外，利用峰和康普頓邊緣形狀的不同，也能很好地把它們區(qū)分開。由于上述特點(diǎn)，以及計算方法比較簡單，導(dǎo)數(shù)法是目前使用最廣泛的一種尋峰方法。圖5-1-8 譜曲線及其各階

41、導(dǎo)數(shù)的圖形 圖5-1-9 重疊峰及其各階導(dǎo)數(shù)的圖形 一階導(dǎo)數(shù)法導(dǎo)峰7：一階導(dǎo)數(shù)法尋峰的計算步驟如下。首先用公式（）計算出各道的平滑后的一階導(dǎo)數(shù)值。按道址順序檢索出一階導(dǎo)數(shù)值并求出由正到負(fù)的零點(diǎn)。如果第m道的一階導(dǎo)數(shù)值為正，第m1道的一階導(dǎo)數(shù)值為負(fù)，則有一個峰的峰位在m，m1道之間。用一階導(dǎo)數(shù)值的線性插值求出峰位的精確數(shù)值。當(dāng)找到一個峰位時，在此峰位的左側(cè)求出一階導(dǎo)數(shù)值的局部正極值的道址m1，在峰位的右側(cè)求出一階導(dǎo)數(shù)負(fù)極小值對應(yīng)的道址m2。令W2m2m1，W2代表了峰的寬度。只有滿足峰寬判定條件。0.8 FWHMW23FWHM()時，才認(rèn)為該峰可能是一個真峰。W2太小時，該峰可能是由統(tǒng)計漲落造成

42、的假峰。W2太大時，該峰可能是一個康普頓邊緣。在這兩種情況下，找到的峰都不是真峰，必需剔除。在實際應(yīng)用中公式（）中的系數(shù)3和0.8需要適當(dāng)?shù)剡M(jìn)行調(diào)整，以適合該測量系統(tǒng)的具體情況。第二個判別條件是峰高判定。峰高判定條件為()式中，mp是一階導(dǎo)數(shù)過零點(diǎn)的道址；ymp是峰位處的譜數(shù)據(jù)；Bmp是峰位處的本底值；TRH為預(yù)先設(shè)定的尋峰閾值，一般取24之間。在實際計算中為了避免本底計算的誤差，對公式（）進(jìn)行如下的簡化：設(shè)譜曲線為連續(xù)曲線，峰函數(shù)為高斯函數(shù)，連續(xù)譜函數(shù)為()式中x為道址，X為峰位，A為峰的凈高度，B為本底函數(shù)。由于本底函數(shù)的斜率很小，譜的一階導(dǎo)數(shù)近似等于而()由(5.1.23)、(5.1.2

43、4)、5)式可以把峰高判定條件改寫為()這個公式就是在一階導(dǎo)數(shù)法尋峰程序中實際應(yīng)用的峰高判定條件。圖5-1-10 一階導(dǎo)數(shù)法峰判定公式中的幾個參數(shù)的圖示在上面討論的一階導(dǎo)數(shù)法尋峰中，使用了式（）作為峰寬判定條件。實際應(yīng)用中發(fā)現(xiàn)這個判定條件在由于測量條件的變化而導(dǎo)致峰的形狀改變時常常不易滿足，從而造成假峰或漏失真峰。為此文獻(xiàn)3中采用了一種先由譜中找出峰區(qū)的左、右邊界，然后再在該峰區(qū)內(nèi)用一階導(dǎo)數(shù)法導(dǎo)峰的方法。其具體步驟如下：首先確定峰區(qū)的左、右邊界首址。沿道址增加方向檢索譜數(shù)據(jù)及其一階導(dǎo)數(shù)值，當(dāng)譜的一階導(dǎo)數(shù)值為正，其值大于某個數(shù)值，且隨著道址增加，譜的一階導(dǎo)數(shù)值也增大時，該道址mL就是峰區(qū)的左邊界

44、道址。在mL道應(yīng)滿足下述條件()(8)(9)中期，ymL是mL道的譜數(shù)據(jù)。是mL道平滑譜的一階導(dǎo)數(shù)值（參見式5-1-10和表5.2）。K為常數(shù)，是對應(yīng)于某一置信度時正態(tài)分布的橫坐標(biāo)。K值與置信度的關(guān)系如表5.5所示。用于尋峰程序中的置信度稱為峰置信度，可以預(yù)先指定，一般在0.600至0.999之間。當(dāng)峰置信度選為0.90時，能夠剔除假峰，而且能夠找到大多數(shù)真峰。中mL道一階導(dǎo)數(shù)值的標(biāo)準(zhǔn)偏差。當(dāng)變換窗口為五點(diǎn)時()表5.5 K值與置信度的關(guān)系置信度K0.9993.0900.9952.5760.9902.3260.9751.9600.9501.6450.9001.2820.8000.8420.75

45、00.6750.7000.5250.6000.2540.5000在譜區(qū)中沿道址增加的方向檢索，在某一道址中，公式（）、（8）、（9）滿足時，此道址就是峰區(qū)的左邊界道址mL，當(dāng)?shù)乐防^續(xù)增大，在某一道址中滿足公式（）、（2）、（3）時，認(rèn)為該道址是峰區(qū)的右邊界道址mR()(8)(9)在峰區(qū)的左、右邊界確定之后，需要檢驗峰區(qū)的寬度。當(dāng)下述公式()滿足時，認(rèn)為是一個正常峰區(qū)，否則認(rèn)為峰區(qū)太窄。而式（）中W是預(yù)先估計的最小峰區(qū)寬度。峰區(qū)確定之后，在mL至mR之間對譜數(shù)據(jù)的平滑一階導(dǎo)數(shù)值進(jìn)行檢索。當(dāng)一階導(dǎo)數(shù)值由正變負(fù)時，則在該道附近存在一個峰。在該道附近對一階導(dǎo)數(shù)值進(jìn)行線性內(nèi)插，求出一階導(dǎo)數(shù)值過零點(diǎn)的精確

46、位置就可以計算出精確峰位。下面給出在峰區(qū)內(nèi)用一階導(dǎo)數(shù)法尋峰的子程序PKFD1。在程序中，實數(shù)組SP中存放原始譜數(shù)據(jù)，實數(shù)組DIF1中存放譜數(shù)據(jù)各點(diǎn)的平滑一階導(dǎo)數(shù)值。為了節(jié)省內(nèi)存，另外設(shè)置了一個長度為256的數(shù)組SPA。首先調(diào)用子程序GETSP把譜數(shù)據(jù)由SP分段地傳送到SPA中。再調(diào)用子程序DIFER1，計算這一段譜的平滑一階導(dǎo)數(shù)值，計算結(jié)果存放在DEF1中。由SPA和DIF1計算峰區(qū)的左、右邊界（IPKSTA、IPKEND），在峰區(qū)內(nèi)判定一階導(dǎo)數(shù)值由正到負(fù)的零點(diǎn)，并用線性內(nèi)插法計算精確峰位，精確峰位的計算結(jié)果存放在PKL中。程序中使用了兩個尋峰閾值：TRH1用于決定峰區(qū)的峰置信度，由操作員根據(jù)

47、譜中峰的形狀指定。SUBROUTINE PKFD1（TRH1，TRH2，TRH3，MB，MF，NP3）LOGICAL IWTESTDIMENSION SPA（256），DIF1（256）COMMON/CURRP/FW，SG，CDL，CDR，STEP，NOSP，NPK，NPS，NPF，NN，NMCOMMON/FPDATA/IP（15），RP（10）COMMON/PSCB/NCSH，CHPK（15），CSIGMA(15)，CCDL(15)，CCDR(15)，CSTEP(15)COMMON/ARRAY/SP（2100），PKL（160），T（128）FSIG（M1）·1*SQPT(4

48、83;*(SPA(M1+2)SPA(M1-2)SPA(M1+1)SPA(M1-1)IGMBNSP25610IF(IF#GE#MF) GOTO 700CALL GETSP（SPA，IG，NSP，MF）CALL DIFER1（SPA，DIF1，NSP）NSP1NSP-3IWTRH2I3IWTEST1W.TL.1050DO 300L1I,NSP1IF(DIF1(L1).LE.0) GOTO 300IF(SPA(L1).LE.0.001) GOTO 300IF(SPA(L1+1).LT.SPA(L1) GOTO 300IF(DIF1(L1+1).LE.DIF1(L1)-FSIG(L1)GOTO 300

49、DTTRH1*SQRT(SPA(L1)IF(DIF1(L1).LE.DT) GOTO 300DT1TRH1*SQRT(SPA(L1+1)IF(DIEQ(L1+1).LE.DT1) GOTO 300GOTO 310300CONTINUEGOTO 31010 IPKSTAL1ISTOPNSPIF(.NOT.IWTEST)ISTOPNSP-1320IF(L1+2.GT.ISTOP) GOTO 600M1L1+2330IF(SPA(M1).LE.0.001) GOTO 200IF(DIF1(M1).GE.0.0) GOTO 390SIGF0.0IF(M1.LT.NSP-2)SIGFFSIG(M1)IF

50、(DIF1 (M11).GT.DIF1 (M1)SIGF) GOTO 390IF(IWTEST)GOTO 340IF(DIF1(M1)-SIGF.LE.SQRT(SPR(M1) GOTO 390340IF(-DIF1(M1)-SIGF.LE.SQRT(SPR(M1)GOTO 200390IF(M1.EQ.ISTOP) GOTO 600M1M11GOTO 330200IPKENDM1WIPKEND-IPKSTA2) GOTO 230JIPKSTAGOTO 250230IIPKENDGOTO 50250JJ1IF(J.GE.IPKEND) GOTO 230IF(DIF1 (J).GE.0.0) G

51、OTO 250NP3NP31IF(NP3.GT.158) GOTO 700PKL(NP3)J2IGDIF1(J1)/(DIF1(J1)-DIF1(J)IF(PKL(NP3).GT.MF) GOTO 700400JJ1IF(J.GE.IPKEND) GOTO 230IF(DIF1 (J).kE.0.0) GOTO 400GOTO 250500LGNSPIG-10GOTO 10600IF(IPKSTA.LE.100) GOTO 650IGIPKSTALG10，GOTO 10650IGIG100GOTO 10700 RETURNENDSUBROUTINE DIFER1(SPA,DIF1,NSP)DI

52、MENSION SPA(1)，DIF1(1)COMMON/CURRP/FW,SG,CDL,CDR,STEP,NOSP,NPK,NPS,NPF,NN,NMDO 10 I1,NSP10DIF1(I)0.0DIF1(1)SPA(2)SPA(1)DIF1(2)0.5*(SPA(3)SPA(1)DIF1(NSP)SPA(NSP)SPA(NSP1)DIF1(NSP1)0.5*(SPA(NSP)SPA(NSP2)NSP2NSP2IF(FW.LT.4) TOGO 30DO 20 I3,NSP220DIF1(I)0.1*(2.0*SPA(I2)SPA(I1)SPA(I1)SPA(I1)2.0*SPA(I2)RE

53、TURN30DO 40 I3,NSP240DIF1(I)0.5(SPA(I1)SPA(I1)RETURNEND一階導(dǎo)數(shù)尋峰方法是最常用的一種尋峰方法。用這種方法能找出大部分單峰。計算方法較簡單，運(yùn)算速度比較快。其主要缺點(diǎn)是不能分辨相距很近的重疊峰。由圖5-1-9可以看出，當(dāng)兩個峰相距很近時，譜曲線只有一個極值，一階導(dǎo)數(shù)中只出現(xiàn)一個由正到負(fù)的過零點(diǎn)。在這種情況下，一階導(dǎo)數(shù)尋峰方法把由兩個組分峰重疊而成的重峰當(dāng)成了一個單峰，從而漏失了一個有意義的峰。另一方面，當(dāng)譜數(shù)據(jù)中的統(tǒng)計漲落比較大時，使用一階導(dǎo)數(shù)尋峰方法可能找到很多個由統(tǒng)計漲落造成的假峰。在用較小的尋峰閾值尋找弱峰時，出現(xiàn)假峰的現(xiàn)象更為嚴(yán)重。對找到的峰進(jìn)行凈面積判定是降低假峰出現(xiàn)幾率的有效方法。當(dāng)峰的凈面積比峰的總面積（峰的凈面積和本底面積之和）的標(biāo)準(zhǔn)偏差大若干倍時，凡確認(rèn)該峰是一個真峰，否則認(rèn)為它是假峰，予以剔除。二階導(dǎo)數(shù)法尋峰：由圖5-1-8中可以看出，在峰頂處由于譜曲線有最大的負(fù)曲率，譜曲線的二階導(dǎo)數(shù)在峰位處出現(xiàn)負(fù)的局部極值。譜曲線在峰位附近的這個特征可以用來確定峰位。由圖5-1-9也可以看出，對于一個重疊峰，雖然譜曲線只有一個局部極大值，識別