本發(fā)明涉及數(shù)字信號處理領(lǐng)域，尤其是科里奧利質(zhì)量流量計(jì)的頻率跟蹤及估計(jì)方法。

背景技術(shù)：

科里奧利質(zhì)量流量計(jì)(Coriolis Mass Flowmeter，CMF)是一種基于諧振原理和科氏效應(yīng)(流體在直線運(yùn)動(dòng)時(shí)產(chǎn)生與質(zhì)量流量成正比的科里奧利力原理)，并能以極高的測量精度和可靠性實(shí)現(xiàn)精確測量流體流量質(zhì)量的高端儀器儀表。目前，國內(nèi)市場上銷售的CMF大多采用模擬驅(qū)動(dòng)方式和模擬信號處理方法，易受外界震動(dòng)和流體流速不穩(wěn)定的影響，使得工業(yè)現(xiàn)場測量質(zhì)量流量精度低于實(shí)驗(yàn)時(shí)標(biāo)定的精度。模擬電路驅(qū)動(dòng)和處理方法是對來自CMF輸出的兩路頻率相同、相位不同的正弦信號進(jìn)行濾波和過零點(diǎn)檢測，將其轉(zhuǎn)換成兩路方波信號，再利用方波信號的上升沿或下降沿來計(jì)數(shù)解算信號的頻率。模擬驅(qū)動(dòng)及信號處理的缺點(diǎn)在于驅(qū)動(dòng)增益有限、頻率跟蹤緩慢，不能解決兩相流、批料流發(fā)生時(shí)流量管固有頻率及阻尼比劇烈變化情況下的流量管停振而無法跟蹤到頻率的問題。

《面向時(shí)變的科里奧利質(zhì)量流量計(jì)信號的處理方法研究與實(shí)現(xiàn)》公開了一種科里奧利質(zhì)量流量計(jì)的信號處理方法，提出將多抽一濾波器、自適應(yīng)格型陷波濾波器和負(fù)頻率修正的滑動(dòng)DTFT遞推算法組合起來作為一套科里奧利流量計(jì)信號處理方法，該方法采用兩級多抽一濾波器，同時(shí)采用30階FIR低通濾波器，雖然能較好的保證線性相位特性，但是其數(shù)據(jù)運(yùn)算量大、處理實(shí)時(shí)性不強(qiáng)，只能作為一種計(jì)算機(jī)仿真驗(yàn)證方法來研究，很難運(yùn)用在存儲(chǔ)容量有限的嵌入式數(shù)字信號處理器上。針對CMF輸出兩路信號的高精度頻率算法估計(jì)和跟蹤表述，該方法對于其采用的自適應(yīng)格型陷波濾波器的參數(shù)描述模糊，雖然提到使用Burg算法進(jìn)行自適應(yīng)調(diào)整，卻沒有詳細(xì)說明參數(shù)整定的具體方法。

合肥工業(yè)大學(xué)徐科軍等的中國專利“基于AFF和SGA的科氏質(zhì)量流量計(jì)數(shù)字信號處理系統(tǒng)，CN1467485”采用自適應(yīng)陷波器(Adaptive Funnel Shaped Filter，AFF)，以自適應(yīng)線性增強(qiáng)信號、跟蹤和測量信號頻率，能兼顧頻率跟蹤能力及跟蹤精度，但問題在于該方案是基于定點(diǎn)型DSP(Digital Signal Processor，DSP)開發(fā)的算法。參考德州儀器TMS320系列DSP數(shù)據(jù)手冊以及文獻(xiàn)《DSP基本體系結(jié)構(gòu)和特點(diǎn)》可知，定點(diǎn)型DSP相比于浮點(diǎn)型DSP，定點(diǎn)DSP的字長每增加1bit，動(dòng)態(tài)范圍擴(kuò)大6dB，在處理低信噪比信號的場合，需要對數(shù)據(jù)不斷地移位定標(biāo)或者截尾。但是頻繁地移位定標(biāo)要耗費(fèi)大量的程序空間和執(zhí)行時(shí)間，截尾處理則會(huì)使計(jì)算的精度劣化。32bit浮點(diǎn)運(yùn)算DSP的動(dòng)態(tài)范圍可以達(dá)到1536dB，這不僅擴(kuò)大了數(shù)據(jù)動(dòng)態(tài)范圍，還因?yàn)榇蟠鬁p少了定標(biāo)，移位和溢出檢查，從而提高了運(yùn)算精度，可節(jié)省大量運(yùn)算時(shí)間和存儲(chǔ)空間。由于浮點(diǎn)DSP的浮點(diǎn)運(yùn)算用硬件來實(shí)現(xiàn)，可以在單周期內(nèi)完成，因而其處理速度遠(yuǎn)高于定點(diǎn)DSP，這一優(yōu)點(diǎn)在實(shí)現(xiàn)高精度復(fù)雜算法時(shí)尤為突出，為復(fù)雜算法的實(shí)時(shí)處理提供了保證。

