CN101713795A - 一种在采样率限制下数字测量频率的方法 - Google Patents

Abstract

本发明涉及一种在采样率限制下数字测量频率的方法，技术特征在于：待测频率信号进行AD转换，转换的待测信号数据按照M个为一组；将相连两秒的两组数据做互相关，得到互相关函数值；计算两组数据的自相关函数值，根据互相关量计算频偏，则待测信号的实际频率值为标称值与频偏量之和。本发明解决频率数字化精密测量的理论和技术问题，为高精度、低成本的基于虚拟仪器技术的频率数字化测量仪器的设计和实现提供一整套完整解决方案。

Description

一种在采样率限制下数字测量频率的方法

技术领域

本发明涉及一种在采样率限制下数字测量频率的方法，能够解决采样率限制的数字精密测频方法，可以解决现有模数转换设备难以达到采样率要求，不能实现中、高频段频率的数字测量。

背景技术

目前，在高精度的频率测量方法和设备方面，主要的设备包括：TimeTech(德国)的比相仪和JPL(美国)的振荡器稳定度分析仪。这两套系统的主要测频原理是双混频时差测量，其主要结构如图1所示。双混频时差测量法具有平衡的结构，能够抵消公共噪声，但是因为系统结构复杂，以及测量能力受公共振荡器性能影响，测量设备造价高达几十万元。

随着微处理器和大规模集成电路的迅速发展，在频率测试系统中，越来越多的传统的测量方法被数字化测量方法所取代。基于数字化技术的电网频率测量技术已经得到广泛的应用。现今的模数转换器件能够容易的实现对电网50Hz左右频率信号的数字化测量转换需求。奈奎斯特采样定理要求若把模拟信号转换成数字信号，采样频率通常至少为其最大频率的两倍，当模拟信号的频率非常高的时候，按此采样频率采样将在器件价格和系统复杂度方面付出高昂的代价。实际上一方面由于受到器件工作速度的限制，当前ADC的采样率仍不能做得很高，另一方面，即使信号采样率达到很高的程度，采样后处理器对信号的处理速度也可能达不到要求。因此为了能够在一定程度上克服了这两方面的不足，许多折衷方案被提出，如各种基于谱估计的欠采样方法，由于欠采样使得信号频谱混叠，必须增加额外的信息解模糊才能得到信号频率的无模糊估计。基于谱估计算法的欠采样频率估计仍是基于恢复原采样信号的思想，目前来看，这类方法总的来说设计复杂度高，并且精度提高潜力并不明显，进一步提高精度对设备的要求也比较高，因此很难满足日益增长的时频测量精度需求。这样的一些原因也直接导致了目前的数字测频技术的发展受到限制，只在电力线(50Hz)等低频段得到了较好的开发利用。

鉴于以上两方面的原因，现有方法中必须解决高精度频率测量系统的结构精简和中低频段频率信号的数字化测量实现两方面的问题。具体要解决三个方面的关键问题：

1)频率信号特征信息的提取，主要是基于提取频率信号细节信息的欠采样方法，且能够实现对频率源瞬时变化的检测；

2)基于数字信号处理技术的频率估计，在欠采样方法的前提下对频率的估计需要设计的数字信号处理方法，要满足在有限信息条件下实现对精密频率源的高精度测量，同时为了较好的实时性需求，算法不宜过于繁复；

3)量化误差、随机干扰噪声的抑制方法，对于任何精密测量方案来说，噪声的分析评估都是至关重要的，由于进行模数转换引入的误差，电磁环境、电力电压因素等等对测量结果的影响以及补偿方法都是必须的。

发明内容

要解决的技术问题

本发明的目的是克服现有技术中的不足之处，提供了一种在采样率限制下数字测量频率的方法，可以实现对5MHz、10MHz等多种频率值信号的高精度数字测量。

技术方案

