CN112595889A - 非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法 - Google Patents

非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法 Download PDF

Info

Publication number
CN112595889A
CN112595889A CN202011238816.XA CN202011238816A CN112595889A CN 112595889 A CN112595889 A CN 112595889A CN 202011238816 A CN202011238816 A CN 202011238816A CN 112595889 A CN112595889 A CN 112595889A
Authority
CN
China
Prior art keywords
sampling
signal
frequency
channel
damping factor
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.)
Granted
Application number
CN202011238816.XA
Other languages
English (en)
Other versions
CN112595889B (zh
Inventor
黄国兴
倪安
卢为党
彭宏
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhejiang University of Technology ZJUT
Original Assignee
Zhejiang University of Technology ZJUT
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Zhejiang University of Technology ZJUT filed Critical Zhejiang University of Technology ZJUT
Priority to CN202011238816.XA priority Critical patent/CN112595889B/zh
Publication of CN112595889A publication Critical patent/CN112595889A/zh
Application granted granted Critical
Publication of CN112595889B publication Critical patent/CN112595889B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
    • G01R23/00—Arrangements for measuring frequencies; Arrangements for analysing frequency spectra
    • G01R23/02—Arrangements for measuring frequency, e.g. pulse repetition rate; Arrangements for measuring period of current or voltage
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Data Mining & Analysis (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Software Systems (AREA)
  • Databases & Information Systems (AREA)
  • Algebra (AREA)
  • Operations Research (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • Computing Systems (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Evolutionary Biology (AREA)
  • Probability & Statistics with Applications (AREA)
  • Measuring Frequencies, Analyzing Spectra (AREA)
  • Complex Calculations (AREA)

Abstract

一种非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法，首先将生成的MEDS信号分流至两条并行通道中，在第一条通道中对信号进行采样，并运行振幅和分量个数最优估计算法，得到相应的估计值，与此同时，在第二条通道中也对信号进行采样，但是该通道中的采样时钟的上升沿与第一条通道中的采样时钟的上升沿相差，Te≤1/fmax，fmax是信号的最大频率，然后运行复频率联合估计算法，得到相应的估计值，从而完成MEDS信号的参数测量。本发明能够适用于频率分量个数未知和模型不匹配的非理想信号环境，且信号重构精度优于现有方法。

Description

非理想多指数衰减正弦信号的欠Nyquist采样与参数测量 方法
技术领域
本发明属于信号处理技术领域，具体涉及一种非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法。
背景技术
在语音分析、核磁共振成像等诸多应用中，多指数衰减正弦信号(MultipleExponentially Damped Sinusoids，MEDS)的参数测量是一个普遍存在的问题。因此，有不少研究者为了能够解决该问题，提出了很多用于解决MEDS信号的参数测量方法。Ruiz等学者提出了一种估计方法，该方法能够减少不完全采样数据产生的确定性误差。Chen等学者提出了一种基于子空间的参数估计方法。然而，上述方法都是基于著名的Nyquist采样理论，该理论要求采样率是信号带宽的两倍。因此，对于宽频带信号，这些方法需要非常高的采样率，导致采样量大，参数测量处理繁重。
考虑具有如下形式的连续时间MEDS信号x(t)
Figure BDA0002767677880000011
其中，ck,k＝1,2,…K是复振幅，并且分量个数K未知，sk＝rk+j2πfk是复频率，rk是阻尼因子，fk是极点频率，σ(t)是模型匹配误差信号。由于信号发生器不能生成理想的MEDS信号，我们设置模型匹配误差信号σ(t)≠0。现对MEDS信号x(t)进行采样，假设采样速率为fs，则采样值可表示为：
Figure BDA0002767677880000012
其中，
Figure BDA0002767677880000021
n∈Z+。要想利用2K个x[n]的样本对参数ck和uk进行估算，我们可以采用频谱估计方法。假设
Figure BDA0002767677880000022
为估计值，由于三角函数的周期性，欠采样会导致频谱混叠现象。
近年来，不少学者对采用欠Nyquist采样法进行MEDS信号参数测量十分关注。但是，使用欠Nyquist采样法有一个弊端，即该方法会导致频谱混叠现象。因此，如何解决频谱混叠现象是一个关键问题。R.Venkataramani等学者提出了一种多速度欠Nyquist采样方案来解决这种频谱混叠问题，但是，该方法需要对频率分量的数量施加一定的限制。R.G.Baraniuk等学者提出了基于压缩感知(Compressed Sensing，CS)理论的一些方法，比如随机解调器和随机采样等信息转换技术也是一种潜在的频率估计方法，但是，硬件难以实现CS所要求的的对信号进行随机采样。S.Qin等人提出了一种基于互质速率的双通道欠奈奎斯特采样和参数测量方案，然而，S.Huang等学者证明了在某些条件下，S.Qin的方法不能对频率进行唯一估计。S.Huang等学者提出了一种基于成对互质速率的三通道欠Nyquist采样和参数测量方案，该方法能够成功地解决频谱混叠问题，然而，该方法需要大量的样本才能获得唯一解。G.Huang等学者针对多个正弦信号提出了两种多通道采样方案，该方案同样能够解决频谱混叠问题。但是，这些系统需要知道频率分量的个数，而且并不是专门为MEDS信号设计的。然而，很多现有的方法都是基于Nyquist采样理论的，需要大量的采样值。因此，如何在通过较少的采样值测量MEDS信号的参数的同时解决欠Nyquist采样带来的频谱混叠问题仍然是一个关键难题。
发明内容
为了克服已有技术的不足，针对非理想MEDS信号的欠Nyquist采样与参数测量方法，本发明提出一种非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法，由两个欠Nyquist采样通道组成，这两个通道各有一个采样时钟，分别为CLK1和CLK2，并且两个通道具有相同的采样速率，CLK1和CLK2的上升沿之间的时间延迟Te小于或等于Nyquist采样间隔，即Te≤1/fmax，这里fmax表示信号的最大频率，从而形成了采样时钟交错；本发明能够适用于频率分量个数未知和模型不匹配的非理想信号环境，且信号重构精度优于现有方法。
本发明解决其技术问题所采用的技术方案是：
一种非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法，包括以下步骤：
步骤一，初始化生成MEDS信号，信号形式如下：
Figure BDA0002767677880000031
其中，ck,k＝1,2,…K是复振幅，并且分量个数K未知，sk＝rk+j2πfk是复频率，rk是阻尼因子，fk是极点频率，σ(t)是模型匹配误差信号，然后将该信号分流到2条并行的采样通道中，分流所到达的两个通道各有一个采样时钟CLK1和CLK2，这两个通道的采样速率相同，CLK1和CLK2的上升沿之间的时间延迟Te小于或等于Nyquist采样间隔，即Te≤1/fmax，fmax表示信号的最大频率，从而形成了采样时钟交错；
步骤二，第一条通道采样：在第一条采样通道中，信号x(t)根据时钟CLK1被均匀采样，并且采样由下式给出：
Figure BDA0002767677880000032
其中，K为未知的分量个数，ck表示复振幅，
Figure BDA0002767677880000033
sk＝rk+j2πfk表示复频率，rk表示阻尼因子，fk表示极点频率，m∈Z+，σ[n]表示σ(t)被时钟CLK1采样后的样本值；
步骤三，运行振幅和分量个数最优估计算法，由于分量个数K未知且模型匹配误差σ(t)≠0，通过最小化σ(t)的能量来找到最优的分量K和最优的参数
Figure BDA0002767677880000041
步骤四，第二条通道采样：在第二条采样通道中，待测MEDS信号
Figure BDA0002767677880000042
由时钟CLK2均匀采样，并且采样公式如下：
Figure BDA0002767677880000043
其中，
Figure BDA0002767677880000044
为复振幅ck的估计值，sk＝rk+j2πfk为复频率，rk为阻尼因子，fk为极点频率，Ts为采样周期，Te为CLK1和CLK2的上升沿之间的时间延迟，n′∈Z+；
步骤五，运行复频率的联合估计算法：输入为参数
Figure BDA0002767677880000045
采样速率fs、交错时间Te以及N′≥K个连续采样值xe[n′]，输出为阻尼因子rk和极点频率fk的估计值
Figure BDA0002767677880000046
进一步，所述步骤三中，算法过程的数学模型描述如下：
3.1，建立优化模型，优化模型表示如下：
Figure BDA0002767677880000047
由于
Figure BDA0002767677880000048
公式(3)转换为：
Figure BDA0002767677880000049
3.2，优化模型离散化，令t＝n/fs，
Figure BDA00027676778800000410
公式(4)离散为如下形式：
Figure BDA00027676778800000411
其中，
Figure BDA0002767677880000051
是指优化对象函数；
3.3，用频谱估计法求解
Figure BDA0002767677880000052
对
Figure BDA0002767677880000053
进行初始化，通过求解公式三将P0优化至Popt，然后通过最小化公式(6)所表示的优化对象函数就得到最优的分量个数K和最优参数
Figure BDA0002767677880000054
再进一步，所述步骤五中，算法过程如下：
5.1，估计阻尼因子rk，阻尼因子rk可以很容易地从
Figure BDA0002767677880000055
中估算出来，估算过程如下：由
Figure BDA0002767677880000056
与sk＝rk+j2πfk得
Figure BDA0002767677880000057
则得阻尼因子rk的估计值如下式：
Figure BDA0002767677880000058
同时得到极点频率fk的估计值如下式：
Figure BDA0002767677880000059
其中，∠(·)表示(·)的幅角主值，且0≤∠(·)＜2π；
5.2，由公式(7)和公式(8)得
Figure BDA00027676778800000510
其中，
Figure BDA00027676778800000511
是已知的，
Figure BDA00027676778800000512
是未知的；
5.3，求解bk的估计值，利用N′≥K个连续采样值xe[n′](1≤n′≤K)构建一个N′×1的向量，该向量形式如下：
xe＝[xe[0],xe[1],…,xe[N′-1]]T (12)
利用
Figure BDA0002767677880000061
的幅角主值创建一个K×N′的范德蒙矩阵，该矩阵形式如下：
Figure BDA0002767677880000062
利用
Figure BDA0002767677880000063
和
Figure BDA0002767677880000064
的幅角主值创建一个N′×K的对角矩阵
Figure BDA0002767677880000065
通过求解
Figure BDA0002767677880000066
得到系数矩阵
Figure BDA0002767677880000067
5.4，求解极点频率fk的估计值，由于三角函数的周期性，估计的归一化频率
Figure BDA0002767677880000068
与
Figure BDA0002767677880000069
的参数相差2πmk，即：
Figure BDA00027676778800000610
其中，∠(·)是(·)的幅角主值，0≤∠(·)＜2π，mk∈Z,k∈{1,2,…,K}，因为
Figure BDA00027676778800000611
fmax表示信号的最大频率，得到如下不等式：
Figure BDA00027676778800000612
进一步化简公式(16)得如下不等式：
Figure BDA00027676778800000613
记
Figure BDA00027676778800000614
公式(17)转换为：
0≤mk≤(Q-1) (18)
因为
Figure BDA00027676778800000615
Q∈Z+，则有
Figure BDA00027676778800000616
如果我设置主、次通道采样时钟上升沿相差的时间间隔
Figure BDA0002767677880000071
则有：
Figure BDA0002767677880000072
与公式(18)相结合，则有：
0≤mkfsTe＜1 (20)
所以对于
Figure BDA0002767677880000073
mk有唯一解，其解如下式：
Figure BDA0002767677880000074
其中，
Figure BDA0002767677880000075
极点频率fk的估计值由
Figure BDA0002767677880000076
唯一确定。
本发明首先将生成的MEDS信号分流至两条并行通道中，在第一条通道中对信号进行采样，并运行振幅和分量个数最优估计算法，得到相应的估计值，与此同时，在第二条通道中也对信号进行采样，但是该通道中的采样时钟的上升沿与第一条通道中的采样时钟的上升沿相差Te(Te≤1/fmax，fmax是信号的最大频率)，然后运行复频率联合估计算法，得到相应的估计值，从而完成MEDS信号的参数测量。
本发明的有益效果主要表现在：能够适用于频率分量个数未知和模型不匹配的非理想信号环境，且信号重构精度优于现有方法。
附图说明
图1是基于欠Nyquist采样的MEDS信号的参数测量系统的原理图。
图2是本发明方法在噪声情况下的性能与其他系统的比较结果，其中，(a)不同系统的均方误差(b)不同系统的平均方差。
具体实施方式
下面结合附图对本发明作进一步描述。
参照图1，一种非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法，包括以下步骤：
步骤一，初始化生成MEDS信号，信号形式如下：
Figure BDA0002767677880000081
其中，ck,k＝1,2,…K是复振幅，并且分量个数K未知，sk＝rk+j2πfk是复频率，rk是阻尼因子，fk是极点频率，σ(t)是模型匹配误差信号，然后将该信号分流到2条并行的采样通道中，分流所到达的两个通道各有一个采样时钟CLK1和CLK2，这两个通道的采样速率相同，CLK1和CLK2的上升沿之间的时间延迟Te小于或等于Nyquist采样间隔，即Te≤1/fmax，fmax表示信号的最大频率，从而形成了采样时钟交错；
步骤二，第一条通道采样：在第一条采样通道中，信号x(t)根据时钟CLK1被均匀采样，并且采样由下式给出：
Figure BDA0002767677880000082
其中，K为未知的分量个数，ck表示复振幅，
Figure BDA0002767677880000083
sk＝rk+j2πfk表示复频率，rk表示阻尼因子，fk表示极点频率，m∈Z+，σ[n]表示σ(t)被时钟CLK1采样后的样本值；
步骤三，运行振幅和分量个数最优估计算法，由于分量个数K未知且模型匹配误差σ(t)≠0，通过最小化σ(t)的能量来找到最优的分量K和最优的参数
Figure BDA0002767677880000084
算法过程的数学模型描述如下：
3.1，建立优化模型，优化模型表示如下：
Figure BDA0002767677880000085
由于
Figure BDA0002767677880000086
公式(3)转换为：
Figure BDA0002767677880000087
3.2，优化模型离散化，令t＝n/fs，
Figure BDA0002767677880000091
公式(4)离散为如下形式：
Figure BDA0002767677880000092
其中，
Figure BDA0002767677880000093
是指优化对象函数；
3.3，用频谱估计法求解
Figure BDA0002767677880000094
对
Figure BDA0002767677880000095
进行初始化，通过求解公式三将P0优化至Popt，然后通过最小化公式(6)所表示的优化对象函数就得到最优的分量个数K和最优参数
Figure BDA0002767677880000096
步骤四，第二条通道采样：在第二条采样通道中，待测MEDS信号
Figure BDA0002767677880000097
由时钟CLK2均匀采样，并且采样公式如下：
Figure BDA0002767677880000098
其中，
Figure BDA0002767677880000099
为复振幅ck的估计值，sk＝rk+j2πfk为复频率，rk为阻尼因子，fk为极点频率，Ts为采样周期，Te为CLK1和CLK2的上升沿之间的时间延迟，n′∈Z+；
步骤五，运行复频率的联合估计算法：输入为参数
Figure BDA00027676778800000910
采样速率fs、交错时间Te以及N′≥K个连续采样值xe[n′]，输出为阻尼因子rk和极点频率fk的估计值
Figure BDA00027676778800000911
其算法过程如下：
5.1，估计阻尼因子rk，阻尼因子rk可以很容易地从
Figure BDA00027676778800000912
中估算出来，估算过程如下：由
Figure BDA00027676778800000913
与sk＝rk+j2πfk得
Figure BDA00027676778800000914
则得阻尼因子rk的估计值如下式：
Figure BDA0002767677880000101
同时得到极点频率fk的估计值如下式：
Figure BDA0002767677880000102
其中，∠(·)表示(·)的幅角主值，且0≤∠(·)＜2π；
5.2，由公式(7)和公式(8)得
Figure BDA0002767677880000103
其中，
Figure BDA0002767677880000104
是已知的，
Figure BDA0002767677880000105
是未知的；
5.3，求解bk的估计值，利用N′≥K个连续采样值xe[n′](1≤n′≤K)构建一个N′×1的向量，该向量形式如下：
xe＝[xe[0],xe[1],…,xe[N′-1]]T (12)
利用
Figure BDA0002767677880000106
的幅角主值创建一个K×N′的范德蒙矩阵，该矩阵形式如下：
Figure BDA0002767677880000107
利用
Figure BDA0002767677880000108
和
Figure BDA0002767677880000109
的幅角主值创建一个N′×K的对角矩阵
Figure BDA00027676778800001010
通过求解
Figure BDA00027676778800001011
得到系数矩阵
Figure BDA00027676778800001012
5.4，求解极点频率fk的估计值，由于三角函数的周期性，估计的归一化频率
Figure BDA00027676778800001013
与
Figure BDA00027676778800001014
的参数相差2πmk，即：
Figure BDA0002767677880000111
其中，∠(·)是(·)的幅角主值，0≤∠(·)＜2π，mk∈Z,k∈{1,2,…,K}，因为
Figure BDA0002767677880000112
fmax表示信号的最大频率，得到如下不等式：
Figure BDA0002767677880000113
进一步化简公式(16)得如下不等式：
Figure BDA0002767677880000114
记
Figure BDA0002767677880000115
公式(17)转换为：
0≤mk≤(Q-1) (18)
因为
Figure BDA0002767677880000116
Q∈Z+，则有
Figure BDA0002767677880000117
如果我设置主、次通道采样时钟上升沿相差的时间间隔
Figure BDA0002767677880000118
则有：
Figure BDA0002767677880000119
与公式(18)相结合，则有：
0≤mkfsTe＜1 (20)
所以对于
Figure BDA00027676778800001110
mk有唯一解，其解如下式：
Figure BDA00027676778800001111
其中，
Figure BDA00027676778800001112
极点频率fk的估计值由
Figure BDA00027676778800001113
唯一确定。
实验对比：为了验证本发明方法的优势，我们将本发明方法与S.Qin提出的双通道延时系统和S.Huang提出的三通道采样系统进行对比，其比较结果如表1所示。
Figure BDA00027676778800001114
Figure BDA0002767677880000121
表1
从表1可以看出，本发明方法仅以3K个样本就能解决频谱混叠问题，这比其他两个采样系统要少得多。此外，当频率分量个数K未知时，本发明方法是适用的，而其他两个系统则不适用。
本发明为了检测系统的有效性，还进行了几个相关的实验。对于公式一所定义的MEDS信号
Figure BDA0002767677880000122
本发明仍然能够成功地完成参数测量。MEDS信号x(t)含有K个频率分量，其最大频率为fmax＝10KHz，复振幅ck随机选取范围为(0,1+i)，阻尼因子rk在(0,1)上随机分布，并且极点频率fk在(0,fmax)上分布。为了方便比较，我们尽可能的设置相似的系统参数：双通道延时系统和本发明方法的信道采样率设置为fs＝fmax/4＝2.5KHz；三通道采样系统的三条信道采样率分别设置为fs＝2,2.5,3.3KHz。
首先，在无噪声的情况下，我们用少量样本来验证本发明所用采样系统的有效性。设置待测MEDS信号x(t)的分量个数K＝5。假设K是已知的，并且从双通道延时系统中选取4K＝20个采样值以及从三通道采样系统中选取188个采样值恢复信号。在本发明方法中，只要3K＝15个采样值(并且K是未知的)就能恢复信号。无噪声情况下的恢复结果对比如表2所示。
Figure BDA0002767677880000123
Figure BDA0002767677880000131
表2
从表2可以看出，尽管Qin和Huang的系统能够估计所有的极点频率，但是不能估计复振幅和阻尼因子。本发明方法可以准确恢复待测MEDS信号的各个参数。更值得注意的是，本发明方法所需的采样值个数比Qin和Huang的系统更少。
然后，我们研究了本发明方法在噪声情况下的性能，并与其他系统进行了比较。比较结果如图2所示。
设置待测MEDS信号x(t)的分量个数K＝3，并且在(0,10)KHz的区间上随机均匀绘制频谱图。用高斯白噪声对信号进行干扰。这些系统设置的目的是获取相同数量的采样值，采样值个数为Num＝188。实验测试进行了1000次，输入SNR电平从0dB增加到50dB。在图2中，我们针对输入SNR绘制了均方误差(MSE)和方差(VAR)曲线图。从图2可以看出，可以从以下三种情况分析估计性能：
1)低SNR(SNR<5dB)。对于低SNR，样本数量在参数估计中起主要作用。由于三个系统的样本数是固定的，因此估计性能保持相似。
2)提高SNR(5dB≤SNR≤25dB)。对于提高SNR，恢复算法在参数估计中更为重要。由于在恢复过程中我们采用了一些优化和去噪技术，显然，本发明方法的MSE和VAR曲线的性能要优于其他系统。
3)高SNR(SNR>25dB)。对于高SNR，使用足够多的采样值数据可以获得稳定的参数估计性能。尽管MSE继续下降，但对于高SNR值，VAR保持恒定。
该实验表明，在噪声条件下，本发明方法的性能优于双通道延时系统和三通道采样系统。

Claims (3)

1.一种非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法，其特征在于，所述方法包括以下步骤：
步骤一，初始化生成MEDS信号，信号形式如下：
Figure FDA0002767677870000011
其中，ck,k＝1,2,…K是复振幅，并且分量个数K未知，sk＝rk+j2πfk是复频率，rk是阻尼因子，fk是极点频率，σ(t)是模型匹配误差信号，然后将该信号分流到2条并行的采样通道中，分流所到达的两个通道各有一个采样时钟CLK1和CLK2，这两个通道的采样速率相同，CLK1和CLK2的上升沿之间的时间延迟Te小于或等于Nyquist采样间隔，即Te≤1/fmax，fmax表示信号的最大频率，从而形成了采样时钟交错；
步骤二，第一条通道采样：在第一条采样通道中，信号x(t)根据时钟CLK1被均匀采样，并且采样由下式给出：
Figure FDA0002767677870000012
其中，K为未知的分量个数，ck表示复振幅，
Figure FDA0002767677870000013
sk＝rk+j2πfk表示复频率，rk表示阻尼因子，fk表示极点频率，m∈Z+，σ[n]表示σ(t)被时钟CLK1采样后的样本值；
步骤三，运行振幅和分量个数最优估计算法，由于分量个数K未知且模型匹配误差σ(t)≠0，通过最小化σ(t)的能量来找到最优的分量K和最优的参数
Figure FDA0002767677870000014
步骤四，第二条通道采样：在第二条采样通道中，待测MEDS信号
Figure FDA0002767677870000015
由时钟CLK2均匀采样，并且采样公式如下：
Figure FDA0002767677870000016
其中，
Figure FDA0002767677870000017
为复振幅ck的估计值，sk＝rk+j2πfk为复频率，rk为阻尼因子，fk为极点频率，Ts为采样周期，Te为CLK1和CLK2的上升沿之间的时间延迟，n′∈Z+；
步骤五，运行复频率的联合估计算法：输入为参数
Figure FDA0002767677870000018
采样速率fs、交错时间Te以及N′≥K个连续采样值xe[n′]，输出为阻尼因子rk和极点频率fk的估计值
Figure FDA0002767677870000019
2.如权利要求1所述的非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法，其特征在于，所述步骤三中，算法过程的数学模型描述如下：
3.1，建立优化模型，优化模型表示如下：
Figure FDA0002767677870000021
由于
Figure FDA0002767677870000022
公式(3)转换为：
Figure FDA0002767677870000023
3.2，优化模型离散化，令t＝n/fs，
Figure FDA0002767677870000024
公式(4)离散为如下形式：
Figure FDA0002767677870000025
其中，
Figure FDA0002767677870000026
是指优化对象函数；
3.3，用频谱估计法求解
Figure FDA0002767677870000027
对
Figure FDA0002767677870000028
进行初始化，通过求解公式三将P0优化至Popt，然后通过最小化公式(6)所表示的优化对象函数就得到最优的分量个数K和最优参数
Figure FDA0002767677870000029
3.如权利要求1或2所述的非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法，其特征在于，所述步骤五中，算法过程如下：
5.1，估计阻尼因子rk，阻尼因子rk可以很容易地从
Figure FDA00027676778700000210
中估算出来，估算过程如下：由
Figure FDA00027676778700000211
与sk＝rk+j2πfk得
Figure FDA00027676778700000212
则得阻尼因子rk的估计值如下式：
Figure FDA00027676778700000213
同时得到极点频率fk的估计值如下式：
Figure FDA00027676778700000214
其中，∠(·)表示(·)的幅角主值，且0≤∠(·)＜2π；
5.2，由公式(7)和公式(8)得
Figure FDA00027676778700000215
其中，
Figure FDA0002767677870000031
是已知的，
Figure FDA0002767677870000032
是未知的；
5.3，求解bk的估计值，利用N′≥K个连续采样值xe[n′](1≤n′≤K)构建一个N′×1的向量，该向量形式如下：
xe＝[xe[0],xe[1],…,xe[N′-1]]T (12)
利用
Figure FDA0002767677870000033
的幅角主值创建一个K×N′的范德蒙矩阵，该矩阵形式如下：
Figure FDA0002767677870000034
利用
Figure FDA0002767677870000035
和
Figure FDA0002767677870000036
的幅角主值创建一个N′×K的对角矩阵
Figure FDA0002767677870000037
通过求解
Figure FDA0002767677870000038
得到系数矩阵
Figure FDA0002767677870000039
5.4，求解极点频率fk的估计值，由于三角函数的周期性，估计的归一化频率
Figure FDA00027676778700000310
与
Figure FDA00027676778700000311
的参数相差2πmk，即：
Figure FDA00027676778700000312
其中，∠(·)是(·)的幅角主值，0≤∠(·)＜2π，mk∈Z,k∈{1,2,…,K}，因为
Figure FDA00027676778700000313
fmax表示信号的最大频率，得到如下不等式：
Figure FDA00027676778700000314
进一步化简公式(16)得如下不等式：
Figure FDA00027676778700000315
记
Figure FDA00027676778700000316
公式(17)转换为：
0≤mk≤(Q-1) (18)
因为
Figure FDA00027676778700000317
Q∈Z+，则有
Figure FDA00027676778700000318
如果我设置主、次通道采样时钟上升沿相差的时间间隔
Figure FDA00027676778700000319
则有：
Figure FDA00027676778700000320
与公式(18)相结合，则有：
0≤mkfsTe＜1 (20)
所以对于
Figure FDA0002767677870000041
mk有唯一解，其解如下式：
Figure FDA0002767677870000042
其中，
Figure FDA0002767677870000043
极点频率fk的估计值由
Figure FDA0002767677870000044
唯一确定。
CN202011238816.XA 2020-11-09 2020-11-09 非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法 Active CN112595889B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011238816.XA CN112595889B (zh) 2020-11-09 2020-11-09 非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011238816.XA CN112595889B (zh) 2020-11-09 2020-11-09 非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法

Publications (2)

Publication Number Publication Date
CN112595889A true CN112595889A (zh) 2021-04-02
CN112595889B CN112595889B (zh) 2022-07-15

Family

ID=75183411

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011238816.XA Active CN112595889B (zh) 2020-11-09 2020-11-09 非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法

Country Status (1)

Country Link
CN (1) CN112595889B (zh)

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101251556A (zh) * 2008-03-04 2008-08-27 北京航空航天大学 正弦信号四参数检测方法和虚拟仪器信号检测装置
KR20150047109A (ko) * 2013-10-23 2015-05-04 삼성전자주식회사 자기 공명 영상 장치 및 그에 따른 자기 공명 영상의 이미징 방법
EP3301461A1 (en) * 2016-09-28 2018-04-04 Siemens Aktiengesellschaft Method for detection of harmonics of a univariate signal
CN108132383A (zh) * 2017-12-26 2018-06-08 哈尔滨工业大学 一种时间交错的多谐波信号欠采样方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101251556A (zh) * 2008-03-04 2008-08-27 北京航空航天大学 正弦信号四参数检测方法和虚拟仪器信号检测装置
KR20150047109A (ko) * 2013-10-23 2015-05-04 삼성전자주식회사 자기 공명 영상 장치 및 그에 따른 자기 공명 영상의 이미징 방법
EP3301461A1 (en) * 2016-09-28 2018-04-04 Siemens Aktiengesellschaft Method for detection of harmonics of a univariate signal
CN108132383A (zh) * 2017-12-26 2018-06-08 哈尔滨工业大学 一种时间交错的多谐波信号欠采样方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
董宁斐等: "基于多路欠采样的多分量LFM信号参数估计", 《数据采集与处理》, no. 01, 15 January 2018 (2018-01-15), pages 36 - 44 *

Also Published As

Publication number Publication date
CN112595889B (zh) 2022-07-15

Similar Documents

Publication Publication Date Title
CN110632387B (zh) 一种基于交流量子电压的谐波电压测量方法
CN108414833B (zh) 一种信号分量频率的精确估计方法
CN102539944B (zh) 一种基于鉴相法测量相位噪声的方法
CN101852826A (zh) 一种电力系统的谐波分析方法及其装置
CN108318737A (zh) 一种高精度数字鉴相方法
Shaw et al. A Kalman-filter spectral envelope preprocessor
CN109889231A (zh) 基于随机解调和有限新息率的脉冲串信号欠采样方法
CN112051583A (zh) Fmcw距离测量系统中拍频信号非线性校正方法
CN114460527B (zh) 关联度延拓Hilbert移相电子式互感器校验仪溯源方法和系统
CN110837003B (zh) 一种基于三角窗的双窗全相位dft同步相量测量方法及系统
CN110008434B (zh) 一种高精度的简谐信号参数估计方法
CN108710029B (zh) 一种信号谐波分量初相位的精确估计方法
CN112595889B (zh) 非理想多指数衰减正弦信号的欠Nyquist采样与参数测量方法
JPWO2003036313A1 (ja) クロック・スキュー測定装置、クロック・スキュー測定方法
CN112964931B (zh) 基于双通道欠采样的非理想多阻尼谐波信号参数测量方法
Zhang Variance analysis of noisy LFM signal in linear canonical Cohen's class
Shao et al. A frequency correcting method combining bilateral correction with weighted average for vortex flow sensor signal
CN105282086A (zh) 变频系统群时延的双音正交测量方法
CN111551785B (zh) 基于无迹卡尔曼滤波的频率与谐波检测方法
CN109030942B (zh) 谐相角分析方法
Liu et al. Digital instantaneous frequency measurement with wide bandwidth for real-valued waveforms using multiple sub-Nyquist channels
CN110297199B (zh) 一种基于全相位fft的铯光泵磁力仪频率测量方法及系统
CN112485518B (zh) 一种双通道多指数衰减正弦信号欠采样方法
CN114301552A (zh) 一种数字调制信号测试方法及系统
CN110133738B (zh) 基于IpDFT的质子磁力仪自由感应衰减信号的频率估计方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant