摘 要：提出了一種基于重疊短漢寧窗離散時(shí)間傅里葉變換（DTFT）算法的科氏流量計(jì)信號(hào)處理方法。首先采用新式自適應(yīng)陷波器對(duì)科氏流量計(jì)時(shí)變信號(hào)進(jìn)行濾波求其頻率，然后采用重疊短漢寧窗DTFT算法實(shí)時(shí)計(jì)算兩路信號(hào)之間的相位差，再根據(jù)頻率和相位差求得時(shí)間差，最終測(cè)得質(zhì)量流量。仿真和實(shí)測(cè)結(jié)果表明本文方法實(shí)用有效。

關(guān)鍵字：科氏流量計(jì) 自適應(yīng)陷波器 相位差 負(fù)頻率

    引言

    科氏流量計(jì)（Coriolis mass flowmeter，CMF）利用科里奧利效應(yīng)，通過測(cè)量傳感器輸出的兩路同頻正弦信號(hào)之間的相位差來計(jì)算流體質(zhì)量流量，是當(dāng)前研究最多、發(fā)展最為迅速、最具代表性的質(zhì)量流量計(jì)之一，具有廣闊的應(yīng)用前景。對(duì)科氏流量計(jì)的信號(hào)進(jìn)行處理，關(guān)鍵在于準(zhǔn)確地測(cè)量出來自兩個(gè)傳感器輸出信號(hào)的頻率和相位。近年來，國內(nèi)外的相關(guān)研究機(jī)構(gòu)積極研究科氏流量計(jì)數(shù)字信號(hào)處理方法，以克服傳統(tǒng)信號(hào)處理方法存在的不足，滿足各種實(shí)際應(yīng)用對(duì)科氏流量計(jì)性能越來越高的要求。

    首先采用自適應(yīng)格型陷波器對(duì)科氏流量計(jì)信號(hào)進(jìn)行濾波，以求其頻率，然后采用滑動(dòng)Goertzel算法來計(jì)算兩路信號(hào)之間的相位差。然而，滑動(dòng)Goertzel算法存在一個(gè)較長時(shí)間的收斂過程，只有在某一點(diǎn)后計(jì)算所得的結(jié)果才有效，對(duì)于時(shí)變信號(hào)不能進(jìn)行長時(shí)跟蹤。采用重疊窗的滑動(dòng)Goertzel算法實(shí)時(shí)計(jì)算兩路信號(hào)每個(gè)采樣點(diǎn)之間的相位差。具有重疊矩形窗的滑動(dòng)Goertzel算法的相位差實(shí)時(shí)性與計(jì)算精度均有所提高，但矩形窗和重疊區(qū)寬度的選取，必須考慮計(jì)算量與計(jì)算精度的均衡；同時(shí)，采用重疊矩形窗而引入的冗余計(jì)算也使得計(jì)算量明顯增大。提出了計(jì)及負(fù)頻率影響的離散時(shí)間傅里葉變換（Discrete time Fourier transform，DTFT）算法，縮短了收斂過程，提高了計(jì)算精度，但該算法不適用于流量變化的時(shí)變信號(hào)。

    本文提出一種基于重疊短漢寧窗DTFT算法的科氏流量計(jì)信號(hào)處理方法。首先，采用新式自適應(yīng)陷波器對(duì)頻率、相位和幅值均按隨機(jī)游動(dòng)模型變化的科氏流量計(jì)時(shí)變信號(hào)進(jìn)行濾波，并求其頻率（即頻率跟蹤），然后采用重疊短漢寧窗DTFT算法計(jì)算兩路信號(hào)的相位差和時(shí)間差。為檢驗(yàn)方法的有效性，本文在Matlab仿真驗(yàn)證的基礎(chǔ)上，用實(shí)采數(shù)據(jù)進(jìn)行了實(shí)驗(yàn)驗(yàn)證。

    1 科氏流量計(jì)信號(hào)處理方法

    1.1 基本思想

    在頻率跟蹤階段，采用的自適應(yīng)格型陷波器具有收斂速度較快、短時(shí)跟蹤信號(hào)頻率隨機(jī)緩慢變化精度較高的特點(diǎn)，但其長時(shí)間持續(xù)跟蹤信號(hào)頻率變化的能力較差，初始參數(shù)值的選取受隨機(jī)噪聲的影響較大。為了提高去噪效果、獲得更快的收斂特性和更高的跟蹤精度，本文采用一種基于新式自適應(yīng)陷波器的頻率跟蹤方法來持續(xù)跟蹤信號(hào)頻率的變化。在相位差計(jì)算階段，本文采用重疊短漢寧窗DTFT算法來實(shí)時(shí)計(jì)算信號(hào)的相位差，考慮了其中負(fù)頻率成分的影響，縮短了收斂過程。同時(shí)，采用重疊短窗截取，提高了實(shí)時(shí)性和計(jì)算精度，也適用于流量變化的實(shí)際情況。集成兩者的優(yōu)點(diǎn)，形成了本文所提出的基于重疊短漢寧窗DTFT算法的科氏流量計(jì)信號(hào)處理方法。

    根據(jù)上述方法基本思想，實(shí)現(xiàn)方法的技術(shù)路線如圖1所示。