本发明的技术特征在于将频率标称值已知的待测频率信号送入数据采集卡的输入通道，同时将具有待测信号标称频率值的参考频率信号送入分频设备的输入端口，采用以下方法实现待测频率信号的分析处理，具体步骤如下：

步骤1：设定分频参数为待测信号标称值除以一个整倍数后有一个固定偏差，偏差可以取1～100赫兹之间的任意值，并将该分频参数作为数据采集卡中ADC转换的采样时钟触发信号；

步骤2：待测频率信号在采样时钟触发信号的控制下进行AD转换，转换的待测信号数据按照M个为一组；所述的M等于每一秒间隔内在采样时钟触发下待测信号经AD转换得到的数据个数；

步骤3：将相连两秒的两组数据做互相关，得到互相关函数值：

$R_{y_i{y}_{i+1}}\left(0\right)=\frac{1}{f_x}\underset{n=0}{\overset{f_s-1}{\mathrm{\Σ}}}{y}_i\left(n\right)y_{i+1}\left(n\right),$ 其中y_i+1(n)＝y_i(n+f_s)；

其中f_s为采样时钟触发信号的频率值；y_i为相连两秒数据的前一秒M个数据；y_i+1为相连两秒数据的后一秒M个数据；n的定义域是[0～(f_s-1)]的整数值。

步骤4：然后与步骤3相同的两组数据分别用 $R_{y_i}\left(0\right)=\frac{1}{f_s}\underset{n=0}{\overset{f_s-1}{\mathrm{\Σ}}}{y}_i\left(n\right)y_i\left(n\right)$ 和 $R_{y_{i+1}}\left(0\right)=\frac{1}{f_s}\underset{n=0}{\overset{f_s-1}{\mathrm{\Σ}}}{y}_{i+1}\left(n\right)y_{i+1}\left(n\right)$ 公式计算出y_i和y_i+1的自相关函数值和

A_i表示采样后待测信号的第i秒的幅度值；

$A_i=\sqrt{2R_{y_i}\left(0\right)}$ $A_{i+1}=\sqrt{2R_{y_{i+1}}\left(0\right)}$

步骤5：根据步骤3得到的互相关量和步骤4得到的A_i和A_i+1的值计算频偏Δf：

$\mathrm{\Δf}=\frac{1}{2\mathrm{\π}}\mathrm{ar}\cos \frac{R_{y_i{y}_{i+1}}\left(0\right)}{\sqrt{R_{y_i}\left(0\right)R_{y_{i+1}}\left(0\right)}};$

步骤6：则待测信号的实际频率值为标称值f₀与频偏量Δf之和。

有益效果

本发明提出在采样率限制下数字测量频率的方法，以频率的数字化处理技术为基础，通过欠采样手段采集包含频率信号特征信息的数据，研究在目前模数转换速率水平条件下频率数字精密测量的方法，推进精密频率测量仪器的虚拟化进程。

本发明能够实现频率的数字化测量，基于欠采样技术、相关处理算法、曲线拟合等方法结合虚拟仪器技术，共同实现高精度频率数字化测量，该发明能对10MHz频率1秒取样的Allan方差稳定度优于10^-13量级，远高于现今数字测频的技术指标。

本发明解决频率数字化精密测量的理论和技术问题，为高精度、低成本的基于虚拟仪器技术的频率数字化测量仪器的设计和实现提供一整套完整解决方案。

附图说明

图1：现有技术中双混频时差测量方法原理图

图2：本发明方法的原理图

图3：一种解决采样率限制的数字精密测频方法的硬件结构图

具体实施方式

现结合实施例、附图对本发明作进一步描述：

本实施例的测频方法的实现包括硬件：数据采集卡、分频设备和计算机三部分，软件实施测频算法安装在计算机上运行。

数据采集卡包括模数转化器件和PCI接口、通信部分，美国国家仪器公司(NI)的PCI-6122型数据采集卡，具有500KHz采样率和16位分辨率。