西安東風(fēng)機(jī)電有限公司的中國專利“一種基于諧波分解的科里奧利質(zhì)量流量計(jì)數(shù)字信號解算裝置及方法，CN1837758”提到其在DSP中采用特殊設(shè)計(jì)的梳狀濾波器實(shí)現(xiàn)帶通濾波來抑制直流和工頻干擾，且其所述的濾波器可以覆蓋其使用的CMF裝置的全部工作頻率。該專利采用FIR濾波器來實(shí)現(xiàn)其所述功能，且要求中心頻率是固定的，若采用該專利所述方法實(shí)現(xiàn)CMF高精度頻率估計(jì)，所要面對的問題在于解決其所述的覆蓋所有CMF裝置的全部工作頻率，此算法設(shè)計(jì)復(fù)雜，雖然可以起到穩(wěn)定濾波器特性的效果，但需要以付出大量的運(yùn)算空間和運(yùn)算時(shí)間為代價(jià)，其使用的最小二乘法對信號進(jìn)行諧波分解，雖然不受非整周期截?cái)嗟挠绊?，但這種方法計(jì)算量大，實(shí)時(shí)跟蹤頻率的效果不佳。

綜上所述，對于高精度估計(jì)科里奧利質(zhì)量流量計(jì)的輸出信號頻率，且保持實(shí)時(shí)跟蹤，需要采用合適的頻率估計(jì)跟蹤算法，盡量優(yōu)化算法流程，在充分發(fā)揮硬件性能的基礎(chǔ)上，保證頻率估計(jì)的精度和算法的快速收斂性，并且不增加算法的時(shí)間復(fù)雜度和空間復(fù)雜度。

技術(shù)實(shí)現(xiàn)要素：

為了克服現(xiàn)有技術(shù)的不足，本發(fā)明為提高科里奧利質(zhì)量流量計(jì)輸出信號頻率估計(jì)的精度，保證跟蹤頻率的穩(wěn)定性，提出采用改進(jìn)的分段相位差頻率估計(jì)方法對科里奧利質(zhì)量流量計(jì)輸出信號進(jìn)行解算，以得到與流體流量質(zhì)量和密度相關(guān)的瞬時(shí)頻率。

本發(fā)明解決其技術(shù)問題所采用的技術(shù)方案包括下述步驟：

步驟一、根據(jù)CMF一次儀表有效輸出信號確定貝塞爾濾波器的通帶截止頻率、通帶最大衰減、阻帶截止頻率、阻帶最小衰減四個(gè)技術(shù)指標(biāo)，將四個(gè)技術(shù)指標(biāo)輸入電路設(shè)計(jì)軟件，得到用于科里奧利質(zhì)量流量計(jì)二次儀表的濾波器電路圖；

步驟二、按照步驟一生成的濾波器電路圖，在CMF二次儀表的模-數(shù)轉(zhuǎn)換器

(Analog-to-Digital Converter，AD)前端搭建有源模擬貝塞爾低通濾波器；

步驟三、用一個(gè)具有4通道同步采樣同步轉(zhuǎn)換功能的AD轉(zhuǎn)換器將CMF一次儀表的輸出信號轉(zhuǎn)換為數(shù)字信號；

步驟四、對于步驟三中AD轉(zhuǎn)換器轉(zhuǎn)換后輸出的數(shù)字信號，用離散正弦信號x[n]表示：

式(1)中A為信號幅度，ω₀為信號數(shù)字角頻率，φ為信號初始相位，n＝0,1,2…N-1，取x[n]的前N/2個(gè)點(diǎn)離散傅里葉變換(Discrete Fourier Transform，DFT)，變換后的頻域譜為X[k]，則

其中k＝0，1，2…N/2-1，f_s為信號采樣頻率，設(shè)k_max為離散信號快速傅里葉變換(Fast Fourier Transformation，F(xiàn)FT)處理后的峰值譜位置，式(2)中的表示為(k_max+δ)，則X[k]寫為：

則該正弦信號的頻率ω₀表示為如下形式：

提取步驟三AD轉(zhuǎn)換器轉(zhuǎn)換的數(shù)字信號的4通道中的一路，截取一組長度為2N-1的輸出數(shù)組，以標(biāo)志位x[i]分成前后兩個(gè)子段，即{x(n)，-N+1≤n≤i}與{x′(n)，i≤n≤N-1}，得到該同一路序列前后兩子段的頻偏估計(jì)：

其中為前向FFT序列的峰值譜的相位譜，為后向FFT序列的峰值譜的相位譜；

步驟五、使用解卷繞補(bǔ)償2π角度使回歸線性相位，綜合(4)、(5)式得到：

公式(6)中為信號的無偏估計(jì)頻率。

本發(fā)明的有益效果是采用了數(shù)字信號處理方法估計(jì)CMF諧振頻率，不受模擬電路時(shí)延和RLC參數(shù)影響，具有靈活、高速與高精度的優(yōu)點(diǎn)，本發(fā)明涉及的CMF高精度頻率跟蹤方法無需使用自適應(yīng)線性增強(qiáng)信號等復(fù)雜算法，在處理數(shù)據(jù)的復(fù)雜度上相比于基于自適應(yīng)漏斗形濾波器(Adaptive Funnel Shaped Filter，AFF)和滑動(dòng)Goertzel算法(Sliding Goertzel Algorithm，SGA)的科氏質(zhì)量流量計(jì)數(shù)字信號處理系統(tǒng)有明顯改善，正確穩(wěn)定跟蹤到CMF輸出頻率精度達(dá)到99.999％。

附圖說明

圖1是本發(fā)明涉及的科里奧利質(zhì)量流量計(jì)高精度頻率跟蹤方法流程圖。

圖2是本發(fā)明的低通濾波器的幅頻特性曲線圖。

圖3是本發(fā)明的低通濾波器的群延遲圖。

具體實(shí)施方式

下面結(jié)合附圖和實(shí)施例對本發(fā)明進(jìn)一步說明。

如圖1所示，本發(fā)明涉及的CMF高精度頻率跟蹤方法需要經(jīng)過確定濾波器性能指標(biāo)、設(shè)計(jì)濾波器、AD同步采集轉(zhuǎn)換信號、解算CMF輸出頻率和補(bǔ)償CMF輸出頻率精度等五大部分組成，其具體實(shí)施步驟如下：

步驟一、根據(jù)CMF一次儀表有效輸出信號確定貝塞爾濾波器的通帶截止頻率、通帶最大衰減、阻帶截止頻率、阻帶最小衰減四個(gè)技術(shù)指標(biāo)，將四個(gè)技術(shù)指標(biāo)輸入德州儀器FilterPro Desktop濾波器設(shè)計(jì)軟件，設(shè)計(jì)一個(gè)貝塞爾濾波器電路，用于CMF原始輸出信號濾波，得到用于科里奧利質(zhì)量流量計(jì)二次儀表的濾波器電路圖；

步驟二、按照步驟一生成的濾波器電路圖，在CMF二次儀表的模-數(shù)轉(zhuǎn)換器(Analog-to-Digital Converter，AD)前端搭建有源模擬貝塞爾低通濾波器，用于濾除疊加在CMF一次儀表輸出信號上的高頻噪聲，同時(shí)利用貝塞爾濾波器優(yōu)良的線性相位特性保證信號相位失真程度最低；

為衰減高頻信號，在模-數(shù)轉(zhuǎn)換器前端設(shè)計(jì)一個(gè)低通濾波器，考慮到高精度頻率估計(jì)算法處理的CMF輸出信號要求相位不失真，故所采用的濾波器必須保持最佳的線性相位特性，本發(fā)明采用一個(gè)有源模擬貝塞爾低通濾波器，可以在CMF輸出信號到達(dá)ADC之前，消除疊加在信號上的高頻噪聲，且保證信號相位失真程度最低。

步驟三、用一個(gè)具有4通道同步采樣同步轉(zhuǎn)換功能的AD轉(zhuǎn)換器將CMF一次儀表的輸出信號轉(zhuǎn)換為數(shù)字信號；

將CMF一次儀表輸出信號送入DSP處理之前需要轉(zhuǎn)換成數(shù)字信號，本發(fā)明采用的具有4通道同步采樣同步轉(zhuǎn)換功能的AD轉(zhuǎn)換器，實(shí)現(xiàn)了對同步性的要求。

步驟四、對于AD輸出的信號，本發(fā)明用離散正弦信號表示為：

取x[n]譜的前N/2點(diǎn)離散傅里葉變換為X[k]：

其中k＝0，1，2…N/2-1，f_s為信號采樣頻率，設(shè)k_max為FFT的峰值譜位置，式(2)中的可表示為(k_max+δ)，則X[k]為：

該正弦信號的頻率ω₀為：

