CN101825660A - 欠采样下的正弦信号频率的高效测量方法及实施装置 - Google Patents

Abstract

本发明属于数字信号处理技术领域。为提供一种欠采样条件下频率等参数估计，还提供能完成精确频率测量的欠采样下的正弦信号频率的高效测量方法和装置，本发明采用的技术方案是，分别做传统FFT谱分析和全相位FFT谱分析，得到峰值谱G(q)、峰值谱G_a(q)，从峰值谱上直接读出相位值

，取G(q)模的平方得到功率谱值P_g(q)，将P_g(q)与G_a(q)模除后，得到幅值估计

；取G(q)的相位值

与G_a(q)的相位值

的差值除以τ＝(N-1)/2后得到频偏估计Δk，最后将ΔkΔω与qΔω叠加即得数字角频率估计

，

和

直接作为所测信号的相位和幅值估计。本发明主要应用于数字信号处理中的欠采样测量。

Description

欠采样下的正弦信号频率的高效测量方法及实施装置

技术领域

本发明属于数字信号处理技术领域。具体涉及在采样速率不足(相比于信号频率)的情况下，一种高精度地估计正弦波频率的新方法及其相应的实验装置，即欠采样下的正弦信号频率的高效测量方法及实施装置。

背景技术

在通信、仪表、电力、光学应用、故障诊断等工程领域，存在大量对高频正弦波的频率和相位进行高精度、高效、快速的估计问题。例如相干解调是通信中最常见的解调方式[1]，相干解调需要接收机从接收信号中提取信息，并产生与发送端同频、同相的本地载波。然而到达接收机的载频信号的频率和相位有可能发生变化，如何用数字化方法高效、精确的测出其频率和相位的变化值，一直是无线通信中的一个难题，这是因为无线载波频率通常很高，为实现不失真的信号波形采样，按照香农定理，就必须至少以两倍信号频率的采样速率进行采样，这样对模数转换器(Analog to Digital Converter，ADC)的要求很高。常用的办法是将信号逐级进行下变频到中频，然后再采样。但是这样做一来需引入变频措施，二来变频过程要依赖非线性元件来实现，这不可避免地会引入非线性失真。另外，目前市场上常用的频率计、相位计[2]和功率计等仪表设备大多数都已经数字化，但是任何频率计都有一个测量范围，当所测量的信号频率很高时，超出了其测量范围会导致仪器失效。目前的频率计的测量方法，还主要是以计数器计数为主，这需要频率计内部产生一个振荡频率非常高的时钟信号，这对硬件和工艺要求非常高。事实上，将高效的数字信号处理算法(如基于传统FFT的频率估计内插方法[3])嵌入到频率计后，可以大大降低仪器成本、提高仪器测量精度和扩大仪器的测量范围。由于信号的频率、幅值和相位的估计问题普遍存在于与国民经济紧密相关的通信、测量仪表[4]、故障诊断[5]等多个生产环节中，因此研发出高性能的相位估计法具有很高的工程价值和经济价值。而研发欠采样情况下的正弦波参数测量方案，更可大大提高现有仪器精度水平的潜力。

正弦波频率参数测量与估计的方法很多，传统依靠模拟器件的方法，如：矢量法、二极管鉴相法、脉冲计数法等，其测量系统复杂，需专用器件，硬件成本高。近年来，频率估计逐渐向数字化方向发展，其优点在于硬件成本低、适应性强，只需单片机、DSP(Digital SignalProcessor，数字信号处理器)、FPGA(Field Programming Gate Array，现场可编程门阵列)等通用器件就可完成，对不同的测量对象仅需改变程序算法即可，且其精度一股高于模拟式测量。因此，选定一套精确的频率估计算法是关键。

然而，为研发出高精度的频率测量算法(即参数估计算法)，仅在工程领域中去考虑问题是远远不够的。因为无论是在哪种应用场合，也无论是测量哪种物理量，采样后信号的表现形式都是离散观测数据，若要精确、快速、有效地从观测数据中提取出参数信息，这就涉及很多的理论问题，所涵盖的知识可延伸到数字信号处理、信号检测与估计理论、信息论、概率论与数理统计、随机过程等多个学科领域，只有加深对这些领域基础知识的理解，从新的角度提出一些优化参数估计性能的措施，才有可能研发出性能更为优良的算法。现有的数字化频率估计主要包括以下几种方法：

(1)脉冲计数法

这是最常用的相位测量法，其频率和相位原理如图1所示，频率测量的主要过程为：用已知频率的计数脉冲对所测信号进行计数，通过读取在正弦波半个周期内(两次过零点的间隔时间，可用触发的办法记录过零点时刻)计数脉冲的个数来估计信号的频率值。而相位测量则需要在本地端产生与所测正弦信号同频且已知初相的参考正弦波，然后分别对这两路正弦波进行脉冲计数(假设其基本计数周期为T₀)，根据其脉冲计数差值(假设n个脉冲差)算出其延时nT₀，进而测算出其相位差2πfnT₀。若由参考相位加上此相位差，还可得到相位估计。