实施例结合国家授时中心的10MHz时频基准信号，通过测量该信号来验证本发明的测频能力。

国家授时中心的10MHz时频基准信号作为分频设备的参考输入信号，待测频率是由铷原子钟输出的10MHz频率信号，幅度范围约在±1V，需要测量该铷原子钟的实际输出频率值。

步骤1：根据待测频率信号的标称频率值和数据采集卡的采样率上限，确定分频倍数N＝100，即分频设备输出的采样时钟频率为100KHz-1Hz，将该采样时钟信号送入数据采集卡的PFI端口，作数据采集卡中ADC转换的采样时钟用；

步骤2：根据步骤1得到的采样时钟在数据采集卡的参数配置表中设置数据采集卡的采样时钟参数为100KHz-1Hz；

步骤3：将待测频率信号接入数据采集卡的四个输入通道之一，如通道0，信号在该通道按照采样时钟触发信号进行AD转换，经过转换后的待测信号数据按照每100K-1个一组储存在磁盘指定的区域，；

步骤4：将相连的两组数据做互相关运算，有， $R_{y_i{y}_{i+1}}\left(0\right)=\frac{A_i{A}_{i+1}}{2}\cos \left(2\mathrm{\π\Δf}\right)$ 本实施例的某一组互相关函数计算得 $R_{y_i{y}_{i+1}}\left(0\right)=\frac{1}{f_s}\underset{n=0}{\overset{f_s-1}{\mathrm{\Σ}}}{y}_i\left(n\right)y_{i+1}\left(n\right)=0.50000000308438;$

步骤5：另外用与步骤4相同的两组数据分别计算出y_i和y_i+1的自相关函数 $R_{y_i}\left(0\right)=0.50000000517618$ 和 $R_{y_{i+1}}\left(0\right)=0.50000000101425$ 的值，从而根据 $A_i=\sqrt{2R_{y_i}\left(0\right)}$ 和 $A_{i+1}=\sqrt{2R_{y_{i+1}}\left(0\right)}$ 得到A_i和A_i+1的值；

步骤6：根据步骤4得到的互相关量和步骤5得到的幅度值求频偏Δf，根据方程 $\mathrm{\Δf}=\frac{1}{2\mathrm{\π}}\mathrm{ar}\cos \frac{R_{y_i{y}_{i+1}}\left(0\right)}{\sqrt{R_{y_i}\left(0\right)R_{y_{i+1}}\left(0\right)}}=1.047755927109550e-006$ 得到Δf值。因为待测频率信号标称值已知为10MHz，因此待测信号的实际频率值为1e7+1.047755927109550e-006Hz。

由上述实施例可以看出，本方法的主要特点是通过设置采样时钟频率，提取待测频率信号的频率偏差信息，降低测量中高频信号的频率对数据采集设备AD转换速率的要求，利用较低采样率实现中高频频率信号的频率测量，该方法采样时尽管不需要满足奈奎斯特定理的采样率采样待测频率信号，但是对于待测信号的频率偏差信号(频率值表达式为N+ΔfHz，N为分频器的分频参数，是整数，Δf是待测信号的频偏量)需采用高于该信号2倍频率的采样率进行采样，该采样率的设置根据用户的测量要求确定，一般情况是高采样率对应高分辨率测量结果，但是超过100倍以上的采样率的性价比较低，建议是10倍～100倍之间。

对于采样噪声，每一点的采样噪声包含两个方面，随机噪声和采样时基偏差，采样后的信号通过相关估计估算频率值，相关处理对信号噪声有抑制作用，相关方法基于整个周期间相似度的比较，以及其求和平均的算法能有效的平滑外部噪声和系统随机噪声对采样数据测量结果的影响。根据随机噪声之间不相关性质，对随机噪声乘积求均值，噪声对结果的影响与采样率有关，不会超过采样点噪声总和的

