段金松， 李 婷， 王 瑋， 麻金龍， 武旭東

（核工業(yè)北京地質(zhì)研究院， 北京 100029）

地面γ 能譜測量作為目標(biāo)區(qū)域放射性礦產(chǎn)勘查的一項重要技術(shù)手段， 其測量穩(wěn)定性與準(zhǔn)確度受到譜線漂移的影響， 而引起譜漂的主要因素為溫度變化［1］。 由于便攜式多道γ能譜儀的探測器設(shè)計多采用閃爍體與光電倍增管組合的形式， 該組合測量方式的溫度效應(yīng)明顯［2］。 閃爍體發(fā)光效率、 發(fā)光衰減時間等指標(biāo)隨溫度變化而變化； 光電倍增管暗電流脈沖、 靈敏度等基本特性受溫度影響； 此外，其外圍電路系統(tǒng)也具有一定的溫度系數(shù)。 這些因素共同導(dǎo)致探測器輸出信號幅度改變、能量峰位漂移， 對能譜數(shù)據(jù)的解析造成困難。因此， 自動穩(wěn)譜算法或穩(wěn)譜系統(tǒng)的設(shè)計與研究是十分必要的。

常見的穩(wěn)譜措施為采用放射源作為參考源的穩(wěn)譜方法， 即在探測器內(nèi)嵌一個低能γ射線參考源或α 射線參考源， 通過穩(wěn)定該參考源的等效峰位達(dá)到穩(wěn)定待測能區(qū)譜峰的目的［3-4］。 由于參考源形成的參考峰在低能區(qū)，當(dāng)環(huán)境溫度差異較大時， 譜漂的非線性導(dǎo)致穩(wěn)譜效果并不理想。 采用單色LED 作為外部光源形成參考峰成為代替上述內(nèi)嵌放射源穩(wěn)譜的一個新的方案， 系統(tǒng)根據(jù)測量得到的峰位道址變化信息， 通過調(diào)整能譜儀的程控放大器及數(shù)控高壓進(jìn)行反饋控制， 使LED 參考峰位穩(wěn)定不變， 從而達(dá)到穩(wěn)定被測特征峰的目的［5-6］。 但該方法只能校準(zhǔn)光電倍增管及外圍電路的溫漂影響， 缺少對閃爍體溫度效應(yīng)的有效補償。 目前關(guān)于使用測量能譜數(shù)據(jù)擬合尋峰， 識別放射性核素并進(jìn)行穩(wěn)譜的方法已有研究［7-8］。 但對于放射性物質(zhì)含量較低的測量對象， 穩(wěn)譜時間可能過長， 對測量結(jié)果造成誤差。 筆者使用Gaussian-LM 算法進(jìn)行擬合尋峰， 通過PID 動態(tài)調(diào)節(jié)的方式實現(xiàn)系統(tǒng)快速自動穩(wěn)譜， 為野外地面γ 能譜準(zhǔn)確測量提供技術(shù)支撐。

1 特征峰穩(wěn)譜原理

圖1 特征峰自動穩(wěn)譜流程圖Fig. 1 Flow chart of spectral auto stabilization

天然放射性環(huán)境測量過程中， 在待測能譜中總能找到K、 U、 Th 中的一個特征峰用作系統(tǒng)自動穩(wěn)譜， 其工作原理如圖1 所示。 系統(tǒng)通過多道數(shù)據(jù)采集卡對全譜數(shù)據(jù)進(jìn)行采集，依據(jù)行業(yè)標(biāo)準(zhǔn)EJ/T 584—2014 《勘查用便攜式γ 輻射儀和四道γ 能譜儀》 設(shè)置能量開窗范圍分別為 K： 1.37～1.56 MeV， U： 1.66～1.88 MeV， Th： 2.42～2.79 MeV。 當(dāng)任一核素能窗范圍內(nèi)累計計數(shù)滿足既定閾值時， 進(jìn)行數(shù)據(jù)擬合并尋峰求得當(dāng)前核素的特征峰位。 比較實測峰位與理論峰位的偏差計算當(dāng)前譜漂移量， 并設(shè)置相應(yīng)的高壓控制值。 系統(tǒng)通過PID控制算法實時調(diào)節(jié)高壓模塊增益， 實現(xiàn)測量周期內(nèi)實測特征峰在標(biāo)準(zhǔn)道址附近很小的區(qū)間內(nèi)波動， 達(dá)到穩(wěn)譜的目的。

2 擬合尋峰方法

2.1 Gaussian 函數(shù)

