王懷平 周建斌 庹先國 王 明 汪雪元,2 馮 林

1（成都理工大學(xué) 核技術(shù)與自動化工程學(xué)院 成都 610059）

2（東華理工大學(xué) 核技術(shù)應(yīng)用教育部工程研究中心 南昌 330013）

3（東華理工大學(xué) 機械與電子工程學(xué)院 南昌 330013）

隨著高性能模數(shù)轉(zhuǎn)換（Analog to Digital Converter，ADC）芯片和可編程硬件邏輯器件的快速發(fā)展及應(yīng)用，輻射探測與能譜測量中數(shù)字核脈沖實時處理方法研究取得顯著進展。成形濾波器是用于處理探測器-前置放大器系統(tǒng)輸出核脈沖信號的核電子學(xué)模塊，它通過將輸入的核脈沖信號成形為某種特殊形狀（如類高斯、梯形等）脈沖，達到抑制電子學(xué)噪聲干擾、消除彈道虧損以及減少脈沖堆積等目的，從而提高輻射探測系統(tǒng)的能量分辨率或脈沖計數(shù)率［1］。相較于傳統(tǒng)模擬成形濾波電路，數(shù)字核脈沖信號處理具有算法設(shè)計靈活、參數(shù)調(diào)節(jié)方便、精度高、穩(wěn)定性好等特點，當(dāng)輻射探測系統(tǒng)采用高速ADC芯片獲取數(shù)字化核脈沖信號時，采用高性能可編程硬件邏輯器件（Field Programmable Gate Array，F(xiàn)PGA）實現(xiàn)數(shù)字脈沖處理算法，可以滿足核脈沖或者核數(shù)據(jù)的實時處理需要［2］。因此，處理核脈沖的數(shù)字成形濾波算法既要考慮系統(tǒng)功能需求，也應(yīng)考慮算法結(jié)構(gòu)是否易于可編程硬件邏輯器件實現(xiàn)。Jordanov等［3-4］利用卷積方法得到快速數(shù)字合成梯形（三角形）成形遞推算法以及cusp-like 成形算法，算法結(jié)構(gòu)簡單、參數(shù)容易調(diào)整、適于實時處理。Imperiale 等［5］通過對成形濾波器傳遞函數(shù)進行Z變換分別得到階躍脈沖、單指數(shù)以及雙指數(shù)脈沖信號梯形成形數(shù)字遞推公式。Wang 等［6］在分析傳統(tǒng)模擬極零相消電路功能基礎(chǔ)上，構(gòu)建具有等效功能的數(shù)字極零相消模型，并設(shè)計出兩種梯形（三角形）成形濾波算法，既能滿足輻射測量系統(tǒng)極零相消需要，又能實現(xiàn)核脈沖信號成形濾波功能。Nakhostin［7］基于模擬CR-（RC）n準(zhǔn)高斯濾波器傳遞函數(shù)的Z變換推導(dǎo)出1～4階CR-（RC）n數(shù)字遞推算法，相較于梯形成形濾波器該算法在抑制并行噪聲方面性能更優(yōu)。Song 等［8］對模擬Sallen-Key 濾波器進行分析，建立了該電路的全參數(shù)（Kirchhoff's Current Law）KCL方程，并通過對所得電路微分方程求解，推導(dǎo)出一種最優(yōu)準(zhǔn)高斯脈沖成形算法，該算法所需參數(shù)容易計算，成形脈沖沒有反沖。此外，Zhao 等［9］還研究了Sine成形濾波算法。陳世國［10］通過分析計算三角形、梯形和高斯脈沖成形方法的性能參數(shù)，得出三角形成形噪聲抑制性能最好，梯形成形彈道虧損抑制性能最好，而高斯成形綜合性能最好。

由于高斯信號的對稱性和完備性，高斯濾波方法在核信號及核數(shù)據(jù)的分析與處理中被大量使用，然而，高斯函數(shù)的復(fù)雜性使得構(gòu)造真高斯成形數(shù)字遞推算法相對困難［11］。因此，數(shù)字高斯成形算法常采用對模擬高斯成形濾波電路（如Sallen-Key、CR-（RC）n）建立的微分方程數(shù)值模型推導(dǎo)而得到，電路等效模型局限和元器件參數(shù)約束等導(dǎo)致成形脈沖并不具備對稱性，Sallen-Key 數(shù)字高斯成形后脈沖存在比較嚴(yán)重的拖尾，單級使用時還會出現(xiàn)反沖等問題［12-13］。本文提出一種類高斯成形數(shù)字卷積實現(xiàn)方法，它以梯形脈沖信號為基礎(chǔ)，通過第一次卷積運算實現(xiàn)雙極性脈沖成形，再次卷積進行累積求和從而得到左右對稱的類高斯成形脈沖，采用Z變換方法推導(dǎo)出該類高斯成形算法的數(shù)字遞推公式，利用MATLAB仿真模擬驗證了算法的有效性和可靠性，并對其幅頻特性進行了詳細(xì)分析，最后將該算法應(yīng)用于搭建的X 射線熒光測量系統(tǒng)實驗平臺，進一步驗證了算法的主要性能。該類高斯成形算法結(jié)構(gòu)相對簡單，易于可編程硬件邏輯器件實現(xiàn)。

