基于图及持续性熵特征的CR信号调制识别方法及系统

基于图及持续性熵特征的cr信号调制识别方法及系统
技术领域
1.本发明属于信号识别与处理领域，具体涉及基于图及持续性熵特征的cr信号调制识别方法及系统。


背景技术：

2.认知无线电(cognitive radio，cr)信号的调制识别是频谱感知中重要任务之一，同时也是解调等后续信号处理的前提，其性能的好坏将影响整个认知信号处理的效果。目前，主要有两类方法：一类是基于似然函数(likelihood-based，lb)的调制识别方法，另一类是基于特征(feature-based，fb)的调制识别方法。前者虽然能得到贝叶斯最优解，但也存在依赖信号先验信息、无法保证封闭解存在、高计算复杂度和概率失配等问题；后者虽然无需先验信息，但特征提取需要人为定义或大量数据训练，对经验依赖性较大。对fb算法而言，如何提取区分性强、鲁棒性强且计算复杂度低的特征，是实现调制识别的关键。
3.近年来，拓扑域信号处理技术，特别是全球定位系统(global positioning system，gps)研究的兴起，为调制识别的特征工程提供新的途径。gsp技术将信号从传统的时域、频域转换到图域进行调制识别，在一定程度上改善算法在低信噪比时的性能，已得到了相关领域研究者的关注。目前，gsp已在信号调制识别方面得到初步应用，其基本思想是将信号及其派生形式(如功率谱、时域样本等)通过一定的规则映射到图域，而后基于所生成的图拓扑结构进行特征分析，获取信号的特征用于识别。但有关gsp在cr调制识别当中应用的文献不多。


技术实现要素：