圖1 研究技術(shù)路線圖

    1.1.1 頻率跟蹤方法

    圖2為新式自適應(yīng)陷波器的結(jié)構(gòu)圖，其傳遞函數(shù)為

        （1）

圖2 新式自適應(yīng)陷波器結(jié)構(gòu)圖

    式中：為n時(shí)刻對(duì)應(yīng)陷阱的陷波頻率，即對(duì)應(yīng)正弦波y（n）的估計(jì)頻率；ρk決定對(duì)應(yīng)陷阱的帶寬；圖2中Δ的合理取值可提高陷波器抑制噪聲影響的能力。

    1.1.2 相位差計(jì)算方法

    圖3～4所示分別為單一頻率正弦信號(hào)加矩形窗和加漢寧窗后得到的幅度頻譜。其中，區(qū)間（0，π）稱為正頻率區(qū)間，反映了信號(hào)的真實(shí)頻率；區(qū)間（-π，0）和（π，2π）稱為負(fù)頻率區(qū)間，分別是正頻率區(qū)間（0，π）關(guān)于ω＝0和ω＝π對(duì)稱的頻譜鏡像。區(qū)間（0，2π）為該信號(hào)的真實(shí)頻譜；區(qū)間（-2π，0）為真實(shí)頻譜的周期延拓。

圖3 單一頻率正弦信號(hào)加矩形窗的幅度頻譜