1 卷積型類高斯成形濾波器實現(xiàn)原理

首先，構(gòu)造如下單位沖激響應(yīng)傳遞函數(shù)h1（t），以實現(xiàn)脈沖信號的雙極性變換：

式中：δ(t)為單位沖激函數(shù)；tc為脈沖信號截斷寬度。

輻射探測器-前置放大器系統(tǒng)的輸出信號通常為具有快速上升沿和緩慢下降沿的雙指數(shù)衰減脈沖信號，其理想化的數(shù)學(xué)模型可采用式（2）進行描述：

式中：A為雙指數(shù)衰減信號的幅值；τ1和τ2為雙指數(shù)衰減信號的兩個時間常數(shù)；u（t）為單位階躍函數(shù)。當(dāng)滿足τ1?τ2條件時，該模型可進一步簡化為單指數(shù)衰減信號模型，即：

以式（3）所描述的單指數(shù)衰減脈沖作為輸入信號，經(jīng)式（1）所描述的系統(tǒng)處理后，可得到雙極性變換輸出信號y1（t），即：

式（4）所描述的單指數(shù)衰減脈沖卷積雙極性成形仿真模擬如圖1所示。

圖1 單指數(shù)衰減脈沖信號卷積雙極性成形仿真模擬Fig.1 Convolutional bipolar pulse shaping simulation for uniexponential decay pulse signal

其次，信號的積分運算，可采用信號與單位階躍函數(shù)的卷積來實現(xiàn)［14］，其方法如下：

式中：參數(shù)k為比例系數(shù)。選取k=1/tc，式（6）所描述的系統(tǒng)輸出信號y2（t）的仿真模擬如圖1所示。當(dāng)輸入信號x（t）為具有左右對稱的梯形脈沖時，合理設(shè)置參數(shù)tc取值，經(jīng)式（1）所描述的系統(tǒng)處理后得到具有對稱零面積的雙極性梯形脈沖信號輸出y1（t），然后利用式（6）所描述的系統(tǒng)對y1（t）進行卷積運算（即累加求和），可成形輸出為左右對稱的準(zhǔn)高斯脈沖信號y2（t）。

2 卷積型類高斯成形數(shù)字遞推算法

基于梯形脈沖變換的數(shù)字類高斯成形算法實現(xiàn)原理如圖2所示。首先，輻射探測器-前置放大器系統(tǒng)輸出的核脈沖信號經(jīng)由高速ADC 實現(xiàn)快速采樣及數(shù)字化處理后得到核脈沖序列Vi（n）；其次，核脈沖序列Vi（n）與數(shù)字系統(tǒng)hT（n）卷積成形為梯形脈沖序列Vo1（n）；然后，梯形脈沖序列Vo1（n）經(jīng)過數(shù)字系統(tǒng)h1（n）處理后成形為具有對稱零面積的雙極性脈沖序列Vo2（n）；最后是脈沖序列Vo2（n）與數(shù)字系統(tǒng)h2（n）卷積實現(xiàn)累加求和，輸出信號Vo（n）為具有左右對稱的類高斯成形脈沖序列。

圖2 卷積型類高斯成形數(shù)字遞推算法實現(xiàn)原理框圖Fig.2 Digital recursive algorithm implementation principle block diagram for convolutional quasi-Gaussian pulse shaping filter

圖2所描述的類高斯成形數(shù)字遞推算法推導(dǎo)如下。輸入核脈沖信號采用理想的單指數(shù)衰減函數(shù)表示：

式中：Vmax為核脈沖幅度：τ為衰減時間常數(shù)。經(jīng)高速ADC數(shù)字化處理得到其脈沖序列Vi（n），即：

式中：Ts為ADC 采樣周期。式（7）的Z變換表達式如下：

其中：d=。

在連續(xù)時域內(nèi)，如圖3 所示梯形脈沖信號的數(shù)學(xué)函數(shù)表達式可寫為［5］：