提取步驟三AD轉(zhuǎn)換器轉(zhuǎn)換的數(shù)字信號的4通道中的一路，截取一組長度為2N-1的輸出數(shù)組，以標(biāo)志位x[i]分成前后兩個(gè)子段，即{x(n)，-N+1≤n≤i}與{x′(n)，i≤n≤N-1}，得到則該同一路序列前后兩子段的頻偏估計(jì)：

其中為前向FFT序列的峰值譜的相位譜，為后向FFT序列的峰值譜的相位譜。

步驟五、反正切計(jì)算得到的相位需要回歸到線性相位特性，使用解卷繞補(bǔ)償2π角度使得回歸線性相位，綜合(4)、(5)式得到：

公式(6)中為信號的無偏估計(jì)頻率。

實(shí)施例：

步驟一、針對國內(nèi)某廠家生產(chǎn)的四個(gè)通用系列型號的CMF流量管，其CMF信號輸出頻率的范圍介于40-240Hz，CMF一次儀表輸出信號頻率除共振頻率之外，還包含有各種高頻、低頻噪聲信號；此外，CMF輸出的信號共振頻率還與流體的性質(zhì)、流速和密度狀態(tài)有關(guān)，當(dāng)流量管處于兩相流或者批料流狀態(tài)時(shí)，流量管的固有頻率及阻尼比都會(huì)發(fā)生劇烈變化，在處理CMF一次儀表輸出的信號時(shí)需對信號進(jìn)行預(yù)處理，本發(fā)明使用一個(gè)低通濾波器濾除CMF一次儀表輸出信號的高頻噪聲。本發(fā)明使用的模擬低通濾波器技術(shù)指標(biāo)確定為：通帶截止頻率Ω_p＝400Hz，通帶最大衰減α_p＝3dB，阻帶截止頻率Ω_s＝600Hz，阻帶最小衰減α_s＝-30dB。

步驟二、本發(fā)明采用兩級多路反饋的3階貝塞爾低通濾波器來實(shí)現(xiàn)CMF原始輸出信號的原始濾波。本發(fā)明設(shè)計(jì)的低通濾波器的通帶增益為1V，濾波器由RC有源電路(Resistence-Capacitance電路，以下簡稱RC電路)組成，RC電路前一部分為一個(gè)1階低通貝塞爾濾波器，RC電路后一部分為一個(gè)2階低通貝塞爾濾波器，兩個(gè)貝塞爾低通濾波器級聯(lián)成一個(gè)3階貝塞爾低通濾波器。如圖2所示，是本發(fā)明得到的三階模擬貝塞爾低通濾波器的幅頻特性曲線，在信號頻率為255.44Hz處，信號已衰減-3.03dB，可見濾波器對240Hz以上的頻率能較好的抑制，圖3是濾波器的群延遲圖，頻率小于257.62Hz時(shí)，濾波器能保持較好的線性相位特性，群延遲時(shí)間為690μs。

步驟三、本發(fā)明基于CMF輸出兩路信號同步性的要求采用了一片具有4通道同步采樣同步轉(zhuǎn)換功能的AD轉(zhuǎn)換器，實(shí)現(xiàn)對CMF一次儀表輸出的A、B兩路信號的同步采樣。根據(jù)奈奎斯特采樣定律，AD采樣頻率必須設(shè)置為阻帶截止頻率的2倍以上，故本發(fā)明使用采樣頻率為4687.5Hz，以此來滿足抗混疊濾波以及FFT運(yùn)算要求。

步驟四、離散時(shí)間情況下的正弦信號一般表達(dá)式為：

取x[n]譜的前N/2點(diǎn)離散時(shí)間信號的DFT觀測為X[k]，A為CMF輸出信號的幅度，工程上要求信號的幅度滿足1Hz/1mV，為方便推理論證，本發(fā)明以幅度A表示，則X[k]可表示為：

本發(fā)明把該正弦信號的頻率ω₀可表示為如下形式：

k_max為FFT的峰值譜位置，本發(fā)明通過提取一路轉(zhuǎn)換后得到的長度為2N-1的輸入數(shù)組，再截取一個(gè)隨機(jī)標(biāo)志位x[i]的前后兩個(gè)子段，即{x(n)，-N+1≤n≤i}與{x′(n)，i≤n≤N-1}，該同一路序列前后兩子段的頻偏估計(jì)為：

其中為前向FFT序列的峰值譜的相位譜，為后向FFT序列的峰值譜的相位譜。

步驟五、本發(fā)明在求解到兩子段的相位差之后，計(jì)算得到的相位需要回歸到線性相位特性，本發(fā)明使用解卷繞補(bǔ)償2π角度使得回歸線性相位，再綜合(4)、(5)式，得到：

公式(6)中為信號的無偏估計(jì)頻率。