CN116559786A - 一种混叠雷达信号脉内特征分析方法及系统 - Google Patents
一种混叠雷达信号脉内特征分析方法及系统 Download PDFInfo
- Publication number
- CN116559786A CN116559786A CN202310371898.2A CN202310371898A CN116559786A CN 116559786 A CN116559786 A CN 116559786A CN 202310371898 A CN202310371898 A CN 202310371898A CN 116559786 A CN116559786 A CN 116559786A
- Authority
- CN
- China
- Prior art keywords
- signal
- frequency
- modulation
- calculating
- value
- 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.)
- Pending
Links
- 238000004458 analytical method Methods 0.000 title claims abstract description 36
- 238000000034 method Methods 0.000 claims abstract description 91
- 238000001228 spectrum Methods 0.000 claims description 81
- 230000008569 process Effects 0.000 claims description 25
- 238000000926 separation method Methods 0.000 claims description 19
- 230000003595 spectral effect Effects 0.000 claims description 14
- 238000004590 computer program Methods 0.000 claims description 11
- 238000004364 calculation method Methods 0.000 claims description 7
- 230000000630 rising effect Effects 0.000 claims description 6
- 238000001914 filtration Methods 0.000 claims description 4
- 238000009499 grossing Methods 0.000 claims description 4
- 238000005516 engineering process Methods 0.000 abstract description 4
- 238000007667 floating Methods 0.000 description 22
- 238000005070 sampling Methods 0.000 description 11
- 230000009466 transformation Effects 0.000 description 10
- 238000004422 calculation algorithm Methods 0.000 description 8
- 239000002131 composite material Substances 0.000 description 8
- 238000001514 detection method Methods 0.000 description 8
- 238000012545 processing Methods 0.000 description 8
- 239000013598 vector Substances 0.000 description 8
- 108091026890 Coding region Proteins 0.000 description 7
- 238000010586 diagram Methods 0.000 description 6
- 238000013461 design Methods 0.000 description 5
- 230000006870 function Effects 0.000 description 4
- 238000004088 simulation Methods 0.000 description 4
- 238000004891 communication Methods 0.000 description 2
- 238000011161 development Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000010606 normalization Methods 0.000 description 2
- 238000010587 phase diagram Methods 0.000 description 2
- 230000004044 response Effects 0.000 description 2
- 238000012935 Averaging Methods 0.000 description 1
- 108010076504 Protein Sorting Signals Proteins 0.000 description 1
- 102100026758 Serine/threonine-protein kinase 16 Human genes 0.000 description 1
- 101710184778 Serine/threonine-protein kinase 16 Proteins 0.000 description 1
- 230000000739 chaotic effect Effects 0.000 description 1
- 230000008878 coupling Effects 0.000 description 1
- 238000010168 coupling process Methods 0.000 description 1
- 238000005859 coupling reaction Methods 0.000 description 1
- 125000004122 cyclic group Chemical group 0.000 description 1
- 238000000354 decomposition reaction Methods 0.000 description 1
- 230000005484 gravity Effects 0.000 description 1
- 238000012423 maintenance Methods 0.000 description 1
- 239000011159 matrix material Substances 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000004806 packaging method and process Methods 0.000 description 1
- 238000012567 pattern recognition method Methods 0.000 description 1
- 238000000819 phase cycle Methods 0.000 description 1
- 238000003825 pressing Methods 0.000 description 1
- 230000005855 radiation Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 239000007787 solid Substances 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
- 238000012800 visualization Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
- G01S7/021—Auxiliary means for detecting or identifying radar signals or the like, e.g. radar jamming signals
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
- G01S7/28—Details of pulse systems
- G01S7/285—Receivers
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
- G01S7/38—Jamming means, e.g. producing false echoes
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
- G01S7/40—Means for monitoring or calibrating
Landscapes
- Engineering & Computer Science (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Computer Networks & Wireless Communication (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Radar Systems Or Details Thereof (AREA)
Abstract
本发明公开了一种混叠雷达信号脉内特征分析方法及系统,涉及雷达对抗侦查领域,该方法包括:步骤一、把时频域交叠的混合信号的分离成单一的信号;步骤二、分别对分离后的每个信号检测信号起止时间;步骤三、分别对分离后的每个信号识别信号调制类型;步骤四、分别对分离后的每个信号估计调制参数。本申请还公开了一种基于混叠雷达信号脉内特征分析技术的计算机系统。本申请具有自动、准确和快速分析混叠雷达信号脉内特征的能力。
Description
技术领域
本发明属于雷达对抗侦察技术领域,具体涉及一种混叠雷达信号脉内特征分析方法及系统。
背景技术
如今,新体制雷达发展迅速并逐渐占居主导地位。新体制雷达中常规脉冲雷达信号比例已减少,线性调频、非线性调频、相位编码等雷达信号逐渐增多。这就对雷达信号的分选和识别提出了更高的要求,需要对雷达信号的脉内调制特征进行分析。通过分析,可以了解该雷达的效能及该国雷达技术的水平;研究该雷达的弱点以便在战时对其进行有效干扰。
雷达信号脉内特征分析是20世纪80年代中期开始研究的新技术,主要根据脉内的细微特征来分辨和识别复杂的调制信号,引导我方作出及时、准确的反应,对雷达对抗起着至关重要的作用。雷达信号脉内特征是雷达信号细微特征的最集中、最重要的表现,相对于雷达信号基本特征而言,是对信号的更深刻、更细微的数学物理描述。目前,雷达信号调制方式识别的方法主要有两大类,一是一是基于假设检验理论的判决理论方法;二是不依赖于假设条件的统计模式识别方法。雷达调制参数估计主要包括:高阶统计量,循环谱相关,时频分布,小波变换,分数阶Fourier变换以及混沌信号处理、Chirplet变换、现代谱估计等方法。
复杂战场电磁环境下,辐射源数量日趋密集,脉冲密度高达每秒数百万级,从通信信号频段到雷达信号频段,各个频段的信号交错纵横。对于中频信号,电子侦察接收机在截获过程中会有多个信号同时到达,特别是低频段,宽带接收时会混杂通信信号等我们不关注的信号,增加了关注信号识别的难度,给雷达信号脉内特征分析带来很大的困扰,导致调制方式识别和调制参数估计结果出现错误。
因此,在雷达信号脉内特征分析时,必须先对脉内时频域交叠的信号进行有效分离,此外,如何能够正确识别出分离信号的调制类型并估计其调制参数以方便对信号的处理应用,也是急需解决的技术问题。
发明内容
为解决上述技术问题,本发明提出一种混叠雷达信号脉内特征分析方法及系统,用于多信号混叠所引起的脉内特征有效分析的技术方案。
本发明第一方面公开了一种混叠雷达信号脉内特征分析方法,所述混叠雷达信号脉内特征分析方法包括如下步骤:
步骤一、将时频域交叠的混合信号分离成独立的信号并存储;
步骤二、根据分离后各个独立信号的包络,检测分离后的每一个所述信号的起止时间;
步骤三、识别分离后的每一个所述信号的调制类型;
步骤四、基于获得的每一个所述信号的调制类型,估计分离后的每一个所述信号的调制参数。
根据本发明第一方面的方法,步骤一中,所述混叠信号分离的过程包括:
步骤11)、对原始I/Q通道信号作快速傅里叶变换,得到待处理信号频谱;
步骤12)、找到所述待处理信号频谱幅值最大值所在位置;
步骤13)、判断所述待处理信号频谱幅值最大值所在位置是否为三阶多项式相位信号;如果是,计算所述幅值最大值附近应该保留的频谱宽度,进入步骤16);如果否,进入步骤14);
步骤14)、计算所述待处理信号频谱幅值最大值所在位置的峰度值和偏度值,依据所述峰度值和偏度值判断所述频谱幅值最大值所在位置是信号或噪声;如果是信号,则进入步骤15),否则退出混叠信号分离过程;
步骤15)、计算所述频谱幅值最大值附近应该保留的频谱宽度;
步骤16)、保留所述频谱宽度的信号,并进行存储,将所述待处理信号频谱的其他部分设为噪声,返回到步骤12)。
根据本发明第一方面的方法,步骤15)中,计算频谱幅值最大值附近应该保留的频谱宽度具体包括:
判断所述频谱幅值最大值信号是否为正弦调频信号,如果是,则计算该正弦调频信号最大频点左右要保留的带宽点数k,左右两边点数设置为相同,确定应该保留的频谱宽度;如果否,计算信号频谱与矩形频谱和三角形频谱的相像系数,根据获得的信号频谱与矩形频谱和三角形频谱的相像系数确定应该保留的频谱宽度。
根据本发明第一方面的方法,步骤二中,所述根据包络检测分离后的每一个信号的起止时间包括:
步骤21:对每一信号的包络作归一化处理,并进行平滑处理;
步骤22:计算平滑后的归一化包络并计算其均值和标准差;
步骤23:确定截取信号时上升沿和下降沿之间的最优比率;如果标准差小于0.1,最优比率取均值和标准差之和,进入步骤24;如果标准差大于或等于0.25,最优比率取0.4,进入步骤24;如果标准差在[0.1,0.25]区间内,则比较均值大小,均值大于0.7,把小值点处信号去掉,均值小于0.25则把大值点处信号去掉,均值在[0.25,0.7]区间内,保留全部信号,根据保留信号下标序号计算出信号的起止时间,跳过步骤24;
步骤24:根据确定的最优比率计算信号的起止时间,确定平滑后的归一化信号包络中大于最优比率的点,计算出信号包络的上升沿和下降沿的下标序号,根据下标序号计算出信号的起止时间。
根据本发明第一方面的方法,步骤三中,所述识别分离后的每一个信号的调制类型包括:用二阶四阶矩方法估计信噪比;对信号进行短时滤波处理;计算信号去卷叠后的瞬时相位;采用20重相位差分方法计算信号的瞬时频率;对时频曲线作一阶线性拟合,计算得到回归系数B,残差r和复相关系数R;根据获得的回归系数B,残差r和复相关系数R来判断每一个信号的调制类型。
根据本发明第一方面的方法,步骤四中,所述估计分离后的每一个信号的调制参数包括:把调制参数分为公有参数和特有参数进行估计,公有参数包括中心频率、带宽、幅度和脉宽,特有参数和调制类型有关,不同的调制类型,设置不同的估计方式来获得不同的特有参数。
本发明第二方面公开了一种混叠雷达信号脉内特征分析系统,所述系统存储有计算机程序指令,通过执行所述计算机程序指令,实现混叠雷达信号脉内特征分析方法,能够按缺省参数运行,也允许用户修改部分参数,并将运行结果可视化显示。
本发明第三方面公开了一种电子设备。电子设备包括存储器和处理器,存储器存储有计算机程序,处理器执行计算机程序时,实现本公开第一方面中任一项的混叠雷达信号脉内特征分析方法中的步骤。
本发明第四方面公开了一种计算机可读存储介质。计算机可读存储介质上存储有计算机程序,计算机程序被处理器执行时,实现本公开第一方面中任一项的混叠雷达信号脉内特征分析方法中的步骤。
综上,本发明提出的方案具备如下技术效果:
1、本申请提出了一种时频域交叠的混合信号的分离方法,能够自动、准确、快速分离混叠信号,为混叠信号的脉内特征分析奠定了坚实的基础。
2、本申请能够有效识别常见的信号调制方式,准确度高、适应性强。
3、本申请能够快速准确估计常见信号的调制参数,参数估计准确度高。
4、本申请既适应于混叠信号脉内特征分析,又适应于非混叠信号脉内特征分析。
5、本申请既可以用缺省参数自动运行,大大降低操作人员的工作量以及专业要求,也可以人工修改部分关键参数,实现重点关注信号的精细分析。
6、本申请提供的计算机软件具备操作简单、运行速度快、分析结果可视化的特点,便于操作人员使用和分析。
附图说明
为了更清楚地说明本发明具体实施方式或现有技术中的技术方案,下面将对具体实施方式或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图是本发明的一些实施方式,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。
图1为本发明的混叠雷达信号脉内特征分析方法信号处理流程示意图。
图2为本发明的混叠雷达信号脉内特征分析方法软件业务流程图。
图3为本发明的混叠雷达信号脉内特征分析方法中信号分解业务流程图。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚,下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例只是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
参见图1,混叠雷达信号脉内特征分析方法包括如下步骤:
步骤一、获得混叠信号后,将时频域交叠的混合信号分离成独立的信号并存储;
步骤二、根据包络检测分离后的每一个信号的起止时间;
步骤三、识别分离后的每一个信号的调制类型;
步骤四、估计分离后的每一个信号的调制参数;
步骤五、输出结果。
具体地,步骤一、所述混叠雷达信号分离包括如下步骤:
步骤11)、对原始I/Q通道信号作快速傅里叶变换,得到待处理信号频谱;
步骤12)、找到所述待处理信号频谱幅值最大值所在位置;
步骤13)、因为三阶多项式相位信号的频谱宽度较宽,且峰度和偏度与正态分布比较接近,容易被判断为噪声,因此先判断是否为三阶多项式相位信号。判断所述待处理信号频谱幅值最大值所在位置是否为三阶多项式相位信号;如果是,只计算最大值附近应该保留的频谱宽度(左k1点,右k2点),进入步骤16);如果否,进入步骤14);
步骤14)、计算所述待处理信号频谱幅值最大值所在位置的峰度值和偏度值,依据所述峰度值和偏度值判断所述频谱幅值最大值所在位置是信号或噪声;如果是信号,则进入步骤15),否则退出混叠信号分离过程;
步骤15)、计算所述频谱幅值最大值附近应该保留的频谱宽度(左k1点,右k2点或左右各k点);
步骤16)、保留最大值频点附近一定带宽(左k1点右k2点或左右各k点)的频谱,其它设为噪声,有待进一步进行调制类型识别及调制参数估计;并把原信号中被提取出的部分设为噪声,返回到步骤12)。
进一步的,步骤12)中,所述找到频谱最大值所在位置的过程包括:
每次只找频谱最大值点,在后续步骤中,如果不能达到退出混叠雷达信号分离的条件,则最多只执行步骤二15次后,自动退出。
进一步的,步骤13)中,所述判断是否为三阶多项式相位信号的过程包括:
频域信号的模经过移动平均法平滑后,比较平滑后的模的最小值、最大值和中位数的关系。如果中位数和最小值之差的两倍小于最大值与中位数之差,则计算信号瞬时频率,并对瞬时频率作二阶线性拟合,如果复相关系数大于0.85同时残差标准差小于3,则认为是三阶多项式相位信号,否则认为不是三阶多项式相位信号。
进一步的,步骤13)中,所述计算最大值附近应该保留的频谱宽度的过程包括:
令med为平滑后的模的中位数,令th=0.5×med。从最大值开始向左边以20为步长移动,计算出这20个点的平均值,如果连续两次平均值与med之差小于th,则停止移动,并得到左边要保留的点数k1;并用同样的方法计算出右边要保留的点数k2。
进一步的,步骤15)中,所述计算最大值附近应该保留的频谱宽度(左k1点,右k2点或左右各k点)的过程包括:
1.首先判断信号是否为正弦调频信号,如果是,则计算正弦调频信号最大频点左右要保留的带宽点数k(左右两边点数一样)。判断和计算方法为:
1)对信号频谱作归一化处理之后,计算归一化处理后的频谱的模F;
2)设置大值比例门限th=max(0.22,mean(F)+3×std(F)),选出大于门限的连续序号;
3)每一段连续序号理论上应该是同一个大值频点,如果大值频点数小于等于2,则进一步判断是不是调制频率比较小的正弦调频信号,判断方法如下,否则进入第4)步;
a.计算信号的瞬时频率;
b.计算瞬时频率的极大值和极小值点;
c.由极大值和极小值点的位置计算瞬时频率的周期;
d.按周期把瞬时频率分段,分别对每一段做正弦拟合。正弦拟合的思路是基于最小二乘法的原理,得到关于待拟合参数的非线性方程组,再用高斯-牛顿迭代法求解该非线性方程组,最终可得到拟合系数、拟合残差、残差标准差;
e.如果所有残差标准差均小于0.6,则认为是正弦调频信号,否则认为不是。
4)把每一个连续序号集合中,连续序号数量超过5个的剔除,再次进行上一步判断;
5)判断大值点是否在大值点序号集合里,如果不在,说明不是正弦调频信号,如果在,进入下一步;
6)以每段序号中最大值对应的序号为准,计算大值频点间的间距,设置容差为3,找到间距小于容差的大值频点(间距小于容差认为间距相等)。如果间距相等的大值谱线数量大于20,则认为不是正弦调频信号,否则进入下一步;
7)如果间距相等的大值频点占比超过80%,则认为是正弦调频信号,并删除间距不相等的点;否则去掉相邻大值差距特别大的,判断依据是相邻大值之差大于5倍大值频点间距的中值;
8)如果最大值频点的位置位于大值频点两端1/3,则认为不是SFM信号,否则继续判断;
9)由于有可能会漏掉大值点,所以需要考虑间距成倍数的点。如果间距和最小间距或最小间距的两倍之差在容差(设为3)范围内的大值频点超过80%,则判断为是正弦调频信号;
10)如果是正弦调频信号,重新计算大值点的自差分,并将大值点的边界左右各移动5个间距,寻找是否有遗漏的大值点,有遗漏则将遗漏的大值点加入,得到要保留的点数k1和k2。最后把计算得到的k1和k2各加上0.5MHz带宽。如果不是正弦调频信号,则进入步骤2;
2.如果不是正弦调频信号,计算信号频谱与矩形和三角形的相像系数,相像系数的计算步骤如下:
求出信号频谱的中心频率和有效带宽并对带宽进行归一化处理,只保留有效带宽内的频谱部分,这样可排除带外噪声的影响。预处理后的频域信号为{G(j),j=1,2,…,N},N为预处理后的信号长度;
1)构造矩形序列:mx是G(j)的最大值;
2)计算G(j)和U(k)的相像系数,计算公式如下:
3)构造三角形序列:
式中,mx是G(j)的最大值;
4)计算{G(j)}和{T(k)}的相像系数,计算公式如下:
5)计算相像系数时,除了得到信号频谱分别与矩形和三角形的相像系数以外,还得到大于设定带宽门限的点数B和大值点的带宽比例门限th<1;
6)如果B≤5,左右两边保留的点数均为k,取值是0.1N/fs的整数部分,其中,fs是信号的采样频率,N是信号点数;为防止k过小或过大,限定k的最小值为7,最大值为11,当B=1时,k的最小值限定为5;
7)如果B>5,且信号频谱和矩形的相像系数大于和三角形的相像系数,则认为是调频类信号,最大值频点附近保留的带宽点数计算方法为:先固定最大频点左边点数为10,在最大频点右边依次以2MHz步长选取频域信号Sf,按上面方法计算Sf与矩形序列的相像系数,直到相像系数不再变化,把当前的点数作为最大频点右边要保留的点数p2;按p2固定最大频点右边部分的信号,在最大频点左边依次以2MHz步长选取频域信号Sf,计算Sf与矩形序列的相像系数,直到相像系数不再变化,把当前的点数作为最大频点左边要保留的点数p1;在得到p1和p2的基础上,缩小步长进行更精细的计算,以求右边点数为例,从粗略值p2开始以0.5MHz的小步长逐渐减少保留的点数,计算与矩形序列的相像系数,直到相像系数发生变化,把当前的点数作为最大频点右边要保留的点数k2。同样的方法计算最大频点左边要保留的点数k1。
8)如果B>5,且信号频谱和矩形的相像系数小于和三角形的相像系数,则认为是编码类信号。由于其频谱是对称的,左右两边保留的带宽一样。最大值频点附近保留的带宽点数计算方法为:在最大频点两边依次以2MHz步长选取频域信号Sf,计算Sf与三角形序列的相像系数,直到相像系数不再变化,把当前的点数作为最大频点两边要保留的点数粗略值p,从p开始以0.5MHz的小步长逐渐减少保留的点数,计算与三角形序列的相像系数,直到相像系数发生变化,把当前的点数作为最大频点右边要保留的点数k;
进一步的,步骤二中,所述对分离出的信号分别根据包络检测起止时间的过程包括:
1)计算信号的包络作归一化处理,并进行平滑;
2)计算平滑后的归一化包络并计算它的均值和标准差;
3)确定截取信号时上升沿和下降沿的最优比率。如果标准差小于0.1,最优比率取均值和标准差之和,进入第4)步;如果标准差大于或等于0.25,最优比率取0.4,进入第4)步;如果标准差在[0.1,0.25]区间,则看均值大小,均值大于0.7,把小值点去掉,均值小于0.25则把大值点去掉,均值在[0.25,0.7]区间,保留全部信号,根据保留信号下标序号计算出起止时间,跳过第4)步;
4)根据确定的最优比率计算信号的起止时间并取出信号。找平滑后的归一化包络中大于最优比率的点,如果有多段,要分段处理,计算出包络的上升沿和下降沿的序号,并除以采样频率得到信号的起止时间。
对于步骤三中,多种脉内调制类型的识别方法,具体包括如下步骤:
步骤31)、首先判断是否为三阶多项式相位信号或正弦调频信号,这一步工作在混叠信号分离与检测中已经做过。如果是,对分离出的下一个信号重新进入步骤31),否则对本信号进入步骤32);
步骤32)、用二阶四阶矩方法估计信噪比;
步骤33)、对信号进行短时滤波处理;
步骤34)、计算信号去卷叠后的瞬时相位;
步骤35)、采用20重相位差分方法计算信号的瞬时频率;
步骤36)、对时频曲线作一阶线性拟合,计算得到回归系数B,残差r和复相关系数R;
步骤37)、判断是否为常规信号,如果是,进入步骤315),否则进入步骤38);
步骤38)、判断是否为线性调频信号,如果是,进入步骤315),否则进入步骤39);
步骤39)、判断是否为双线性调频信号,如果是,进入步骤315),否则对进入步骤310);
步骤310)、检测频率跳变点,找到跳变点的序号,如果跳变点个数小于或等于3个,则认为是未知调制类型信号,进入步骤315),否则进入步骤311);
步骤311)、计算瞬时频率的峰度,如果峰度大于4,进入步骤312),否则进入步骤313);
步骤312)、如果复相关系数大于0.45,判断为线性调频-二相编码复合调制信号,否则判断为相位编码信号,具体编码类型需要进一步判断,判断完成后进入步骤315);
步骤313)、把步骤十检测到的跳变点去除,重新检测跳变点,如果第二次检测到的跳变点个数超过第一次检测到跳变点个数的60%,判定为频率编码信号,进入步骤315),否则进入步骤314);
步骤314)、对瞬时频率作一阶线性拟合,如果复相关系数大于0.85且残差标准差小于6,判定为线性调频-二相编码复合调制信号,否则判定为频率编码-二相编码复合调制信号;
步骤315)、是否分离出的所有信号都已经完成调制类型识别,如果是,则进入步骤316),否则,对没有完成调制类型识别的信号进入步骤31)。
步骤316)、在混叠信号分离时,可能把频率编码信号拆分成常规信号,也有可能把频率编码-二相编码复合调制信号拆分成多个二相编码信号,因此要判断和处理这两种情况,处理完成后进入调制参数估计。
进一步的,步骤37)中,所述判断是否常规信号的过程包括:
1)如果残差r的标准差大于5,则认为不是常规信号,否则作进一步判断;
2)如果残差r的标准差小于0.8且R同时小于0.5,则认为是常规信号,否则作进一步判断;
3)用3dB带宽估计信号中心频率,如果中心频率小于0.05,认为脉内无调制;否则计算信号的归一化频谱,找到频谱的大值点,如果大值点只有一个,则认为是常规信号;考虑到噪声影响,如果大值点多于1个,进一步判断大值点个数和连续性,设定门限L表示信号长度,fs表示采样频率,round表示四舍五入取整,限定3≤p≤5,如果大值点个数小于p同时大值点序号间距大于3的在大值点总数中占比小于10%,则认为是常规信号。
进一步的,步骤38)中,所述判断是否线性调频信号的过程包括:
1)如果残差r的标准差大于6或者复相关系数R小于0.2,则认为不是线性调频信号,否则作进一步判断;
2)用3dB带宽估计信号中心频率;
3)估计信号去载频后的瞬时相位;
4)对相位作二阶线性拟合,计算得到回归系数B2,残差r2和复相关系数R2;
5)如果同时满足r的标准差小于2.5,r2的标准差小于0.45和R2大于0.95,则认为是线性调频信号;
6)如果R2小于0.5或者r、r2的标准差都大于5,判断为不是线性调频信号,需进一步判断;
7)如果不满足以上两条条件,根据频谱进一步区分。计算信号频谱的幅值,以超过最大幅值一半的为大值频点,计算大值频点之间的差值δf,如果同时满足R大于0.45,大值点频谱宽度大于3,所有的δf小于0.3,则认为是线性调频信号;如果R大于0.7同时R2大于0.95,则计算频点相差超过门限(超过门限的大值频点带宽大于20时,门限取1.5,其它情况取1)的点数,如果大值频点个数大于3同时频点相差超过门限的点数占比低于10%,则认为是线性调频信号,否则做进一步判断。
进一步的,步骤39)中,所述判断是否双线性调频信号的过程包括:
1)找到瞬时频率数据的尖峰点,从尖峰点把数据分为两段;
2)对两段数据分别作一阶线性拟合;
3)如果两段的相关系数中有小于0.5的,判定为不是DLFM信号,需要进一步判断;
4)如果两端的一阶系数(斜率)互为相反数同时任意一个相关系数大于0.85,判断是DLFM信号,否则需要进一步判断。
进一步的,步骤310)中,所述检测频率跳变点的方法的过程包括:
1)瞬时频率phiN每相隔40个点作差,即δ=phiNi+40-phiNi;
2)按下式计算跳变点检测门限pjump;
其中,fs是采样频率,N=20是多重相位差分重数,p=0.45;
3)如果δ的绝对值大于pjump,则判断为跳变点。
进一步的,步骤312)中,所述具体编码类型进一步判断方法的过程包括:
1)计算时域信号的4次方后,做快速傅里叶变换,计算频域幅度的平方,并将零频点移到频域的中心,得到xn;
2)设置门限mx=0.5×max(xn),找到大于门限的点,称为大值点;
3)如果大值点多于1个且大值点序号不连续,则判断为是多相编码信号,否则进一步判断;
4)用3dB带宽进行载频粗估计;
5)基于多尺度Haar小波变换估计信号的码元周期;
6)如果码元周期大于1us,提取一个码元周期的数据再次估计载频;
7)把信号变换到基带后,作自相关运算,得到其虚部并规范化;
8)如果自相关虚部全为0,判定为二相编码信号,否则需要进一步判断;
9)对时域信号的平方作快速傅里叶变换,计算频谱的幅度并归一化处理,找到最大值,用xfm表示最大值前后5MHz的频谱。给定系数th,系数计算方法为xfm的均值加xfm的两倍标准差再加0.1,并限制th<0.5,令频谱的最大值乘以系数th为门限,寻找频谱中大于门限的点,如果数量小于3个,判定为二相编码信号,否则判定为四相编码信号。
进一步的,步骤316)中,所述判断和处理频率编码信号或频率编码-二相编码复合调制信号被拆分的过程包括:
1.频率编码信号被拆分成常规信号
1)记录所有识别为常规信号的起止时间和时域信号;
2)按照频率编码信号的每个片段时间不重叠,且相邻两个片段时间间隔在0.15us之内的原则,找到所有可能是频率编码信号的片段组合;
3)把剩下可能是频率编码拆开的常规信号片段拼接起来,按前面叙述的方法计算信噪比,识别调制类型,如果识别结果是频率编码信号,则存储识别结果。
2.频率编码-二相编码复合调制信号被拆分成二相编码信号
1)记录所有识别为二相编码信号的起止时间和时域信号;
2)按照信号的每个片段时间不重叠,且相邻两个片段时间间隔在0.15us之内的原则,找到所有可能是频率编码-二相编码复合调制信号信号的片段组合;
3)分别对每一个片段组合进行分析,把这些片段拼接起来,按前面叙述的方法计算信噪比,识别调制类型,如果识别结果是频率编码-二相编码复合调制信号,则存储识别结果。
所述步骤四中,对识别出调制类型的信号,估计其调制参数的方法具体包括如下步骤:
步骤41)、公有参数估计。时域信号sn作快速傅里叶变换到频域,计算频域幅度的平方,记为xn,找到xn中大于0.5×max(xn)的值所在的位置,记为ff,信号长度记为L,采样频率记为fs,令由如下公式计算中心频率fcenter、带宽B、幅度A和脉宽T。
步骤42)、特有参数估计。
具体的估计方法根据信号类型不同,估计方式不同,具体设置如下:
1.常规信号。
常规信号需要估计的参数有载频fc、初始相位对于fc的估计主要有修正的频率插值算法(M-Rife算法),其原理是:对信号进行离散傅里叶变换后,利用信号频谱的最大谱线与最近邻的较大谱线进行频域插值,从而得到真实频率的精确估计。得到fc后,利用去载频后的相位序列求均值即可得到/>。
2.线性调频信号。
线性调频信号需要估计的参数有起始频率f0、调频斜率k和初始相位。
估计过程包括:
1)线性调频信号序列可表达为:
s(n)=A(n)exp{j(2πf0nΔt+πkn2Δt2+φ0)}
其中A(n)是信号幅度,Δt是采样间隔,k、f0和φ0分别是信号的调频斜率、起始频率和初始相位,定义中心频率为fM=f0+kNΔt/2,则信号序列可改写为
其中φ1=φ0+πfMNΔt-πkN2Δt2/4。
2)将线性调频信号延迟相关后得到正弦波序列:
y(n)=s(n+N/2)·s*(n)=A(n)2exp{j(πNknΔt2+φ2)},n=0,…,N/2
其中φ2=πfMNΔt-πkN2Δt2/4。y(n)是一个载频为kNΔt/2的常规信号
序列,用M-Rife算法对该序列做常规信号频率估计,求得估计频率fk,并计算出调频斜率的估计值k=fk/(NΔt/2);
3)构造去调频斜率序列:
z(n)=exp{-j(πkΔt2(n-N/2)2)}
将s(n)与z(n)相乘得到一个调频斜率很小并且中心频率fM的线性调频序列。由于调频斜率非常小,该序列可看作是一个载频为fM的常规信号,因此用M-Rife算法对这个序列做频率估计就可得到fM的估计值,并得到起始频率的估计值f0=fM-kNΔt/2;
4)计算去载频后的相位,对相位作二阶线性拟合,拟合系数的常数项即初始相位
3.双线性调频信号。
双线性调频信号的待估计的参数有载频fc、调频斜率k、初始相位和调频中点fm。只需取频率尖峰点以前的部分信号,按线性调频信号的估计方法即可得到载频fc、调频斜率k和初始相位/>调频中点可由下式计算:
其中,L是信号点数,fs是采样频率,round()表示四舍五入取整。
4.正弦调频信号。
正弦调频信号待估计的参数有载频fc、调制系数mf、调制频率fm和初始相位
首先用快速傅里叶变换求得信号频谱,找到所有大于最大值0.9倍的频率点,取均值即为信号载频fc;对信号去载频后经过20重相位差分求得瞬时频率,和瞬时频率频谱最大值对应的且小于采样频率一半的频率即为调制频率fm;调制系数mf为2倍的瞬时频率频谱最大幅度与fm的比值。相位去载频后,作一阶线性拟合,拟合系数的常数项即初始相位
5.三阶多项式相位信号。
三阶多项式信号待估计的参数有调频系数ai,i=1,2,3,载频f0及初始相位
下面给出三阶延迟差分:
DP1[s(n),τ]=s(n)
DP2[s(n),τ]=s(n)×s(n-τ)
DP3[s(n),τ]=s(n)×[s(n-τ)]2s(n-2τ)
式中τ为延时长度,0≤τ≤N-1。
m阶离散多项式相位变换DPTm定义为DPm[s(n),τ]的离散傅里叶变换,参数估计过程包括:
1)初始化:设m=3;
2)选择延时τm=N/m,估计am;
3)s(m-1)(n)=s(m)(n)exp{-jam(nΔt)m};
4)m自减1,即m=m-1,如果m≥1,返回步骤2),否则至步骤5);
5)信号自相关运算之后,第一个出现的最大频谱对应的频率就是载频f0;
6)计算的相位即/>
6.相位编码信号
二相编码和四相编码信号待估计调制参数有载频f0、码元周期Tb、编码序列{αj},初始相位多相编码信号主要调制参数有载频f0、码元周期Tb、码元个数N、多相编码种类type,对相位编码信号进行参数估计的过程为:
1)首先利用傅里叶变换得到雷达信号的频谱分布,对频谱进行平滑,将3dB带宽内的功率谱重心作为载频的粗估计;
2)在上一步得到载频f0的估计值基础上,对原始信号去载频后利用Haar小波变换得到相位编码信号的码元周期估计;如果是二相编码、四相编码信号,执行第3)、4)步,如果是多相编码信号,执行第5)、6)步;其中,用Haar小波变换估计相位码元周期过程如下:
a.根据载频的估计值对信号去载频,下变频到基带;
b.估计基带信号的3dB带宽,选取三个尺度因子a1=0.4/Bω,a2=0.6/Bω,a3=0.75/Bω;
c.分别计算基带信号在三个尺度因子下的Haar小波变换幅度的平方,并把三个结果叠加;
d.对叠加后的结果作离散傅里叶变换,并滤除直流分量,估计峰值谱线对应的频率fb;
e.得到相位编码信号的码元周期估计值Tb=1/fb:
3)对基带调制信号还原相位,第一个相位值即为初始相位/>
4)基于瞬时自相关极性判断的方法估计编码序列{αj};
5)码元个数round(·)表示四舍五入取整;
6)利用最大信号检测法识别多相编码种类,识别过程为:
a.当码元个数N不是完全平方数时,信号只能是P3码或P4码;当是奇数时,信号是Frank码、P1码、P3码或P4码中的一种;
当是偶数时,信号可能是多相码任意一种;
b.根据码元个数N,构造对应的本地多相码,将接收到的复信号与构造的多相码共轭相乘;
c.乘积转换到频域并将零频点移到中心;
d.以频谱最大值的0.5倍为门限,找到频谱大于门限的点的序号,序号集合的最后一个减去第一个所得之差,记为Δ;
e.Δ最小者对应的构造的本地多相码即为多相码种类;
f.如果最小值不止一个,则将接收到的复信号转置后与构造的多相码信号矩阵相乘,并计算绝对值,最大绝对值所对应的多相码信号种类即为多相编码种类检测结果。
7.频率编码信号
频率编码信号的待估参数包括码元周期Tc、子码频率fc和编码宽度w,对频率编码信号进行参数估计的过程为:
1)对复信号作20重相位差分,得到瞬时频率pf,并将pf按40的步长做自差分,得到Δpf;
2)对Δpf按3σ准则检测跳变点;
3)根据跳变点位置将复信号分段,利用3dB带宽对每段的载频进行估计;
4)利用STFT变换(短时傅里叶变换)结合小波变换估计频率编码信号的码元周期,具体过程为:
a.对信号做STFT变换,参数设置为窗宽256,每一段的重叠样本数为窗宽的87.5%并四舍五入取整,STFT变换中的傅里叶变换长度为512。STFT的返回值有频率向量F,频谱图计算的时刻点T和能量谱密度P;
b.求出每一个滑动窗内的中心频率。假设每一段能量谱密度的最大值所在位置为tm,那每一个滑动窗内的中心频率fcenter的计算方法为:
c.伸缩因子设为8,对中心频率做离散小波变换,得到序列,并对幅度平方得CWT;变换公式为(下式把中心频率简写为):
yi=yi-1+2fi+3-fi-1-fi+7 i≠1。
d.设检测门限为0.01263fs,对CWT检测跳变幅度超过门限对应的点,对这些点差分,选取第一段连续点的中间值与最后一段连续点的中间值之间的点,与这些点对应的CWT值作希尔伯特变换再作快速傅里叶变换;
e.去除直流分量后,找到快速傅里叶变换结果的最大峰值点,取倒数求得码元周期Tc。
5)跳变点序号作自差分除以采样频率fs得到每段跳频序列的宽度,并将宽度小于0.3Tc的频率值去掉;
6)检测相邻两个频率差值小于0.5的,如果有则合并,它们的算术平均作为合并后的新频率,并得到新的编码宽度w;
7)根据编码宽度修正码元周期,其中n1是频点个数,由编码宽度w除以最小编码宽度之商四舍五入后求和计算所得;
8)将删除合并后的频率值扩充到所有频点上,即得到每个子码频率的估计值。
8.频率编码-二相编码复合调制信号
频率编码-二相编码复合调制信号的待估参数包括频率编码码元周期Tf、频率编码个数Nf、编码频率ff、相位编码码元周期Tp、相位编码码元个数Np、相位编码序列{αj}、初始相位参数估计过程如下:
1)将信号平方,去掉相位编码调制;
2)利用前面介绍的STFT变换结合小波变换估计频率编码信号的码元周期Tf,频率编码个数L是信号长度,fs是采样频率,round(·)表示四舍五入取整;
3)把原信号等分为Nf段,用3dB带宽对每段的载频进行估计得到编码频率ff的估计值;
4)把原信号等分为Nf段,用前面介绍的二相编码信号的参数识别方法对每段信号进行识别;
5)出现次数最多的码元周期值被认为是相位编码码元周期Tp的估计值;
6)对所有码元周期值值为Tp的段落,计算它们所对应的编码序列的长度,选择出现次数最多的长度即为相位编码码元个数Np的估计值,Np所对应的编码序列的第一个作为相位编码序列{αj}的估计值,每段信号估计出的初始相位的平均值是最终初始相位的估计值。
9.线性调频-二相编码复合调制信号
线性调频-二相编码的待估参数包括载频fc、调频斜率k0、相位编码码元周期Tp、相位编码码元个数Np、相位编码序列{αj}、初始相位参数估计过程如下:
1)信号平方,用识别线性调频信号的方法得到的结果除以2得到载频fc、调频斜率k0;
2)信号去调频,用估计二相编码信号参数的方法估计出相位编码信号的码元周期和初始相位,即为相位编码码元周期Tp和初始相位的估计值;
3)用前面介绍的方法对二相编码信号的每个采样点,估计其编码,将编码中连续的1表示为1个1,连续的0表示为1个0,记为向量v,并计算出每段连续出现1或0的个数,表示为l。计算向量 其中fs表示采样频率,round(·)表示四舍五入取整,ceil(·)表示向+∞方向取整;/>
4)将v中的1或0按num中对应的数扩充,得到最终的相位编码序列{αj}的估计值。
此外,本申请还提供一种混叠雷达信号脉内特征分析系统,所述软件系统业务实现包括如下步骤:
步骤一、软件系统技术选型。由于QT是基于C++的二次开发平台,其本身具备良好的封装机制使得QT的模块化程度非常高,可重用性较好。此外,QT提供了一种称为signals/slots的安全类型来替代callback,这使得各个元件之间的协同工作变得十分简单。因此,我们采用最通用的跨平台、开源编程工具QT,支持软件在最流行的两大平台Window和Linux上运行,覆盖了90%以上的用户。在脉内信号图形显示上,我们采用了基于QT的开源图形显示技术QWT,它具备良好的执行效率,可以提供毫秒级的显示响应。
步骤二、软件系统架构设计。本申请在架构设计中遵循了前后端分离、模块化的设计原则。
前后端分离的概念源于解耦,是降低代码复杂度,是开发便于管理和维护。在多年来软件开发实践中有过非常良好的实践。本软件也采用了这个理念,将图形、数据的显示和操作,与后台的算法分析进行了良好的隔离,中间用一套稳定的调用接口来实现连接。
(1)后端算法被封装成动态链接库,依照设计规范提供稳定的调用接口;
(2)前端利用QT中的view概念,基于Qmainwindow、Qwidget等图形类,形成一套与数据解耦的显示界面;
(3)前端通过预先约定的接口调用后端算法获取显示数据,由显示数据驱动显示的效果,但界面的业务逻辑不会受到计算结果的影响。
在设计中,本申请将软件前端分解成了控制整体窗口的主窗口模块、显示图形的图形显示模块、显示参数的树形图模块,显示文件的文件列表模块。基本保障了良好的内聚度,又将耦合度降低到了便于维护的程度。
步骤三、软件系统功能实现。
参见图2,为本发明的混叠雷达信号脉内特征分析方法软件业务流程图。在位置0处为部件为由某实侦设备采集待分析的中频数据文件,在位置1处为所需要的配置参数,由于不同采集设备可能有不用的参数特征,故此提供此配置表,供用户根据实际情况修改。
通过软件打开执行位置0处的数据文件,并读取执行位置1处的配置文件,两者获取成功后,由执行位置2处的信号分离模块进行信号的拆解,形成一个信号数据列表。信号完成分离后,数据进入执行位置3处的进行调制类型识别。
在识别算法对信号列表进行逐一分辨。识别结果为下述列表中的一个:
0 未知
5 脉内无调制
12 判断为常规(NS)信号
1 判断为线性调频(LFM)信号
2 判断为多项式相位(PPS)信号
3 判断为二相编码(BPSK)信号
4 判断为四相编码(QPSK)信号
13 判断为多相编码(MPSK)信号
9 判断为频率编码(FSK)信号
11 判断为双线性(DLFM)调频
21 判断为正弦调频(SFM)信号
81 判断为FSK_BPSK复合调制信号
82 判断为LFM_BPSK复合调制信号
信号调制类型识别完成后,进入执行位置4处进行调制参数估算,根据识别的类型进行参数估算,估算时,根据调制类型的不同有不同的参数,调制类型与参数对照表如下所列:
NS:包含有载频f0(浮点数)、初始相位phi0(浮点数))
LFM:包含有载频f0(浮点数)、初始相位phi0(浮点数)、调频斜率k0(浮点数))
DLFM:包含有载频f0(浮点数)、初始相位phi0(浮点数)、调频斜率k0(浮点数)、调频中点fM(浮点数))
BPSK:包含有载频f0(浮点数)、初始相位phi0(浮点数)、码元周期Tc(浮点数)、编码序列cn(整数,行向量))
QPSK:包含有载频f0(浮点数)、初始相位phi0(浮点数)、码元周期Tc(浮点数)、编码序列cn(整数,行向量))
FSK:包含有码元周期Tc(浮点数)、编码频率fvalue(浮点数,行向量)、编码宽度width(浮点数,行向量))
PPS:包含有载频fc(浮点数)、初始相位phi0(浮点数)、调频系数k(浮点数,列向量)、信号幅度估计b(浮点数))
参数估算完成后,转换为信号波形图、频谱图、时频图和时相图的横纵坐标组。将此坐标组返回到图形显示界面。流程进入执行位置5处进行图表显示模块。
图表显示模块即为本软件的界面部分,在本发明中,根据模块化原则,将图表显示模块分解为4个模块。
(1)主窗口子模块:用于形成窗口主体,提供其它模块运行的容器。
(2)文件列表子模块:用于显示用于分析的数据列表,点击按钮“添加中频文件”,即可打开文件选择窗口,供用户从磁盘存储器添加要分析的文件;点击“清空”按钮则可清除列表中的所有文件。
(3)图形显示子模块:用于显示信号数据在时间、频率坐标上的分布。包括波形图(以到达时间为横坐标,以幅度为纵坐标)、频谱图(以频率为横坐标,以幅度为纵坐标)、时频图(以到达时间为横坐标,以频率为纵坐标)、时相图(以到达时间为横坐标,以相位值为纵坐标)。并为图形提供了缩放、局部数据选择、截屏、重置等操作。
(4)参数显示子模块:以属性表的格式显示中频信号的参数。根据调制类型的不同,有不同的显示列表。本软件将公共参数提出分为一栏,包含中心频率、原始信噪比、信噪比、幅度、脉宽、带宽等;将各个调制类型特有的参数放在另一个栏动态展示。
信号分离过程说明参见图3,信号分离模块读取到中频数据和配置参数后,开启分析步骤,参见图3,步骤202,从中频数据中获取信号频率起始点,然后沿着这段数据获取信号结束点。之后执行步骤203,根据起止时间计算信号的长度,如果持续时间在0.5us之内,则放弃此段数据,不做分析;如果大于0.5us,则进入步骤204,将数据存入一个列表,供信号调制类型识别。
本实例还给出了一种混叠雷达信号脉内特征分析方法,分析仿真数据,该混叠雷达信号脉内特征分析方法包括如下步骤:
步骤一、仿真生成脉宽为10us的线性调频和脉宽为5us的常规信号以2us的时延混叠的信号作为输入信号,信噪比5dB。
线性调频信号的调制参数为起始频率f0=400MHZ、调频斜率k=5MHZ/us和初始相位常规信号的调制参数为载频fc=350MHZ,初始相位/>
步骤二、将输入信号按照本申请的所述方法进行混叠分离,得到两个信号;
步骤三、分别对分离出的两个信号检测起止时间,得到第一个信号的起止时间为0.03us-9.85us,第二个信号的起止时间为2.00us-7.09us;
步骤四、分别对两个信号按照本申请所述方法进行调制类型识别,识别结果为第一个信号是线性调频信号。第二个信号是常规信号;
步骤五、分别对两个信号按照本申请所述方法进行调制参数估计,估计结果为线性调频信号:起始频率f0=400.159MHZ、调频斜率k=4.999MHZ/us和初始相位常规信号:载频fc=350.003MHZ,初始相位/>
本实例处理的是一个混叠信号,使用缺省参数运行,在处理器为AMD Ryzen 74700U的家用笔记本电脑上运行时间为0.275s,混叠信号分离正确,调制类型识别正确,起止时间检测精度高,参数估计结果精度高。
另一实例给出了一种混叠雷达信号脉内特征分析技术,参考图1,分析仿真数据,该混叠雷达信号脉内特征分析方法包括如下步骤:
步骤一、仿真生成脉宽为3us的常规信号和脉宽为6us的二相编码信号以3us的时延混叠的信号作为输入信号,信噪比0dB。
常规信号的调制参数为fc=300MHZ,初始相位二相编码信号的调制参数为载频f0=300MHZ,码元周期Tb=1us、编码序列{1,0,1,0,1,1},初始相位/>
步骤二、将输入信号按照本申请的所述方法进行混叠分离,得到两个信号;
步骤三、分别对分离出的两个信号检测起止时间,得到第一个信号的起止时间为0.00us-3.05us,第二个信号的起止时间为3.91us-9.97us;
步骤四、分别对两个信号按照本申请所述方法进行调制类型识别,识别结果为第一个信号是常规信号。第二个信号是二相编码信号;
步骤五、分别对两个信号按照本申请所述方法进行调制参数估计,估计结果为常规信号:fc=299.992HZ,初始相位二相编码信号的调制参数为载频f0=300.047MHZ,码元周期Tb=1.010us、编码序列{1,0,1,0,1,1},初始相位/>
本实例处理的是一个频谱重叠,时间不重叠的混叠信号,使用缺省参数运行,在处理器为AMD Ryzen 7 4700U的家用笔记本电脑上运行时间为0.164s,混叠信号分离正确,调制类型识别正确,起止时间检测精度高,参数估计结果精度高。
本实例给出了一次中频信号分析,实现了本申请声明的设计流程,完成了对某多相编码数据的分析,信号的载频为691MHZ,码元周期为11.5us,识别结果与图中可视化信号特征吻合。
请注意,以上实施例的各技术特征可以进行任意的组合,为使描述简洁,未对上述实施例中的各个技术特征所有可能的组合都进行描述,然而,只要这些技术特征的组合不存在矛盾,都应当认为是本说明书记载的范围。以上实施例仅表达了本申请的几种实施方式,其描述较为具体和详细,但并不能因此而理解为对发明专利范围的限制。应当指出的是,对于本领域的普通技术人员来说,在不脱离本申请构思的前提下,还可以做出若干变形和改进,这些都属于本申请的保护范围。因此,本申请专利的保护范围应以所附权利要求为准。
Claims (9)
1.一种混叠雷达信号脉内特征分析方法,其特征在于,所述混叠雷达信号脉内特征分析方法包括如下步骤:
步骤一、将时频域交叠的混合信号分离成独立的信号并存储;
步骤二、根据分离后各个独立信号的包络,检测分离后的每一个所述信号的起止时间;
步骤三、识别分离后的每一个所述信号的调制类型;
步骤四、基于获得的每一个所述信号的调制类型,估计分离后的每一个所述信号的调制参数。
2.根据权利要求1所述的混叠雷达信号脉内特征分析方法,其特征在于,步骤一中,所述将时频域交叠的混合信号分离的过程包括:
步骤11)、对原始I/Q通道信号作快速傅里叶变换,得到待处理信号频谱;
步骤12)、找到所述待处理信号频谱幅值最大值所在位置;
步骤13)、判断所述待处理信号频谱幅值最大值所在位置是否为三阶多项式相位信号;如果是,计算所述幅值最大值附近应该保留的频谱宽度,进入步骤16);如果否,进入步骤14);
步骤14)、计算所述待处理信号频谱幅值最大值所在位置的峰度值和偏度值,依据所述峰度值和偏度值判断所述频谱幅值最大值所在位置是信号或噪声;如果是信号,则进入步骤15),否则退出混叠信号分离过程;
步骤15)、计算所述频谱幅值最大值附近应该保留的频谱宽度;
步骤16)、保留所述频谱宽度的信号,并进行存储,将所述待处理信号频谱的其他部分设为噪声,返回到步骤12)。
3.根据权利要求2所述的混叠雷达信号脉内特征分析方法,其特征在于,步骤15)中,计算所述频谱幅值最大值附近应该保留的频谱宽度具体包括:
判断所述频谱幅值最大值信号是否为正弦调频信号,如果是,则计算该正弦调频信号最大频点左右要保留的带宽点数k,左右两边点数设置为相同,确定应该保留的频谱宽度;如果否,计算信号频谱与矩形频谱和三角形频谱的相像系数,根据获得的信号频谱与矩形频谱和三角形频谱的相像系数确定应该保留的频谱宽度。
4.根据权利要求1所述的混叠雷达信号脉内特征分析方法,其特征在于,步骤二中,所述根据信号包络检测分离后的每一个所述信号的起止时间包括:
步骤21:对每一信号的包络作归一化处理,并进行平滑处理;
步骤22:计算平滑后的归一化包络并计算其均值和标准差;
步骤23:按照以下方式确定截取信号时上升沿和下降沿之间的最优比率:如果标准差小于0.1,最优比率取均值和标准差之和,进入步骤24;如果标准差大于或等于0.25,最优比率取0.4,进入步骤24;如果标准差在[0.1,0.25]区间内,则比较均值大小,均值大于0.7,把小值点处信号去掉进行计算,均值小于0.25则把大值点处信号去掉进行计算,均值在[0.25,0.7]区间内,保留全部信号,根据保留信号下标序号计算出信号的起止时间,跳过步骤24;
步骤24:根据确定的最优比率计算信号的起止时间,确定平滑后的归一化信号包络中大于最优比率的点,计算出信号包络的上升沿和下降沿的下标序号,根据下标序号计算出信号的起止时间。
5.根据权利要求1所述的混叠雷达信号脉内特征分析方法,其特征在于,步骤三中,所述识别分离后的每一个信号的调制类型包括:用二阶四阶矩方法估计信噪比;对信号进行短时滤波处理;计算信号去卷叠后的瞬时相位;采用20重相位差分方法计算信号的瞬时频率;对时频曲线作一阶线性拟合,判断每一个信号的调制类型。
6.根据权利要求1所述的混叠雷达信号脉内特征分析方法,其特征在于,步骤四中,所述估计分离后的每一个信号的调制参数包括:把调制参数分为公有参数和特有参数进行估计,公有参数包括中心频率、带宽、幅度和脉宽,特有参数和调制类型有关,不同的调制类型,设置不同的估计方式来获得不同的特有参数。
7.一种混叠雷达信号脉内特征分析系统,其特征在于,所述系统存储有计算机程序指令,通过执行所述计算机程序指令,实现权利要求1~6任意一项所述的混叠雷达信号脉内特征分析方法,能够按缺省参数运行,允许用户修改部分参数,并将运行结果可视化显示。
8.一种电子设备,其特征在于,所述电子设备包括存储器和处理器,存储器存储有计算机程序,处理器执行计算机程序时,实现权利要求1~6任意一项所述的混叠雷达信号脉内特征分析方法。
9.一种计算机可读存储介质,其特征在于,所述计算机可读存储介质上存储有计算机程序,所述计算机程序被处理器执行时,实现权利要求1~6任意一项所述的混叠雷达信号脉内特征分析方法。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202310371898.2A CN116559786A (zh) | 2023-04-10 | 2023-04-10 | 一种混叠雷达信号脉内特征分析方法及系统 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202310371898.2A CN116559786A (zh) | 2023-04-10 | 2023-04-10 | 一种混叠雷达信号脉内特征分析方法及系统 |
Publications (1)
Publication Number | Publication Date |
---|---|
CN116559786A true CN116559786A (zh) | 2023-08-08 |
Family
ID=87490592
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202310371898.2A Pending CN116559786A (zh) | 2023-04-10 | 2023-04-10 | 一种混叠雷达信号脉内特征分析方法及系统 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN116559786A (zh) |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN117111016A (zh) * | 2023-10-20 | 2023-11-24 | 南京航天工业科技有限公司 | 复杂电磁环境下基于信道化的实时脉内分析方法及系统 |
-
2023
- 2023-04-10 CN CN202310371898.2A patent/CN116559786A/zh active Pending
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN117111016A (zh) * | 2023-10-20 | 2023-11-24 | 南京航天工业科技有限公司 | 复杂电磁环境下基于信道化的实时脉内分析方法及系统 |
CN117111016B (zh) * | 2023-10-20 | 2023-12-19 | 南京航天工业科技有限公司 | 复杂电磁环境下基于信道化的实时脉内分析方法及系统 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Yang et al. | Multicomponent signal analysis based on polynomial chirplet transform | |
JP6026531B2 (ja) | レーダー用デジタル受信機を用いるレーダーパルス検出 | |
JP5965587B2 (ja) | ジッタ及びノイズ分析方法並びに試験測定機器 | |
US6269135B1 (en) | Digital phase discriminations based on frequency sampling | |
EP1267172B1 (en) | Apparatus and method for spectrum analysis-based serial data jitter measurement | |
JP4774543B2 (ja) | ジッタ推定装置及び推定方法 | |
US6594595B2 (en) | Apparatus for and method of measuring cross-correlation coefficient between signals | |
CN116559786A (zh) | 一种混叠雷达信号脉内特征分析方法及系统 | |
US20020176525A1 (en) | Apparatus for and method of measuring clock skew | |
JP3271504B2 (ja) | 周波数推定回路およびそれを用いたafc回路 | |
US20050185708A1 (en) | Apparatus for measuring jitter, method of measuring jitter and computer-readable medium storing a program thereof | |
Grajal et al. | Real time FPGA implementation of an automatic modulation classifier for electronic warfare applications | |
KR101280512B1 (ko) | 실시간 pmop/fmop 신호 처리를 위한 디지털 수신기 및 그 방법 | |
CN110708267B (zh) | 频偏信息估计值确定方法 | |
Zeng et al. | An approach to intra-pulse modulation recognition based on the ambiguity function | |
US8433739B2 (en) | Methods and systems for detecting repetitive synchronized signal events | |
CN115529213B (zh) | 一种用于分离lfm脉冲重叠信号的方法及装置 | |
CN115902528B (zh) | 一种直流牵引网振荡与短路故障辨识方法 | |
US8438203B2 (en) | Methods and systems for processing and displaying data | |
Paisana et al. | An alternative implementation of a cyclostationary detector | |
Djukanović et al. | A parametric method for non-stationary interference suppression in direct sequence spread-spectrum systems | |
Chan et al. | Evaluation of various FFT methods for single tone detection and frequency estimation | |
CN112116917B (zh) | 基于相位跃变度的电抗器本体与风机声信号分离方法 | |
KR100911689B1 (ko) | 실시간 음악 비트 주기 추출 방법 및 실시간 음악 비트주기 추출 장치 | |
RU2365052C1 (ru) | Адаптивный классификатор сложных широкополосных импульсных сигналов |
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 |