圖3 梯形脈沖信號Fig.3 Trapezoidal pulse signal

式（14）所描述的梯形脈沖信號類高斯成形過程可表示為：

式（15）離散化后進行Z變換得：

式中：na=ta/Ts，nb=tb/Ts，nc=tc/Ts，Vmax為梯形脈沖幅度。當(dāng)輸入信號Vi（t）為式（7）表示的單指數(shù)衰減脈沖時，數(shù)字類高斯成形系統(tǒng)傳遞函數(shù)的Z變換表達式如下：

為了降低系統(tǒng)設(shè)計的復(fù)雜度，因此選擇級聯(lián)結(jié)構(gòu)來實現(xiàn)，式（17）可表示為：

H1（z）、H5（z）及H6（z）相當(dāng)于微分運算，H2（z）、H3（z）及H4（z）相當(dāng)于積分運算。對于線性移不變系統(tǒng)，交換子系統(tǒng)的級聯(lián)順序并不會影響最終的輸出，因此將微分H1放在第一級可防止溢出，微分H6放在最后一級有利于減小誤差［15-16］。

對式（19）進行逆Z變換，可得到單指數(shù)衰減脈沖信號成形為類高斯脈沖信號的數(shù)字遞推公式如下：

該卷積型數(shù)字類高斯成形濾波器的實現(xiàn)邏輯框圖如圖4所示。

圖4 卷積型數(shù)字類高斯成形濾波器邏輯框圖Fig.4 Logical block diagram for convolutional digital quasi-Gaussian pulse shaping filter

利用MATLAB 仿真模擬來驗證該算法的有效性和可靠性。設(shè)置輸入信號為單指數(shù)衰減脈沖，即V（it）=Vmax·exp（-t/τ）·u（t），其中Vmax=1 000 mV、τ=3.2 μs、Ts=50 ns，分別利用該類高斯成形算法和梯形成形算法進行數(shù)字核脈沖處理，選擇相同成形時間（即達峰時間）：類高斯成形算法的成形參數(shù)設(shè)置為na=16、nb=16、nc=32，成形后輸出類高斯脈沖寬度為64 個點；梯形成形算法成形參數(shù)設(shè)置為na=20、nb=44，成形后輸出梯形脈沖寬度也同樣為64 個點，其成形仿真模擬如圖5 所示，圖中黑色實線為輸入單指數(shù)衰減脈沖、藍色實線為梯形成形脈沖、紅色實線為類高斯成形脈沖。仿真模擬研究發(fā)現(xiàn)：式（20）數(shù)字遞推算法如果運算順序安排不當(dāng)，將會導(dǎo)致輸出信號發(fā)散溢出，因此需將微分運算放在第一級；梯形成形算法參數(shù)設(shè)置要求nc=na+nb，且nb≥na；該類高斯成形算法是以梯形脈沖為基礎(chǔ)的，成形參數(shù)設(shè)置要求nc≥na+nb，且nb≥na，當(dāng)nc=na+nb時，類高斯成形脈沖沒有平頂，而當(dāng)nc＞na+nb時，梯形脈沖經(jīng)過第一次卷積成形得到的雙極性脈沖，在正負(fù)極性兩個脈沖中間會輸出一個零幅值間隙（gap），其輸出為有平頂?shù)念惛咚姑}沖信號。另外，式（20）成形輸出類高斯脈沖信號的幅度除以成形參數(shù)nb可得到與輸入信號的幅度嚴(yán)格相等。

圖5 脈沖成形寬度相同時類高斯成形算法與梯形成形算法仿真模擬Fig.5 Quasi-Gaussian and trapezoidal pulse shaping algorithms simulation with identical shaping pulse width

圖5的仿真模擬還給出了在相同成形脈沖寬度條件下兩種成形算法對堆積脈沖的分離能力。堆積脈沖①兩個脈沖的峰間距tb大于或等于成形脈沖寬度tw（即tb≥tw），梯形成形或者類高斯成形后可以得到兩個完全獨立成形脈沖，并能恢復(fù)堆積脈沖的初始幅度；堆積脈沖②為兩個脈沖峰間距tw/2＜tb＜tw，類高斯成形算法在接近tw/2時依然具有較好的堆積脈沖分離能力，且能恢復(fù)堆積脈沖初始幅度，而梯形成形算法已無法實現(xiàn)堆積脈沖分離；堆積脈沖③為兩個脈沖峰間距tb≤tw/2，此時類高斯成形算法和梯形成形算法均已無法實現(xiàn)堆積脈沖分離，輸出為大于成形脈沖寬度的畸形脈沖。

3 卷積型類高斯成形濾波器的幅頻特性