这种方法的缺陷在于：(1)需专门产生脉冲计数信号，且测量精度很大程度上取决于基本脉冲宽度T₀，只有减小T₀才能提高测量精度，扩大仪器测量范围，这对硬件设备的要求很高；(2)仅能对单频信号进行测量估计；(3)需产生参考正弦波才能完成相位计数，且要求参考信号与所测正弦信号的频率完全一致，若略有偏离，则会给测量结果带来很严重的偏差。

(2)希尔伯特变换法[6]

对于余弦类信号x(t)＝cos(2πft+θ₀)，为估计t＝0时的相位θ₀，以采样频率f_s对x(t)进行采样后得到序列x(t)＝cos(2πft+θ₀)，n＝0，1，...，N-1；若对x(n)作希尔伯特变换可得其解析信号

n＝0，1，...，N-1，再根据式

即可得到θ₀的估计。

然而希尔伯特变换法对采样频率f_s要求很高，当信号频率f＝kf_s/N，k∈z⁺时，上式测得的相位是准确的；一旦采样频率f_s稍微发生偏离，使得这种整数倍的关系不成立，其估计就不正确了。此方法的抗噪性能差；另外，此方法仍只能对单频信号进行估计。

对两个存在延时关系的序列分别采用希尔伯特变换法测出其相位，由其相位差即可得到信号频率的估计。

(3)正弦曲线拟合法[7-9]

John Kuffel提出了正弦曲线拟合法[7，8]，该方法可获得很高的估计精度。它分为四参数法(频率、幅度、相位和直流分量均未知)和三参数法(除频率已知外，其它3个均未知)两种情况，文献[9]指出，四参数正弦曲线拟合过程并不是闭合的线性过程，尚无确切的数学公式可直接计算出拟合参数，若拟合初始条件选择不当，易使得迭代过程发散或收敛到局部最优，且拟合需大量运算时间。为提高效率，文献[9]提出先估算频率，再进行三参数的正弦曲线拟合(为闭合线性过程，且绝对收敛)。然而，此方法要求频率估计非常精确，否则一个直接后果就是会把频率估计误差带入到相位估计中，另外，此方法仍只能对单频信号进行相位估计。

(4)各类基于FFT频谱校正的频率估计法

基于FFT的相位估计法是近年来工程界和学术界研究得最热的、也是实际应用最广的参数估计法。原因是：其一，由于FFT(Fast Fourier Transform，快速傅立叶变换)是DFT(DiscreteFourier Transform，离散余弦变换)的快速算法，具有其它算法无法比拟的方便、简洁的计算优势；其二，在理论上，当N足够大时，若对离散采样值进行FFT可获得最大似然解的近似逼近；其三，各种丰富的信息可在FFT的结果中直接得到反映，如FFT谱线位置的分布反映了信号的频率大小，而FFT结果为复数，因而隐含了信号的相位信息等，这是其它正交变换所不具备的，如DCT(离散余弦变换，Discrete Cosine Transform)、沃尔什变换、K-L变换(Karhunen-Loeve Transform，卡洛变换)、Wigner-ville(维格纳-威利)变换等。虽然FFT的栅栏效应限制了其频率分辨率，但其栅栏效应可通过对FFT谱线进行内插的方法加以解决，这样就衍生出了多种离散频谱校正法，文献[10]总结了4类频谱校正法，包括能量重心法[11]、比值法[12-13]、DFT+FFT谱细化法[14]、相位差法[15-17]，这些方法都可对存在多种频率成份的信号进行较精确的相位估计，因而广泛应用于振动分析及故障检测、电力系统谐波分析、电介质损耗角测量、雷达测速等领域。

然而通过实验和理论证明，这些基于FFT频谱校正的相位估计法存在如下缺陷：(1)FFT存在谱泄漏现象，当信号包含多种频率成份时，各频率成份间会产生谱间干扰，这会很大程度地影响相位的测量精度；(2)由于存在谱泄漏，这些校正法尤其不适合于存在密集频谱的相位估计场合[18]；(3)经调研发现，这些FFT频谱校正法在估计相位时，都是按照先估计频率，得到频率偏离值后再进行相位估计的步骤进行的，因此这样做的一个直接影响是会把频率估计的误差带入到相位估计中。

发明内容

本发明的目的在于克服现有技术的不足，在采样速率明显跟不上甚至远远低于正弦波信号频率(即欠采样)的情况下，高效、高精度地完成对正弦波信号频率等参数估计(如幅值和相位测量也能由本方案实现)，并且还提供能完成精确频率测量的硬件实施装置。

为达到上述目的，本发明采用的技术方案如下：欠采样下的正弦信号频率的高效测量方法，包括下列步骤，设要测出频率、幅值和相位值的正弦波信号为x(t)＝Acos(2πf₀t+θ₀)，

首先以低于真实频率值f₀的采样速率f_s对信号x(t)进行模数转换而得到欠采样后的离散信号g(n)＝x(nT_s)，T_s＝1/f_s；

