唐 琳,趙衛(wèi)東，余松科,劉 澤,余小東，孟 源,黃興祿

1. 成都大學(xué)電子信息與電氣工程學(xué)院，四川 成都 610106 2. 數(shù)學(xué)地質(zhì)四川省重點(diǎn)實(shí)驗(yàn)室(成都理工大學(xué))，四川 成都 610059 3. 成都大學(xué)科研處，四川 成都 610106 4. 成都大學(xué)計(jì)算機(jī)學(xué)院，四川 成都 610106

引 言

輻射測(cè)量中最重要的兩個(gè)指標(biāo)，一個(gè)是能量分辨率，另一個(gè)則是計(jì)數(shù)率。但能量分辨率和計(jì)數(shù)率是相互矛盾的，改善能量分辨率將降低計(jì)數(shù)率，相反提高計(jì)數(shù)率將損失能量分辨率[1-3]。探測(cè)器的能量分辨率由多個(gè)因素決定，比如探測(cè)信號(hào)的統(tǒng)計(jì)漲落、信號(hào)處理器的噪聲、外部干擾以及溫度等。統(tǒng)計(jì)漲落決定了給定探測(cè)器能量分辨率的理論極限，而其他因素的影響則可以通過適當(dāng)?shù)脑肼暈V除和電子技術(shù)來降低[4-6]。

以往關(guān)于能量分辨率的研究大多利用譜反卷積對(duì)獲取到的能譜進(jìn)行后處理，從而降低特征峰的半高寬。這些后處理方法是通?；趯@取到的能譜建模為輸入能譜和探測(cè)器響應(yīng)函數(shù)這兩個(gè)隨機(jī)變量的函數(shù)，而探測(cè)器的響應(yīng)函數(shù)又可以建模為輸入脈沖高度與輸出脈沖高度的聯(lián)合概率分布。目前至少存在三種對(duì)譜線進(jìn)行后處理的方法可以提高能量分辨率，包括使用正則化、極大似然和最大熵。然而，這三種后處理方法都涉及了為光譜反褶積建模探測(cè)器的能量響應(yīng)函數(shù)，因此這些方法將是一項(xiàng)計(jì)算量非常大的任務(wù)[7]。

傳統(tǒng)的X射線光譜數(shù)據(jù)處理平臺(tái)采用的是多道譜的成譜方法，其主要特征在于譜圖中每一個(gè)道址上計(jì)數(shù)的增量是基于單個(gè)脈沖幅度的，也就是說，每一個(gè)脈沖幅度被接收，相應(yīng)的道址上的計(jì)數(shù)值就會(huì)加一。早在2005年，國(guó)外就有研究人員初步提出了多脈沖成譜的相關(guān)方法[8]，但國(guó)內(nèi)尚無相關(guān)成果公布。筆者在前期研究工作[9-10]中已經(jīng)提出了脈沖剔除法和脈沖修復(fù)法處理突變脈沖得到真實(shí)可靠的脈沖幅度信息，并因此降低了統(tǒng)計(jì)漲落。本文在上述研究基礎(chǔ)上提出一種X射線光譜數(shù)據(jù)處理平臺(tái)的優(yōu)化設(shè)計(jì)，采用可在線實(shí)施的多脈沖局部平均(MPLA)成譜技術(shù)，選擇一個(gè)恰當(dāng)范圍的平均窗口對(duì)多個(gè)脈沖幅度進(jìn)行平均，平均后的幅度值被接收，再將相應(yīng)的道址上的計(jì)數(shù)值加一，該技術(shù)的實(shí)施使得X射線光譜數(shù)據(jù)處理平臺(tái)得以優(yōu)化。本文通過理論推導(dǎo)揭示了MPLA概率密度變換的基本理論特征，也通過實(shí)驗(yàn)驗(yàn)證了MPLA算法的特征，證明了在具有正態(tài)分布PDF的頻譜峰值的典型情況下，即使僅對(duì)兩個(gè)脈沖高度進(jìn)行平均，變換后峰的FWHM也會(huì)變窄。

1 原 理

MPLA算法與傳統(tǒng)多道譜的根本區(qū)別在于前者使用多個(gè)脈沖高度的平均值來確定需要增加計(jì)數(shù)的道址，而后者則是每讀取一個(gè)脈沖幅度值，就給對(duì)應(yīng)道址上的計(jì)數(shù)加一。MPLA算法的核心在于多個(gè)脈沖幅度的平均，但這里的平均并不是針對(duì)所有的脈沖幅度。如果對(duì)所有測(cè)量到的脈沖幅度都進(jìn)行平均，那么最終得到的譜圖就會(huì)完全失去其原有的統(tǒng)計(jì)特征，也會(huì)破壞不同元素特征峰之間的計(jì)數(shù)差異，因此MPLA算法中的平均必須是針對(duì)選定窗口內(nèi)一定脈沖幅度數(shù)量進(jìn)行平均，再通過平均得到的脈沖幅度值來確定增加計(jì)數(shù)的道址。