類高斯成形算法在設(shè)計時注重于脈沖信號的時域特征，但對信號成形處理從頻域來看還起到了濾波作用［17-18］，根據(jù)類高斯成形濾波器和梯形成形濾波器的傳遞函數(shù)，設(shè)置與圖5相同的成形參數(shù)，兩者幅頻特性對比如圖6（a）所示，類高斯成形濾波器具有更快的阻帶衰減速度。

圖6 類高斯成形算法的幅頻特性(a) 達峰時間相同時類高斯成形與梯形成形的幅頻特性對比，(b) 不同成形參數(shù)類高斯成形的幅頻特性Fig.6 Amplitude-frequency characteristics of quasi-Gaussian pulse shaping algorithm (a) Amplitude-frequency characteristics comparison between quasi-Gaussian and trapezoidal pulse shaping algorithms with identical peaking times,(b) Amplitude-frequency characteristics of quasi-Gaussian pulse shaping algorithms with different shaping parameters

成形參數(shù)取值不同時類高斯成形濾波器幅頻特性如圖6（b）所示，na取值增大時，通頻帶減小，低頻成分幅度相對增加，高頻噪聲抑制作用增強；nc取值如果大于na和nb之和，輸出為有平頂?shù)念惛咚姑}沖，從幅頻特性來看，nc取值增大，通頻帶減小，低頻成分幅度相對增加，但也增加了高頻阻帶振蕩衰減時間。值得注意的是成形參數(shù)na和nc取值增大時，將導(dǎo)致類高斯脈沖變寬，增加了脈沖堆積概率。

4 實驗應(yīng)用及測試

實驗選用科頤維KYW2000A型X光管（額定管壓和額定管流分別為50 kV、1 mA，風(fēng)冷制冷方式）和Ag靶搭建的X射線熒光測量系統(tǒng)，設(shè)置X光管管壓和管流分別為18.8 kV、78.4 μA，X 光管出射的X射線經(jīng)過12.5 μm的Be窗濾片后照射標(biāo)準(zhǔn)Mn樣品，輻射探測器選用Amptek 的Fast SDD 探測器（型號為XR-100SDD），數(shù)據(jù)獲取采用自制20 MSPS 14 bit ADC 數(shù)據(jù)采集板。利用MATLAB 對獲取到的核脈沖數(shù)據(jù)進行離線處理，實測核脈沖信號類高斯成形（成形參數(shù)設(shè)置na=16、nb=16、nc=32，成形脈沖寬度為3.2 μs）、梯形成形（成形參數(shù)設(shè)置na=20、nb=44，成形脈沖寬度為3.2 μs），結(jié)果如圖7 所示。由圖7 可見，成形得到的獨立類高斯脈沖信號具有非常好的左右對稱性；在相同達峰時間條件下，類高斯成形算法具有更好的堆積脈沖分離與初始幅度恢復(fù)能力。

圖7 實測核脈沖信號及其類高斯、梯形成形Fig.7 Measured nuclear pulse signal and quasi-Gaussian and trapezoidal pulse shaping

在上述實驗條件下，設(shè)置600 s 測量時間，探測器輸出脈沖的計數(shù)率為20.39 k·s-1，所獲取的原始脈沖數(shù)據(jù)經(jīng)MATLAB 離線處理后，原始脈沖幅度譜、原始脈沖采用梯形成形算法以及類高斯成形算法處理后得到的能譜如圖8 所示，其能量分辨率半高寬（Full Width at Half Maximum，F(xiàn)WHM）分別為184 eV@5.89 keV、132 eV@5.89 keV 和130 eV@5.89 keV，導(dǎo)致能量分辨率差異的主要原因是：原始脈沖幅度譜在獲取脈沖幅度時沒有消除彈道虧損以及進行脈沖堆積判棄，因此造成特征峰（5.89 keV）處的脈沖幅度統(tǒng)計漲落增大，而且堆積脈沖導(dǎo)致本底增加，其能量分辨率較差；梯形成形和類高斯成形都具有恢復(fù)脈沖幅值的能力，在一定條件下還可實現(xiàn)脈沖堆積的分離，所以它們的能量分辨率要優(yōu)于原始脈沖幅度譜，兩者能量分辨率相當(dāng)，相較而言，類高斯成形所獲能譜的特征峰具有更高的峰面積和總計數(shù)，其原因在于類高斯成形具有更好的脈沖堆積分離能力。