然后对g(n)分别做传统FFT谱分析和全相位FFT谱分析，而得到传统FFT谱G(k)k＝0，1，...，N-1和全相位FFT谱G_a(k)，k＝0，1，...，N-1，取G(k)的模值即得传统振幅谱|G(k)|，取G(k)的相角即得传统FFT相位谱

同样地，取G_a(k)的模值即得全相位FFT振幅谱|G_a(k)|，取G_a(k)的相角即得全相位FFT相位谱

再对全相位FFT振幅谱|G_a(k)|进行谱峰搜索，并记录谱峰所在的k＝q的位置，从该峰值位置直接读出全相位FFT相位谱值即可作为相位测量值

另外，在谱峰k＝q位置很容易计算出传统FFT的功率谱值P_g(q)＝|G_a(q)|²，将P_g(q)除以谱峰位置的全相位振幅谱值|G_a(q)|后，得到幅值估计

取谱峰k＝q位置的FFT相位值与全相位FFT相位值

的差值除以τ＝(N-1)/2后得到频偏估计Δk，最后将ΔkΔω与qΔω叠加即得数字角频率估计

和

直接作为所测信号的相位和幅值估计；

确定出f₀和f_s的整数倍值p，将p和f_s值与

叠加后，即可作为最终的频率估计值

全相位FFT谱分析是，用长为(2N-1)的卷积窗w_c对中心样点x(0)前、后(2N-1)个数据进行加权，然后将间隔为N的数据两两进行重叠相加，再对重叠相加后的数据进行DFT即得全相位谱分析结果X_a(k)，k＝0，1，...，N-1，再取其模值后，即得全相位FFT振幅谱|X_a(k)|，取其相角后，即得全相位FFT相位谱

一种欠采样下的正弦信号频率的高效测量实施装置，包括：

信号调理电路，用于对输入信号进行模拟预处理，对信号幅度范围进行必要调整，并去除外干扰噪声；

模数转化器A/D，用于采样得到样本序列x(n)，以并行数字输入的形式进入数字信号处理器DSP；

数字信号处理器DSP，用于对g(n)分别做传统FFT谱分析和全相位FFT谱分析，由传统FFT得到g(n)的峰值谱G(q)，由全相位FFT谱分析得到g(n)的峰值谱G_a(q)，从全相位FFT得到的峰值谱上直接读出相位值

取G(q)模的平方得到功率谱值P_g(q)，将P_g(q)与G_a(q)模除后，得到幅值估计

取FFT峰值谱G(q)的相位值

与全相位FFT峰值谱G_a(q)的相位值

的差值除以τ＝(N-1)/2后得到频偏估计Δk，最后将ΔkΔω与qΔω叠加即得数字角频率估计 和

直接作为所测信号的相位和幅值估计；

确定出f₀和f_s的整数倍值p，将p和fs值与

叠加后，即可作为最终的频率估计值

输出驱动及显示模块，用于显示出信号估计值。

附图说明

图1基于时钟脉冲计数的频率和相位估计。

图2欠采样高频正弦波参数估计流程。

图3全相位FFT振幅谱和相位谱分析流程(阶数N＝4)。

图4满足香农定理采样时的传统FFT谱和全相位FFT谱。

图5原波形、欠采样波形及其全相位FFT谱。

图6数字角频率分布区域。

图7两种峰值振幅谱分布。

图8全相位FFT/传统FFT综合插值法参数估计流程。

图9严重失真下的欠采样波形及其全相位FFT谱。

图10本发明的欠采样频率测量硬件实现装置

图11与图10对应的DSP内部处理流程。

具体实施方式

本发明提出基于全相位FFT方法进行相位估计，全相位FFT是专利申请人王兆华在专利文献[19]中提出的一种新型谱分析方法。谱分析包括振幅谱分析和相位谱分析两个方面，专利文献[19]只研究了全相位FFT的振幅谱特性，没有涉及相位谱的研究，也没有研究如何根据谱分析的结果进行参数估计。

本专利所解决的就是在用欠采样速率对信号进行采样的情况下，综合本专利所发现的欠采样情况下全相位FFT和传统FFT的振幅谱和相位谱特征，完成对信号频率的高精度、高效估计。

1.欠采样情况下正弦波信号的频率、相位与幅值的测量原理

本方案的基于欠采样情况下正弦波信号的频率、相位与幅值的测量流程如图2所示。图2中，为测出未知正弦波信号x(t)＝Acos(2πf₀t+θ₀)的频率、幅值和相位值，首先以低于真实频率值f₀的采样速率f_s对信号x(t)进行模数转换而得到离散信号g(n)＝x(nT_s)，T_s＝1/f_s，然后对g(n)分别做FFT谱分析和全相位FFT谱分析，其中从全相位FFT峰值谱上直接读出相位值

再依据基于全相位FFT与传统FFT的综合插值法(后面详细介绍)计算得到欠采样信号g(n)的幅值估计和数字角频率估计

其中

和

可直接作为所测信号的相位和幅值估计；而频率估计还需要确定出f₀和f_s的整数倍值p，将p f_s值与

叠加后，即可作为最终的频率估计值