MPLA算法涉及兩項(xiàng)可變參數(shù)，一是平均窗口的大小，本文以參數(shù)r表示，另一項(xiàng)參數(shù)則是每一次平均的脈沖幅度數(shù)量，本文以參數(shù)n表示。在下文對(duì)MPLA算法原理的描述中，我們?nèi)∑骄翱趓等于2，每一次參與平均的脈沖幅度數(shù)量為2。

在MPLA算法的執(zhí)行流程包含以下幾個(gè)步驟，如圖1所示。首先讀取第一個(gè)脈沖幅度時(shí)定位一個(gè)平均窗口，當(dāng)r為2時(shí)窗口內(nèi)包含的道址數(shù)為5，脈沖幅度讀取成功后更新當(dāng)前平均窗口的幅度和脈沖個(gè)數(shù)，每次更新后即對(duì)平均窗口內(nèi)的脈沖個(gè)數(shù)進(jìn)行判斷，當(dāng)其小于預(yù)設(shè)的參數(shù)n時(shí)，則繼續(xù)讀取下一個(gè)脈沖幅度，當(dāng)平均窗口內(nèi)的脈沖幅度數(shù)量等于參數(shù)N時(shí)則對(duì)相應(yīng)平均窗口內(nèi)的脈沖幅度進(jìn)行平均，得出的平均數(shù)即為需要更新計(jì)數(shù)的道址，然后再對(duì)取平均值的窗口內(nèi)脈沖幅度和脈沖個(gè)數(shù)進(jìn)行清零。

為了詳細(xì)闡述MPLA算法的實(shí)現(xiàn)過程，下文將以圖文結(jié)合的方式呈現(xiàn)一個(gè)MPLA實(shí)現(xiàn)示例。在當(dāng)前示例中，MPLA算法的平均窗口r設(shè)置為2，每一次參與平均的脈沖幅度數(shù)量設(shè)置為2。在圖2所示的MPLA算法步驟一中，讀取到的第一個(gè)脈沖幅度值為660，平均窗口的中心則為660，而窗口大小為2則代表窗口的范圍被限定在658～622之間。

此時(shí)活動(dòng)平均窗口僅有一個(gè)，并且該窗口內(nèi)可用于平均的脈沖個(gè)數(shù)小于2，開始執(zhí)行MPLA算法示例的步驟二，如圖3所示，本次讀取的脈沖幅度為650，平均窗口中心移到650，平均窗口范圍為648～652，該范圍與前一個(gè)脈沖所處的范圍并無重疊，這就產(chǎn)生了一個(gè)新的平均窗口，對(duì)該窗口內(nèi)的脈沖幅度和以及脈沖個(gè)數(shù)進(jìn)行更新后再次判斷是否有達(dá)到平均條件的窗口，若沒有繼續(xù)執(zhí)行步驟三。

圖1 MPLA算法實(shí)現(xiàn)流程Fig.1 MPLA algorithm implementation flow

圖2 MPLA算法實(shí)現(xiàn)示例步驟一Fig.2 Step 1 of MPLA algorithm implementation example

在執(zhí)行完步驟一和步驟二后可以發(fā)現(xiàn)，在執(zhí)行多脈沖平均的過程中，可能會(huì)存在多個(gè)活動(dòng)的平均窗口，如果下一個(gè)讀取的脈沖高度正好在包含平均和的兩個(gè)通道之間，則會(huì)隨機(jī)將脈沖高度添加到兩個(gè)通道之一； 如果下一個(gè)讀取的脈沖幅度所產(chǎn)生的平均窗口覆蓋了該區(qū)域內(nèi)已經(jīng)存在的某個(gè)脈沖幅度和，那么則用新的脈沖幅度與原有脈沖幅度相加，并將該平均窗口內(nèi)的脈沖個(gè)數(shù)加一。

如圖4所示，在讀幅度為659的脈沖之前，所示脈沖幅度和的區(qū)域內(nèi)已經(jīng)有650和660兩個(gè)脈沖幅度和，每一個(gè)脈沖幅度和對(duì)應(yīng)的脈沖數(shù)量都為一，還未達(dá)到窗口內(nèi)取平均的條件。下一個(gè)讀取的脈沖幅度659產(chǎn)生的平均窗口范圍為657～661，剛好包含了原本存在660，因此更新后的脈沖幅度和就為1 319，并且該窗口內(nèi)脈沖個(gè)數(shù)為2，達(dá)到了取平均的條件。通過脈沖幅度和除以脈沖個(gè)數(shù)計(jì)算出平均脈沖幅度為660，已經(jīng)被平均的脈沖幅度和以及脈沖個(gè)數(shù)清零，其余窗口內(nèi)的計(jì)數(shù)保持不變，同時(shí)將取平均值得到的脈沖幅度所對(duì)應(yīng)的道址上的計(jì)數(shù)加一，如圖5所示。

