CN116484198A - 无人机遥控图传信号的参数盲估计方法 - Google Patents
无人机遥控图传信号的参数盲估计方法 Download PDFInfo
- Publication number
- CN116484198A CN116484198A CN202310416626.XA CN202310416626A CN116484198A CN 116484198 A CN116484198 A CN 116484198A CN 202310416626 A CN202310416626 A CN 202310416626A CN 116484198 A CN116484198 A CN 116484198A
- Authority
- CN
- China
- Prior art keywords
- time
- frequency
- remote control
- matrix
- frequency matrix
- 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
- 238000000034 method Methods 0.000 title claims abstract description 72
- 230000005540 biological transmission Effects 0.000 title claims abstract description 16
- 239000011159 matrix material Substances 0.000 claims abstract description 112
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 37
- 230000009466 transformation Effects 0.000 claims abstract description 13
- 238000003064 k means clustering Methods 0.000 claims abstract description 9
- 230000002068 genetic effect Effects 0.000 claims description 10
- 238000005192 partition Methods 0.000 claims description 6
- 238000012804 iterative process Methods 0.000 claims description 3
- 238000001514 detection method Methods 0.000 description 9
- 238000010586 diagram Methods 0.000 description 6
- 238000007430 reference method Methods 0.000 description 5
- 238000004088 simulation Methods 0.000 description 5
- 238000004364 calculation method Methods 0.000 description 3
- 238000002474 experimental method Methods 0.000 description 3
- 238000001228 spectrum Methods 0.000 description 3
- 230000007423 decrease Effects 0.000 description 2
- 238000006073 displacement reaction Methods 0.000 description 2
- 238000000605 extraction Methods 0.000 description 2
- 239000011521 glass Substances 0.000 description 2
- 230000010354 integration Effects 0.000 description 2
- 238000012360 testing method Methods 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 239000006185 dispersion Substances 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000002401 inhibitory effect Effects 0.000 description 1
- 230000002452 interceptive effect Effects 0.000 description 1
- 230000035772 mutation Effects 0.000 description 1
- 230000000737 periodic effect Effects 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 238000012552 review Methods 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 230000008054 signal transmission Effects 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
- G06F18/211—Selection of the most significant subset of features
- G06F18/2111—Selection of the most significant subset of features by using evolutionary computational techniques, e.g. genetic algorithms
-
- 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/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/23—Clustering techniques
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
- G06N3/00—Computing arrangements based on biological models
- G06N3/12—Computing arrangements based on biological models using genetic models
- G06N3/126—Evolutionary algorithms, e.g. genetic algorithms or genetic programming
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02D—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN INFORMATION AND COMMUNICATION TECHNOLOGIES [ICT], I.E. INFORMATION AND COMMUNICATION TECHNOLOGIES AIMING AT THE REDUCTION OF THEIR OWN ENERGY USE
- Y02D30/00—Reducing energy consumption in communication networks
- Y02D30/70—Reducing energy consumption in communication networks in wireless communication networks
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Data Mining & Analysis (AREA)
- General Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Mathematical Physics (AREA)
- Health & Medical Sciences (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- General Engineering & Computer Science (AREA)
- Evolutionary Biology (AREA)
- Artificial Intelligence (AREA)
- Evolutionary Computation (AREA)
- Biophysics (AREA)
- Computational Mathematics (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Physiology (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Computing Systems (AREA)
- Pure & Applied Mathematics (AREA)
- General Health & Medical Sciences (AREA)
- Software Systems (AREA)
- Databases & Information Systems (AREA)
- Algebra (AREA)
- Genetics & Genomics (AREA)
- Biomedical Technology (AREA)
- Computational Linguistics (AREA)
- Molecular Biology (AREA)
- Selective Calling Equipment (AREA)
Abstract
本发明提供了一种无人机遥控图传信号的参数盲估计方法，通过计算接收信号的短时傅里叶变换得到时频矩阵；对时频矩阵进行去噪处理，基于驻留时时长的k‑means算法或基于时频方差的k‑means聚类算法对去噪之后的时频矩阵进行聚类以去除定频干扰，并重构得到重构时频矩阵；提取重构时频矩阵的时频脊线；对时频脊线作Haar小波变换，并取Haar小波变换后的模值；用逐差法求模值各峰之间的时间间隔以及模值各峰对应的时间；根据时间间隔以及模值各峰对应的时间，估计遥控信号的起跳时间；基于起跳时间，估计遥控信号的实际频率。相比于现有技术，本发明可以去除掉定频干扰，并实现对遥控信号的参数进行估计。
Description
技术领域
本发明属于无人机技术领域，具体涉及一种无人机遥控图传信号的参数盲估计方法。
背景技术
随着无人机的普遍，对无人机违规行为处理的前提是探测到违规的具体无人机。目前已有无人机探测手段包括光电探测、雷达探测以及无线电频谱探测。为了提高无人机的适用性，商业无人机公司一直在小型化、操纵性等方面提高技术。
现有无人机尺寸越来越小机动性越来越强，采用探测低小慢这类无人机的方案，如采用光电探测或者雷达探测，识别率将大打折扣。目前市场上主流的消费级无人机都会携带运动相机以及相应的图像信号传输系统，在无人机飞行过程中将一直向遥控地面站发射信号，因此可以根据无人机图传信号对其进行检测与识别，跳频图案反映了无人机遥控图传信号的重要时频特征，根据跳频图案可以得到无人机遥控图传信号的跳频频率、跳速、跳频周期等重要参数，为防御方对违规的无人机实施反制提供了重要信息。因此技术人员通过对无人机遥控信号的参数进行估计可以实现违规无人机的检测。
现有技术中，无人机遥控信号参数估计方法是：首先对接收信号进行短时间傅里叶变换(Short-Time Fourier Transform,STFT)，将接收信号变换至时频域，然后提取时频分布的时频脊线。遥控信号的时频脊线表现为时频面上的不断发生突变的折线段，横向表示时间维度，纵向表示频率维度，通过对时频脊线的分析可以获得信号的瞬时频率，进而得到其随时间的变化规律，进而得以估计遥控信号的各项参数。
现有技术方案中由于定频干扰的存在，时频脊线会集中分布在该定频干扰所在的频率附近，破坏原本的时频脊线分布特征，当干扰功率进一步增强时，其时频脊线近似为一条直线，从而使得基于时频脊线的参数估计方法失效。
发明内容
为了解决现有技术中存在的上述问题，本发明提供了一种无人机遥控图传信号的参数盲估计方法。本发明要解决的技术问题通过以下技术方案实现：
本发明提供了一种无人机遥控图传信号的参数盲估计方法包括：
步骤1，接收信号，并计算所述接收信号的短时傅里叶变换得到时频矩阵；
步骤2，对所述时频矩阵进行去噪处理，得到去噪之后的时频矩阵；
步骤3，基于驻留时时长的k-means算法或基于时频方差的k-means聚类算法对去噪之后的时频矩阵进行聚类以去除定频干扰，并重构得到重构时频矩阵；
步骤4，提取所述重构时频矩阵的时频脊线；
步骤5，对所述时频脊线作Haar小波变换，并取Haar小波变换后的模值；
步骤6，用逐差法求模值各峰之间的时间间隔以及模值各峰对应的时间；
步骤7，根据所述时间间隔以及模值各峰对应的时间，估计遥控信号的起跳时间；
步骤8，基于所述起跳时间，估计遥控信号的实际频率。
本发明的有益效果：
本发明提供了一种无人机遥控图传信号的参数盲估计方法，通过计算接收信号的短时傅里叶变换得到时频矩阵；对时频矩阵进行去噪处理，基于驻留时时长的k-means算法或基于时频方差的k-means聚类算法对去噪之后的时频矩阵进行聚类以去除定频干扰，并重构得到重构时频矩阵；提取重构时频矩阵的时频脊线；对时频脊线作Haar小波变换，并取Haar小波变换后的模值；用逐差法求模值各峰之间的时间间隔以及模值各峰对应的时间；根据时间间隔以及模值各峰对应的时间，估计遥控信号的起跳时间；基于起跳时间，估计遥控信号的实际频率。相比于现有技术，本发明可以去除掉定频干扰，并实现对遥控信号的参数进行估计。
以下将结合附图及实施例对本发明做进一步详细说明。
附图说明
图1是本发明提供的一种无人机遥控图传信号的参数盲估计方法的示意图；
图2a是特定频率处残余噪声、遥控信号以及定频干扰在时间维上分布情况图；
图2b是本发明的基于驻留时长的聚类结果；
图2c是本发明的基于时频方差的聚类结果；
图3a是本发明中改进时频脊线的遥控信号参数盲估计方法重构的时频图及时频脊线图；
图3b是本发明中基于时频方差聚类的遥控信号参数盲估计方法重构的时频图及时频脊线图；
图4a是本发明的跳周期的估计结果图；
图4b是本发明的跳频频率的估计结果图；
图5a是本发明的不同干扰功率下周期估计性能比较图；
图5b是本发明的不同干扰功率下频率估计性能比较图。
具体实施方式
下面结合具体实施例对本发明做进一步详细的描述，但本发明的实施方式不限于此。
本发明提供了两种方案，分别为改进时频脊线的遥控信号参数盲估计方法和基于时频方差聚类的遥控信号参数盲估计方法。
1、改进时频脊线的遥控信号参数盲估计方法的技术原理是：首先在接收到信号后，计算接收信号的STFT，得到对应的时频矩阵；接着进行时频矩阵的预处理，得到去噪后的时频矩阵；对去噪后的矩阵按照基于驻留时时长的k-means算法进行聚类，去除定频干扰，对时频矩阵进行重构，得到重构时频矩阵；之后提取重构时频矩阵的时频脊线，并对时频脊线作Haar小波变换，取其模值|W(a,t)|；用逐差法求|W(a,t)|各峰之间的时间间隔为Th，其倒数即为遥控信号的跳速；求|W(a,t)|各峰对应的时间，估计遥控信号的起跳时间；取各跳周期的中点，估计遥控信号的频率。
2、基于时频方差聚类的遥控信号参数盲估计方法的技术原理是：首先在接收到信号后，计算接收信号的STFT，得到对应的时频矩阵；之后利用遗传算法求解最佳的截断门限，得到新的时频矩阵；计算时频矩阵的时频方差S，通过k-means算法对进行聚类(聚类数目为2)，去除定频干扰并重构得到重构时频矩阵；提取重构时频矩阵的时频脊线，利用Haar小波对其进行检测；按公式估计出跳频信号的周期和跳速；对重构时频矩阵进行ISTFT，并转换到频域求出遥控信号的频率，完成遥控信号的参数的估计。
如图1所示，本发明提供了一种无人机遥控图传信号的参数盲估计方法包括：
步骤1，接收信号，并计算所述接收信号的短时傅里叶变换得到时频矩阵；
步骤2，对所述时频矩阵进行去噪处理，得到去噪之后的时频矩阵；
步骤3，基于驻留时时长的k-means算法或基于时频方差的k-means聚类算法对去噪之后的时频矩阵进行聚类以去除定频干扰，并重构得到重构时频矩阵；
步骤4，提取所述重构时频矩阵的时频脊线；
步骤5，对所述时频脊线作Haar小波变换，并取Haar小波变换后的模值；
步骤6，用逐差法求模值各峰之间的时间间隔以及模值各峰对应的时间；
步骤7，根据所述时间间隔以及模值各峰对应的时间，估计遥控信号的起跳时间；
步骤8，基于所述起跳时间，估计遥控信号的实际频率。
下面对两种方案的不同以及详细过程做详细介绍。
(1)改进时频脊线的遥控信号参数盲估计方法，包括以下步骤：
S1，计算信号的STFT，得到对应的时频矩阵STFTx(,f)；
S2，进行时频矩阵的去噪处理，得到去噪后的时频矩阵STFT(,f)；
在信噪比较高的情况下，利用遗传算法对时频矩阵中的噪声元素进行去除；在信噪比较低的情况下，利用基于迭代法的时频矩阵去噪法对时频矩阵中的噪声元素进行去除。
本发明可以利用遗传算法求解最佳的截断门限ε，并利用截断门限ε对所述时频矩阵STFTx(t,f)进行去噪，得到去噪之后的时频矩阵STFT(t,f)。
利用遗传算法求解最佳的截断门限ε，并利用截断门限ε对所述时频矩阵STFTx(t,f)进行去噪，得到去噪之后的时频矩阵STFT(t,f)包括：
步骤2.1，取时频矩阵STFTx(t,f)的最大值M和最小值m，按照式ε1＝(M+m)/2求出截断门限ε的初始阈值ε1；
步骤2.2，以初始值ε1为界将时频矩阵STFTx(t,f)分为两个部分Xs和Xn，Xs表示信号包含定频干扰的时频成分，Xn表示噪声的时频成分；
步骤2.3，统计Xs和Xn各部分的数目Ns和Nn，然后计算Xs和Xn这两部分的平均值，分别表示为Xs1和Xn1；
步骤2.4，将Xs1和Xn1这两部分作为新的信号成分，求Xs1和Xn1的平均值，从而得到新的阈值ε2，使用新的阈值ε2替换所述初始阈值ε1；
步骤2.5，重复步骤2.2至步骤2.5不断更新得到新的阈值，直至εk+1＝εk时，结束迭代过程，得到最终的阈值ε＝εk+1；
步骤2.6，将得到最终阈值作为截断门限ε，将时频矩阵中低于截断门限ε的元素去除，得到去噪后的时频矩阵STFT(t,f)。
S3，对去噪后的时频矩阵STFT(t,f)按照基于驻留时时长的k-means算法进行聚类，去除定频干扰，以对时频矩阵进行重构，得到重构时频矩阵STFTR(t,f)；
如图2a所示，信号在时频图中的驻留时间长度表征了信号的驻留时间，对于噪声来说，其在时间维度上离散分布，驻留时间非常短；遥控信号是一个跳频信号，其在时频面上是一段一段存在的，驻留时间就是一个跳频周期；定频干扰在整个观测时间内都存在，是一条横贯时间轴的长线段。因此，可以首先对时频矩阵各频率分量进行驻留时间的统计，对于一个M×N的时频矩阵，首先要判断矩阵中的元素是否为0，若为0则判断下一个元素，若不为0则记录驻留时间长度为1，继续判断下一个元素，若为0则返回驻留时长为1，若不为0则驻留时长加1，直至输出驻留时长。根据驻留时间的不同，利用k-means聚类算法将噪声、遥控信号与定频干扰分开。
基于驻留时长的聚类算法以去除定频干扰，并重构得到重构时频矩阵的具体过程如下：
步骤3.1a，统计去噪后的时频矩阵STFT(t,f)中各频率分量的驻留时间，设为{len1,len2,len3....lenn}，随机选取3个中心点{μ1,μ2,μ3}将数据分为3个簇划分{λ1,λ2,λ3}；
步骤3.2a，计算所有数据距离聚类中心点的欧氏距离dij＝|leni-μj|，并将leni分入对应的类别λi；
步骤3.3a，计算每个簇划分λi中所有样本点的平均值作为新的中心点，重复步骤3.1a、3.2a的过程，直至所有中心点都不再发生变化，完成聚类，如图2b所示。
步骤3.4a，将聚类为遥控信号的频率分量保留，将聚类为定频干扰与残余噪声的频率分量置为0，更新原来的时频矩阵得到重构时频矩阵STFTR(t,f)，如图3a所示。
S4，提取重构时频矩阵的时频脊线如图3a所示。
对于STFT来说，其时频脊线是指它的时频分布中每一时刻的峰值频率，即X(t,f)表示信号的时频分布。
步骤5，对fi(t)作Haar小波变换，并取其模值；
相邻的频率分量之间存在着跳变，由此产生了一定的类似于方波的阶跃现象。由于小波变换具有类似于“放大镜”的作用，能够对信号细节部分进行一定程度的放大，因此选择合适的小波函数能够对信号中的一些奇异点进行检测。小波函数有很多，针对不同的信号特征可以选用不同的小波函数，由于时频脊线近似为方波，而Haar小波对这类跳变具有很好的检测能力，可以在局部检测出跳变的时刻点，因此可以利用Haar小波检测出跳变时刻，进而精确估计出信号的跳频周期等参数。小波变换的公式为
式中，a表示尺度参数(通常取a＝1)，fi(t)表示提取的时频脊线，*表示共轭，τ表示位移，ψ(t)表示Haar小波的母函数，其表达式为
若在同一周期内，由于遥控信号的频率不发生改变，设此时频率为fk，则对应的小波变换结果为：
若在不同周期内，由于积分区间跨越不同的跳周期，遥控信号的频率发生了突变，假设在τ0处，频率由fk变为fk+1。
当时，小波变换的结果为：
当时，小波变换的结果为：
由上述计算可知，在同一周期内，|W(a,τ)|恒为0，而在不同的跳周期内，|W(a,τ)|的值不为0，当处在频率发生跳变的位置时，取得最大值。随着τ在整个区间的滑动，|W(a,τ)|会在0与峰值之间不断变化。|W(a,τ)|的变化反映了遥控信号频率的变化，因为在同一个跳频周期内，遥控信号的频率不发生变化，只有在不同跳频周期内，遥控信号的频率才会不同，因此两个突变时刻的间隔也就是信号的跳频周期。
S6，用逐差法求|W(a,t)|各峰之间的时间间隔为Th，其倒数即为遥控信号的跳速；
各谱峰的时间间隔就是该跳频信号跳周期估计Th，其倒数就是跳速，为了获得更为准确的周期估计值，采用奇数项和偶数项对消的逐差法计算，即
其中，表示峰值数目，tpp(k)表示第k个峰值出现的时间。
步骤7，求|W(a,t)|各峰对应的时间，估计遥控信号的起跳时间；
忽略接收机的时延，第一个峰值(第二个跳频)所对应的时间即为遥控信号的第一跳(起跳)时间。由于噪声的存在，时频脊线的提取会有微小的误差，因此各个峰值出现的时间并不是完全准确的估计。可以根据多个峰值取时间平均以抵消噪声带来的误差，即
S8，取各跳周期的中点，估计遥控信号的频率。
对于遥控信号频率的估计，可以根据提取到的时频脊线直接估计得到，考虑到噪声的影响，通常会取各跳周期的中间点对应的频率作为实际频率的估计值，即
(2)基于时频方差聚类的遥控信号参数盲估计方法包括以下步骤：
S1，计算接收信号的STFT，得到对应的时频矩阵STFTx(t,f)。
S2，利用遗传算法求解最佳的截断门限ε，得到时频矩阵STFT(t,f)。
本发明可以利用遗传算法求解最佳的截断门限ε，并利用截断门限ε对所述时频矩阵STFTx(t,f)进行去噪，得到去噪之后的时频矩阵STFT(t,f)。
利用遗传算法求解最佳的截断门限ε，并利用截断门限ε对所述时频矩阵STFTx(t,f)进行去噪，得到去噪之后的时频矩阵STFT(t,f)包括：
步骤2.1，取时频矩阵STFTx(t,f)的最大值M和最小值m，按照式ε1＝(M+m)/2求出截断门限ε的初始阈值ε1；
步骤2.2，以初始值ε1为界将时频矩阵STFTx(t,f)分为两个部分Xs和Xn，Xs表示信号包含定频干扰的时频成分，Xn表示噪声的时频成分；
步骤2.3，统计Xs和Xn各部分的数目Ns和Nn，然后计算Xs和Xn这两部分的平均值，分别表示为Xs1和Xn1；
步骤2.4，将Xs1和Xn1这两部分作为新的信号成分，求Xs1和Xn1的平均值，从而得到新的阈值ε2，使用新的阈值ε2替换所述初始阈值ε1；
步骤2.5，重复步骤2.2至步骤2.5不断更新得到新的阈值，直至εk+1＝k时，结束迭代过程，得到最终的阈值ε＝εk+1；
步骤2.6，将得到最终阈值作为截断门限ε，将时频矩阵中低于截断门限ε的元素去除，得到去噪后的时频矩阵STFT(,f)。
S3，计算STFT(,f)的时频方差S，通过k-means算法对S进行聚类(聚类数目为2)，去除定频干扰并重构得到矩阵STFTR(,f)。
假设遥控信号去噪之后如下式所示：
定义时频矩阵STFT(,f)的时频方差为S，则有
其中，j＝1,2,…,n。进而可以得到全部的时频方差的表达式：
其中，J(,f)表示干扰信号。
时频方差表征了在某一频率上，信号沿着时间维分布的离散程度。对于跳频信号来说，其分布集中在在固定的时频域内，时频方差较大；而对于定频干扰和噪声来说，由于其较均匀的分布于整个时间维上，因而时频方差较小，借助于k-means聚类算法，可以将两者很好的区分出来。
基于时频方差的k-means聚类算法对去噪之后的时频矩阵进行聚类以去除定频干扰，并重构得到重构时频矩阵具体过程为：
步骤3.1b，计算去噪之后的时频矩阵STFT(t,f)中每一行的时频方差si(i∈[1,N])，时频方差的集合为设为{s1,s2,s3....sN}，随机选取2个中心点{α1,α2}将数据分为2个簇划分{β1,β2}；
步骤3.2b，计算所有据距离聚类中心点的欧氏距离Dij＝|si-αj|，(j∈1,2)，并将si分入对应的类别βj；
步骤3.3b，计算每个簇划分βj中所有样本点的平均值作为新的中心点；
步骤3.4b，重复步骤3.2b、3.3b的过程，直至所有中心点都不再发生变化，完成聚类，如图2c所示。
步骤3.5b，将聚类为遥控信号的频率分量保留，将聚类为定频干扰与残余噪声的频率分量置为0，更新原来的时频矩阵得到重构时频矩阵STFTR(t,f)，如图3b所示。
S4，提取时频脊线如图3b所示。
对于STFT来说，其时频脊线是指它的时频分布中每一时刻的峰值频率，即其中，X(t,f)表示信号的时频分布。
S5，对fi(t)作Haar小波变换，并取其模值|W(a,t)|；
相邻的频率分量之间存在着跳变，由此产生了一定的类似于方波的阶跃现象。由于小波变换具有类似于“放大镜”的作用，能够对信号细节部分进行一定程度的放大，因此选择合适的小波函数能够对信号中的一些奇异点进行检测。小波函数有很多，针对不同的信号特征可以选用不同的小波函数，由于时频脊线近似为方波，而Haar小波对这类跳变具有很好的检测能力，可以在局部检测出跳变的时刻点，因此可以利用Haar小波检测出跳变时刻，进而精确估计出信号的跳频周期等参数。小波变换的公式为
式中，a表示尺度参数(通常取a＝1)，fi()表示提取的时频脊线，*表示共轭，τ表示位移，ψ(t)表示Haar小波的母函数，其表达式为：
若在同一周期内，由于遥控信号的频率不发生改变，设此时频率为fk，则对应的小波变换结果为：
若在不同周期内，由于积分区间跨越不同的跳周期，遥控信号的频率发生了突变，假设在τ0处，频率由fk变为fk+1。
当时，小波变换的结果为：
当时，小波变换的结果为：
由上述计算可知，在同一周期内，|W(a,τ)|恒为0，而在不同的跳周期内，|W(a,τ)|的值不为0，当处在频率发生跳变的位置时，取得最大值。随着τ在整个区间的滑动，|W(a,τ)|会在0与峰值之间不断变化。|W(a,τ)|的变化反映了遥控信号频率的变化，因为在同一个跳频周期内，遥控信号的频率不发生变化，只有在不同跳频周期内，遥控信号的频率才会不同，因此两个突变时刻的间隔也就是信号的跳频周期。
S6，用逐差法求|W(a,t)|各峰之间的时间间隔为Th，其倒数即为遥控信号的跳速；
各谱峰的时间间隔就是该跳频信号跳周期估计Th，其倒数就是跳速，为了获得更为准确的周期估计值，采用奇数项和偶数项对消的逐差法计算，即
其中，表示峰值数目，tpp()表示第k个峰值出现的时间。
S7，求|W(a,t)|各峰对应的时间，估计遥控信号的起跳时间；
忽略接收机的时延，第一个峰值(第二个跳频)所对应的时间即为遥控信号的第一跳(起跳)时间。由于噪声的存在，时频脊线的提取会有微小的误差，因此各个峰值出现的时间并不是完全准确的估计。可以根据多个峰值取时间平均以抵消噪声带来的误差，即
S8，截取某一跳的重构时频矩阵STFTR(,f)，并进行ISTFT求出遥控信号的频率，完成遥控信号的参数的估计。
ISTFT的公式如下，
本发明根据时间间隔Th确定跳频周期。截取某一跳的重构时频矩阵，对某一跳的重构时频矩阵求短时傅里叶逆变换，得到遥控信号的实际频率。
在已经求得跳频周期的基础上，截取某一跳信号重构的时频矩阵，然后对重构的时频矩阵求短时傅里叶逆变换(inverse short time Fourier transform,ISTFT)，将其转化到频域上，通过快速傅里叶变换得到其频率值。因为快速傅里叶变换的分辨率可以按照需求进行调整，因此采用这种方法得到的频率的估计值比利用时频脊线进行频率估计的精度更高。
下面通过仿真实验对本发明的效果做进一步说明。
1.仿真条件
实验数据取经过下变频后的中频信号，遥控信号的一化频率集为[0.025,0.05,0.075,0.10,0.125,0.15,0.175,0.20]，跳周期为5ms，定频干扰的归一化频率为0.25，采样点数为40000，STFT采用长度为507的Hamming窗。定义信噪比为信号与噪声功率之比的对数形式，干信比(interference-to-signal ratio，ISR)为干扰与信号功率之比的对数形式，周期估计相对误差为 频率估计相对误差为/>估计结果的归一化均方误差(NRMSE)定义为如下形式，
其中，表示参数的第i次估计结果，Ns为仿真次数。
2.仿真内容
实验1定频干扰功率设为与信号功率相同，在信噪比[-15,5]dB范围内，分别利用三种算法对每组数据进行100次Monte Carlo试验，得到跳周期的归一化均方误差和跳频频率的归一化均方误差，仿真结果如图4a和图4b所示。其中，参考方法1表示参数估计方法，参考方法2表示二次差分法，本发明方法1表示改进时频脊线的参数估计方法，本发明方法2表示基于时频方差聚类的参数估计方法。
实验2考虑不同的定频干扰功率对估计结果的影响，定义估计相对误差E1、小于1％时为一次正确估计。在[-15,5]dB信噪比范围内，分别仿真了干信比在[0，10]dB时的估计情况。利用方法3对每组数据进行100次Monte Carlo试验，得到各个ISR条件下的正确估计概率，结果如图5a和图5b所示。
1.仿真结果分析
由图4a可知，随着信噪比的提高，参考方法1对周期的估计NEMSE逐渐下降，但整体上方法1的NEMSE比较大，这说明当定频干扰存在时，参考方法1无法对周期进行正确估计；参考方法2虽然可以对周期进行估计，但由于二次差分毛刺较多，对噪声仍旧比较敏感，本发明方法1和本发明方法2这两种算法的NEMSE随着信噪比的增大逐渐下降，并且这两种算法对噪声有较好的抑制效果，在信噪比较低时仍具有较高的估计精度，当信噪比大于-5dB时，NRMSE小于0.1，当信噪比大于0dB时，NRMSE迅速收敛至0，当信噪比低于-8dB时，由于所引入的噪声过大，无法很好的提取到有效的时频区间，导致后续基于驻留时长或基于时频方差的聚类方法得到的结果发生恶化，从而使NRMSE变大，整体上看，本发明方法2相比于本发明方法1精度略高，收敛速度略快，但其运算量较大，基于驻留时长的聚类过程需要更多的时间完成。
由图4b可知，对遥控信号的频率的估计结果与对周期的估计结果类似。参考方法1、参考方法2和本发明方法1的NRMSE最终收敛至0.1。这是因为STFT无法使时域和频域的分辨率同时达到最大，其时域分辨率与频域分辨的乘积恒大于一个不为零的数。因此，利用上述方法得到的频率估计值的NRMSE存在不为零的下限。而本发明方法2选择在频域对跳频频率进行估计，不存在分辨率不够的问题，当信噪比大于-8dB时，NRMSE能够迅速收敛至0。
通过分析由图4a和图4b可知，方法3在运算量、处理时间及对跳频频率的估计误差方面的性能优于方法2，因此在实验2中，本发明重点分析方法3的抗干扰性能。
由图5a和图5b可知，随着干扰功率的增强，需要在更高的信噪比条件下才能得到跳周期的正确估计结果，跳周期的正确估计概率曲线发生右移。当ISR大于6dB时，由于干扰的功率过大，遗传算法不能正确寻找到合适阈值，无法提取有效的时频区间，导致算法发生恶化。当ISR等于0dB时，算法在信噪比高于-5dB时快速收敛。当ISR小于6dB时，在信噪比低于5dB的情况下，算法正确估计概率可达到90％以上。对跳频频率的估计结果与上述结果类似，当干扰功率过大时，由于信号较弱而被当作噪声去除，使得ISTFT得到的时域信号发生了畸变和丢失，对跳频频率的估计发生错误。
综上所述，基于时频方差聚类的参数估计算法在ISR小于5dB时能够同时消除噪声和干扰对估计结果带来的影响，且具有较高的估计精度和较强的抗干扰能力。
此外，术语“第一”、“第二”仅用于描述目的，而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此，限定有“第一”、“第二”的特征可以明示或者隐含地包括一个或者更多个该特征。在本发明的描述中，“多个”的含义是两个或两个以上，除非另有明确具体的限定。
尽管在此结合各实施例对本申请进行了描述，然而，在实施所要求保护的本申请过程中，本领域技术人员通过查看所述附图、公开内容、以及所附权利要求书，可理解并实现所述公开实施例的其他变化。在权利要求中，“包括”(comprising)一词不排除其他组成部分或步骤，“一”或“一个”不排除多个的情况。
以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明，不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干简单推演或替换，都应当视为属于本发明的保护范围。
Claims (10)
1.一种无人机遥控图传信号的参数盲估计方法，其特征在于，包括：
步骤1，接收信号，并计算所述接收信号的短时傅里叶变换得到时频矩阵；
步骤2，对所述时频矩阵进行去噪处理，得到去噪之后的时频矩阵；
步骤3，基于驻留时时长的k-means算法或基于时频方差的k-means聚类算法对去噪之后的时频矩阵进行聚类以去除定频干扰，并重构得到重构时频矩阵；
步骤4，提取所述重构时频矩阵的时频脊线；
步骤5，对所述时频脊线作Haar小波变换，并取Haar小波变换后的模值；
步骤6，用逐差法求模值各峰之间的时间间隔以及模值各峰对应的时间；
步骤7，根据所述时间间隔以及模值各峰对应的时间，估计遥控信号的起跳时间；
步骤8，基于所述起跳时间，估计遥控信号的实际频率。
2.根据权利要求1所述的无人机遥控图传信号的参数盲估计方法，其特征在于，步骤2包括：
利用遗传算法求解最佳的截断门限ε，并利用截断门限ε对所述时频矩阵STFTx(t，f)进行去噪，得到去噪之后的时频矩阵STFT(t，f)。
3.根据权利要求2所述的无人机遥控图传信号的参数盲估计方法，其特征在于，所述利用遗传算法求解最佳的截断门限ε，并利用截断门限ε对所述时频矩阵STFTx(t，f)进行去噪，得到去噪之后的时频矩阵STFT(t，f)包括：
步骤2.1，取时频矩阵STFTx(t，f)的最大值M和最小值m，按照式ε1＝(M+m)/2求出截断门限ε的初始阈值ε1；
步骤2.2，以初始值ε1为界将时频矩阵STFTx(t，f)分为两个部分Xs和Xn，Xs表示信号包含定频干扰的时频成分，Xn表示噪声的时频成分；
步骤2.3，统计Xs和Xn各部分的数目Ns和Nn，然后计算Xs和Xn这两部分的平均值，分别表示为Xs1和Xn1；
步骤2.4，将Xs1和Xn1这两部分作为新的信号成分，求Xs1和Xn1的平均值，从而得到新的阈值ε2，使用新的阈值ε2替换所述初始阈值ε1；
步骤2.5，重复步骤2.2至步骤2.5不断更新得到新的阈值，直至εk+1＝εk时，结束迭代过程，得到最终的阈值ε＝εk+1；
步骤2.6，将得到最终阈值作为截断门限ε，将时频矩阵中低于截断门限ε的元素去除，得到去噪后的时频矩阵STFT(t，f)。
4.根据权利要求3所述的无人机遥控图传信号的参数盲估计方法，其特征在于，步骤3中基于驻留时时长的k-means算法对去噪之后的时频矩阵进行聚类以去除定频干扰，并重构得到重构时频矩阵包括：
步骤3.1a，统计去噪后的时频矩阵STFT(t，f)中各频率分量的驻留时间，设为{len1，len2，len3....lenn}，随机选取3个中心点{μ1，μ2，μ3}将数据分为3个簇划分{λ1，λ2，λ3}；
步骤3.2a，计算所有数据距离聚类中心点的欧氏距离dij＝|leni-μj|，并将leni分入对应的类别λi；
步骤3.3a，计算每个簇划分λi中所有样本点的平均值作为新的中心点，重复步骤3.1a、3.2a的过程，直至所有中心点都不再发生变化，完成聚类；
步骤3.4a，将聚类为遥控信号的频率分量保留，将聚类为定频干扰与残余噪声的频率分量置为O，更新原来的时频矩阵得到重构时频矩阵STFTR(t，f)。
5.根据权利要求3所述的无人机遥控图传信号的参数盲估计方法，其特征在于，步骤3中基于时频方差的k-means聚类算法对去噪之后的时频矩阵进行聚类以去除定频干扰，并重构得到重构时频矩阵包括：
步骤3.1b，计算去噪之后的时频矩阵STFT(t，f)中每一行的时频方差si(i∈[1，N])，时频方差的集合为设为{s1，s2，s3....sN}，随机选取2个中心点{α1，α2}将数据分为2个簇划分{β1，β2}；
步骤3.2b，计算所有据距离聚类中心点的欧氏距离Dij＝|si-αj|，(j∈1，2)，并将si分入对应的类别βj；
步骤3.3b，计算每个簇划分βj中所有样本点的平均值作为新的中心点；
步骤3.4b，重复步骤3.2b、3.3b的过程，直至所有中心点都不再发生变化，完成聚类；
步骤3.5b，将聚类为遥控信号的频率分量保留，将聚类为定频干扰与残余噪声的频率分量置为O，更新原来的时频矩阵得到重构时频矩阵srFTR(t，f)。
6.根据权利要求5所述的无人机遥控图传信号的参数盲估计方法，其特征在于，步骤3.1b中去噪之后的时频矩阵如下式：
时频矩阵STFT(t，f)的时频方差为S，则有
其中，j＝1，2，...，n，全部的时频方差的表达式为：
其中，J(t，f)表示干扰信号。
7.根据权利要求4或6所述的无人机遥控图传信号的参数盲估计方法，其特征在于，步骤4包括：
提取到重构时频矩阵STFTR(t，f)的时频脊线
对于STFT来说，其时频脊线是指它的时频分布中每一时刻的峰值频率，即X(t，f)表示信号的时频分布。
8.根据权利要求7所述的无人机遥控图传信号的参数盲估计方法，其特征在于，步骤6包括：
用逐差法求|W(a，t)|各峰之间的时间间隔为Th；
求取Th的倒数，获得遥控信号的跳速。
9.根据权利要求8所述的无人机遥控图传信号的参数盲估计方法，其特征在于，步骤8包括：
当采用取各跳周期的中点对应的频率，作为实际频率的估计值。
10.根据权利要求8所述的无人机遥控图传信号的参数盲估计方法，其特征在于，
根据时间间隔Th确定跳频周期；
截取某一跳的重构时频矩阵，对某一跳的重构时频矩阵求短时傅里叶逆变换，得到遥控信号的实际频率。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310416626.XA CN116484198A (zh) | 2023-04-18 | 2023-04-18 | 无人机遥控图传信号的参数盲估计方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310416626.XA CN116484198A (zh) | 2023-04-18 | 2023-04-18 | 无人机遥控图传信号的参数盲估计方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN116484198A true CN116484198A (zh) | 2023-07-25 |
Family
ID=87214974
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202310416626.XA Pending CN116484198A (zh) | 2023-04-18 | 2023-04-18 | 无人机遥控图传信号的参数盲估计方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN116484198A (zh) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN118399835A (zh) * | 2024-07-01 | 2024-07-26 | 太原理工大学 | 电机的转速估计方法 |
-
2023
- 2023-04-18 CN CN202310416626.XA patent/CN116484198A/zh active Pending
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN118399835A (zh) * | 2024-07-01 | 2024-07-26 | 太原理工大学 | 电机的转速估计方法 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110865357B (zh) | 一种基于参数优化vmd的激光雷达回波信号降噪方法 | |
| US20220349986A1 (en) | Radar communication with interference suppression | |
| CN111934711B (zh) | 一种时频混叠跳频信号的参数估计方法 | |
| CN107153178B (zh) | 外辐射源雷达参考信号含有多径干扰时的目标检测方法 | |
| CN102510363A (zh) | 一种强干扰源环境下的lfm信号检测方法 | |
| CN116484198A (zh) | 无人机遥控图传信号的参数盲估计方法 | |
| CN111507047A (zh) | 一种基于SP-CUnet的逆散射成像方法 | |
| CN107843875A (zh) | 基于奇异值分解降噪的贝叶斯压缩感知雷达数据融合方法 | |
| CN111881858A (zh) | 一种微震信号多尺度去噪方法、装置及可读存储介质 | |
| CN112505675A (zh) | 目标角度和距离定位方法、装置、雷达和存储介质 | |
| CN114362851B (zh) | 一种基于机器学习的无线信道数据去噪方法 | |
| CN116699526A (zh) | 一种基于稀疏与低秩模型的车载毫米波雷达干扰抑制方法 | |
| CN116756491A (zh) | 一种基于蜣螂优化算法优化小波阈值的阀门信号降噪方法 | |
| CN110336585A (zh) | 一种基于mwc的跳频信号参数估计方法 | |
| CN116687392A (zh) | 基于时频信息矩阵的毫米波雷达跌倒检测的杂波消除方法 | |
| CN117574187A (zh) | 一种多辐射源信号分选时延估计方法及系统 | |
| CN105656577A (zh) | 面向信道冲激响应的分簇方法和装置 | |
| CN109752633B (zh) | 一种对变电站局部放电信号进行定位的方法及系统 | |
| CN116299176A (zh) | 一种基于霍夫变换的目标空间特征提取及融合定位方法 | |
| CN112327260B (zh) | 一种sar回波数据中脉冲式干扰信号的抑制方法和装置 | |
| CN116155319A (zh) | 一种电磁态势监视与分析系统及方法 | |
| CN112711001B (zh) | 一种精细去噪辅助的激光雷达波形分解方法 | |
| CN109270509B (zh) | 数据损失情况下基于矩阵填充的doa估计方法及系统 | |
| CN111914802B (zh) | 一种基于迁移学习的电离层返回散射传播模式识别方法 | |
| CN116500553B (zh) | 雷达互相干扰抑制方法、装置、计算机设备及存储介质 |
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 |