图2中的全相位FFT是申请人在专利文献[19]中提出的，其处理流程图见图3所示，但文献[19]只实现了幅值谱估计，事实上由全相位FFT还能得到相位谱，其完整的振幅谱和相位谱的估计过程如图3所示。

只需用长为(2N-1)的卷积窗w_c对中心样点x(0)前、后(2N-1)个数据进行加权，然后将间隔为N的数据两两进行重叠相加，再对重叠相加后的数据进行DFT即得全相位谱分析结果X_a(k)，k＝0，1，...，N-1，再取其模值后，即得全相位FFT振幅谱|X_a(k)|，取其相角后，即得全相位FFT相位谱

2.欠采样测量原理

(1)满足香农定理采样情况下的信号离散谱分布

对于单频正弦波信号x(t)＝Acos(2πf₀t+θ₀)来说，理想情况下信号的傅立叶谱是正负频率轴上在f＝f₀处和f＝-f₀处的冲击函数。当以满足香农定理的采样速率f_s(f_s＞2f₀)对x(t)x(t)采样后，再对采样序列x(n)做阶数为N的FFT可得到离散谱X(k)，由于信号频率f₀与FFT的频率分辨率Δf＝f_s/N之间通常不能保持整数倍关系，这使得单频信号的离散傅立叶谱通常不会仅仅是两根谱线，这时将会在序号为q＝[Nf₀/f_s]的谱位置生成峰值谱线(‘[·]’表示取整运算)，而在k＝q附近通常会泄漏出多根谱线。显然，因谱泄漏而产生的多根谱线降低了谱分析的直观性和可读性。另外，这时FFT相位谱也变得很紊乱，从峰值相位谱线上无法直接得出相位信息。例如以采样速率f_s＝16Hz对频率f₀为3.4Hz、幅值A为2、初相值为40°的正弦波进行采样，再做阶数N＝16的FFT谱分析，则可得到如图4(a)所示的振幅谱和如图4(b)所示的相位谱。很显然所生成的峰值谱线落在q＝3处，若直接以该处的对应频率值qΔf＝3Hz去估算信号频率则会产生0.4Hz的频率估计误差，而读出的振幅值约为2|X(q)|＝2×0.7＝1.4，与真实振幅值相差0.6，而图4(b)中的相位值

超出100°，远远偏离了实际值40°。申请人曾在文献[19]中提出全相位FFT谱分析方法，其全相位FFT振幅谱|X_a(k)|如图4(c)所示，相位谱

如图4(d)所示。从图4(c)、图4(d)可明显看出，|X_a(k)|相比于传统FFT的|X(k)|，其谱泄漏程度大为降低，并且其峰值相位谱

值正好等于实际值40°，直接实现了初相的精确估计，同时从图4(d)中我们还发现在谱序号为k＝N-3＝13处，测出的相位值为-40°，即为真实相位值的相反数，因此这两根对称的谱线都包含有真实的相位信息。

从图4可看出，即使是在满足香农定理采样情况下，传统FFT仍无法直接得到信号的频率、幅值和相位值；而全相位FFT谱泄漏程度大为降低，能直接测出信号的相位值，但也不能直接测到信号的频率和幅值，因此需借助本专利提出的比值法插值来得到其精确估计值。若在不满足香农定理情况(即欠采样情况)下对高频信号进行采样，要测出高频信号的参数值，其情况就更为复杂。因为除了保留图4情况下由于离散谱分析所带来的估计误差外，欠采样情况下还必然会因采样频率太低而引起波形失真，这种波形失真必然会进一步增大原信号的频率和幅值的估计难度。

本专利提出如图2的处理流程来测量高频信号的频率、幅值和相位值，总的策略是先采用全相位FFT/传统FFT综合插值法去估计欠采样失真信号的频率、幅值和相位值，基于此再把这些估计值转化为原高频信号的频率、幅值和相位估计。之所以可以完成这种测量值转化，这是由欠采样信号与原始信号频谱间的内在联系而决定的。

(2)不满足香农定理欠采样情况下的信号离散谱分布

由于相比于传统FFT谱分析，全相位FFT谱分析可以大大改善频谱泄漏效应，故本专利采用全相位FFT对欠采样信号进行谱分析。

图5(a)给出了单频余弦信号x(t)＝Acos(2π100t+40°)的波形，同时给出了以采样速率f_s＝53对该信号进行采样的离散样点，图5(b)将这些离散样点连接起来形成波形g(n)。由于采样频率f_s＝53Hz远小于2倍的信号频率(2f₀＝200Hz)，故属于欠采样情况，采样后的波形频率偏低，出现了明显失真。图5(c)给出了采样序列的全相位FFT振幅谱(谱分析阶数N＝16)，从中可看出其谱泄漏情况并不突出。图5(d)给出其全相位FFT相位谱。