圖3 MPLA算法實(shí)現(xiàn)示例步驟二Fig.3 Step 2 of MPLA algorithm implementation example

圖4 MPLA算法實(shí)現(xiàn)示例步驟三Fig.4 Step 3 of MPLA algorithm implementation example

圖5 MPLA算法實(shí)現(xiàn)示例步驟四Fig.5 Step 4 of MPLA algorithm implementation example

2 理論推導(dǎo)與分析

理論上來說，傳統(tǒng)的多道譜(MCA)體現(xiàn)的是能量上的概率密度(pdf)，每一個(gè)脈沖幅度都對(duì)最終得到的譜圖的計(jì)數(shù)率有所貢獻(xiàn)，而MPLA將這種概率密度作為輸入，并通過在指定范圍內(nèi)的平均將其轉(zhuǎn)換為相同能量上新的概率密度。

設(shè)r為窗口平均參數(shù)，n為平均參數(shù)的脈沖高度數(shù)。X0是一個(gè)通過多脈沖幅度和本地取平均值得到脈沖幅度，且X0的范圍為[X0-r,X0+r]。令概率密度函數(shù)為f(x)，累積分布函數(shù)為F(x)。若X1是一個(gè)隨機(jī)變量描述的落在數(shù)理統(tǒng)計(jì)范圍內(nèi)的第二脈沖幅度，則X1與X0有相同的概率密度函數(shù)，但區(qū)間為[X1-r,X1+r]。道址上的計(jì)數(shù)將在采集到2個(gè)數(shù)理統(tǒng)計(jì)范圍內(nèi)的脈沖幅度時(shí)進(jìn)行遞增，相應(yīng)的隨機(jī)變量P0=(X0+X1)/2。輸入概率密度函數(shù)f(x)將完全指定X0和X1的分布，以及MPLA轉(zhuǎn)換產(chǎn)生的P0的概率密度函數(shù)。

因?yàn)閒(x)是關(guān)于μ對(duì)稱的，所以對(duì)于任意值b可以得出f(μ-b)=f(μ+b)。同理，如果f(x)是關(guān)于μ對(duì)稱的，那么fp(y)也關(guān)于μ對(duì)稱。如果f(x)是對(duì)稱的，并且不斷的趨于平均值μ，那么P0的方差小于X0。推導(dǎo)過程如下：

在不損失一般性的前提下，我們假設(shè)為μ=0。該理論簡(jiǎn)化后為式(1)所示。

(1)

表達(dá)式(1)的左側(cè)可進(jìn)一步展開為式(2)所示。

(2)

(3)

由于f(x)是對(duì)稱的，并且一直朝著均值遞增(此處我們假定了均值為0)，我們可以得出假如x0=μ=0那么E[X1|X0=x0]=x0，同樣的，假如x0≠μ那么|E[X1|X0=x0]|<|x0|，在此基礎(chǔ)上，式(3)可進(jìn)一步推導(dǎo)出式(4)。

(4)

上述推論表明，不管MPLA的參數(shù)值如何，MPLA轉(zhuǎn)換都不會(huì)改變分布，并保持對(duì)稱性。如果原始分布是單個(gè)高斯峰，則所得分布將具有相同的平均值，維持峰值位置，并且也將是對(duì)稱的。理想的MPLA轉(zhuǎn)換降低了原始分布的半高寬，銳化了峰形。對(duì)于單個(gè)高斯峰，半高寬與分布的方差有關(guān)。

3 實(shí)驗(yàn)部分

MPLA算法是一種在線實(shí)時(shí)進(jìn)行的脈沖處理方法，其實(shí)施平臺(tái)在硬件上包括前端的探測(cè)器、后端的PC處理軟件以及核心部分的數(shù)字信號(hào)處理器。如前文所述，MPLA算法本身并沒有使用范圍的限制，它是將傳統(tǒng)的概率密度作為輸入，并通過在指定范圍內(nèi)的平均將其轉(zhuǎn)換為相同能量上新的概率密度，其本質(zhì)上就是一種概率密度的轉(zhuǎn)換，并沒有對(duì)探測(cè)器有任何的限定。本文所采用的探測(cè)器是AMPTEK生產(chǎn)的高性能FAST-SDD探測(cè)器。數(shù)字信號(hào)處理器包括信號(hào)轉(zhuǎn)換、多級(jí)放大、模數(shù)轉(zhuǎn)換等前端處理單元、脈沖幅度分析單元、MPLA單元以及脈沖成譜單元，其結(jié)構(gòu)框圖如圖6所示。