圖8 標(biāo)準(zhǔn)Mn樣品X熒光分析實驗所得原始脈沖幅度譜(a)、梯形成形(b)及類高斯成形(c)處理能譜Fig.8 Primitive energy spectrum (a),energy spectrum with trapezoidal pulse shaping (b) and quasi-Gaussian pulse shaping (c) from X-ray fluorescence analysis experiments of standard Mn sample

設(shè)置X 射線光管管壓和管流分別為18.8 kV、78.4 μA，獲取的原始脈沖離線數(shù)據(jù)采用梯形成形和類高斯成形處理后所得能譜5.89 keV特征峰能量分辨率（FWHM）與達峰時間關(guān)系如圖9（a）所示，隨著達峰時間增大，兩種成形算法所得能譜的能量分辨率先提高后損失，達峰時間為6.4 μs時，獲得最佳能量分辨率。

圖9 能量分辨率與達峰時間的關(guān)系(a)和能量分辨率與X射線光管電流的關(guān)系(b)Fig.9 Relationship between energy resolution and peaking time (a) and the relationship between energy resolution and X-ray tube current (b)

設(shè)置X 射線光管管壓為18.8 kV，達峰時間為3.2 μs，測量時間為20 s，采用梯形成形和類高斯成形處理后所得能譜5.89 keV特征峰能量分辨率與X射線光管管流的關(guān)系如圖9（b）所示。在給定的X射線熒光分析測量系統(tǒng)條件下，設(shè)置相同的X 射線光管管壓，可認(rèn)為X 射線光管電流與探測器的脈沖計數(shù)率之間存在正相關(guān)性。因此，圖9（b）也可看作能量分辨率與探測器入射計數(shù)率之間的關(guān)系。可見，隨著X射線光管電流的增大（探測器對應(yīng)入射計數(shù)率增大），將導(dǎo)致兩者能量分辨率相應(yīng)變差，其原因主要是X射線光管電流增大引起核脈沖信號堆積概率提高。

5 結(jié)語

本文提出一種卷積型類高斯濾波成形算法，它以梯形脈沖信號為基礎(chǔ)，通過第一次卷積實現(xiàn)雙極性脈沖成形，再次卷積累加求和得到左右對稱類高斯脈沖信號輸出，利用Z變換方法推導(dǎo)了該算法的數(shù)字遞推公式，算法模型易于可編程硬件邏輯器件實現(xiàn)，適用于輻射測量實時核脈沖信號處理。采用MATLAB 仿真模擬驗證了該算法的有效性和可靠性，并研究了該類高斯成形濾波器的幅頻特性：

1）在相同達峰時間條件下，該類高斯成形算法比梯形成形算法具有更快的阻帶衰減速度；

2）類高斯成形算法參數(shù)取值要求nb≥na且nc≥na+nb，當(dāng)nc=na+nb時，輸出為無平頂?shù)念惛咚姑}沖，當(dāng)nc＞na+nb時，輸出為有平頂?shù)念惛咚姑}沖，參數(shù)na和nc增大時，通頻帶變小，低頻成分幅度相對增加，成形脈沖變寬，增加了脈沖堆積概率。分別將該類高斯成形算法和梯形成形算法應(yīng)用于X射線熒光分析測量實驗，對獲取到的核脈沖信號進行離線處理，進一步驗證了算法的有效性和堆積脈沖分離能力，該類高斯成形算法輸出類高斯脈沖信號具有很好的左右對稱性而且無下沖現(xiàn)象，達峰時間相同時，比梯形成形算法具有更好的堆積脈沖分離能力；

3）在給定X 光管管壓和管流實驗條件下，原始核脈沖信號存在較嚴(yán)重脈沖堆積，設(shè)置相同達峰時間，梯形成形和類高斯成形處理后所得能譜能量分辨率相當(dāng)，但后者具有更高的峰面積和總計數(shù)；

4）X 射線光管管壓和管流一定時，兩種成形算法所得能譜的能量分辨率均呈現(xiàn)先提高后損失趨勢；

5）X 射線光管管壓和達峰時間一定時，兩種成形算法所得能譜的能量分辨率隨著X光管電流增大（對應(yīng)探測器入射計數(shù)率）而相應(yīng)變差。

作者貢獻聲明王懷平、周建斌負(fù)責(zé)文章初稿寫作、算法設(shè)計及有效性實驗驗證；周建斌負(fù)責(zé)文章的審閱與修訂；庹先國負(fù)責(zé)文章的審閱與修訂；王明負(fù)責(zé)測量系統(tǒng)實驗平臺的搭建、實驗測試、實驗數(shù)據(jù)的收集和整理；汪雪元參與實驗應(yīng)用及測試工作的討論，提出建議；馮林主要負(fù)責(zé)文章的公式校對和輔助編輯工作。