从图5中可发现，即使是在欠采样情况下，全相位FFT仍可直接估计出原高频信号的相位信息(图5d的相位谱图上，在峰值谱k＝2处，其相位值等于实际值的相反数-40°，在峰值谱k＝N-2＝14处，其相位值为40°)。但频率信息还不能直接提取出来，因为当谱分析阶数N＝16时，从图5(c)的谱峰位置k＝2处，可大致估计出信号的频率是2Δ＝2×f_s/N＝6.6250Hz，从图5(c)的谱峰位置k＝14处，可大致估计出信号的频率是14Δ＝14×f_s/N＝46.3750Hz，它们与真实值100Hz均相差甚远。这可从图5(a)、5(b)的波形看出，因为若以满足香农定理进行采样，则要求在一个周期内至少能采到2个样点，而图5(a)中的原信号在一个周期内平均采不到1个样点，故图5(b)采样后的波形严重失真，失真波形在图5(d)离散谱分析中的反映是无法标识真实谱应在的位置。

(3)欠采样情况下的信号参数信息的隐蔽性分析

虽然欠采样后信号波形出现失真，但原高频信号的参数信息却没有丢失，而是隐藏在欠采样后的波形里，因而需对欠采样后的波形的参数隐蔽性进行深入分析。

以采样速率fs对频率为f0的余弦信号进行采样后，其数字化角频率为ω＝2πf₀/f_s，因而在满足香农定理情况下进行采样时(f_s＞2f₀)，其数字角频率ω∈[0，2π)，我们称该区间为满足香农定理的主值区间(如图6所示)；当在不满足香农定理情况下进行欠采样时(f_s＜2f₀)，则其数字角频率必然超出ω∈[0，2π)的主区间，而跨越到如图6所示的欠采样区间。

另一方面，信号离散化后，数字角频率ω也必然相应的周期化，其周期为2π。这样对于欠采样情况，图6欠采样区间中的任一数字角频率ω模除周期2π后，总可以映射到图6中主值区间[0，2π)内的某个值ω^*上，即满足ω＝p2π+ω^*(p∈z⁺)。而在数字频率与模拟频率的对应关系上，ω与信号真实频率f₀相对应，2π与采样频率f_s相对应，故正整数p实际上就是真实频率f₀与采样频率f_s的整数倍值(p是粗略估计的整数)，在实际应用中往往容易事先得知。因而要想精确估计出信号的真实频率，需要精确地去估计主值区间内的ω^*值。而ω^*正是欠采样后的波形g(n)的数字角频率值。

借助离散傅立叶变换(Discrete Fourier Transform，DFT，其快速算法为FFT)，便可进一步近似测出在主值区间[0，2π)内的ω^*值：对采样的信号做阶数为N的FFT，假定所生成的FFT峰值谱的谱线序号为‘q’，则可以粗略测算出ω^*的近似数值为qΔω(q∈z⁺，Δω＝2π/N)，受频率分辨率所限，qΔω与ω^*的误差将会在一个频率间隔Δω内。因而需设计算法去尽量减小这一误差，这就需要借助插值算法，估算出ω^*的真实位置。而频谱泄漏是影响该插值算法精度的重要因素，由于全相位FFT相比于传统FFT，其谱泄漏程度大为降低，故基于全相位FFT谱分析的频率插值算法的估计精度肯定比传统FFT高。故本专利先将问题转化为用基于全相位FFT的谱插值算法去估计ω^*，粗略得到整数‘p’估计后，从而得到原信号的角频率ω估计，再依据f₀＝ω/2π·f_s估算频率。

而幅值和相位信息无需转化，可直接从欠采样后的波形得到。这是因为欠采样后的信号g(n)与原信号x(t)＝Acos(2πf₀t+θ₀)的关系为

$g\left(n\right)=x\left(t\right){|}_{t=n{T}_s}=A\cos (2\mathrm{\π}{f}_{0}n{T}_s+{\mathrm{\θ}}_0)---\left(1\right)$

故虽然欠采样后数字角频率ω^*＝2πf₀T_s值需间接求出，但从式(1)可看出g(n)仍是余弦波表达式，x(t)的幅值即为g(n)的幅值，t＝0时的x(t)的初相位θ₀即为g(n)在n＝0时的初相位。因而，直接估计正弦波序列g(n)的幅值和初相即估计出了x(t)的幅值和初相。

进一步分解式(1)，有

$g\left(n\right)=\frac{A}{2}{e}^{j2\mathrm{\π}{f}_{0}n{T}_s}{e}^{j{\mathrm{\θ}}_0}+\frac{A}{2}{e}^{-j2\mathrm{\π}{f}_{0}n{T}_s}{e}^{-j{\mathrm{\θ}}_0}=\frac{A}{2}{e}^{j{\mathrm{\ω}}^{*}n}{e}^{j{\mathrm{\θ}}_0}+\frac{A}{2}{e}^{-j{\mathrm{\ω}}^{*}n}{e}^{-j{\mathrm{\θ}}_0}---\left(2\right)$

从式(2)可看出，单频信号x(t)＝Acos(2πf₀t+θ₀)实际包含频率为ω^*和-ω^*的两个复指数成分(这与图4c和图5c中的左、右两个对称的谱峰相对应)，前者的相位值为θ₀，则后者的相位值为-θ₀。但当采样频率取值不同时，由全相位FFT分析得到的这两个频率成分对应的谱线位置会有所不同。