4.本发明针对现有技术中的不足，提供一种基于图及持续性熵特征的cr信号调制识别方法及系统。
5.为实现上述目的，本发明采用以下技术方案：
6.基于图及持续性熵特征的cr信号调制识别方法，其特征在于，包括以下步骤：
7.步骤1：在进行调制识别前先对观测信号做短时滤波的预处理；
8.步骤2：计算观测信号的平方谱，加窗后进行图域转换；根据图的gini系数判断观测信号为i类信号还是ii类信号；
9.步骤3：对i类信号和ii类信号使用不同的识别特征和识别规则进行识别：
10.对i类信号，计算观测信号的四次方谱，加窗后进行图域转换；根据图的度之和判断观测信号为qpsk信号还是16qam信号；
11.对ii类信号，计算观测信号的时频曲线，并求时频曲线的持续性图和持续性熵，根据持续性熵，判断观测信号为2ask信号还是bpsk信号。
12.为优化上述技术方案，采取的具体措施还包括：
13.进一步地，所述步骤1的预处理过程包括以下子步骤：
14.步骤1.1：将离散采样后的原始观测信号进行分段，其形式表示为：
15.ri(n)＝si(n)+ωi(n)，i(n
s-1)≤n≤(i+1)(n
s-1)
16.其中，si(n)是发送端发送的调制信号，ωi(n)是加性高斯白噪声，ns为信号段的长度，i表示段数；
17.步骤1.2：计算ri(n)的ns点的离散傅里叶变换，得到ri(k)＝dft[ri(n)]；
[0018]
步骤1.3：设计一个传输特性如下的带通滤波器：
[0019][0020]
其中，ks为|ri(k)|的最大谱线位置，|
·
|表示取模，ds为滤波的点数；将ri(k)输入设计的滤波器，其输出表示为：
[0021]ri
′
(k)＝hs(k)ri(k)
[0022]
步骤1.4：计算ri′
(k)的ns点逆离散傅里叶变换，得到ri′
(n)＝idft(ri′
(k))；
[0023]
步骤1.5：将每个经分段重构后的时域信号按原先的分段顺序重新组合成新的观测信号
[0024]
进一步地，所述步骤2包括以下子步骤：
[0025]
步骤2.1：对经过步骤1处理后的观测信号做平方运算，进而通过离散傅里叶变换得到观测信号的平方谱；
[0026]
步骤2.2：以谱线最大值为中心位置加矩形窗，取出矩形窗内的谱线，窗口的大小为取出的谱线数量；对平方谱χ(τ)，其经过加窗后的输出序列表示为：
[0027]
χw(τ)＝χ(τ)
·
w(τ)
[0028]
式中，
[0029][0030]
其中，τw为|χ(τ)|的最大值位置，dw为窗口宽度；
[0031]
步骤2.3：对加窗处理后的平方谱进行图域转换，并计算图的gini系数gc，根据经验设置判决门限ηb，0.1＜ηb＜0.5；若gc＜ηb则将观测信号判为i类信号，否则将观测信号判为ii类信号；其中，i类信号包括qpsk和16qam，ii类信号包括2ask和bpsk。
[0032]
进一步地，所述步骤3中，对i类信号的识别过程包括以下步骤：
[0033]
步骤3.1.1：对经过步骤1处理后的观测信号做四次方运算，并对其做傅里叶变换后取模，得到其四次方谱；以谱线最大值为中心位置加矩形窗，取出矩形窗内的谱线，窗口的大小为取出的谱线数量；
[0034]
步骤3.1.2：将加窗后的信号进行图域转换，计算图的度之和d
sum
：
[0035][0036]
其中，表示由图上各个顶点的度的值构成的向量，d
δ
为第δ个顶点的度，n0表示图的顶点数；
[0037]
步骤3.1.3：根据经验设定门限35＜ηi＜40，若d
sum
＜ηi，则判为qpsk信号，否则判为16qam信号。
[0038]
进一步地，所述步骤3中，对ii类信号的识别过程包括以下步骤：
[0039]
步骤3.2.1：使用短时傅里叶变换获取经过步骤1处理后的观测信号的时频曲线，将信号分成n
short
个分块数据，n
short
为短时窗个数，对每个分块的信号段做频率估计，得到瞬时频率的估计值进而得到时频曲线样本
[0040]
步骤3.2.2：对得到的时频曲线样本，采用子集滤流的方法生成时频曲线的持续性图；对于两个配对的点(l，m)，m，l＝1，...，n
short
，绝对值为其持续性值；将每个配对点映射成并绘制在二维直角坐标系中，从而得到持续性图；
[0041]
步骤3.2.3：对于得到的持续性图，其有n
pd
个配对点，映射到持续性图上的点分别为为和分别是l和m中的取值，表示第i个配对成功的点，i＝1，...，n
pd
；计算每个配对点的持续性值pi：
[0042][0043]
并对得到的持续性值归一化得到归一化后的持续性值npi：
[0044][0045]
则持续性熵h
pd
为：
[0046][0047]
步骤3.2.4：根据经验设定门限η
ii
＝3，识别规则为：若h
pd
＜η
ii
，则判为2ask信号，否则判为bpsk信号。
[0048]
本发明还提出了基于图及持续性熵特征的cr信号调制识别系统，其特征在于，包括：
[0049]
预处理模块，用于在进行调制识别前先对观测信号做短时滤波的预处理；
[0050]
判断模块，用于计算观测信号的平方谱，加窗后进行图域转换；根据图的gini系数判断观测信号为i类信号还是ii类信号；
[0051]
识别模块，用于对i类信号和ii类信号使用不同的识别特征和识别规则进行识别：
[0052]
对i类信号，计算观测信号的四次方谱，加窗后进行图域转换；根据图的度之和判断观测信号为qpsk信号还是16qam信号；
[0053]
对ii类信号，计算观测信号的时频曲线，并求时频曲线的持续性图和持续性熵，根据持续性熵，判断观测信号为2ask信号还是bpsk信号。
[0054]
本发明还提出了一种计算机可读存储介质，存储有计算机程序，其特征在于，所述计算机程序使计算机执行如上所述的基于图及持续性熵特征的cr信号调制识别方法。
[0055]
本发明还提出了一种电子设备，其特征在于，包括：存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行计算机程序时，实现如上所述的基于图及持续性熵特征的cr信号调制识别方法。
[0056]
本发明的有益效果是：本发明针对2ask、bpsk、qpsk和16qam四种调制信号，提出了一种基于图及持续性熵特征的识别方法及系统。首先通过短时滤波技术对观测信号进行预处理，用以补偿因非线性变换带来的信噪比损失，以增加识别特征的可区分性。而后分类间识别和类内识别两个步骤，通过提取经处理后的观测信号生成的图的gini系数、度之和以及持续性熵作为识别统计量，实现调制识别。仿真实验表明，在低信噪比时，本发明仍具较好的识别能力。
附图说明
[0057]
图1是本发明的识别方法流程图。
[0058]
图2为不同信噪比下四种调制信号经加窗后平方谱转换成的图的gini系数对比图。
[0059]
图3为qpsk及16qam信号四次方谱生成图的度之和对比示意图。
[0060]
图4为ii类信号时频曲线的持续性熵对比示意图。
[0061]
图5为所提出方法在不同信噪比条件下的识别正确率统计示意图。
具体实施方式
[0062]
现在结合附图对本发明作进一步详细的说明。
[0063]
如图1所示的基于图及持续性熵特征的cr信号调制识别方法，包括如下步骤：
[0064]
步骤1：滤波预处理：在进行调制识别前先对观测信号做短时滤波的预处理，用以补偿非线性变换带来的信噪比损失。短时滤波的点数根据经验一般选择10～30。
[0065]
假设经过了awgn信道的观测信号r(t)可以表示为：
[0066]
r(t)＝s(t)+ω(t)
[0067]
其中，s(t)为发送端发送的已调制信号，ω(t)为经过了awgn信道之后叠加的零均值实高斯白噪声，其方差为σ2。
[0068]
信噪比较低时，信号将出现淹没在噪声中的情况，识别特征可区分性下降；并且在提取特征的过程中往往需要进行非线性运算，这将增加噪声功率，降低信号处理的信噪比，因此在对cr信号调制识别前，先进行滤波处理，以减少低信噪比环境带来的影响。具体如下：
[0069]
步骤1.1：将离散采样后的原始观测信号按适当的长度进行分段，其形式可以表示为：
[0070]ri
(n)＝si(n)+ωi(n)，i(n
s-1)≤n≤(i+1)(n
s-1)
[0071]
其中，si(n)是发送端发送的调制信号，ωi(n)是加性高斯白噪声，ns为信号段的长度，i表示段数。
[0072]
步骤1.2：计算ri(n)的ns点的离散傅里叶变换(discrete fourier transform，dft)，得到ri(k)＝dft[ri(n)]。
[0073]
步骤1.3：设计一个传输特性如下的带通滤波器：
[0074][0075]
其中，ks为|ri(k)|的最大谱线位置(|
·
|表示取模)，ds为滤波的点数。将ri(k)输
入设计的滤波器，其输出可表示为：
[0076]ri
′
(k)＝hs(k)ri(k)
[0077]
步骤1.4：计算ri′
(k)的ns点逆离散傅里叶变换(inverse discrete fourier transform，idft)，得到ri′
(n)＝idft(ri′
(k))。
[0078]
步骤1.5：将每个经分段重构后的时域信号按原先的分段顺序重新组合成新的观测信号
[0079]
步骤2：类间识别：计算观测信号平方谱，加窗后进行图域转换。若得到的图的gini系数gc＜ηb，则将观测信号判为i类信号，否则将观测信号判为ii类信号。具体如下：
[0080]
步骤2.1：对重新组合的新观测信号做平方运算，进而通过离散傅里叶变换得到观测信号平方谱。
[0081]
步骤2.2：以谱线最大值为中心位置加矩形窗，取出最大值及其左右一定数量的谱线，窗口的大小为取出的谱线数量。对实离散序列χ(τ)，其经过加窗后的输出序列表示为：
[0082]
χw(τ)＝χ(τ)
·
w(v)
[0083]
式中，
[0084][0085]
其中，τw为|χ(τ)|的最大值位置，dw为窗口宽度。根据经验，类间识别的窗口宽度为20。
[0086]
步骤2.3：对加窗处理后的平方谱进行图域转换，并计算图的gini系数gc，根据经验可在(0.1～0.5)之间设置判决门限ηb。若gc＜ηb则将观测信号判为i类信号，否则将观测信号判为ii类信号。其中，i类信号包括qpsk和16qam，ii类信号包括2ask和bpsk。
[0087]
步骤3：类内识别：对i、ii两类信号使用不同的识别特征与规则，具体如下：
[0088]
i类信号类内识别：计算观测信号四次方谱，加窗后进行图域转换。若得到的图的度之和d
sum
＜ηi，则将观测信号判为qpsk信号，否则将观测信号判为16qam信号。
[0089]
ii类信号类内识别：计算观测信号的时频曲线，并求时频曲线的pd和持续性熵，若持续性熵h
pd
＜η
ii
，则将观测信号判为2ask信号，否则将观测信号判为bpsk信号。
[0090]
若步骤2中判为i类信号：
[0091]
步骤3.1.1：计算四次方谱，并加窗。先对经过了短时滤波的观测信号做四次方运算，并对其作傅里叶变换后取模，得到其四次方谱。以谱线最大值为中心位置加矩形窗，取出最大值及其左右一定数量的谱线，窗口的大小为取出的谱线数量，此时所用窗口宽度为25。
[0092]
步骤3.1.2：计算图的度之和。将加窗后的信号进行图域转换，计算图的度之和：
[0093][0094]
步骤3.1.3：调制识别。根据经验可以在(35～40)范围内设定门限ηi，若d
sum
＜ηi，则判为qpsk信号，否则判为16qam信号。
[0095]
若步骤2中判为ii类信号：
[0096]
步骤3.2.1：提取时频曲线。计算经过了短时滤波预处理的ii类信号的时频曲线，本方法使用短时傅里叶变换获取信号的时频曲线。对于观测信号r(t)，其短时傅里叶变换定义为：
[0097][0098]
式中，hs(u)为短时窗，*为取共轭运算。本质上，短时傅里叶变换，就是对短时窗内的分块信号做傅里叶变换，这些分块信号可以分开也可以重叠。假设短时窗有n
short
个，将信号分成了n
short
个分块数据，对每个分块的信号段做频率估计，得到瞬时频率的估计值如此就得到了时频曲线样本
[0099]
步骤3.2.2：持续性图(pd)获取。对得到的时频曲线样本，采用子集滤流的方法生成时频曲线的pd，关键点配对规则如下：
[0100]
(1)当引入一个新的连通分量时，称产生这个分量的局部极小值可代表该分量；
[0101]
(2)当遇到一个局部极大值，同时有两个连通分量待配对时，将极大值与代表这两个分量的两个局部极小值中的较高者配对，另一个极小值代表合并产生的连通分量。对于两个配对的点(l，m)，m，l＝1，...，n
short
，函数的差值的绝对值即为其持续性值。将每个配对点映射成并绘制在二维直角坐标系中，从而得到持续性图。横坐标为birth，代表出现时的函数值，纵坐标为death，代表消失时的函数值。
[0102]
步骤3.2.3：计算持续性熵。对于得到的pd，假设其有n
pd
个配对点，其函数值映射到pd上的点分别为计算每个配对点持续性值：
[0103][0104]
并对得到的持续性值归一化得到归一化后的持续性值：
[0105][0106]
则持续性熵为：
[0107][0108]
步骤3.2.4：调制识别。根据经验设定门限η
ii
＝3，识别规则为：若h
pd
＜η
ii
，则判为2ask信号，否则判为bpsk信号。
[0109]
图2～图4为类间及类内识别时不同调制类型检验统计量对比。
[0110]
类间识别中用于检测生成的图是否为完全图的检验统计量为图的gini系数gc，识别规则为：若gc＜ηb，则判定观测信号为i类信号，否则判定观测信号为ii类信号。其中ηb为判决门限，其取值范围根据经验选择。
[0111]
图2所示为信号长度为1024点，码元速率为200kbaud，载波频率为600khz，采样频率为4mhz，信噪比为-8～-3db，窗口宽度为20，量化级数为6的四种调制信号的gini系数gc
的对比示意图，图中每种信噪比下的gini系数都是1000次仿真的均值。由图中可见：i类信号(qpsk和16qam)和ii类信号(2ask和bpsk)在不同信噪比下的gini系数都有较大的差异，因此可以用gini系数作为特征来区分两类信号。
[0112]
图3所示为信号长度为1024点，码元速率为200kbaud，载波频率为600khz，采样频率为4mhz，信噪比为-8～-3db，窗口宽度为25，量化级数为10的两种i类信号的度之和对比图，图中每种信噪比下的度之和都是1000次仿真的均值。由图中可见，在不同信噪比时qpsk信号和16qam信号度之和都有一定的差异，并且随着信噪比的增大，度之和的差异也在变大。因此，在适用snr下，可以用图的度之和作为区分两种信号类型的特征。
[0113]
图4所示为信号长度为1024点，码元速率为200kbaud，载波频率为600khz，采样频率为4mhz，信噪比为-8～-3db条件下，两种ii类信号的持续性熵h
pd
随信噪比的变化对比示意图，图中每种信噪比下的持续性熵都是1000次仿真的均值。从图中可以看出在不同的信噪比下，bpsk信号和2ask信号的持续性熵都有较明显的差异，可以作为识别统计量用于调制识别。
[0114]
图5所示为所提出方法在不同信噪比条件下的平均识别正确率。可见，图及持续性熵特征法对调制信号的平均识别正确率随着信噪比的增加而增大。当信噪比为-5db时，信号的平均识别正确率接近80％；当信噪比为-3db时，信号的平均识别正确率可达99％以上。
[0115]
在另一实施例中，本发明还提出了基于图及持续性熵特征的cr信号调制识别系统，包括：
[0116]
预处理模块，用于在进行调制识别前先对观测信号做短时滤波的预处理；
[0117]
判断模块，用于计算观测信号的平方谱，加窗后进行图域转换；根据图的gini系数判断观测信号为i类信号还是ii类信号；
[0118]
识别模块，用于对i类信号和ii类信号使用不同的识别特征和识别规则进行识别：
[0119]
对i类信号，计算观测信号的四次方谱，加窗后进行图域转换；根据图的度之和判断观测信号为qpsk信号还是16qam信号；
[0120]
对ii类信号，计算观测信号的时频曲线，并求时频曲线的持续性图和持续性熵，根据持续性熵，判断观测信号为2ask信号还是bpsk信号。
[0121]
在另一实施例中，本发明提出了一种计算机可读存储介质，存储有计算机程序，计算机程序使计算机执行如上所述的基于图及持续性熵特征的cr信号调制识别方法。
[0122]
在另一实施例中，本发明提出了一种电子设备，包括：存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，处理器执行计算机程序时，实现如上所述的基于图及持续性熵特征的cr信号调制识别方法。
[0123]
在本技术所公开的实施例中，计算机存储介质可以是有形的介质，其可以包含或存储以供指令执行系统、装置或设备使用或与指令执行系统、装置或设备结合使用的程序。计算机存储介质可以包括但不限于电子的、磁性的、光学的、电磁的、红外的、或半导体系统、装置或设备，或者上述内容的任何合适组合。计算机存储介质的更具体示例会包括基于一个或多个线的电气连接、便携式计算机盘、硬盘、随机存取存储器(ram)、只读存储器(rom)、可擦除可编程只读存储器(eprom或快闪存储器)、光纤、便捷式紧凑盘只读存储器(cd-rom)、光学储存设备、磁储存设备、或上述内容的任何合适组合。
[0124]
以上仅是本发明的优选实施方式，本发明的保护范围并不仅局限于上述实施例，
凡属于本发明思路下的技术方案均属于本发明的保护范围。应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理前提下的若干改进和润饰，应视为本发明的保护范围。