由于放射性統(tǒng)計漲落、 閃爍體發(fā)光效率的不均勻性及儀器固有的漲落特性， 使閃爍探測器的輸出脈沖幅度圍繞某均值小幅漲落，故 γ 能譜單能峰近似服從高斯分布［9-10］， 如式（1）所示。

式中： 變量x—γ 射線能量或能譜曲線道址；h—峰高； μ—標(biāo)準(zhǔn)峰位對應(yīng)的能量或道址； σ與半峰寬的關(guān)系—FWHM＝2.355 σ。 通過求h、 μ 和 σ 可有效確定特征峰面積、 峰位及其邊界值。 經(jīng)試驗驗證， 通常選擇 μ±2.58 σ 的置信區(qū)間， 即置信概率大于99.0%來確定特征峰邊界較為合適。

2.2 LM 算法對參數(shù)尋優(yōu)

Levenberg-Marquart 算法 （簡稱 LM 算法）是基于高斯－牛頓法和梯度下降法提出的參數(shù)向量最優(yōu)化算法［11］， 本質(zhì)是用模型函數(shù)對待估參數(shù)在其臨域內(nèi)做線性近似。 它在高斯－牛頓法基礎(chǔ)上引入阻尼因子控制每一步迭代的步長以及方向， 通過自適應(yīng)調(diào)整阻尼因子實現(xiàn)收斂， 使模型未知參數(shù)不斷逼近真值， 在進(jìn)行能譜數(shù)據(jù)擬合優(yōu)化過程中得到穩(wěn)定的結(jié)果。 迭代初期， 首先需要根據(jù)能窗范圍內(nèi)譜數(shù)據(jù)預(yù)估參數(shù)初始值。 該初始值的選取是決定算法迭代速度和計算精度的關(guān)鍵因素， 選取的初值越接近真實值， 收斂和求解全局最優(yōu)解的速度越快。

LM 算法具體實現(xiàn)流程如圖2 所示， 初始化參數(shù)（h、 μ、 σ）及阻尼因子 λ， 計算當(dāng)前參數(shù)值處的殘差向量與Jacobian 矩陣， 根據(jù)阻尼因子及 Jacobian 矩陣獲取 Hessian 矩陣及迭代步長， 由迭代步長計算新的可能的待估參數(shù)， 進(jìn)而獲取新的殘差向量， 對比殘差決定參數(shù)及阻尼因子的更新方向。 在曲線擬合實踐中， 依條件選用不同的調(diào)整系數(shù)α、 β 實現(xiàn)參數(shù)的快速收斂， α 稱為正向調(diào)整系數(shù)通常選取0.1， β 稱為逆向調(diào)整系數(shù)通常選取 10。 λ′＝α λ 或＝β λ， 其中 λ′表示調(diào)整后的阻尼因子。

3 系統(tǒng)實時穩(wěn)譜

為了實現(xiàn)系統(tǒng)快速穩(wěn)譜， 防止在校正譜漂的過程中出現(xiàn)高壓增益調(diào)節(jié)不足或過調(diào)節(jié)的現(xiàn)象， 根據(jù)上述尋峰結(jié)果， 結(jié)合系統(tǒng)控制原理， 采用軟件數(shù)字PID 控制算法進(jìn)行穩(wěn)譜處理， 通過對高壓增益調(diào)節(jié)曲線分析， 制定優(yōu)化模糊規(guī)則， 并將模糊規(guī)則參數(shù)代入PID實際計算模型中， 增強系統(tǒng)魯棒性， 其工作原理如圖3 所示。

系統(tǒng)以目標(biāo)高壓控制量D 作為PID 輸出控制量， 選取擬合峰位對應(yīng)的道址P 作為自變量， 對于單一控制變量， PID 控制算法為：

圖2 LM 算法實現(xiàn)流程Fig. 2 Levenberg-Marquart algorithm implementation process

圖3 PID 控制算法Fig. 3 PID control algorithm

式中： Pref—理論峰位對應(yīng)道址； Pn、 Pn-1、 Pn-2—第 n 次、 第 n-1 次、 第 n-2 次得到的擬合峰位道址； KP、 KI、 KD—比例單元、 積分單元和微分單元的控制因子； 其中KI的主要作用為增益控制， 增加KI可降低系統(tǒng)的穩(wěn)態(tài)誤差，提高系統(tǒng)響應(yīng)速度。 KP和KD與峰位道址P 的變化率相關(guān)， 能夠預(yù)測偏差， 產(chǎn)生超前的校正作用， 較好地改善動態(tài)性能［12］。

4 試驗結(jié)果及討論

4.1 擬合尋峰結(jié)果

選用鉀、 鈾、 釷單一核素模型針對能窗劃分范圍內(nèi)的特征峰進(jìn)行擬合， 采用LM 算法對Gaussian 擬合函數(shù)進(jìn)行參數(shù)尋優(yōu)， 結(jié)果如圖4 所示， 分別求得不同核素模型的擬合峰位、 半峰全寬（FWHM）如表 1 所示。

圖4 Gaussian 擬合及LM 算法尋優(yōu)結(jié)果Fig. 4 Gaussian fitting result and Levenberg-Marquart algorithm optimization

由表1 可知， 采用LM 算法對特征峰進(jìn)行擬合取得較好效果，R2（確定系數(shù)）均大于0.95，其擬合峰位道址與理論峰位較接近， 針對1024 道便攜式多道 γ 能譜儀， 誤差在±1 道。因此， 在儀器自動穩(wěn)譜過程中可以使用該算法進(jìn)行特征峰尋峰操作。

4.2 穩(wěn)譜效果

將便攜式多道γ 能譜儀置于可編程溫濕度環(huán)境實驗箱 （型號： CTH-WK4503-02FO）中， 開展一系列環(huán)境溫度漸變試驗， 考察基于Gaussian 與LM 算法的特征峰自動穩(wěn)譜效果。 試驗溫度范圍設(shè)置為-20～50 ℃， 實驗箱以 3 ℃/min 的速率按照 20→-20→50 ℃的順序進(jìn)行降溫和升溫操作。 分別選用鉀、 鈾、 釷單一核素模型， 以及環(huán)境本底模擬多種γ 能譜儀自動穩(wěn)譜環(huán)境， 設(shè)置單次譜數(shù)據(jù)采集時間為60 s， 實時穩(wěn)譜效果如圖5 所示。

圖5 中紅色線為特征峰位， 黑色線為高壓值。圖5a 中， 環(huán)境本底情況下， K 的特征峰（能量為 1.46 MeV） 較明顯， 故采用 K 峰進(jìn)行穩(wěn)譜。 由圖5 可知， 隨著環(huán)境溫度變化， 儀器可以實現(xiàn)自動高壓補償并確保特征峰在很小的范圍內(nèi)波動。

表1 擬合尋峰結(jié)果分析Table 1 Finding result Analysis of peak fitting

圖5 便攜式多道γ 能譜儀穩(wěn)譜驗證試驗Fig. 5 Verification experiment of spectral stabilization for portable multi-channel gamma spectrometer

基于上述儀器實時穩(wěn)譜功能， 選擇混合模型， 從-20～50 ℃， 以 10 ℃的溫度間隔分別設(shè)置8 個溫度梯度， 到溫后保持2 h 進(jìn)行譜數(shù)據(jù)采集， 驗證不同溫度梯度下擬合尋峰，自動穩(wěn)譜精度， 結(jié)果如表2 所示。

由表2 可知， 不同溫度梯度下， 開啟自動穩(wěn)譜功能后的實測峰位道址最大譜漂不超過±3 道， 該結(jié)果表明了Gaussian 擬合函數(shù)與LM 算法參數(shù)尋優(yōu)在自動穩(wěn)譜過程中的可靠性， 以及上述通過PID 進(jìn)行高壓模塊增益調(diào)節(jié)在實時穩(wěn)譜中的可行性。 但本文介紹的自動穩(wěn)譜算法對于重合峰分辨能力較弱， 下階段可采用多階Gaussian 函數(shù)擬合， 增強重峰分辨能力。

表2 不同溫度梯度下的特征峰穩(wěn)譜效果Table 2 Effect of spectral stabilization under different temperature gradients

5 結(jié)語

針對便攜式多道γ 能譜儀測量過程中因溫度影響造成譜線漂移的問題， 提出采用LM算法利用逐次逼近原理對Gaussian 擬合函數(shù)進(jìn)行參數(shù)尋優(yōu)， 獲取當(dāng)前特征峰峰位道址，通過峰位偏差計算當(dāng)前譜漂移量， 再依據(jù)譜漂移量利用PID 控制算法調(diào)節(jié)高壓控制字繼而進(jìn)行相應(yīng)的高壓補償， 實現(xiàn)系統(tǒng)自動穩(wěn)譜功能， 該方法具有硬件資源少、 響應(yīng)速度快、數(shù)據(jù)運算量小等優(yōu)點， 很適合便攜式儀器。經(jīng)試驗驗證， 在-20～50 ℃范圍內(nèi)實測峰位道址最大譜漂為±3 道。 目前該方法已應(yīng)用于γ能譜儀實際測量中， 運行穩(wěn)定， 結(jié)果可靠，達(dá)到預(yù)期要求。