对于香农采样f_s＞2f₀情况，则有p＝0、ω^*＝2πf₀/f_s＜π，除以频率分辨率Δω＝2π/N后，则有：k＝[ω^*/Δω]＜N/2，结合式(2)，这样必然在全相位FFT谱分析相位图的左半轴峰值谱上出现初相信息，而在右半轴峰值谱上出现负初相信息(如图4(d)所示)。

而对于欠采样情况，其相位谱正、负相位分布则要进一步细化考虑。由于ω＝2pπ+ω^*，当ω∈[2pπ，2pπ+π]时，则必然有ω^*∈[0，π)，这时由全相位FFT分析得到的左半轴峰值相位谱值是真实值θ₀；当ω∈[2pπ+π，2pπ+2π)时，则必然有ω^*∈[π，2π)，这时由全相位FFT分析得到的左半轴峰值相位谱值是真实值的相反数-θ₀(图5(d)即符合这种情况)。因而事先要有ω大致是落在[2pπ，2pπ+π)内还是[2pπ+π，2pπ+2π)内的信息，这就对应着真实频率f₀是落在[pf_s，pf_s+f_s/2)内还是[pf_s+f_s/2，pf_s+f_s)内的粗略信息，这在实际中是不难做到的。

因而关键在于精确地而不是大致地估计出ω^*值，本专利提出采用基于全相位FFT谱分析/传统FFT谱分析的综合估计法来完成，基于此从而得到幅值A的估计。

(4)基于全相位FFT谱/传统FFT谱分析的频率估计与幅值估计

由于全相位FFT同传统FFT一样，都是离散谱分析，所以必然存在栅栏效应，这使得实际观测到的谱峰位置只能是位于频率分辨率Δω的整数倍位置上(如图7所示，假定其谱序号为q)，因而频谱分析的结果必然是在峰值谱线周围等间隔地泄漏出旁谱线，次高谱线既可能出现在峰值谱线的右边(如图7a所示)，也可能出现在峰值谱线的左边(如图7(b)所示)，而信号实际谱峰位置却不能直接观测得到，它存在于峰值谱线和次高谱线之间(如图7(a)、7(b)虚线所示)，两者存在Δk＜0.5个频率间隔，因而需要估测Δk值，这实际上是个插值问题。

可综合利用信号FFT谱的相位信息和全相位FFT谱的相位信息与幅值信息来估计信号的数字角频率和幅值，如图8所示，对欠采样后的信号g(n)进行传统FFT得到其峰值谱G(q)，进行全相位FFT得到其峰值谱G_a(q)，取G(q)模的平方得到功率谱值P_g(q)，将P_g(q)与G_a(q)模除后，即可得到幅值估计取FFT峰值谱G(q)的相位值与全相位FFT峰值谱G_a(q)的相位值

的插值除以τ＝(N-1)/2后得到频偏估计Δk，最后将ΔkΔω与qΔω叠加即得数字角频率估计

图8中，假设FFT谱阶数为N，则FFT需N个样点，全相位FFT需2N-1个样点，图2其他过程并不耗费多余样点，故整个参数估计过程仅需2N-1个样点。

(5)实验效果

图5是针对信号频率还不够高，欠采样后失真不够大时的波形参数估计情况。图9同时给出了以采样速率f_s＝23kHz对频率为200kHz、初相角为40°、幅值为2的高频余弦信号进行欠采样的波形及其全相位FFT谱分析(其谱分析阶数N＝32)结果。从图9(b)可以看出，这时欠采样波形已经严重失真，已经完全偏离了余弦信号的形状。而从图9(d)的峰值相位谱图中，仍可以很准确地读出真实相位值40°。按本发明提出的参数估计流程，可估算出其频率为200000.026Hz，误差仅为0.00013％；而其幅值估计为1.9888374，误差仅为0.5581％。

表1给出了对频率为200kHz的高频信号，分别在19.1KHz～209.6KHz范围内以不同的欠采样频率f_s(均低于香农不失真频率400kHz)对该信号进行采样，再按照图2本专利方法对频率、幅值和相位进行估计，其估计结果如表1所示(谱分析阶数N＝16)。

表1不同采样速率下对高频正弦信号频率、幅值和相位的估计结果

采样速率fs(×103Hz)	频率估计(×103Hz)(真实值：200×103Hz)	幅值估计值(真实值A＝2)	相位估计值(真实值40°)	整数p估计	左半谱峰位置q
						19.1	200.07898486	2.11885088	38.216839	10	8

采样速率fs(×103Hz)	频率估计(×103Hz)(真实值：200×103Hz)	幅值估计值(真实值A＝2)	相位估计值(真实值40°)	整数p估计	左半谱峰位置q
						39.2	200.00228700	1.95308060	39.999575	5	3
59.3	200.00380478	1.99935268	39.999398	3	7
						68.4	200.02125990	1.97919965	39.957112	2	2
89.5	199.99986548	1.97989198	39.999999	2	5
						109.9	199.99977610	1.99540015	39.999965	1	4
138.5	200.00445032	1.99663685	39.999713	1	8
						168.3	200.00098598	1.99994489	39.999950	1	4