圖4 單一頻率正弦信號(hào)加漢寧窗的幅度頻譜

    正頻率區(qū)間（0，π）中的譜峰A反映了該信號(hào)的真實(shí)頻率，由于加窗截短造成了頻譜泄漏的影響，使得負(fù)頻率區(qū)間中譜峰A′的旁瓣可能會(huì)延伸并迭加到正頻率區(qū)間中譜峰A的旁瓣或主瓣上，從而造成旁瓣或主瓣干涉的影響。當(dāng)譜峰A靠近正頻率區(qū)間（0，π）的中點(diǎn)時(shí)，即當(dāng)信號(hào)頻率接近或等于Nyquist頻率的一半，亦即等于采樣頻率的1/4時(shí)，旁瓣干涉的影響最小，如圖3（a）所示；當(dāng)譜峰A靠近正頻率區(qū)間（0，π）的兩端時(shí)，即當(dāng)信號(hào)頻率較低或接近Nyquist頻率時(shí)，旁瓣干涉的影響會(huì)迅速增大；當(dāng)譜峰A的位置過于靠近ω＝0或ω＝π時(shí)，會(huì)導(dǎo)致出現(xiàn)譜峰A與A′的主瓣部分重疊的情況，從而造成旁瓣和主瓣干涉的影響，如圖3（b）和圖3（c）所示。此時(shí)，若仍按常規(guī)的方法來進(jìn)行頻譜分析，將會(huì)造成較大的偏差，甚至出現(xiàn)無法分析的情況。

    由圖3～4可以看出：在進(jìn)行頻譜分析時(shí)，當(dāng)信號(hào)頻率接近或等于Nyquist頻率的一半（即采樣頻率的1/4）時(shí)，可以忽略負(fù)頻率成分的影響；當(dāng)信號(hào)頻率較低或接近Nyquist頻率時(shí)，必須計(jì)及負(fù)頻率成分的影響。此外，由圖3與圖4比較可以看出，相同條件下，采用加漢寧窗所取得的效果比加矩形窗的效果要好。本文基礎(chǔ)上進(jìn)行改進(jìn)，提出了一種基于重疊短漢寧窗DTFT算法的科氏流量計(jì)信號(hào)處理方法，以適用于流量變化情況下的時(shí)變信號(hào)。

    設(shè)對(duì)濾波后的增強(qiáng)信號(hào)采用重疊的短漢寧窗進(jìn)行截取，得到兩路信號(hào)，由于選取的漢寧窗寬度值N很小，可近似認(rèn)為兩路信號(hào)在窗內(nèi)是時(shí)不變的。加漢寧窗后的兩路信號(hào)x₁（n）和x₂（n）在（短窗內(nèi)最后一采樣時(shí)刻的頻率估計(jì)值）處的DT-FT分別為

        （2）

    式（2）中，兩路信號(hào)x₁（n）和x₂（n）可表示為

        （3）

    式中：ω＝2πf₀/f_s；n＝0，1，…，N-1；A₁，A₂分別為兩路信號(hào)的幅度；f₀為信號(hào)頻率，f_s為采樣頻率；θ₁，θ₂分別為兩路信號(hào)的初始相位。ω（n）為漢寧窗函數(shù)的時(shí)域表達(dá)式

        （4）

    用φω1，φω2分別表示X_ω1（）和X_ω2（）的相位。式中

    

    則兩路信號(hào)的相位差Δθ可表示為

    

    計(jì)及負(fù)頻率成分的影響，經(jīng)過推導(dǎo)可得到兩路信號(hào)的相位差為

        （5）

    

    1.2 實(shí)現(xiàn)步驟

    針對(duì)時(shí)變信號(hào)，基于重疊短漢寧窗DTFT算法的科氏流量計(jì)信號(hào)處理方法實(shí)現(xiàn)步驟如下：

    （1）采用新式自適應(yīng)陷波器對(duì)科氏流量計(jì)的兩路時(shí)變信號(hào)進(jìn)行濾波，得到兩路增強(qiáng)信號(hào)x₁（n），x₂（n），并實(shí)時(shí)跟蹤求得信號(hào)的頻率估計(jì)值（n）；

    （2）采用重疊短漢寧窗（漢寧窗寬度為N，重疊區(qū)寬度為N-1）對(duì)濾波后的增強(qiáng)信號(hào)進(jìn)行截取，得到對(duì)應(yīng)漢寧窗截取的兩路信號(hào)序列分別為x_ω1（n）和x_ω2（n）；

    （3）采用常規(guī)的DTFT算法分別計(jì)算每個(gè)短漢寧窗內(nèi)x_ω1（n）和x_ω2（n）在處的DTFT，求出tan_φω1和tan_φω2。

    

    （4）由，N計(jì)算出m1～m4，并同tan_φω1和tan_φω2一起代入式（5），從而求得相位差；

    （5）最后再根據(jù)計(jì)算出的和，求每一時(shí)刻的時(shí)間差

    1.3 具體算法

    根據(jù)上述基本思想及實(shí)現(xiàn)步驟，基于重疊短漢寧窗DTFT算法的科氏流量計(jì)信號(hào)處理方法的具體算法如下：

    步驟1：設(shè)定變量及初始值。

    

    步驟2：采用新式ANF（Adaptive notch filter）估計(jì)信號(hào)頻率，并進(jìn)行實(shí)時(shí)跟蹤，得到信號(hào)頻率（n）。

    

    

    步驟3：用基于重疊短漢寧窗DTFT算法來實(shí)時(shí)計(jì)算兩路信號(hào)的相位差，并進(jìn)行實(shí)時(shí)跟蹤，得到相位差Δθ（n）。

    

    2 實(shí)驗(yàn)結(jié)果及分析

    2.1 仿真結(jié)果及分析

    為了檢驗(yàn)本文方法的效果，使用Matlab對(duì)算法進(jìn)行仿真。由于所研究的科氏流量計(jì)信號(hào)頻率通常在（100±4）Hz范圍內(nèi)變化，相位差小于±4°，仿真實(shí)驗(yàn)中選擇20000個(gè)采樣點(diǎn)，相關(guān)參數(shù)取值如下：f_s＝2000Hz，A（0）＝10，ω（0）＝≈0.3142，σ_e＝0.6，σ_A＝10^-3，σ_ω＝10^-5，σ_Φ＝10^-5。

    由于新式自適應(yīng)陷波器的仿真結(jié)果在中已給出，這里給出陷波器一些變量的初始值。

    

    將本文方法與采用的具有重疊短矩形窗的滑動(dòng)Goertzel算法（以下簡(jiǎn)稱“SGA”方法）進(jìn)行仿真比較分析。這里取初始相位差為0.1°的相位差和時(shí)間差估計(jì)結(jié)果。從圖5～6可以看出，SGA方法存在一定的滯后性，而本文方法的相位差和時(shí)間差的估計(jì)曲線較好地與真實(shí)值曲線吻合，說明本文方法具有更好的跟蹤能力。其主要原因是本文方法考慮了負(fù)頻率成分的影響，且所截取窗較短，能及時(shí)地跟蹤相位差和時(shí)間差的實(shí)時(shí)變化。為了表達(dá)方便，圖中的橫坐標(biāo)“時(shí)間（n）”選為采樣點(diǎn)數(shù)，通過信號(hào)頻率和采樣頻率不難計(jì)算出實(shí)際時(shí)間，后續(xù)圖中也使用了此表示方法。

