一種用于核輻射探測的數(shù)字濾波器優(yōu)化設計方法
【技術(shù)領域】
[0001] 本發(fā)明涉及核輻射探測領域，特別涉及一種用于核輻射探測中的數(shù)字濾波器優(yōu)化 設計方法。
【背景技術(shù)】
[0002] 蓋革-米勒計數(shù)管(Geiger-Μ??Ier counter，以下簡稱GM計數(shù)管)是核福射探測領 域中常用的探測器，但其受到固有死時間與恢復時間的影響，線性范圍與測量上限受到限 制。為了提高GM計數(shù)管的量程，近年來提出了一種通過控制GM計數(shù)管工作電壓實現(xiàn)寬量程 測量的方法，文獻《一種新的GM計數(shù)管Time-To-Count測量方法》（核電子學與探測技術(shù)， 2002年第03期）、《Time-t〇-Count技術(shù)測量輻射強度的原理》(核電子學與探測技術(shù)，2006年 第01期）、《Time-t 〇-C〇Unt方法擴展GM計數(shù)管量程的實驗研究》（原子能科學技術(shù)，2013年04 期）等均報道了此方法。通過GM計數(shù)管的工作電壓在正常計數(shù)坪電壓與低于起始電壓下切 換，先計算GM計數(shù)管加上計數(shù)坪電壓時間內(nèi)的平均計數(shù)率（即GM計數(shù)管具有探測能力時間 內(nèi)的平均計數(shù)率），再由平均計數(shù)率來計算劑量率。該方法的關(guān)鍵在于計算平均等待時間 (即GM計數(shù)管加計數(shù)坪電壓的平均時間），但由于測量平均等待時間服從負指數(shù)分布，欲取 得較精確的平均值，需采用大量的測量結(jié)果進行平均或加權(quán)平均，這就要求測量系統(tǒng)必須 具有較大的存儲量與較快的計算速度。文獻《單管寬量程劑量儀關(guān)鍵技術(shù)分析》(第十六屆 全國核電子學與核探測技術(shù)學術(shù)年會論文集(上冊），2012)從數(shù)字信號處理的角度，將平均 過程看作測量信號通過低通濾波器的結(jié)果，采用低通濾波器處理等待時間，使等待時間樣 本通過濾波器后具有較小的方差，并建議通過改變設計指標的方式設計IIR數(shù)字低通濾波 器。雖然這種方法設計的低通濾波器可滿足測量精度要求，但未考慮濾波器響應時間長短， 導致測量結(jié)果不能即時反映輻射強度或劑量率的變化，從而使得測量結(jié)果明顯滯后于輻射 場的變化。

【發(fā)明內(nèi)容】

[0003] 為解決上述技術(shù)問題，本發(fā)明提供一種用于核輻射探測的數(shù)字濾波器優(yōu)化設計方 法，該方法能在滿足測量精度的前提下減少存儲空間、提高計算效率與反應速度，達到減輕 測量裝置的存儲與計算壓力、提升響應速度的目的。
[0004] 本發(fā)明的基本構(gòu)思是:將現(xiàn)有的平均或加權(quán)平均過程看作測量所得等待時間通過 FIR濾波器的結(jié)果，并采用IIR濾波器代替現(xiàn)有的平均或加權(quán)平均處理等待時間，通過粒子 群優(yōu)化算法對IIR濾波器的系數(shù)[ &1，&2，···，&(1] JboAib2,…，!^]進行優(yōu)化設計，使設計的 IIR濾波器在一定測量精度（由濾波器噪聲增益表征)條件下，具有最快的響應速度。
[0005] 本發(fā)明技術(shù)方案具體包括如下步驟：
[0006]步驟1:初始化參數(shù)
[0007]設置濾波器噪聲增益g、濾波器階數(shù)q、修正因子f (0〈f〈l)、粒子群規(guī)模M、最大迭代 次數(shù)K、權(quán)重CQ、C1、C2、粒子最大速度Vmax，隨機構(gòu)造 M個向量Xi= {xil，Xi2,…，Xi(2q-1)}分別代 表M個粒子當前位置，設置迄今為止全局最優(yōu)位置Pg= {pgl，Pg2,…，pg(2q-1)}為任意值，對應 的濾波器系數(shù)ag= {ag〇，agi，ag2，…，agq}、bg= {bg〇，bgi，bg2，…，bgq}為任思值，響應時間Tg為 無窮大，各粒子的最優(yōu)位置Pi= {pil，Pi2,…，Pi(2q-1)}為任意值(其中i為粒子序號），對應的 濾波器系數(shù)ai = {ai〇，aii，ai2，…，aiq}、bi = {bi〇，bii，bi2，…，biq}為任意值，各最優(yōu)位置對應 的響應時間Ti為無窮大；
[0008] 步驟2:采用粒子群優(yōu)化算法進行濾波器優(yōu)化設計；
[0009] 步驟2.1:令當前迭代次數(shù)k = l;
[0010] 步驟2.2:對各粒了進行設計；
[0011] 步驟2.2.1:令當前粒子序號i = l;
[0012] 步驟2.2.2:利用式（1)更新粒子Xi:
[0013] ViH = CdViH+Ciri (DiH-XiHHcm(Pgd-Xid)
[0014](1)
[0015] Xid = Xid+Vid
[0016] 式中，Xid代表粒子Xi的第d維變量，Vid代表粒子Xi的第d維速度；
[0017] 步驟2.2.3:對當前粒子i，令aj = xij( j = l，2，···，q)構(gòu)造濾波器系數(shù)[I，ai，a2，···， aq]，并計算濾波器的極點，若濾波器某極點位于單位圓外，則將該極點的倒數(shù)作為濾波器 的極點，并重新計算濾波器系數(shù)[I，ai，a2，…，aq]，同時令Xij = aj(j = 1，2，…，q)更新粒子的 前q維；
[0018] 步驟2·2·4:對當前粒子i，令bj = Xi(q+j+l)(j = 0，l，2，···，q-2)，bq-l = 0、bq = 0，構(gòu)造 濾波器系數(shù)[bo，bi，b2，…，bq]，根據(jù)系數(shù)[I，ai，a2，…，aq]、[bo，bi，…，bq]由式(2)計算濾波 器的單位響應hi (η):
[0019]
(2)
[0020] 式中，δ(η)為單位函數(shù)，當η = 〇時，其值為1，否則其值為〇;
[0021 ]步驟2.2.5:按式計算濾波器最小噪聲增益gmin:
[0022]
(3)
[0023] 步驟2.2.6:如果gmin小于等于g，轉(zhuǎn)步驟2.2.7,否則按式(4)更新粒子i的后q-1維， 并轉(zhuǎn)步驟2.2.4;
[0024]
(4).;：
[0025] 步驟2.2.7:對當前粒子i，由式(5)構(gòu)建矩陣心、A2，按式(6)計算向量c:
[0026]
(5)
[0027] c = [l 0 0 ··· OKA1+A2)-1 (6)；
[0028] 步驟2.2.8:由式(7)~(9)構(gòu)建矩陣HH3,按式（10)將向量c拆為向量 C1、c2,按 式(11)構(gòu)建向量b1;
[0034] 步驟2.2.9:由式（12)計算參數(shù)8〇、(：111與(3112，由式（13)計算參數(shù)11、8 :
[0029]
[0030]
[0031]
[0032]
[0033]
[0035] (12)
[0036] (13)?
[0037] 步驟2.2.10:由式（14)計算一元二次方程系數(shù)A、B、C:
[0038]
(14),
[0039]步驟2.2.11:計算一元二次方程4#+8丨+0 = 0的根，如果兩根為實根，分別將兩根 作為濾波器系數(shù)bq-i，并按式(15)計算濾波器系數(shù)bq，轉(zhuǎn)步驟2.2.12;否則轉(zhuǎn)步驟2.2.15;
[0040]
(15)；
[00411 步驟2.2.12:根據(jù)步驟2.2.2所得濾波器系數(shù)ai，a2，…，a q、bo，bi，…，bq-2，結(jié)合步驟 2.2.11計算所得^^與^的兩組值構(gòu)建兩個濾波器，按式（16)分別計算兩濾波器的階躍響 應，并以階躍響應s (η)達到0.9時的η值，確定響應時間tm、tn2，選響應時間短的濾波器作為 設計結(jié)果，取得濾波器系數(shù)[I，ai，a2，…，aq]、[bo，bi，…，bq]與響應時間tn;
[0042]
(16)
[0043] 其中，單位階躍函數(shù)定義如下：
[0044]
(17),
[0045] 步驟2.2.13:對于第i個粒子，將響應時間tn與其歷史最優(yōu)位置對應的響應時間Ti 進行比較，如果tn小于Ti，令Ti = tn、pi = Xi、ai = [I，ai，a2,…，aq]、bi = [bo，bi，…，bq];
[0046] 步驟2·2· 14:將Ti 與 Tg 進行比較，如果Ti 小于Tg，令 Tg = Ti、pg = pi、ag = ai、bg = bi;
[0047] 步驟2.2.15:令i = i+1，如果i〈M，那么轉(zhuǎn)步驟2.2.2，否則轉(zhuǎn)步驟2.3;
[0048] 步驟2.3:令k = k+l，如果k〈K，轉(zhuǎn)步驟2.2.1，否則結(jié)束優(yōu)化設計；
[0049] 步驟3:將全局最優(yōu)位置對應的濾波器輸出，作為優(yōu)化設計結(jié)果。
[0050] 本發(fā)明的優(yōu)勢在于：
[0051 ] (1)可根據(jù)需要指定濾波器的噪聲增益，控制測量誤差；
[0052] (2)優(yōu)化設計出的濾波器與傳統(tǒng)平均或加權(quán)平均算法相比，在相同噪聲增益情況 計算量很??；
[0053] (3)優(yōu)化設計出的濾波器與現(xiàn)有濾波器設計工具設計的濾波器相比，在相同噪聲 增益情況下具有更快的響應速度。
【附圖說明】
[0054] 圖1數(shù)字濾波器優(yōu)化設計方法程序流程圖
[0055] 圖2傳統(tǒng)方法設計濾波器與本方法設計濾波器階躍響應比較圖
【具體實施方式】
[0056] 設核輻射測量裝置要求測量所得輻射強度（劑量率）的誤差小于10%，根據(jù)30-法 貝1J，要求平均等待時間測均方根誤差小于1/30,那么平均時所需的測量數(shù)據(jù)應不少于900， 此處選擇1000,即平均處理后，噪聲功率變?yōu)樵瓉淼?/1000。按照文獻《單管寬量程劑量儀 關(guān)鍵技術(shù)分析》(第十六屆全國核電子學與核探測技術(shù)學術(shù)年會論文集(上冊），2012)修改 濾波器指標的方法，用傳統(tǒng)濾波器設計方法，設計橢圓低通濾波器EF。設置濾波器EF的通帶 歸一化寬度為〇. 001，通帶衰減為IdB，阻帶歸一化寬度為0.003，阻帶衰減為40dB，設計出的 濾波器為3階濾波器，該濾波器的噪聲增益為0.00102,與平均方法噪聲增益相當。采用1000 個測量數(shù)據(jù)進行平均、加權(quán)平均及采用濾波器EF處理三者的存儲量、計算量對比見表1所 不。
[0057] 表1平均、加權(quán)平均、濾波器EF處理三者的存儲量、計算量對比
[0059]從表1可以看出，濾波器EF的計算量與存儲量均明顯少于平均或加權(quán)平均方法，也 就是說，在平滑性能相近情況下，核探測中選用IIR濾波器進行平滑處理可降低測量系統(tǒng)存 儲與計算要求。濾波器EF的系數(shù)見表2所示。
[0060] 本發(fā)明以上述設計低通濾波器EF為參照，設置噪聲增益為0.00102(與濾波器EF相 同），濾波器階數(shù)同樣為3階。由于兩濾波器階數(shù)相同，其計算量與存儲量相同，二者噪聲增 益相同，其測量精度也相同，對比二者響應速度便能分析二者的優(yōu)劣。具體步驟如下：
[0061 ]步驟1:初始化參數(shù)
[0062] 設置濾波器噪聲增益g = 0.00102、濾波器階數(shù)q = 3、修正因子f = 0.9、粒子群規(guī)模 M= 100、最大迭代次數(shù)K= 10000、Co = I、C1 = 2、C2 = 2、粒子最大速度Vmax= 100，隨機構(gòu)造 M個 向量11={如412，"_，叫2(1-1)}分別代表1個粒子當前位置，設置迄今為止全局最優(yōu)位置 1^ = (Pgi，Pg2，…，Pg(2q-1)}為任思值，對應的濾波器系數(shù)ag= {ago，agi，ag2，…，agq}、bg = {bg〇，bgi， bg2，…，bgq}為任意值，響應時間Tg為無窮大，各粒子的最優(yōu)位置Pi = {pil，Pi2，···，P