本发明属于通信信号处理领域，尤其涉及一种基于反馈结构的多谐波信号欠采样方法及系统。
背景技术：
多谐波信号是一种在实际生活中被广泛使用的信号，经常用在音频、通讯和电力系统中。根据奈奎斯特采样定理我们可知，为了从采样样本中完全重建模拟信号，采样率必须大于或等于信号带宽的两倍。随着多谐波信号带宽的逐步增加，采样设备的压力也随之增大。因而，奈奎斯特(Nyqiust)采样定理逐渐成为多谐波信号采样系统设计的瓶颈，制约着信号处理技术的发展。
多谐波信号是一种可以由有限数量的参数来表征的信号，即频率和幅值参数。对于一个由K个频率成份组成的多谐波信号，可以表示为：
其中，fk表示信号的第k个频率成份，表示对应的信号复数幅值参数，T表示信号的持续时间，e表示自然对数的底数，j表示虚数单位。假设信号的最大频率为fmax，即0≤fk<fmax,由于多谐波信号y(t)可以由一组频率和幅值参数唯一确定，该信号在时间[0,T)内的自由度为2K。
如果以fs的采样率对信号y(t)进行均匀采样，那么采集到的样本可以表示为：
其中，fs表示主采样率，N表示采集到的样本数，Z表示整数。当采样率满足Nyquist采样定理，即当fs≥fmax时，参数可以利用谱估计方法唯一确定。然而，当采样率不满足Nyquist采样定理，即当fs<fmax时，将会导致频率混叠问题。根据三角函数的周期性，对于每一个估计频率参数能够得到一组可能解：
其中，表示满足公式(a2)的最小可能解，mk表示整数。当信号主采样率fs满足这种特殊情况时，由于频率参数fa和fb难以区分。将这种特殊情况称之为镜像频率混叠。
目前已有许多关于多谐波信号的欠Nyquiat采样方案被提出。国外的Zoltowski等人提出了利用延时采样结合空间谱估计来估计频率的方法，然而在某些条件下，估计的频率并不完全准确。有学者基于压缩感知理论提出模拟-信息转换(AIC)技术，如随机解调器(RD)、随机采样(RS)、随机出发的调制宽带压缩采样(RT-MWCS)等结构都可以实现对多谐波信号的欠采样及信号的频率估计，不过该方法中的随机采样很难在硬件中实现，而且估计的精度取决于频率网格的密度，通常为了考虑计算的复杂度，网格密度都不高，因此通常恢复的精度都不高。有学者提出利用双通道采样实现对多谐波信号的频率估计，但是该方法估计出的频率有时并不唯一，为了改进此方法国内华中科技大学的研究团队提出利用三个采样通道互质采样的方法来解决频率混叠，但是该方法需要大量的样本。迄今为止，针对多谐波信号还没有一个所需样本数少、易于实现的抗混叠欠采样方案。因此设计一种简单有效，且能够解决频率混叠、镜像频率混叠问题的欠采样结构十分必要。
技术实现要素：
本发明针对多谐波信号欠采样出现的频率混叠、镜像频率混叠等问题，旨在提出一种基于反馈结构的多谐波信号欠采样方法及系统。
本发明的一种基于反馈结构的多谐波信号欠采样方法，包括下述步骤：
S1.主采样通道与信号发送单元连接，以对所述信号发送单元发送的多谐波信号进行低速采样获取主采样样本，根据普罗尼算法估计主采样频率参数和主采样幅值参数的主采样可能解组合；
S2.根据估计的主频率参数集生成反馈采样率，及切换控制信号，所述信号发送单元根据所述切换控制信号断开与所述主采样通道的连接，与反馈采样通道连接；
S3.所述反馈采样通道根据所述反馈采样率对所述多谐波信号进行低速采样获取反馈采样样本，采用普罗尼算法估计反馈频率参数和反馈幅值参数的反馈可能解组合；
S4.根据所述主采样可能解组合和所述反馈可能解组合联合获取信号的估计频率和估计幅值参数。
优选的，在所述步骤S2中，所述反馈采样率f′s满足：
其中，rem((·),1)表示(·)除以1后的余数，ε表示由噪声强度决定的阈值，fa和fb均表示频率参数a,b∈{0,1,…,K-1}，K表示信号的谐波个数，fs表示主采样率，m表示整数。
优选的，在所述步骤S4中获取信号的估计频率和估计幅值参数的过程为：
S41.计算所述主采样通道估计的主采样频率参数可能解：
S42.计算所述反馈采样通道估计的反馈采样频率参数可能解：
S43.根据主采样频率参数可能解和反馈采样频率参数可能解确定真实的估计频率参数；
S44.根据真实的频率参数构造复指数序列：
S45.根据复指数序列中的索引集确定真实的估计幅值参数。
本发明还提供了一种基于反馈结构的多谐波信号欠采样系统，包括：
发送单元，用以发送多谐波信号；
主采样通道，用以对所述多谐波信号进行低速采样获取主采样样本；
第一普罗尼算法单元，用以根据所述主采样样本采用普罗尼算法估计主采样频率参数和主采样幅值参数的主采样可能解组合；
反馈采样率生成器，用以根据估计的主频率参数集生成反馈采样率，并生成切换控制信号，所述发送单元根据所述控制信号，断开与所述主采样通道的连接，所述发送单元与反馈采样通道连接；
所述反馈采样通道，用以根据所述反馈采样率对所述多谐波信号进行低速采样获取反馈采样样本；
第二普罗尼算法单元，用以根据所述反馈采样样本采用普罗尼算法估计反馈频率参数和反馈幅值参数的反馈可能解组合；
联合估算单元，用以根据所述主采样可能解组合和所述反馈可能解组合联合获取信号的估计频率和估计幅值参数。
优选的，所述反馈采样率f′s满足：
其中，rem((·),1)表示(·)除以1后的余数，ε表示由噪声强度决定的阈值，fa和fb均表示频率参数a,b∈{0,1,…,K-1}，K表示信号的谐波个数，fs表示主采样率，m表示整数。
优选的，所述联合估算单元获取信号的估计频率和估计幅值参数的过程为：计算所述主采样通道估计的主采样频率参数可能解：计算所述反馈采样通道估计的反馈采样频率参数可能解：根据主采样频率参数可能解和反馈采样频率参数可能解确定真实的估计频率参数；根据真实的频率参数构造复指数序列：根据复指数序列中的索引集确定真实的估计幅值参数。
上述技术特征可以各种适合的方式组合或由等效的技术特征来替代，只要能够达到本发明的目的。
本发明的有益效果在于，本发明可以实现对多谐波信号的欠采样。当信号频率很高时，采用本发明提出的采样方法完成采样和参数估计的速率远小于奈奎斯特采样频率的速率，极大的降低采样设备的压力。多谐波信号欠采样系统可以有效解决频率混叠，以及镜像频率混叠的问题。
附图说明
图1为本发明所述的基于反馈结构的多谐波信号欠采样方法的一种实施例的方法流程图；
图2为本发明所述的基于反馈结构的多谐波信号欠采样方法的原理图；
图3为本发明所述的基于反馈结构的多谐波信号欠采样系统的一种实施例的系统框图；
图4为不同信噪比下各方法的估计效果曲线图；
图5为不同采样样本对估计效果的影响曲线图；
图6为不同频率分量的数量对估计效果的影响曲线图。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。
需要说明的是，在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。
下面结合附图和具体实施例对本发明作进一步说明，但不作为本发明的限定。
如图1所示，一种基于反馈结构的多谐波信号欠采样方法，包括下述步骤：
S1.主采样通道与信号发送单元连接，以对所述信号发送单元发送的多谐波信号进行低速采样获取主采样样本，根据普罗尼算法(Prony)估计主采样频率参数和主采样幅值参数的主采样可能解组合；
S2.根据估计的主频率参数集生成反馈采样率，及切换控制信号，所述信号发送单元根据所述切换控制信号断开与所述主采样通道的连接，与反馈采样通道连接；
S3.所述反馈采样通道根据所述反馈采样率对所述多谐波信号进行低速采样获取反馈采样样本，采用普罗尼算法估计反馈频率参数和反馈幅值参数的反馈可能解组合；
S4.根据所述主采样可能解组合和所述反馈可能解组合联合获取信号的估计频率和估计幅值参数。
在本实施例中，多谐波信号欠采样方法可以实现对多谐波信号的欠采样。当信号频率很高时，采用本发明提出的采样方法完成采样和参数估计的速率远小于奈奎斯特采样频率的速率，极大的降低采样设备的压力；同时，还可以有效解决频率混叠，以及镜像频率混叠的问题。
结合图2所示，以K个频率成份的多谐波信号为例采用基于反馈结构的多谐波信号欠采样方法进行采用的具体过程如下：
步骤一、初始化：
假设待采样的多谐波信号y(t)由K个频率成份组成，可以表示为：
其中，fk表示信号的第k个频率成份，ck(ck≠0,ck∈C)表示对应的信号复数幅值参数，C表示常数，T表示信号的持续时间，e表示自然对数的底数，j表示虚数单位。假设信号的最大频率fmax先验已知，即0≤fk<fmax,
步骤二，主采样通道低速采样：
多谐波信号y(t)在主采样通道以主采样率fs完成低速均匀采样，采样值可以表示为：
步骤三，采用Prony算法估计主采样频率参数和主采样幅值参数的主采样可能解组合：
主采样通道获取的N≥2K个连续的采样值y[n],n＝0,1,…,N-1后，利用Prony算法估计出2P个主采样参数其中为主采样频率参数最小解，P表示主采样通道中镜像频率非混叠的独立谐波个数。当不存在镜像频率混叠时，P＝K，且估计的幅值参数为准确的值，即当存在镜像频率混叠时，P≠K，且估计的幅值参数也存在混叠。
步骤四，生成反馈通道的采样率并切换采样通道：
根据主采样通道估计的主采样频率参数集生成反馈采样通道的反馈采样率f′s。该反馈通道的采样率应满足：
其中，fa和fb表示频率参数，m的取值范围为：
步骤五，反馈采样通道低速采样：
将采样通道切换到反馈采样通道之后，多谐波信号y(t)在以反馈采样率f′s低速均匀采样，采样值可以表示为：
其中，T1<T表示主采样通道的观测时长，n′＝0,1,…,N′-1，N′表示反馈采样通道采集的样本个数，T′s表示反馈采样通道的采样间隔；
步骤六，采用Prony算法估计反馈频率参数和反馈幅值参数的反馈可能解组合：
反馈采样通道获取的N′≥2K个连续的采样值y′[n′],n′＝0,1,…,N′-1后，利用Prony算法估计出2P′个参数其中，P′表示反馈采样通道中镜像频率非混叠的独立谐波个数，为反馈频率参数最小解。当不存在镜像频率混叠时，P′＝K，且估计的反馈幅值参数为准确的值，即当存在镜像频率混叠时，P′≠K，且估计的反馈幅值参数也存在混叠。
步骤七，采用参数联合估计算法确定估计频率参数和估计幅值参数的准确值。
进一步地，在步骤四中噪声情况下，为了提高采样系统的稳定性，反馈通道的反馈采样率f′s需要满足：
其中，rem((·),1)表示(·)除以1后的余数，ε表示由噪声强度决定的阈值，fa和fb均表示频率参数a,b∈{0,1,…,K-1}，K表示信号的谐波个数，fs表示主采样率，m表示整数。然后利用Cadzow迭代算法处理采样数据，提高估计的鲁棒性。
在优选的实施例中，在所述步骤S4中获取信号的估计频率和估计幅值参数的过程为：
S41.计算所述主采样通道估计的主采样频率参数可能解：
根据三角函数的周期性，所有可能解集合为：
真实的主频率参数
其中，表示满足公式(2)的最小可能解，mk表示整数，fmax表示信号的最大频率，fs表示主采样率；
S42.计算所述反馈采样通道估计的反馈采样频率参数可能解：
根据三角函数的周期性，所有可能解集合为：
真实的频率参数
其中，表示满足公式(5)的最小可能解，m′k表示整数，f′s表示反馈采样率；
S43.根据主采样频率参数可能解和反馈采样频率参数可能解确定真实的估计频率参数；
S44.根据真实的频率参数构造复指数序列：
其中，表示主采样估计的频率参数复指数形式，e表示自然对数的底数，j表示虚数单位，，表示反馈采样估计的频率参数复指数形式；
获得和中不同元素的索引集κ和κ′，其中，κ∪κ′＝{0,1,…,K-1}；
S45.根据复指数序列中的索引集确定真实的估计幅值参数；
其中，表示主采样估计的幅值参数，表示表示反馈采样估计的幅值参数。
如图3所示，本发明还提供了一种基于反馈结构的多谐波信号欠采样系统，包括：
发送单元1，用以发送多谐波信号；
主采样通道2，用以对所述多谐波信号进行低速采样获取主采样样本；
第一普罗尼算法单元5，用以根据所述主采样样本采用普罗尼算法估计主采样频率参数和主采样幅值参数的主采样可能解组合；
反馈采样率生成器3，用以根据估计的主频率参数集生成反馈采样率，并生成切换控制信号，所述发送单元1根据所述控制信号，断开与所述主采样通道2的连接，所述发送单元1与反馈采样通道4连接；
所述反馈采样通道4，用以根据所述反馈采样率对所述多谐波信号进行低速采样获取反馈采样样本；
第二普罗尼算法单元6，用以根据所述反馈采样样本采用普罗尼算法估计反馈频率参数和反馈幅值参数的反馈可能解组合；
联合估算单元7，用以根据所述主采样可能解组合和所述反馈可能解组合联合获取信号的估计频率和估计幅值参数。
在本实施例中，基于反馈结构的多谐波信号欠采样系统实现了对多谐波信号的欠采样。当信号频率很高时，完成采样和参数估计的速率远小于奈奎斯特采样频率的速率，极大的降低采样设备的压力，并且有效解决了频率混叠，以及镜像频率混叠的问题。在优选的实施例中，所述反馈采样率f′s满足：
其中，rem((·),1)表示(·)除以1后的余数，ε表示由噪声强度决定的阈值，fa和fb均表示频率参数a,b∈{0,1,…,K-1}，K表示信号的谐波个数，fs表示主采样率，m表示整数。
在优选的实施例中，所述联合估算单元获取信号的估计频率和估计幅值参数的过程为：计算所述主采样通道估计的主采样频率参数可能解：计算所述反馈采样通道估计的反馈采样频率参数可能解：根据主采样频率参数可能解和反馈采样频率参数可能解确定真实的估计频率参数；根据真实的频率参数构造复指数序列：根据复指数序列中的索引集确定真实的估计幅值参数。
为了验证本发明方法的性能，将本发明中提出的采样系统和现有的延时采样系统、三通道互质采样系统进行对比分析。三种方法的对比如表1所示，从中可以看出，本发明的多谐波信号欠采样方法仅需要4K个样本即可同时解决频率混叠问题和镜像频率混叠问题。
表1
为了更清晰的说明本发明提出的采样系统的优势，进行了仿真实验来验证。信号采用公式(1)中定义的多谐波信号，包含K个频率成份，最大频率为10GHz且最小频率间隔大于1MHz。幅值参数在(0，1+i)内随机选取。信号持续时间T＝10-6sec。各欠采样系统的仿真参数设置如表2所示：
表2
为了定量的描述参数估计的准确性，便于比较。引入归一化均方误差(NMSE)作为评价指标。
其中，fk是真实的频率，是估计频率。因为幅值的误差和频率是成正比的，所以仅用NMSE计算频率来分析方法的有效性。
实验一，考虑有镜像频率混叠、无噪声的情况。频率分量个数K设置为5，分别为[0.4,3.5,5.2,6.8,9.1]GHz。因为所以f2＝3GHz和f3＝5.5GHz的镜像频率会发生混叠。反馈采样系统用4K＝20个样本，时间延时采样系统用4K＝20个样本，三通道互质采样系统的快拍数L取4时，共需要4×(3×4+4×5+3×5)＝188个样本数。时间延时采样和三通道互质采样的参考文献中没有提到幅度的恢复，因此这里仅对其频率进行恢复。通过参数恢复比较表3可以看出来，本发明的多谐波信号欠采样方法可以精确地恢复频率及幅值参数。
表3
实验二，下面分析提出的方法在噪声情况下的表现，频率分量个数K设置为3，并且频率在(0,10)GHz内随机取。在信号上叠加高斯白噪声。实验中，反馈采样系统和时间延迟采样系统的样本数设置为Num＝60K＝180，三通道互质采样的快拍数L设置为4，总的样本数为188。输入信噪比从-20dB变化到50dB，每次实验做100次，获得平均的恢复结果如图4所示，从中可以看出，本发明的多谐波信号欠采样方法具有很好的抗噪性。
实验三，下面这个实验来分析采样样本数对三种方法估计效果的影响。频率分量个数K设置为3，且在在(0,10)GHz内随机选取，并在信号上叠加20dB的高斯白噪声。反馈采样和时间延迟采样的样本数以12的步进从48增加到336，三通道互质采样的快拍数L从1增加到7，总样本数以7的步进从47增加到329。每次实验重复100次，平均的估计结果和样本数的关系如图5所示，本发明多谐波信号欠采样方法要优于其他两种方法，而且随着样本数的增加，优势更加明显。
实验四，本实验来分析频率分量的个数对各方法性能的影响。实验中，频率分量个数K从1变化到25，每个频率都在(0,10)GHz内随机选取。加入信噪比为20dB的高斯白噪声，反馈采样系统和时间延时采样的样本数设置为324，三通道互质采样的样本数设置为329，其中快拍数L设置为7。每次实验重复100次，从图6可以看出在频率分量个数较少的时候，本发明方法的优势非常明显。随着频率分量个数的增加，各采样结构的重构性能趋于稳定。
综上所述，本发明可以实现对多谐波信号的欠采样，最低的等效采样率仅为当信号频率很高时，本发明提出的采样方法可以以远小于奈奎斯特采样频率的速率完成采样和参数估计，可以极大的降低采样设备的压力。同时，本发明方法还可以有效解决频率混叠，以及镜像频率混叠的问题。
虽然在本文中参照了特定的实施方式来描述本发明，但是应该理解的是，这些实施例仅仅是本发明的原理和应用的示例。因此应该理解的是，可以对示例性的实施例进行许多修改，并且可以设计出其他的布置，只要不偏离所附权利要求所限定的本发明的精神和范围。应该理解的是，可以通过不同于原始权利要求所描述的方式来结合不同的从属权利要求和本文中所述的特征。还可以理解的是，结合单独实施例所描述的特征可以使用在其他所述实施例中。