圖5 初始相位差為0.1°時(shí)的相位差估計(jì)曲線

圖6 初始相位差為0.1°時(shí)的時(shí)間差估計(jì)曲線

    表1列出了在不同初始相位差條件下，本文方法與SGA方法相位差和時(shí)間差估計(jì)值的均方誤差值，其計(jì)算公式如式（6，7）所示。MSE=

        （6）

        （7）

    從表1可以看出，不同初始相位差條件下，本文方法得到的相位差和時(shí)間差均方誤差值均小于SGA方法，說明了本文方法的有效性。

    SGA方法所選取的矩形窗長度一般為400，重疊區(qū)寬度為350；而本文方法所選窗長為N，重疊區(qū)寬度為N-1（本文中N＝16），相比而言，本文方法的計(jì)算量大大減少，計(jì)算效率更高，實(shí)時(shí)性更好，跟蹤精度更高，這也說明本文方法更具優(yōu)越性。

表1 不同初始相位差下兩種方法相位差和時(shí)間差估計(jì)均方誤差值的比較

    2.2 實(shí)測(cè)結(jié)果及分析

    實(shí)際采集某型號(hào)科氏流量計(jì)傳感器輸出的兩路信號(hào)，采用SGA方法和本文方法對(duì)實(shí)采數(shù)據(jù)進(jìn)行比較分析。流量計(jì)振動(dòng)信號(hào)的頻率約為146Hz，采樣頻率為10kHz，分別針對(duì)不同流量變化情況下進(jìn)行采樣。由于現(xiàn)有技術(shù)條件的限制，無法求得每點(diǎn)的實(shí)際流量值，圖7僅給出了兩種方法估計(jì)曲線的比較。從圖7可以看出，兩種算法變化趨勢(shì)基本一致，可以較好反映真實(shí)相位差和時(shí)間差的變化情況，但SGA方法存在一定的滯后性，與前述仿真結(jié)果（圖5～6）相符，這也充分說明了仿真結(jié)果的正確性和本文方法的有效性。表2列出了5種不同流量下本文方法得到的頻率、相位差和時(shí)間差估計(jì)均值。圖8所示為時(shí)間差均值和質(zhì)量流量顯示值之間的關(guān)系。根據(jù)科氏流量計(jì)原理可知：在U型管結(jié)構(gòu)與材料確定的情況下，質(zhì)量流量qm與兩路信號(hào)的時(shí)間差Δt成線性關(guān)系，與角頻率ω?zé)o關(guān)。由圖8可見，質(zhì)量流量與時(shí)間差之間存在著較好的線性關(guān)系，與其理論相符合。這也表明，本文方法在工業(yè)生產(chǎn)調(diào)試現(xiàn)場(chǎng)是實(shí)用有效的。

圖7 不同流量下的相位差和時(shí)間差估計(jì)曲線

圖8 時(shí)間差與質(zhì)量流量的關(guān)系

表2 不同流量下的頻率、相位差和時(shí)間差估計(jì)值

    3 結(jié)論

    針對(duì)現(xiàn)有的科氏流量計(jì)信號(hào)處理方法，頻率跟蹤方法存在長時(shí)頻率跟蹤精度和穩(wěn)定性較差的問題，而相位差計(jì)算方法又存在冗余計(jì)算、計(jì)算精度和實(shí)時(shí)性均有待提高的問題，本文提出一種基于重疊短漢寧窗DTFT算法的科氏流量計(jì)信號(hào)處理方法，并進(jìn)行了仿真及實(shí)測(cè)結(jié)果的比較分析，得到如下結(jié)論：（1）整套信號(hào)處理方法有效地提高了科氏流量計(jì)的質(zhì)量測(cè)量精度，仿真及實(shí)測(cè)結(jié)果均表明了本文方法的有效性和實(shí)用性。（2）基于重疊短漢寧窗DTFT算法的相位差計(jì)算方法，計(jì)及了負(fù)頻率成分的影響，能及時(shí)地跟蹤相位差和時(shí)間差的實(shí)時(shí)變化，提高了實(shí)時(shí)測(cè)量精度。（3）本文方法所采用的新式自適應(yīng)陷波器有著更好的收斂特性和長時(shí)持續(xù)跟蹤能力，跟蹤精度更高。