189.7	199.99772890	1.99440317	39.999646	1	2
						209.6	200.04270896	1.98116780	39.925003	0	2

从表1可看出，本专利方法对采样频率的要求很宽松，即使是以低于信号频率10倍的采样频率对信号进行采样，按照本专利提出的估计方法，都能很精确地估计出频率、幅值和相位的参数值。

由于本专利方法每次进行参数估计仅需2N-1个样点，故获得表1实验数据仅仅需要2×16-1＝31个样点。由表1可看出，当f_s＝89.5KHz时，其频率估计误差仅为0.14Hz，相位估计误差不到0.00001度，幅值估计误差不到0.02。故其成本小，而精度较高。

应指出，各种不同采样频率下的测量精度差异，并不是没有规律的。从表1中最后一列谱峰位置可看出，当谱峰位置q＝8(对应f_s＝19.1KHz、138.5KHz)、q＝2(对应f_s＝68.4KHz、189.7KHz、209.6KHz)时，其测量误差相对来说比较大，最高达到几十Hz；而当谱峰位置q＝5(对应f_s＝89.5KHz)、q＝4(对应f_s＝109.9KHz、168.3KHz)时，其测量误差相对比较小，可以控制在1Hz以内。本例中谱分析阶数N＝16，也就是说，当谱峰位置约为q＝N/4附近时，本专利方法的测量精度是很高的。这可从图9的振幅谱图上得到解释，因为频率轴是周期延拓的，而实信号的谱都包含左、右两簇谱峰(即两个边带)，故当谱峰位于直流附近(q＝2)或奈奎斯特频率对应的N/2附近(q＝8)，左、右两边带谱实际上是紧挨的，相互会形成较大的谱间干扰从而降低其测量精度。而当谱峰位置位于q＝N/4附近时情况则不然，这时左、右两边带谱隔得很开，谱间干扰小，故其测量精度高。因而可通过识别谱峰位置是否位于q＝N/4附近来作为采样速率选择的判断依据。

本发明提出的基于欠采样情况下参数估计方法，应用在各工程领域，可以产生如下的有益效果。

第一.对模数转换器的采样速率不做限制，降低了硬件成本。

表1的实验结果说明了，ADC的采样速率可以远远低于信号频率，即使在这种情况下也能高精度地估计出信号的频率值，幅值和相位值也能得出。因此本文方案没有对ADC的采样速率进行限制。事实上，采样速率越高，ADC的工艺就越复杂，市面上的高速ADC价格都很贵。本专利方法对ADC采样速率不做限制，可以大大节省硬件成本。

第二.测量精度高、抗噪能力强。

表1的数据说明了该优点，这也可以从原理上给予解释：原有的基于时钟计数的测量方法，是在信号波形的两个过零点之间，通过对时钟脉冲计数来达到测频的目的，因此当噪声干扰时，会使得过零点产生漂移，计数自然就不准，从而导致测频精度降低。而本专利方法则不一样，因为本专利方法是频域测量方法，噪声是全频带分布，故在频域中，噪声干扰会分摊到所有的N根谱线上，而本专利的测频方法仅用到了峰值谱线和次高旁谱线，故受噪声干扰小、精度高。另外，谱泄漏是影响测量精度的另一个因素，本专利引入的全相位FFT，文献[19]指出其谱泄漏程度比传统FFT小得多，这也是本专利方法精度高的一个主要原因。

第三.可大大扩展频率计仪表的频率测量范围。

由于取消了采样速率和信号频率间在数值上的限制，因而当信号频率很高时，完全可以用本文方法以低速率对信号进行采样，当信号频率较低时，也可以按传统方式，用比信号更高的采样速率进行采样。故无论所测量信号的频率是低还是高，仪器内部的ADC都可以自适应地应付其变化，故可大大扩展频率计仪表的测量范围。

第四.存储量小，仅需内部RAM即可，无需外置存储器。

由于本专利进行频率估计仅需要2N-1个样点，且表1实验数据表明，谱分析阶数N不需要很大。当N＝16只需存储31个样点，便可将频率估计误差限制在几个Hz以内，因而无需外部扩展存储器，只需数字处理器件(DSP)的内部RAM就行，这易于仪表集成化。

下面结合附图和实施例进一步详细说明本发明。

下面首先对实施本发明的硬件予以简单说明。参见图10，为精确估计出输入信号x(t)的频率参数，需借助信号调理电路对输入信号进行模拟预处理，以对信号幅度范围进行必要调整，并去除外干扰噪声等；再经过A/D(模数转化器)采样得到样本序列x(n)，以并行数字输入的形式进入DSP器件(数字信号处理器，Digital Signal Processor)，经过DSP器件的内部算法处理，从而得到信号参数的估计，最后借助输出驱动及其显示模块显示出信号估计值，即图10的整个系统构成了一个“高精度宽范围频率参数估计仪”。

