CN104897962B - 基于互素感知的单频信号短样本高精度测频方法及其装置 - Google Patents
基于互素感知的单频信号短样本高精度测频方法及其装置 Download PDFInfo
- Publication number
- CN104897962B CN104897962B CN201510346171.4A CN201510346171A CN104897962B CN 104897962 B CN104897962 B CN 104897962B CN 201510346171 A CN201510346171 A CN 201510346171A CN 104897962 B CN104897962 B CN 104897962B
- Authority
- CN
- China
- Prior art keywords
- frequency
- fourier transform
- discrete fourier
- spectral
- peak
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Fee Related
Links
- 230000008447 perception Effects 0.000 title claims abstract description 19
- 238000000691 measurement method Methods 0.000 title claims abstract description 8
- 238000001228 spectrum Methods 0.000 claims abstract description 84
- 230000003595 spectral effect Effects 0.000 claims abstract description 79
- 238000000034 method Methods 0.000 claims abstract description 78
- 238000012937 correction Methods 0.000 claims abstract description 32
- 238000005070 sampling Methods 0.000 claims description 51
- 238000012545 processing Methods 0.000 claims description 22
- 238000001914 filtration Methods 0.000 claims description 19
- 238000004364 calculation method Methods 0.000 claims description 13
- 238000005259 measurement Methods 0.000 claims description 12
- 230000001934 delay Effects 0.000 abstract 1
- 230000008569 process Effects 0.000 description 15
- 238000002474 experimental method Methods 0.000 description 9
- 238000010586 diagram Methods 0.000 description 7
- 238000012935 Averaging Methods 0.000 description 5
- 238000013461 design Methods 0.000 description 4
- 238000010183 spectrum analysis Methods 0.000 description 4
- 230000007547 defect Effects 0.000 description 3
- 230000005284 excitation Effects 0.000 description 3
- 230000006870 function Effects 0.000 description 3
- OIGNJSKKLXVSLS-VWUMJDOOSA-N prednisolone Chemical compound O=C1C=C[C@]2(C)[C@H]3[C@@H](O)C[C@](C)([C@@](CC4)(O)C(=O)CO)[C@@H]4[C@@H]3CCC2=C1 OIGNJSKKLXVSLS-VWUMJDOOSA-N 0.000 description 3
- 230000004044 response Effects 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 2
- 230000008859 change Effects 0.000 description 2
- 230000000739 chaotic effect Effects 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 230000008094 contradictory effect Effects 0.000 description 1
- 125000004122 cyclic group Chemical group 0.000 description 1
- 238000012938 design process Methods 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000036039 immunity Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000010606 normalization Methods 0.000 description 1
- 230000002093 peripheral effect Effects 0.000 description 1
- 238000000819 phase cycle Methods 0.000 description 1
- 238000003672 processing method Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 229910052709 silver Inorganic materials 0.000 description 1
- 239000004332 silver Substances 0.000 description 1
- 230000001629 suppression Effects 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 238000010200 validation analysis Methods 0.000 description 1
Landscapes
- Radar Systems Or Details Thereof (AREA)
- Measurement Of Resistance Or Impedance (AREA)
Abstract
本发明公开了一种基于互素感知的单频信号短样本高精度测频方法及其装置，方法包括：对由多通道稀疏信号组成的多通道输出序列进行离散傅里叶变换，得到离散傅里叶变换后的输出序列；对离散傅里叶变换后的输出序列作谱峰扫描，得到谱峰位置的索引，通过索引获取峰值离散傅里叶变换谱线，采用相位差法进行谱校正，获取校正后的频率；将校正后的频率对作为余数，按照闭合解析形式的中国余数定理进行处理，重构出原始信号频率。装置包括：数字信号处理器、输出驱动及其显示电路，数字信号处理器用于上述的计算，重构出原始信号频率；输出驱动及其显示电路对原始信号频率进行显示。本发明降低了频谱估计的时间延迟，提高了样本利用率。
Description
技术领域
本发明涉及数字信号处理领域，尤其涉及信号欠采样、互素感知、多相滤波、频谱校正、基于中国余数定理的参数重构。
背景技术
谱估计是与许多实际物理量测量(如大气动力学研究中的雷达风速测量[1]、振动分析中的转速测量[2]等)紧密相关的信号处理方法[3,4]，在电子学、电气学和自动化等各个学科都有广泛的应用，是数字信号处理的基础问题。传统谱估计只有满足香农采样定理(即要求一个信号周期内至少采集2个样点)，信号才不会发生失真，才有可能对信号进行有效的恢复。然而，在许多工程应用(如雷达[5]、阵列波达方向(DOA)估计[6]等领域)中，受到现场环境、采样设备等条件限制(如最高模数转换速率限制)，难以满足香农采样，并且随着信号频率的提高采样难度和设备成本也会大幅度提高，因而稀疏欠采样样本的频谱估计(采样速率fs远小于奈奎斯特采样速率，即2倍信号频率f0)成为学术界和工程界迫切需要解决的问题。
针对此问题，国内外学者提出了很多稀疏谱估计方法，如文献[7-9]提出将中国余数定理(Chinese Remainder Theorem,CRT)与离散傅里叶变换(Discrete FourierTransform,DFT)相结合的频谱测量方法，该方法通过对高频信号进行不同点数的DFT，求取多路欠采样信号的谱余数，将高频信号的问题转换成余数低频信号的问题，利用得到的谱余数来重构未知的高频频率，降低了信号采样难度和对采样设备的要求，然而该方法仅能实现1级下采样，其低数据率下的处理能力有待提升；文献[10]提出基于样本嵌入(Nestedsample)的测量算法，该方法要求对一个信号进行两种不同速率的采样，将同一信号的高速率样本嵌入到低速率样本中，利用样本位置之间的差分阵列(difference-co-array)增加自由度(freedoms)的性质，通过导出两路样本的互相关特性来获得信号的功率谱估计，该方法在整体上降低了信号采样率，并且对于较密集的频率恢复也有较好的效果，然而该方法在局部无法摆脱高速采样器，仍需进行一定的高速采样，因而它并不是真正意义上的稀疏谱估计方法。因此，找到真正意义上的稀疏谱估计方法，彻底摆脱高速采样器的约束，并且可以应用在短样本情况是一个亟待解决的难题。
为解决稀疏谱估计问题，近年来，一种新型的谱估计方法—互素感知(co-primesensing)理论[11-14]受到越来越多的关注，该方法首先对单个模拟输入信号作两路并行的稀疏采样(要求两路下采样因子M、N数值满足互素关系)，然后对得到的两路稀疏样本分别作多相滤波，将单通道信号转换成多通路子信号，再同时对两路并行的多相滤波输出信号做IDFT(Inverse DFT)，最后对IDFT输出的多通道信号做互相关扫描以估计信号频谱位置，从而获得真正意义上的高分辨率谱。该方法因经历了两级下采样(互素下采样和多相滤波下采样)，故非常适合于低数据率应用场合，且结构简单，算法清晰，一定程度上解决了稀疏谱估计问题。
但文献[11]提到的方法有两个不容忽视的缺陷，一是频率估计精度不够高，这是因为在频谱估计过程中并未加入谱校正，最终的谱线只能大致定位在整数位置，人为地忽略了信号频率的小数偏移，使得频谱估计产生了不可避免的误差；二是耗费样本多，这是由于需要对两路信号做互相关处理，为了得到准确的谱估计，需要耗费大量的快拍(snapshot,1个快拍数量等于两路下采样因子的乘积，即MN)平均，这在很大程度上增加了频谱估计的时间延迟，提高了对硬件设备的要求和成本，不利于大规模编程和实际应用。由于实际信号的复杂性，这些情况很可能大量出现，严重影响频谱估计的性能。以上两个缺陷都不能仅仅通过优化算法中某个环节的性能而消除，必须对算法结构进行改进。
发明内容
本发明提供了一种基于互素感知的单频信号短样本高精度测频方法及其装置，本发明降低了频谱估计的时间延迟，提高了样本利用率，详见下文描述：
一种基于互素感知的单频信号短样本高精度测频方法，所述单频信号短样本高精度测频方法包括以下步骤：
对由多通道稀疏信号组成的多通道输出序列进行离散傅里叶变换，得到离散傅里叶变换后的输出序列；
对离散傅里叶变换后的输出序列作谱峰扫描，得到谱峰位置的索引，通过索引获取峰值离散傅里叶变换谱线，采用相位差法进行谱校正，获取校正后的频率；
将校正后的频率对作为余数，按照闭合解析形式的中国余数定理进行处理，重构出原始信号频率。
其中，所述对由多通道稀疏信号组成的多通道输出序列进行离散傅里叶变换，得到离散傅里叶变换后的输出序列的步骤具体为：
对输入信号x(n)进行两路下采样，下采样因子分别为M和N，得到两路稀疏信号；
对两路稀疏信号分别进行多相滤波，输出两个多通道稀疏信号；
在时刻n，对每个由多通道稀疏信号组成的多通道输出序列分别进行离散傅里叶变换，得到离散傅里叶变换后的输出序列。
其中，所述对离散傅里叶变换后的输出序列作谱峰扫描，得到谱峰位置的索引，通过索引获取峰值离散傅里叶变换谱线，采用相位差法进行谱校正，获取校正后的频率的步骤具体为：
分别对离散傅里叶变换后的输出序列作谱峰扫描，分别得到谱峰位置的索引；分别将搜索到的谱峰位置索引代入到相邻快拍，分别获取两根峰值离散傅里叶变换谱线，用相位差法对两根峰值离散傅里叶变换谱线进行谱校正，分别得到校正后的频率。
一种基于互素感知的单频信号短样本高精度测频装置，所述单频信号短样本高精度测频装置包括：
数字信号处理器，所述数字信号处理器用于对由多通道稀疏信号组成的多通道输出序列进行离散傅里叶变换，得到离散傅里叶变换后的输出序列；对离散傅里叶变换后的输出序列作谱峰扫描，得到谱峰位置的索引，通过索引获取峰值离散傅里叶变换谱线，采用相位差法进行谱校正，获取校正后的频率；将校正后的频率对作为余数，按照闭合解析形式的中国余数定理进行处理，重构出原始信号频率；
输出驱动及其显示电路，所述输出驱动及其显示电路用于对原始信号频率进行输出显示。
其中，所述数字信号处理器包括：
第一计算模块，所述第一计算模块用于对输入信号x(n)进行两路下采样，下采样因子分别为M和N，得到两路稀疏信号；
第二计算模块，所述第二计算模块用于对两路稀疏信号分别进行多相滤波，输出两个多通道稀疏信号；在时刻n，对每个由多通道稀疏信号组成的多通道输出序列分别进行离散傅里叶变换，得到离散傅里叶变换后的输出序列。
其中，所述数字信号处理器包括：
第三计算模块，所述第三计算模块用于分别对离散傅里叶变换后的输出序列作谱峰扫描，分别得到谱峰位置的索引；
第四计算模块，用于分别将搜索到的谱峰位置索引代入到相邻快拍，分别获取两根峰值离散傅里叶变换谱线，用相位差法对两根峰值离散傅里叶变换谱线进行谱校正，分别得到校正后的频率。
其中，所述数字信号处理器为DSP。
本发明提出的基于互素感知的单频信号短样本高精度测频方法，若用于频率估计及实际工程领域，可产生如下有益效果：
第一 保证频率估计精度的前提下大大节约了样本数量。
对于传统互素感知频谱估计方法，由于算法最终使用两路稀疏信号的互相关性求取信号的频谱，而相关处理需要大量的样本进行平均，由文献[11]可知，样本平均过程会产生很大的时间延迟，而若想提高频率估计的精度，样本数量必然不能太少，估计精度与估计效率相互矛盾；而本发明提出的谱估计方法，由于最终直接使用CRT算法进行频率重构，不必进行大量样本的平均，只需保证信号处理过程不发生大的失真的样本数量即可，而且，在提取余数和谱余数校正的处理中，只需选取两个处于滤波稳态的相邻的样本即可，大大节省了样本数量，减少了时间延迟。
例如实验1提到的在样本数很少的情况下互素感知方法谱估计会有很多干扰谱，而本发明提出的方法即使在样本数很少的情况下依然能较准确的估计出给定信号频率。
第二 提高了2个以上数量级的频率估计精度。
传统的互素感知频谱估计方法中，由于没有加入谱校正，使得最终重构的频率只能四舍五入到整数，人为地忽略了信号频率的小数偏移，不可避免地引入了误差，而本发明提出的谱估计方法在提取出谱余数之后，对其进行了频谱校正，将其小数偏移加到整数余数之中，使得最终重构的频率精确到小数，提高了频率估计的精度。
由实验1可以看出，互素感知方法只能将频率精确到整数，而本发明提出的方法成功将真实频率的小数偏移估计出来，估计精度更高，实验已验证，本发明可提高2个以上数量级的频率估计精度。
第三 具有较高的抗噪性能。
传统互素感知谱估计方法，最终需要进行大量的样本平均来达到抑制噪声的效果；而本发明提到的方法利用鲁棒性的闭式CRT对于噪声的容差特性，不需要进行大量样本的平均，依然可以很好的抑制噪声干扰。
例如在实验2中，只需SNR≥3dB，就可以在一个很大的信噪比范围内准确的重构信号频率。
第四 适于在低速率应用场合进行频率场合。
这是由采用互素谱感知的结构决定的，由于图1的互素谱感知采用了两级下采样结构(互素下采样和多相滤波下采样)，且IDFT、谱校正、中国余数定理重构均即是在两级下采样后工作，故非常适合于低数据率应用场合，非常有利于降低成本。
附图说明
图1为短样本稀疏谱估计器设计流程图；
图2为多相滤波结构图；
图3为互素感知处理后的谱余数分布图；
其中，(a)为u路谱余数分布示意图；(b)为v路谱余数分布示意图。
图4为互素感知估计器功率谱分布图；
图5为不同信噪比下的RMSEs示意图；
图6为本发明的硬件实施图；
图7为DSP内部程序流图。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚，下面对本发明实施方式作进一步地详细描述。
针对以上不足，本发明从两个方面予以改进，第一，在提高精度方面，使用文献[15]的方法对两路DFT输出结果进行校正，将DFT的峰值谱提供的整数频点值与校正得到的小数偏移值加起来作为“余数”估计值，从而减小了DFT谱估计偏差；第二，在降低耗费样本数量方面，利用CRT对两路DFT谱校正输出的余数估计值进行综合而直接重构得到信号的频率估计，从而避开了文献[11]的互素谱分析流程的互相关求平均过程，因而节省了大量的样本周期，大大降低了频谱估计的时间延迟，提高了样本利用率。
本发明可以解决现有互素感知理论的不足，并且从根本上对算法结构进行了优化和改进，能更好地应用到稀疏谱估计之中，在雷达、波达方向估计和波束成行等领域有广阔的应用前景。
101：对由多通道稀疏信号组成的多通道输出序列进行离散傅里叶变换，得到离散傅里叶变换后的输出序列；
其中，该步骤具体包括以下步骤：
1)对输入信号x(n)进行两路下采样，下采样因子分别为M和N(M和N为任意的互素的整数，具体的取值根据实际应用中的需要进行设定)，得到两路稀疏信号xu(n)和xv(n)；
2)对两路稀疏信号xu(n)和xv(n)分别进行多相滤波，输出两个多通道稀疏信号：xu,p(n)和xv,q(n)，其中0≤p≤M-1，0≤q≤N-1；
3)在每个时刻n，对多通道输出序列{xu,p(n),p＝0,...,M-1}进行M点IDFT(离散傅里叶变换)，得到离散傅里叶变换后的输出序列{uk(n),k＝0,...,M-1}。
类似地，对多通道输出序列{xv,q(n),q＝0,...,N-1}进行N点IDFT，得到离散傅里叶变换后的输出序列{vl(n),l＝0,...,N-1}。
102：对离散傅里叶变换后的输出序列作谱峰扫描，得到谱峰位置的索引，通过索引获取峰值离散傅里叶变换谱线，采用相位差法进行谱校正，获取校正后的频率；
即，分别对IDFT输出{uk(n),k＝0,...,M-1}和{vl(n),l＝0,...,N-1}作谱峰扫描，得到谱峰位置的索引ku和lv。
对于u路，将搜索到的谱峰位置索引ku代入到相邻快拍(即第n个快拍和第n+1个快拍)的IDFT输出uk(n)和uk(n+1)，获取两根峰值IDFT谱线(即和)，用相位差法进行谱校正，得到校正后的频率
对于v路，作类似处理，得到校正后的频率
103：将校正后的频率对作为余数，按照闭合解析形式的中国余数定理进行处理，重构出原始信号频率
以M、N作为模值，将校正后的频率对作为余数，按照闭合解析形式的中国余数定理进行处理，重构出原始信号频率
本发明的总体信号处理框架如下：
当输入信号x(n)为只有一个频率时，设
其中，a0是信号幅度，f0是归一化频率，f0＝(m+δ)Δf，Δf＝1/(MN)，m为整数且m∈[0，MN-1)，δ是小数频率偏移(-0.5≤δ＜0.5)，θ0是初始相位。
对于u路，经过N倍下采样后的稀疏信号xu(n)的角频率值为(m+δ)2π/M，u路的M点IDFT输出(对应频率分辨率为2π/M)可检测出该角频率的去除未知整模糊度nu(nu∈Z+)后的信息，即
(m+δ)2π/M＝nu2π+(k+δ)2π/M (2)
类似地，经过M倍下采样后的稀疏信号xv(n)的角频率值为(m+δ)2π/N,v路的N点IDFT输出(对应频率分辨率为2π/N)可检测出该角频率的去除未知整数模糊度nv(nv∈Z+)后的信息，即
(m+δ)2π/N＝nv2π+(l+δ)2π/N (3)
联立式(2)、式(3)，有
可发现，式(4)与中国余数定理的形式完全一致。由于互素的模值M、N已知，因而只需在u路和v路IDFT输出端通过谱峰搜索确定ku和lv值，通过谱校正得到的两个估计值频偏由中国余数定理即可完全确定未知的整模糊度nu、nv，进而精确估计f0＝(m+δ)Δf。
需要强调的是，根据图1，整个重构过程需要用到多相滤波，因此，对于每个多相支路而言，需要消耗一定数量的“快拍”才能达到稳态，假定原型滤波器的长度为L，需要消耗的“快拍”数的临界值为C，则C与L需满足以下关系式
通过上述的步骤101-104，本发明降低了频谱估计的时间延迟，提高了样本利用率，满足了实际应用中的需要。
本发明实现了可对时域稀疏欠采样信号(采样速率可远低于奈奎斯特采样速率)做谱估计；可提高原互素谱分析器的频率估计精度；其精度相比于到原互素谱分辨率提高2个以上数量级。可大幅度降低原互素谱分析器耗费的样本，仅需两个快拍即可获得高精度频率估计结果。
实施例2
本发明实施例2对实施例1中的涉及已有的双子段相位差频谱校正法、闭式鲁棒中国余数定理重构法进行详细的介绍，详见下文描述：
201：信号互素感知处理；
对于u路处理，经过N倍下采样后的稀疏信号xu(n)有如下形式：
同理，v路M倍下采样后的稀疏信号xv(n)有如下形式
然后进行多相滤波过程，其中u路多相滤波的结构如图2所示。
其中Ep(z)(0≤p≤M-1)为截止频率为π/M的低通滤波器H(z)的线性多相子滤波器，满足
由于各路多相子滤波器Ep(z)的长度相等，故各相子滤波会产生一近似的群延时τ，该群延时会引起时延相位θp，同时，会产生一幅度失真cp，从而各相子滤波输出信号xu,p(n)为
同理，v路各相子滤波输出信号xv,q(n)为
其中，是由于滤波造成的幅度失真，是相位偏移。θp、cp、和都是关于p的函数。然后对u路输出的M路子信号按列进行IDFT(即对变量p做IDFT)，得到如下形式的信号
同理，对v路输出的N路子信号按列进行IDFT(即对变量q做DFT)，得到如下形式的信号
对于u路而言，由式(11)可以看出，由于IDFT是针对xu,p(n)的变量p而言的，且式(9)中xu,p(n)的频率(m+δ)2π/M很可能超出2π，故结合IDFT的循环周期性，只有当式(11)的谱序号k取值为(m+δ)模除M的四舍五入取整结果时，即满足
k＝[(m+δ)modM]＝m mod M (13)
IDFT谱线才出现峰值。
同理，对于v路而言，只有当式(12)的谱序号l取值为(m+δ)模除N的四舍五入取整结果时，即满足
l＝[(m+δ)modN]＝m mod N (14)
DFT谱线才出现峰值。
以两个单成分的指数信号为例，取M＝11，N＝7，令式(1)中对应的频率f0＝36.3456Δf，Δf＝1/(MN)，幅值取值a＝1，初始相位取值θ＝45π/180。则根据式(13)和式(14)算出的理论上的信号谱余数的取整值分别为k＝4，l＝2。按照互素谱分析流程进行谱分析，得到的IDFT幅值谱分别如图3(a)、图3(b)所示。
图3(a)为u路谱余数分布，(b)为v路谱余数分布。由图3可以看出，两路IDFT谱的峰值谱位置ku、lv与理论谱余数的取整值完全一致，可对其进行谱峰搜索得到。
202：谱余数校正；
如前所述，由于经过IDFT变换提取的谱余数位于整数位置，而真实的余数并非都是整数，如果用这些位置的谱余数作为频率重构的余数，必然会引起估计误差。本发明通过对相邻输入快拍(即第n个快拍和第n+1个快拍)的两根峰值IDFT谱线(即和)用相位差法进行谱校正，来获得小数部分的余数，从而可消除这部分误差。原理及过程如下。
由于信号的互素感知处理不涉及非线性器件，从整体来看，是一个线性时不变(LTI)***，因此，对于u路，输入信号x(n)与峰值ku处的IDFT输出之间的信号处理均是线性处理，可用下式表示
其中，是幅度响应，是相位响应。因为整个***是是不变的，输入信号x(n)幅度恒定，因此，输出的幅度响应不会随时间发生变化(即不会随变量n发生变化)。
对于输入为x(n+MN)的信号(即以相邻的第n+1个“快拍”作为激励，每一个“快拍”包含MN个样本)，将n＝n+MN代入式(1)得
将式(16)作为激励，则IDFT输出有如下形式
由线性***的齐次性，从式(16)可推出
再结合式(15)、式(17)和式(18)，可得出两个相邻“快拍”输出之间的相位差
其中，k*是由于相位的2π周期性造成的相位整周期模糊。为了移除这个整数模糊度，需要对式(19)取2π模，则
由于频率小数偏移δ∈[-0.5,0.5)，需要对得到的相位差进行相位调整，即
最后用2π归一化式(21)得到估计的频率偏移
对于v路，确定IDFT谱峰位置lv后，根据以上步骤，可估计出小数频偏则两个估计的余数可表示为
203：CRT重构频率。
利用闭式中国余数定理可对原始信号各成分的频率值进行重构，其中CRT所需的两个余数值由式(23)给出，两个模值选取为两个互素的下采样因子M和N(最大公约数为1)。为得到频率估计闭合解析形式的中国余数定理的具体重构过程如下：
1)由所给的余数中，计算其中：
其中[·]表示四舍五入取整。
2)计算模除N的余数
其中，是M模除N的模乘积的逆，“≡”表示模，可以提前算出：
3)计算式(4)中的折叠整数(即整数模糊度)
其中，b是1模除N的模乘积的逆。
b≡1mod N (27)
4)根据式(4)，计算折叠整数
最后，将估计出的两个余数以及两个折叠整数代入CRT重构的联立方程式(4)可得最终的频率估计：
为了减小误差，取平均值：
根据以上CRT重构步骤，用三个单频信号(对应的频率分别为f1＝57.3534Hz，f2＝43.0224Hz，f3＝15.5069Hz)进行实验，频率分别为可得到如表1所示的各成分的频率估计结果(包括折叠整数的估计结果和最终的频率估计结果)。
表1基于中国余数定理的频率估计结果
由表1可以看出，经过CRT从重构后，各频率成分估计结果与真实值保持一致，这证明了本发明估计方法的高精度性。
实施例3
(1)谱估计器短样本精度验证
实验1
分别用传统互素谱分析方法和本方法对单频信号x(t)＝exp(j2π107.3456t)进行频率估计，互素的下采样因子M＝17，N＝13，奈奎斯特采样率为Fs＝MN＝221Hz(对应的采样间隔T＝1/Fs＝1/221)，u路和v路的采样速率分别为fs1＝1/NT＝17Hz，fs2＝1/MT＝13Hz，均远低于奈奎斯特采样限制，频率分辨率Δf＝Fs/MN＝1Hz，滤波器长度L＝67。根据式(5)，初始化“快拍”数的临界值为6。用Q表示需要消耗的“快拍”数，每次试验中，选取第Q和第(Q-1)个“快拍”用于频谱校正。表2列出了Q＝2～8的各个情况下，使用本发明的频谱估计方法得到的各个参数值，图4为Q＝8时使用原始互素感知估计器得到的频谱图。
表2基于中国余数定理的频率估计结果
由图4可以看出，传统互素感知谱估计器因受频率分辨率Δf＝Fs/MN＝1Hz的限制，忽略了频率的小数偏移，故重构出的频率为107Hz，误差为0.3456Hz；
而由表2本发明提出的估计方法可以得出，由于加入了谱校正技术，本发明提出的谱估计方法成功将频率的小数偏移提取出来，加到频率余数上，使得最终被重构的频率保留了真实频率的小数偏移，Q＝3,4,5,6,7,8时误差值分别为0.083Hz，0.0055Hz，0.0065Hz，0.0014Hz，0Hz，0Hz，即使在“快拍”数未能达到临界初始化数量的情况下(3≤Q＜6)，依然有较高的精度，而一旦“快拍”数多于临界值，误差即为0。表明本发明提出的频率估计方法相对于原始互素感知谱估计方法在短样本情况下具有较高的估计精度，同时由于所需样本数较少，也在很大程度上降低了时间延迟。
(2)加噪情况下鲁棒性实验
实验2
选择与实验1相同的信号，使用0均值高斯白噪声构造不同信噪比(SNR)的激励，“快拍”数选取临界初始化数Q＝C＝6，每种信噪比情况取1000次实验结果的平均，使用本发明提出的频率估计方法得不同信噪比情况下频率估计的均方根误差(root-mean-squareerrors,RMSEs)，如图5所示。
由图5可以看出，在信噪比SNR≥3dB时，RMSEs小于0.1Hz，误差很小，SNR<3dB时，误差变大，频率估计发生大的失真。说明只要保证信噪比不是太低(本实验中只要SNR≥3dB)，式(29)中的折叠整数估计结果总是正确的，最终的频率估计的误差总能控制在0.02Hz以内，获得比频率分辨率Δf＝Fs/MN＝1Hz高出2个以上数量级的精度。
实施例4
一种基于互素感知的单频信号短样本高精度测频装置，参见图6，,单频信号短样本高精度测频装置包括：
数字信号处理器，数字信号处理器用于对由多通道稀疏信号组成的多通道输出序列进行离散傅里叶变换，得到离散傅里叶变换后的输出序列；对离散傅里叶变换后的输出序列作谱峰扫描，得到谱峰位置的索引，通过索引获取峰值离散傅里叶变换谱线，采用相位差法进行谱校正，获取校正后的频率；将校正后的频率对作为余数，按照闭合解析形式的中国余数定理进行处理，重构出原始信号频率；
输出驱动及其显示电路，所述输出驱动及其显示电路用于对原始信号频率进行输出显示。
其中，数字信号处理器包括：
第一计算模块，所述第一计算模块用于对输入信号x(n)进行两路下采样，下采样因子分别为M和N，得到两路稀疏信号；
第二计算模块，所述第二计算模块用于对两路稀疏信号分别进行多相滤波，输出两个多通道稀疏信号；在时刻n，对每个由多通道稀疏信号组成的多通道输出序列分别进行离散傅里叶变换，得到离散傅里叶变换后的输出序列。
其中，数字信号处理器包括：
第三计算模块，所述第三计算模块用于分别对离散傅里叶变换后的输出序列作谱峰扫描，分别得到谱峰位置的索引；
第四计算模块，用于分别将搜索到的谱峰位置索引代入到相邻快拍，分别获取两根峰值离散傅里叶变换谱线，用相位差法对两根峰值离散傅里叶变换谱线进行谱校正，分别得到校正后的频率。
其中，数字信号处理器为DSP。
实施例5
下面结合图6、图7，以及DSP具体的内部操作步骤对实施例4中的单频信号短样本高精度测频装置进行详细说明，详见下文描述：
在图6中，首先将实际采集的信号以及滤波器系数存入外部RAM中，再将信号以及滤波器系数实时输入到DSP中，经过DSP内部核心算法，对信号进行下采样、多相滤波、IDFT、谱校正以及构建CRT核心算法等，对信号频率进行重构，最后借助输出驱动显示及其显示模块显示出频率值。
其中，图6的DSP(Digital Signal Processor，数字信号处理器)为核心器件，在频率重构的过程中，完成如下主要功能：
调用内部核心算法，完成实际采集信号的下采样、多相滤波、频谱估计和频谱校正等过程；
控制下采样值以及信号样本，实时对其进行调整，使其符合实际需要；
将频率估计结果实时输出至驱动和显示模块。
需指出，由于采用了数字化的估计方法，因而决定图6***的复杂度、实时程度和稳定性的主要因素并不是图6中DSP器件的***连接，而是DSP内部程序存储器所存储的核心估计算法。DSP器件的内部程序流程如图7所示。
本发明将所提出的“交叉项干扰抑制的短样本稀疏谱高精度估计方法及其装置”的核心算法植入DSP器件内，基于此完成高精度、低复杂度、高效的频率估计。
图7流程分为如下几个步骤：
首先根据实际需要，设置信号的下采样值(M和N，互素的整数)，并设计所需的滤波器系数，用于多相滤波过程；
然后，CPU主控器从I/O端口读取参设定的参数，进入内部RAM；
按图1本发明的处理过程进行频率估计器的设计是DSP算法最核心的部分，运行该算法后，即可得到所采集的信号的频率估计；
判断本发明方法是否满足实际需求，若不满足，程序返回，重新根据要求设定信号参数；
直至设计结果符合实际要求，然后通过DSP的输出总线输出至外部显示驱动设备，将频率估计结果进行数码显示。
需指出，由于采用了DSP实现，使得整个频谱估计器设计变得更为灵活快捷，可根据频谱估计器设计过程中的实际需要，灵活变换所需参数，使之最终符合工程需要。
参考文献
[1]青海银,张援农,周晨,等.基于MST雷达垂直风速的大气温度剖面反演[J].物理学报,2014,63(9):94301-094301.
[2]Ding Z,Yao X S,Liu T,et al.Long-range vibration sensor based oncorrelation analysis of optical frequency-domain reflectometry signals[J].Optics express,2012,20(27):28319-28329.
[3]黄翔东,孟天伟,丁道贤,等.前后向子分段相位差频率估计法[J].物理学报,2014,63(21):214304-214304.
[4]丛超,李秀坤,宋扬.一种基于新型间歇混沌振子的舰船线谱检测方法[J].物理学报,2014,63(6):64301-064301.
[5]Wu Q,Liang Q.Coprime sampling for nonstationary signal in radarsignal processing[J].EURASIP Journal on Wireless Communications andNetworking,2013,2013(1):1-11.
[6]Zhang Y D,Qin S,Amin M G.DOA estimation exploiting coprime arrayswith sparse sensor spacing[J].Proc.IEEE ICASSP,Florence,Italy,2014:2267-2271.
[7]Wang W,Xia X G.A closed-form robust Chinese remainder theorem andits performance analysis[J].Signal Processing,IEEE Transactions on,2010,58(11):5655-5666.
[8]Li X,Liang H,Xia X G.A robust Chinese remainder theorem with itsapplications in frequency estimation from undersampled waveforms[J].SignalProcessing,IEEE Transactions on,2009,57(11):4314-4322.
[9]黄翔东,丁道贤,南楠,等.基于中国余数定理的欠采样下余弦信号的频率估计[J].物理学报,2014,63(19):198403-198403.
[10]Pal P,Vaidyanathan P P.Nested arrays:a novel approach to arrayprocessing with enhanced degrees of freedom[J].Signal Processing,IEEETransactions on,2010,58(8):4167-4181.
[11]Vaidyanathan P P,Pal P.Sparse sensing with co-pprime samplers andarrays[J].IEEE Transactions on Signal Processing,2011,59(2):573-586.
[12]Vaidyanathan P P,Pal P.Theory of sparse coprime sensing inmultiple dimensions[J].Signal Processing,IEEE Transactions on,2011,59(8):3592-3608.
[13]Vaidyanathan P P,Pal P.Sparse coprime sensing withmultidimensional lattice arrays[C]//Digital Signal Processing Workshop andIEEE Signal Processing Education Workshop(DSP/SPE),2011IEEE.IEEE,2011:425-430.
[14]Vaidyanathan P P,Pal P.Coprime Sampling and Arrays in One andMultiple Dimensions[M]//Multiscale Signal Analysis and Modeling.Springer NewYork,2013:105-137.
[15]Huang X,Xia X.A Fine Resolution Frequency Estimator Based onDouble Sub-segment Phase Difference[J].IEEE Signal Processing Letters,2015,22(8):1055-1059.
本发明实施例对各器件的型号除做特殊说明的以外，其他器件的型号不做限制，只要能完成上述功能的器件均可。
本领域技术人员可以理解附图只是一个优选实施例的示意图，上述本发明实施例序号仅仅为了描述，不代表实施例的优劣。
以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。
Claims (5)
1.一种基于互素感知的单频信号短样本高精度测频方法，其特征在于，所述单频信号短样本高精度测频方法包括以下步骤：
对由多通道稀疏信号组成的多通道输出序列进行离散傅里叶变换，得到离散傅里叶变换后的输出序列；
对离散傅里叶变换后的输出序列作谱峰扫描，得到谱峰位置的索引，通过索引获取峰值离散傅里叶变换谱线，采用相位差法进行谱校正，获取校正后的频率；
将校正后的频率对作为余数，按照闭合解析形式的中国余数定理进行处理，重构出原始信号频率；
其中，所述对离散傅里叶变换后的输出序列作谱峰扫描，得到谱峰位置的索引，通过索引获取峰值离散傅里叶变换谱线，采用相位差法进行谱校正，获取校正后的频率的步骤具体为：
分别对离散傅里叶变换后的输出序列作谱峰扫描，分别得到谱峰位置的索引；分别将搜索到的谱峰位置索引代入到相邻快拍，分别获取两根峰值离散傅里叶变换谱线，用相位差法对两根峰值离散傅里叶变换谱线进行谱校正，分别得到校正后的频率；
其中，所述用相位差法对两根峰值离散傅里叶变换谱线进行谱校正，分别得到校正后的频率具体为：
将估计出的两个余数以及两个折叠整数代入CRT重构的联立方程下式可得最终的频率估计：
其中，ku、lv分别为两路IDFT谱的峰值谱位置，为频率偏移，为小数频偏，均为折叠整数，M、N均为下采样因子，Δf为频率分辨率；u、v分别为稀疏信号的路号；
为了减小误差，取平均值：
2.根据权利要求1所述的一种基于互素感知的单频信号短样本高精度测频方法，其特征在于，所述对由多通道稀疏信号组成的多通道输出序列进行离散傅里叶变换，得到离散傅里叶变换后的输出序列的步骤具体为：
对输入信号x(n)进行两路下采样，下采样因子分别为M和N，得到两路稀疏信号，M和N为任意的互素的整数；
对两路稀疏信号分别进行多相滤波，输出两个多通道稀疏信号；
在时刻n，对每个由多通道稀疏信号组成的多通道输出序列分别进行离散傅里叶变换，得到离散傅里叶变换后的输出序列。
3.一种用于实施权利要求1-2中任一权利要求所述的基于互素感知的单频信号短样本高精度测频方法的测频装置，其特征在于，所述测频装置包括：
数字信号处理器，所述数字信号处理器用于对由多通道稀疏信号组成的多通道输出序列进行离散傅里叶变换，得到离散傅里叶变换后的输出序列；对离散傅里叶变换后的输出序列作谱峰扫描，得到谱峰位置的索引，通过索引获取峰值离散傅里叶变换谱线，采用相位差法进行谱校正，获取校正后的频率；将校正后的频率对作为余数，按照闭合解析形式的中国余数定理进行处理，重构出原始信号频率；
输出驱动及其显示电路，所述输出驱动及其显示电路用于对原始信号频率进行输出显示。
4.根据权利要求3所述的测频装置，其特征在于，所述数字信号处理器包括：
第一计算模块，所述第一计算模块用于对输入信号x(n)进行两路下采样，下采样因子分别为M和N，得到两路稀疏信号，M和N为任意的互素的整数；
第二计算模块，所述第二计算模块用于对两路稀疏信号分别进行多相滤波，输出两个多通道稀疏信号；在时刻n，对每个由多通道稀疏信号组成的多通道输出序列分别进行离散傅里叶变换，得到离散傅里叶变换后的输出序列。
5.根据权利要求3所述的测频装置，其特征在于，所述数字信号处理器包括：
第三计算模块，所述第三计算模块用于分别对离散傅里叶变换后的输出序列作谱峰扫描，分别得到谱峰位置的索引；
第四计算模块，用于分别将搜索到的谱峰位置索引代入到相邻快拍，分别获取两根峰值离散傅里叶变换谱线，用相位差法对两根峰值离散傅里叶变换谱线进行谱校正，分别得到校正后的频率。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201510346171.4A CN104897962B (zh) | 2015-06-19 | 2015-06-19 | 基于互素感知的单频信号短样本高精度测频方法及其装置 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201510346171.4A CN104897962B (zh) | 2015-06-19 | 2015-06-19 | 基于互素感知的单频信号短样本高精度测频方法及其装置 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN104897962A CN104897962A (zh) | 2015-09-09 |
| CN104897962B true CN104897962B (zh) | 2019-02-22 |
Family
ID=54030745
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201510346171.4A Expired - Fee Related CN104897962B (zh) | 2015-06-19 | 2015-06-19 | 基于互素感知的单频信号短样本高精度测频方法及其装置 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN104897962B (zh) |
Families Citing this family (13)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105259410B (zh) * | 2015-10-26 | 2017-12-05 | 天津大学 | 一种强噪声干扰下的欠采样波形的频率估计方法及其装置 |
| CN105446940B (zh) * | 2015-11-23 | 2018-07-03 | 西安交通大学 | 一种基于同步压缩变换重构的幅值校正方法 |
| CN106027179A (zh) * | 2016-05-12 | 2016-10-12 | 天津大学 | 一种基于综合互素分析的宽带频谱感知方法及其装置 |
| CN106301631B (zh) * | 2016-10-18 | 2019-03-05 | 天津大学 | 一种基于子空间分解的互素欠采样频谱感知方法及其装置 |
| CN106506102B (zh) * | 2016-10-18 | 2019-02-15 | 天津大学 | 一种互素欠采样下高精度、低时延的谱感知方法及其装置 |
| EP3462146A1 (en) * | 2017-07-29 | 2019-04-03 | Tata Consultancy Services Limited | Systems and methods for autonomous stroboscopic machine inspection for multi-point and multi-frequency vibration measurement |
| CN108037481B (zh) * | 2017-12-01 | 2022-02-08 | 天津大学 | 鲁棒性可分级的稀疏阵列频率和doa估计方法及装置 |
| CN109672489A (zh) * | 2018-11-26 | 2019-04-23 | 天津大学青岛海洋技术研究院 | 一种基于互素分析的海上全景频谱感知分析方法及其装置 |
| CN109917360A (zh) * | 2019-03-01 | 2019-06-21 | 吉林大学 | 一种混叠脉冲的参差pri估计方法 |
| CN110212889B (zh) * | 2019-05-29 | 2020-11-13 | 北京机电工程研究所 | 一种数字信号采样装置及方法 |
| CN110896308B (zh) * | 2019-10-31 | 2023-09-12 | 中国工程物理研究院电子工程研究所 | 一种单音信号重构方法 |
| CN112098850B (zh) * | 2020-09-21 | 2024-03-08 | 山东工商学院 | 基于sdo算法的锂离子电池电压故障诊断方法及*** |
| CN112485518B (zh) * | 2020-11-09 | 2022-07-15 | 浙江工业大学 | 一种双通道多指数衰减正弦信号欠采样方法 |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101825660A (zh) * | 2010-05-05 | 2010-09-08 | 天津大学 | 欠采样下的正弦信号频率的高效测量方法及实施装置 |
| CN104007316A (zh) * | 2014-05-29 | 2014-08-27 | 天津大学 | 一种欠采样速率下的高精度频率测量方法及其测量仪 |
| US9037429B2 (en) * | 2011-06-06 | 2015-05-19 | Siemens Industry, Inc. | Methods and apparatus for measuring the fundamental frequency of a line signal |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US8660218B2 (en) * | 2012-01-30 | 2014-02-25 | L-3 Communications Integrated Systems Lp | Systems and methods for signal detection and estimation |
-
2015
- 2015-06-19 CN CN201510346171.4A patent/CN104897962B/zh not_active Expired - Fee Related
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101825660A (zh) * | 2010-05-05 | 2010-09-08 | 天津大学 | 欠采样下的正弦信号频率的高效测量方法及实施装置 |
| US9037429B2 (en) * | 2011-06-06 | 2015-05-19 | Siemens Industry, Inc. | Methods and apparatus for measuring the fundamental frequency of a line signal |
| CN104007316A (zh) * | 2014-05-29 | 2014-08-27 | 天津大学 | 一种欠采样速率下的高精度频率测量方法及其测量仪 |
Non-Patent Citations (4)
| Title |
|---|
| Sparse coprime sensing with multidimensional lattice arrays;Vaidyanathan P P等;《Digital Signal Processing Workshop and IEEE Signal Processing Education Workshop.IEEE》;20111031;全文 |
| 基于中国余数定理的欠采样下余弦信号的频率估计;黄翔东 等;《物理学报》;20141008;第63卷(第19期);第383-390页 |
| 基于稀疏采样的互素谱分析方法;闫子阳 等;《第八届全国信号和智能信息处理与应用学术会议会刊》;20140601;第20-23页 |
| 欠采样稀疏频率估计方法及研究;李莉;《中国优秀硕士学位论文全文数据库 信息科技辑》;20110815;全文 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN104897962A (zh) | 2015-09-09 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN104897962B (zh) | 基于互素感知的单频信号短样本高精度测频方法及其装置 | |
| CN104375111B (zh) | 对密集频谱进行快速高精度细化校正的方法 | |
| Eichstädt et al. | Deconvolution filters for the analysis of dynamic measurement processes: a tutorial | |
| CN104007316B (zh) | 一种欠采样速率下的高精度频率测量方法及其测量仪 | |
| JP6173721B2 (ja) | 周波数解析装置、当該周波数解析装置を用いた信号処理装置、および、当該信号処理装置を用いた高周波測定装置 | |
| CN104991119B (zh) | 一种消除伪峰、谱泄漏效应的互素谱分析方法及其装置 | |
| CN106645780B (zh) | 一种基于dsp的转速检测方法及*** | |
| CN109407501B (zh) | 一种基于相关信号处理的时间间隔测量方法 | |
| CN103983849A (zh) | 一种实时高精度的电力谐波分析方法 | |
| Sudani et al. | A novel robust and accurate spectral testing method for non-coherent sampling | |
| CN109117816A (zh) | 基于六阶样条插值小波的信号奇异点检测方法 | |
| CN103543331B (zh) | 一种计算电信号谐波和间谐波的方法 | |
| CN103969508A (zh) | 一种实时高精密的电力谐波分析方法及装置 | |
| CN103605904B (zh) | 基于误差估算的自补偿电力***幅值算法 | |
| CN105759255B (zh) | 一种cic多相内插滤波超声相控阵波束延时方法 | |
| CN103592513B (zh) | 电力信号谐波分析方法和装置 | |
| CN106324342A (zh) | 一种基于查表的谐波检测方法 | |
| CN110297199B (zh) | 一种基于全相位fft的铯光泵磁力仪频率测量方法及*** | |
| CN105282086A (zh) | 变频***群时延的双音正交测量方法 | |
| JP2000180484A (ja) | 高調波測定装置 | |
| JP5553334B2 (ja) | 正弦波パラメータ推定方法 | |
| CN108957118A (zh) | 一种无功功率计算方法 | |
| CN110285881A (zh) | 一种基于全相位滤波的密集谱频率估计法 | |
| CN109490625B (zh) | 一种基于滑动窗口和半定规划的谐波信号分析方法 | |
| Slosarcik | FFT-based algorithm for metering applications |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant | ||
| CF01 | Termination of patent right due to non-payment of annual fee | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20190222 |