例如，对1Hz正弦波进行10kHz采样，即一个周期内采样10000个点，能平滑噪声至原来的 $1/\sqrt{10000}=0.01,$ 有利于提高测量能力。根据前文所述的实例，本方法测量频率分辨率达到1e-13量级。可以在当前工艺水平条件下大幅度提高测量分辨率，并且测量结果便于使用计算机进行实时显示和保存，后端处理更加灵活。

本发明提出了一种高精度数字测频的解决方案，基于欠采样技术的精密数字测频。主要是在目前模数转换速率和数字信号处理技术水平下，数字高精度测频的实现方法。目前常用的双混频时差测量方法采用模拟结构，测量精度受电路噪声和公共源影响较大，依靠电路工艺改进提高测量精度是个缓慢的过程；而数字测频受限于采样率和数字信号处理能力的限制尚不能用于中高频信号的测量。基于此种现状，通过研究数字采样和精密频率测量的特点，提出了基于提取特征信息的欠采样频率测量方法。以远低于原信号频率的采样率取样，通过研究特别的算法实现频率测量。并充分开发数字信号处理技术对电路噪声的平滑作用，在现有工艺水平下提高测量精度。

多通道测量结构设计和虚拟仪器控制面板设计降低了设备成本。

Claims (1)

1.一种在采样率限制下数字测量频率的方法，其特征在于：将频率标称值已知的待测频率信号送入数据采集卡的输入通道，同时将具有待测信号标称频率值的参考频率信号送入分频设备的输入端口，采用以下方法实现待测频率信号的测量，具体步骤如下：

步骤1：设定分频参数为待测信号标称值除以一个整倍数后有一个固定偏差，偏差可以取1～100赫兹之间的任意值，并将该分频参数作为数据采集卡中ADC转换的采样时钟触发信号；

步骤2：待测频率信号在采样时钟触发信号的控制下进行AD转换，转换的待测信号数据按照M个为一组；所述的M等于每一秒间隔内在采样时钟触发下待测信号经AD转换得到的数据个数；

步骤3：将相连两秒的两组数据做互相关，得到互相关函数值：

$R_{y_i{y}_{i+1}}\left(0\right)=\frac{1}{f_s}\underset{n=0}{\overset{f_s-1}{\mathrm{\Σ}}}{y}_i\left(n\right)y_{i+1}\left(n\right),$ 其中y_i+1(N)＝y_i(n+f_s)；

其中f_s为采样时钟触发信号的频率值；y_i为相连两秒数据的前一秒M个数据；y_i+1为相连两秒数据的后一秒M个数据；n的定义域是[0～(f_s-1)]的整数值。

步骤4：然后与步骤3相同的两组数据分别用 $R_{y_i}\left(0\right)=\frac{1}{f_s}\underset{n=0}{\overset{f_s-1}{\mathrm{\Σ}}}{y}_i\left(n\right)y_i\left(n\right)$ 和 $R_{y_{i+1}}\left(0\right)=\frac{1}{f_s}\underset{n=0}{\overset{f_s-1}{\mathrm{\Σ}}}{y}_{i+1}\left(n\right)y_{i+1}\left(n\right)$ 公式计算出y_i和y_i+1的自相关函数值

和

A_i表示采样后待测信号的第i秒的幅度值；

$A_i=\sqrt{2R_{y_i}\left(0\right)}$ $A_{i+1}=\sqrt{2R_{y_{i+1}}\left(0\right)}$

步骤5：根据步骤3得到的互相关量和步骤4得到的A_i和A_i+1的值计算频偏Δf：

$\mathrm{\Δf}=\frac{1}{2\mathrm{\π}}\mathrm{ar}\cos \frac{R_{y_i{y}_{i+1}}\left(0\right)}{\sqrt{R_{y_i}\left(0\right)R_{y_{i+1}}\left(0\right)}};$

步骤6：则待测信号的实际频率值为标称值f₀与频偏量Δf之和。

Images