圖6 數(shù)字信號(hào)處理器結(jié)構(gòu)圖Fig.6 Structure chart of the digital signal processor

數(shù)字信號(hào)處理器的前端處理部分不是本文講述的重點(diǎn)，因此下文中將重點(diǎn)介紹MPLA單元的硬件實(shí)現(xiàn)過程。MPLA單元在FPGA中實(shí)現(xiàn)，采用的芯片是Xilinx公司Spartan-3家族的XC3S400-4TQG144C。如圖6所示，MPLA單元是整個(gè)數(shù)字信號(hào)處理器的核心部分，其算法的實(shí)現(xiàn)主要包含了算法控制、平均窗口定位、窗口內(nèi)平均以及相應(yīng)道址上的計(jì)數(shù)更新四個(gè)部分。具體的實(shí)現(xiàn)流程在前文中已經(jīng)有詳細(xì)的描述，此處不再贅述。在MPLA單元處理完成后，更新的計(jì)數(shù)則通過SPI總線傳遞給MCU，MCU與上位機(jī)之間通過CAN總線進(jìn)行通信。此外，MCU的另一SPI總線連接到前端處理電路中，用于控制前端電路中的放大參數(shù)和偏移量的調(diào)整。

4 結(jié)果與討論

由圖7可以看出經(jīng)過MPLA算法處理后的譜線半高寬明顯降低，也就是說能量分辨率有了明顯的提高，但與此同時(shí)，多脈沖窄窗口求平均的成譜方法也造成了計(jì)數(shù)率的降低。

由于放射性衰變是隨機(jī)的，所以計(jì)數(shù)率總是在圍繞一個(gè)平均值不停地?cái)[動(dòng)，即放射源在每單位時(shí)間內(nèi)發(fā)生衰變的原子數(shù)目是不同的，時(shí)多時(shí)少，有起有伏，但是它比較集中地在某一范圍內(nèi)波動(dòng)，而這種現(xiàn)象就是放射性衰變的統(tǒng)計(jì)漲落。由概率統(tǒng)計(jì)理論可知，當(dāng)統(tǒng)計(jì)平均值較大時(shí)，核衰變是服從正態(tài)分布的，這也是核衰變的重要統(tǒng)計(jì)特點(diǎn)。在實(shí)驗(yàn)環(huán)節(jié)中，以鐵礦樣品為測(cè)量對(duì)象，采用MPLA算法并進(jìn)行快速多脈沖成譜處理后的結(jié)果如圖7所示。

圖7 采用MPLA算法前后得到的譜圖Fig.7 Spectrum obtained before and after MPLA

圖中黑色譜線是通過傳統(tǒng)成譜方法得到的，紅色譜線是采用MPLA算法處理后得到的，相比之下，經(jīng)MPLA算法處理后的譜線有效地降低了特征峰的半高寬。這就驗(yàn)證了前文理論推導(dǎo)部分得出的結(jié)論，理想的MPLA轉(zhuǎn)換降低了原始分布的半高寬，銳化了峰形。

5 結(jié) 論

給出了多脈沖局部平均變換的理論描述及實(shí)現(xiàn)原理，證明了MPLA算法的一些基本特性，對(duì)具有單個(gè)峰的對(duì)稱分布，MPLA變換總是保持均值不變并減小FWHM。理論推導(dǎo)和實(shí)驗(yàn)結(jié)果表明，MPLA算法具有降低譜峰的FWHM的優(yōu)點(diǎn)。但與此同時(shí)，MPLA的一個(gè)缺點(diǎn)是與傳統(tǒng)成譜方法相比，多脈沖平均的過程中使得計(jì)數(shù)減少。MPLA譜中計(jì)數(shù)減少的總數(shù)受MPLA算法中平均窗口內(nèi)脈沖幅度的個(gè)數(shù)控制。有很多脈沖處理技術(shù)可以彌補(bǔ)這種計(jì)數(shù)損失，例如，可以收集大量的脈沖高度測(cè)量，然后重復(fù)、隨機(jī)地進(jìn)行采樣以應(yīng)用MPLA變換。筆者將在接下來的工作中對(duì)這些脈沖處理技術(shù)進(jìn)行更深入的研究，以期為MPLA算法找到一種非常契合的脈沖處理技術(shù)來補(bǔ)償計(jì)數(shù)率損失。即便存在計(jì)數(shù)損失的缺陷，但不可否認(rèn)，MPLA算法為降低探測(cè)器噪聲提供了一種實(shí)時(shí)、高效、通用的處理方法，并且具有可證明的理論保證，這對(duì)提高探測(cè)器的能量分辨率具有極大的應(yīng)用意義。