需指出的是，(1)为扩大频率计的测量范围，同时为保证高精度的测量结果，图10的采样率可以由DSP控制。这只需要控制DSP的CLKOUT输出时钟CP₂(该时钟小于外部主时钟CP₁)即可实现，这是因为ADC的采样速率完全由外部馈给的时钟CP₁决定。而不同频率的信号，可能需要不同的采样速率，如高频信号需要欠采样速率，而低频信号则需要正常采样速率，不管哪种情况，可以通过判断谱峰位置是否位于q＝N/4附近来判断，当所选的采样速率，使得谱分析得到的谱峰位置位于q＝N/4附近，就可认为精度足够高了。而DSP是可以很方便的通过控制输出时钟来改变外部ADC的采样速率；(2)需要把和各种采样速率对应的正整数p存入到内部寄存器中，这个值是个粗略值，因为用户在进行测量前，肯定可以估算信号频率大概处于哪个范围内，p整数值表征的就是粗略范围估计。

图10内部的DSP参考程序流程如图11所示。

本发明将所提出的“欠采样频率估计”的核心估计算法结合图11的处理流程植入DSP器件内，基于此完成高精度、高效的信号频率参数估计。需指出，由于采用了数字化的估计方法，因而决定图10系统的复杂度、实时程度和稳定度的主要因素并不是图10中的DSP器件的外围连接，而是DSP内部程序存储器所存储的核心估计算法。DSP器件的内部程序流程如图11所示。

在器件选择上，图10的DSP选用通用器件(如TMS5000、6000序列，AD公司的ADSP21XX系列等)即可，因为通用DSP都内置RAM以及具有时钟输出管脚，可以通过编程控制CLKOUT管脚的输出时钟频率。而模数转换器ADC则建议选择支持中速转换的逐次比较型ADC(如AD7679，最高转化速率570Kbps)即可。

Claims (3)

1.一种欠采样下的正弦信号频率的高效测量方法，其特征是，包括下列步骤，设要测出频率、幅值和相位值的正弦波信号为x(t)＝Acos(2πf₀t+θ₀)，

首先以低于真实频率值f₀的采样速率f_s对信号x(t)进行模数转换而得到欠采样后的离散信号g(n)＝x(nT_s)，T_s＝1/f_s；

然后对g(n)分别做传统FFT谱分析和全相位FFT谱分析，而得到传统FFT谱G(k)k＝0，1，...，N-1和全相位FFT谱G_a(k)，k＝0，1，...，N-1，取G(k)的模值即得传统振幅谱|G(k)|，取G(k)的相角即得传统FFT相位谱

同样地，取G_a(k)的模值即得全相位FFT振幅谱|G_a(k)|，取G_a(k)的相角即得全相位FFT相位谱

再对全相位FFT振幅谱|G_a(k)|进行谱峰搜索，并记录谱峰所在的k＝q的位置，从该峰值位置直接读出全相位FFT相位谱值即可作为相位测量值另外，在谱峰k＝q位置很容易计算出传统FFT的功率谱值P_g(q)＝|g_a(q)|²，将P_g(q)除以谱峰位置的全相位振幅谱值|G_a(q)|后，得到幅值估计取谱峰k＝q位置的FFT相位值

与全相位FFT相位值

的差值除以τ＝(N-1)/2后得到频偏估计Δk，最后将ΔkΔω与qΔω叠加即得数字角频率估计

和

直接作为所测信号的相位和幅值估计；

确定出f₀和f_s的整数倍值p，将p和f_s值与

叠加后，即可作为最终的频率估计值

2.根据权利要求1所述的一种欠采样下的正弦信号频率的高效测量方法，其特征是，全相位FFT谱分析是，用长为(2N-1)的卷积窗w_c对中心样点x(0)前、后(2N-1)个数据进行加权，然后将间隔为N的数据两两进行重叠相加，再对重叠相加后的数据进行DFT即得全相位谱分析结果X_a(k)，k＝0，1，...，N-1，再取其模值后，即得全相位FFT振幅谱|X_a(k)|，取其相角后，即得全相位FFT相位谱

3.一种欠采样下的正弦信号频率的高效测量实施装置，其特征是，包括：

信号调理电路，用于对输入信号进行模拟预处理，对信号幅度范围进行必要调整，并去除外干扰噪声；

模数转化器A/D，用于采样得到样本序列x(n)，以并行数字输入的形式进入数字信号处理器DSP；

数字信号处理器DSP，用于对g(n)分别做传统FFT谱分析和全相位FFT谱分析，由传统FFT得到g(n)的峰值谱G(q)，由全相位FFT谱分析得到g(n)的峰值谱G_a(q)，从全相位FFT得到的峰值谱上直接读出相位值θ₀，取G(q)模的平方得到功率谱值P_g(q)，将P_g(q)与G_a(q)模除后，得到幅值估计A；取FFT峰值谱G(q)的相位值与全相位FFT峰值谱G_a(q)的相位值

的差值除以τ＝(N-1)/2后得到频偏估计Δk，最后将ΔkΔω与qΔω叠加即得数字角频率估计

和

直接作为所测信号的相位和幅值估计；

确定出f₀和f_s的整数倍值p，将p和f_s值与

叠加后，即可作为最终的频率估计值

输出驱动及显示模块，用于显示出信号估计值。

Images