基于微多普勒效应的悬停无人机特征检测方法

本发明属于雷达技术领域，更进一步涉及一种悬停无人机特征检测方法。本发明可用于低重频、强地杂波环境下的目标检测。

背景技术：

近年来，随着无人机技术的快速发展，无人机的管控面临着巨大的困难。无人机这种低慢小目标由于其飞行高度“低”、飞行速度“慢”、雷达反射截面积“小”，很容易被淹没在强烈的地杂波中，传统的频域滤波检测方法检测性能严重下降，难以完成对无人机等低慢小目标的检测。

北京航空航天大学在其申请的专利文献“一种基于双发双收相干雷达的无人机旋翼长度和转速的估计方法”(申请号201810779125.7，公开号cn108957443a)中公开了一种对无人机的旋翼长度和转速进行精确估计的方法。该方法首先通过对发射信号进行相位补偿，形成波束对准目标；然后，利用两路接收天线收集旋翼回波信号，并将两路回波信号输入相干器，得到相干信号，通过相干信号的自相关函数提取旋翼的旋转半周期；最后，利用时频谱特定时刻具有幅度极大值的频率的差值，估计无人机旋翼的长度和转速，实现无人机的特征检测和识别。但是，由于无人机旋翼回波的主瓣极窄，当雷达的脉冲重复频率较小时，使用该方法会产生雷达采不到无人机旋翼回波主瓣的问题，导致其无法适用于低重频情况下的无人机参数估计和目标检测。

方鑫在其发表的论文“小型无人机目标雷达探测关键技术研究”(电子科技大学2019博士论文)中提出了一种基于旋翼微多普勒特征的无人机目标检测方法。该方法根据不同叶片微多普勒频率差值与雷达载波频率、叶片长度、叶片旋转速度以及叶片旋转初始角之间的关系，给出了根据微多普勒特征区分不同叶片回波信号的条件，基于时频分析、hrt和二维cfar实现了旋翼叶片的微多普勒特征提取和无人机目标的检测。该方法由于没有考虑杂波的影响，且使用旋翼微多普勒在时频域的正弦曲线对无人机进行检测，而实际数据中无人机旋翼在时频域的曲线仅在高信噪比、高脉冲采样频率的情况下才能被观察到的，因而导致其在低重频、低信杂比情况下的检测效果较差。

技术实现要素：

本发明的目的在于针对上述已有技术的不足，提出一种基于微多普勒效应的悬停无人机特征检测方法，以提取无人机金属转子和旋翼的微多普勒特征，在低重频、低信杂比情况下实现对无人机的有效检测。

为实现上述目的，本发明的技术方案包括如下：

(1)对悬停无人机的雷达回波信号在快时间域进行下变频和脉冲压缩处理，并提取每一个距离单元的慢时间域信号；

(2)将慢时间域数据分为目标单元回波和非目标单元回波两部分，取一半目标单元回波和一半非目标单元回波组成训练数据集，将剩余的回波数据作为测试数据集；

(3)根据训练数据集和测试数据集的慢时间域数据，在不同变换域中计算表征无人机金属转子和旋翼微多普勒特性的特征值：

(3a)在时域计算训练数据的时域赫斯特指数h0和测试数据的时域赫斯特指数h1；

(3b)在频域计算训练数据的频谱熵值ex，0和测试数据的频谱熵值ex，1；

(3c)在时频域计算训练数据时频脊的均值mtfr，0和方差vartfr，0，同时计算测试数据时频脊的均值mtfr,1和方差vartfr,1；

(3d)在倒频域计算训练数据的倒谱熵值eq，0和测试数据的倒谱熵值eq,1；

(4)用上述计算出训练数据集中的所有特征值组成训练特征向量，对支撑向量机进行训练，生成支撑向量机检测器；

(5)使用上述计算出测试数据集的所有特征值组成特征向量，输入到支撑向量机检测器，判断是否有目标存在，实现对悬停无人机的特征检测。

本发明利用时域赫斯特指数、频谱熵值、时频脊的均值和方差、倒谱熵值这5个微多普勒特征对悬停无人机进行特征检测，相比现有技术具有如下优点：

第一，由于这些特征既包含旋翼的微多普勒效应也包含金属转子的微多普勒效应，因此，目标单元回波的特征和非目标单元回波的特征在低重频条件下依旧有明显的区别，可在低重频条件下对悬停无人机进行有效的检测，克服了现有基于旋翼微多普勒的无人机目标检测技术在低重频条件下检测效果差的问题。

第二，由于微多普勒特征受杂波的影响较小，因此，微多普勒特征能够有效区分纯杂波信号和淹没在杂波中的无人机信号，在信杂比较低时也能对悬停无人机进行有效的检测，克服了现有技术难以对淹没在杂波中的无人机进行有效检测的问题，提高了在强杂波环境下的检测效果。

附图说明

图1为本发明的实现流程图；

图2为仿真的目标单元回波；

图3为用本发明对图2进行特征检测的仿真实验结果图。

具体实施方式

下面结合附图，对本发明的实施例和效果作进一步详细描述。

参照图1，对实施的实现步骤如下：

步骤1，对悬停无人机的雷达回波信号在快时间域进行下变频和脉冲压缩处理。

(1.1)获取悬停无人机的雷达回波。

无人机为一悬停的四旋翼无人机，雷达位于一处高楼的楼顶，雷达视线方向正对着无人机进行观测，从雷达的接收机获取悬停无人机的回波信号sr：

sr＝s0+sm+sr+sc+sn

其中，s0为无人机机身的回波，sm为无人机金属转子的回波，sr为无人机旋翼的回波，sc为杂波，sn为噪声；

(1.2)获取悬停无人机的下变频信号：

在快时间域内对悬停无人机的雷达回波信号进行下变频，得到悬停无人机的下变频信号h：

h＝srexp(-j2πztr)，

其中，exp表示以自然常数e为底的指数操作，j表示虚数单位符号，z表示根据工程需求设定的雷达机动目标回波信号的载频，tr表示雷达接收回波信号快时间域的采样时间点，tr的取值范围为fs表示根据工程需求设定的雷达接收机动回波信号的采样频率；

(1.3)获取悬停无人机的脉压信号：

在快时间域内对下变频信号在距离维进行脉冲压缩，得到悬停无人机的脉冲压缩信号g：

其中，表示卷积操作，rect(·)表示矩形窗因子，当时矩形窗因子取值为0，当时矩形窗因子取值为1，tp表示雷达接收的回波信号的脉冲宽度，|·|表示取绝对值操作，k表示雷达接收的回波信号的调频斜率。

步骤2，提取每一个距离单元的慢时间域信号。

在悬停无人机的脉压信号中，提取每一个距离单元的慢时间域信号，包括非目标单元和目标单元的慢时间域回波信号，分别表示如下：

非目标单元的慢时间域回波信号x(n)为：

x(n)＝sc+sn，

其中，n∈[1,n]，n为单帧数据的脉冲个数，sc为杂波，sn为噪声；

目标单元的慢时间域回波信号x(n)为：

其中，a和b为数值不同的两个常数，t＝ntr为慢时间域采样点，tr为脉冲重复时间，s0为机身回波，sr,k(t)为第k个旋翼的回波，sm,k(t)为第k个金属转子的回波。

所述的无人机旋翼回波sr,k(t)和金属转子的回波的公式分别表示如下：

其中，l为旋翼的长度，fd为无人机径向运动的多普勒频率，λ为雷达信号的波长，r0为悬停无人机到雷达的初始距离，nb为单个旋翼的叶片数量,β为无人机相对雷达视线方向的俯仰角，ψt,k＝2πfω+θk为旋翼的旋转角，即雷达与旋翼所在的直线和旋翼叶片所在直线之间的夹角，fω为旋转速度，θk为第k个旋翼的初始旋转角，d为常数，j0(·)为第一类0阶贝塞尔函数，r为金属转子的半径，γt,k为金属转子表面的某散射点p所在切线和雷达与旋翼所在直线的夹角，ψt,k为金属转子表面的散射点p的旋转角，u1(t)和u2(t)为周期矩形脉冲信号：

u2(t)＝1-u1(t)。

步骤3，对慢时间域数据进行分集。

将慢时间域数据分为目标单元回波和非目标单元回波两部分，取一半目标单元回波和一半非目标单元回波组成训练数据集，将剩余的回波数据作为测试数据集。

步骤4，根据训练数据集和测试数据集的慢时间域数据，在不同变换域中计算表征无人机金属转子和旋翼微多普勒特性的特征值。

表征无人机金属转子和旋翼微多普勒特性的特征值，包括时域赫斯特指数h、频谱熵值ef、时频脊的均值mtfr、时频脊的方差vartfr和倒谱熵值eq。这些特征值可以表征悬停无人机在不同变换域的微多普勒特征。

由于地面散射体都是固定不动的，而悬停无人机的旋翼和金属转子都是在快速转动的，所以地杂波的波动性较弱，而悬停无人机的回波的波动性较强。因此，目标单元回波与非目标单元回波的时域波动性有明显的差异，使用时域赫斯特指数h可以计量时域信号的波动性，作为悬停无人机特征检测的一个特征。

由于目标单元回波中有悬停无人机的微多普勒分量，所以其频谱的能量分布较为分散，频谱的混乱度较大。而非目标单元回波中没有微多普勒分量，其频谱的能量主要集中于杂波所在的零频处，混乱度较小。因此，目标单元和非目标单元的频谱混乱度有明显的差异，使用频谱熵值ef可以计量信号频谱的混乱度，作为悬停无人机特征检测的一个特征。

由于目标单元回波中有悬停无人机的微多普勒分量，所以当短时傅里叶的窗长较长时，目标单元回波的时频谱中会出现关于零频对称的相互平行的微多普勒谱线，即金属转子的微多普勒时频脊。对于目标单元回波，除零频外，每一时刻瞬时频率的幅度最大值都在上述时频脊上。对于非目标单元回波，除零频外，每一时刻瞬时频率的幅度最大值都是随机分布的噪声分量。因此，可以提取每一时刻瞬时频率的幅度最大值所在的频率，计算其均值和方差，即时频脊的均值mtfr和时频脊的方差vartfr，作为悬停无人机特征检测的特征。

由于目标单元回波的频谱中有周期性出现的微多普勒谱线，所以其倒谱的能量分布较为分散，倒谱的混乱度较大。而非目标单元回波中没有微多普勒分量，其倒谱的能量主要集中在倒频为零处，混乱度较小。因此，目标单元和非目标单元的倒谱混乱度有明显的差别，使用倒谱熵值eq可以计量信号倒谱的混乱度，作为悬停无人机特征检测的一个特征。

上述特征值的计算分别如下：

(4.1)在时域计算训练数据的时域赫斯特指数h0和测试数据的时域赫斯特指数h1：

(4.1.1)分别从训练数据集中取慢时间域信号x0(n)，计算训练数据的累积离差y0(n)，从测试数据集中取慢时间域信号x1(n)，计算测试数据的累积离差y1(n)：

其中，n∈[1,n]，n为单帧数据的脉冲个数，<·>表示取均值的操作；

(4.1.2)分别从训练数据的累积离差序列y0(n)的开头和末尾开始划分出个连续且长度为m的短序列，构成2nm个短序列，同时，对测试数据的累积离差序列y1(n)执行同样的划分操作，也构成2nm个短序列；

(4.1.3)使用最小二乘法分别对序列y0(n)和序列y1(n)划分出的短序列进行拟合，分别得到训练数据的拟合多项式y0(i)和测试数据的拟合多项式y1(i)，其中，i的取值范围为[1,2nm]；

(4.1.4)分别对训练数据和测试数据进行消除趋势估计，获得训练数据的拟合方差yυ，0(s,m)和测试数据的拟合方差yυ，1(s,m)：

其中，s为拟合多项式的序号，s∈[1,2nm]；

(4.1.5)对(4.1.4)得到的两个拟合方差分别取均值，再开方，得到训练数据的二阶波动函数f0(m)和测试数据的二阶波动函数f1(m)：

(4.1.6)根据训练数据的二阶波动函数f0(m)和测试数据的二阶波动函数f1(m)，使用最小二乘法，分别计算训练数据的时域赫斯特指数h0和测试数据的时域赫斯特指数h1：

(4.2)在频域计算训练数据的频谱熵值ef，0和测试数据的频谱熵值ef，1：

(4.2.1)利用傅里叶变换公式，分别计算训练数据x0(n)的频谱x0(f)和测试数据x1(n)的频谱x1(f)：

其中，f为频率，f∈[-fr/2,fr/2]，fr为雷达的脉冲重复频率；

(4.2.2)分别对训练数据的频谱x0(f)和测试数据的频谱x1(f)进行归一化，得到训练数据的归一化频谱x0,1(f)和测试数据的归一化频谱x1,1(f)：

(4.2.3)分别计算训练数据的频谱熵值ef，0和测试数据的频谱熵值ef，1：

(4.3)在时频域计算训练数据时频脊的均值mtfr，0和方差vartfr，0，同时计算测试数据时频脊的均值mtfr,1和方差vartfr,1：

(4.3.1)将短时傅里叶变换的短时窗设置为汉明窗，为使金属转子的微多普勒特征易于观察和提取，设窗长应大于旋翼旋转周期的2倍，即：

lw＞2t，

其中，lw为短时傅里叶变换的窗长，t为无人机旋翼和金属转子的旋转周期；

(4.3.2)分别对训练数据x0(n)和测试数据x1(n)作短时傅里叶变换，得到训练数据的时频谱stftx0(n,f)和测试数据的时频谱stftx1(n,f)：

其中，n∈[1,n]，f∈[-fr/2,fr/2]，w(·)为短时窗函数；

(4.3.3)分别对训练数据和测试数据的时频谱进行截取，得到去除机身信号的训练数据时频谱stftx0,0(n,f)和去除机身信号的测试数据时频谱stftx1,0(n,f)：

(4.3.4)根据去除机身信号的时频谱，分别计算训练数据中无人机金属转子的时频脊r0(n)和测试数据中无人机金属转子的时频脊r1(n)：

r0(n)＝{|fm||stftx0,0(n,fm)＝max[stftx0,0(n,f)]}，

r1(n)＝{|fm||stftx1,0(n,fm)＝max[stftx1,0(n,f)]}，

(4.3.5)根据(4.3.4)的结果，分别计算训练数据时频脊的均值mtfr，0和测试数据时频脊的均值mtfr,1：

(4.3.6)根据(4.3.4)和(4.3.5)的结果，分别计算训练数据时频脊的方差vartfr，0和测试数据时频脊的方差vartfr,1：

(4.4)在倒频域计算训练数据的倒谱熵值eq，0和测试数据的倒谱熵值eq,1：

(4.4.1)利用倒谱变换公式，分别计算训练数据x0(n)的倒谱c0(q)和测试数据x1(n)的倒谱c1(q)：

c0(q)＝|ifft{log|fft[x0(n)]|}|，

c1(q)＝|ifft{log|fft[x1(n)]|}|，

其中，q为倒频率，fft(·)为傅里叶变换，ifft(·)为逆傅里叶变换；

(4.4.2)分别对训练数据和测试数据的倒谱进行截取，得到去除机身信号的训练数据倒谱c0,0(q)和去除机身信号的测试数据倒谱c1,0(q)：

c0,0(q)＝c0(q)，q∈(t/3,ntr/2)，

c1,0(q)＝c1(q)，q∈(t/3,ntr/2)；

(4.4.3)对截取的倒谱进行归一化，分别计算训练数据的归一化倒谱c0,1(q)和测试数据的归一化倒谱c1,1(q)：

(4.4.4)根据(4.4.3)的结果，分别计算训练数据的倒谱熵值eq,0和测试数据的倒谱熵值eq,1：

步骤5，根据训练数据集中的所有特征值，对支撑向量机进行训练，生成支撑向量机检测器。

(5.1)分别使用训练数据集中每一个慢时间域信号的特征值，组成训练特征向量f0,i：

f0,i＝(h0,ef,0,mtfr,0,vartfr,0,eq,0)，

其中，i∈[1,m0]，m0为训练数据集中慢时间域回波的数量；

(5.2)为训练特征向量f0,i标注类别，得到训练样本集ω0：

ω0＝{(f0,i,yi)|i＝1,2,...,m0}，

其中，yi为慢时间域信号的类别，若慢时间域信号来自目标单元，则yi＝1，若慢时间域信号来自非目标单元，则yi＝-1；

(5.3)根据训练样本集构造一个最优超平面g(f)＝ω^tf+b作为决策曲面，使得目标单元回波和非目标单元回波样本之间的隔离边缘最大化，该最优超平面需满足如下条件：

yi(ω^tf0,i+b)≥1，i＝1,2,...,m0

其中，ω为超平面的法向量，b为超平面的位移项，(·)^t为转置操作；

(5.4)根据(5.3)中的限制条件，将最优超平面参数的求解问题转化为如下的二次规划问题：

再通过拉格朗日乘数法，将上述二次规划问题转化为对偶问题q(α)：

其中，为拉格朗日系数，αi为拉格朗日系数α中的第i个值，αi≥0，i∈[1,m0]；

(5.5)对(5.4)中的对偶问题q(α)进行求解，得到拉格朗日系数α的解为再由拉格朗日系数的解α^*计算超平面法向量的解ω^*和超平面位移项的解b^*：

(5.6)将超平面法向量的解ω^*和超平面位移项的解b^*代入支撑向量机的决策函数j(f)，得到训练好的支撑向量机，即支撑向量机检测器：

其中，sgn(·)为符号函数，f为输入的特征向量，若f是目标单元回波的特征向量，则j(f)＝1，若f是非目标单元回波的特征向量，则j(f)＝-1。

步骤6，利用支撑向量机检测器对悬停无人机进行特征检测。

(6.1)分别使用测试数据集中每一个慢时间域信号的特征值，组成测试特征向量f1,l：

f1,l＝(h1,ef,1,mtfr,1,vartfr,1,eq,1)，

其中，l∈[1,m1]，m1为测试数据集中慢时间域回波的数量；

(6.2)将每一个测试特征向量f1,l分别输入到支撑向量机检测器，得到测试数据所属类别的判决值j(f1,l)：

(6.3)根据j(f1,l)，判断慢时间域信号所在的距离单元是否有目标存在，若j(f1,l)的值为1，则f1,l是目标单元回波的特征向量，即慢时间域信号所在的距离单元有目标存在，若j(f1,l)的值为-1，则f1,l是非目标单元回波的特征向量，即慢时间域信号所在的距离单元无目标存在。

下面结合仿真实验对本发明做进一步的描述：

1.仿真实验条件：

仿真实验的硬件平台为：处理器为intel(r)core(tm)i5-6500cpu，主频为3.2ghz，内存8gb。

仿真实验的软件平台为：windows7操作系统和matlabr2019a。

仿真实验的雷达参数为：雷达发射信号的波长λ＝0.03124m，雷达接收机的采样频率fs＝100mhz，单帧数据的相参处理时间为cpi＝100ms。

仿真实验的目标参数为：目标为一悬停的四旋翼无人机，无人机到雷达的初始距离r0＝500m，目标速度v＝0m/s，方位角α＝0°，俯仰角β＝0°，无人机旋翼转速fω＝100r/s，无人机单个旋翼长度l＝0.12m，无人机金属转子半径r＝0.01m，无人机旋翼和金属转子的初始旋转角[θ0,θ1,θ2,θ3]分别为区间[0，2π]上的四个随机数。

仿真实验的杂波参数为：杂波的分布为韦布尔分布，韦布尔分布的形状参数为p＝5，韦布尔分布的尺度参数为q＝450，杂波的功率谱为立方谱，立方谱的转折频率为fc＝9hz。

2.仿真内容及结果分析：

仿真1，在雷达脉冲重复频率为2khz，信噪比为10db，信杂比为0db的条件下，仿真产生低重频情况下的悬停无人机回波，结果图2所示，其中：

图2(a)为信号时域波形图，

图2(b)为信号频谱图，

图2(c)为信号时频谱图，

图2(d)为信号倒谱图。

由图2可以看出，在脉冲重复频率较低时，悬停无人机回波中仍有明显的微多普勒特征，并且图中的微多普勒特征是由金属转子产生的。因此，可以证明根据无人机金属转子的微多普勒特征能实现低重频情况下的悬停无人机特征检测。

仿真2，在信噪比为10db，雷达脉冲重复频率分别为20khz和2khz的条件下，在区间[-40db,20db]上，每隔5db取一个值作为信杂比，在每一个信杂比下仿真产生1000个目标单元回波和1000个非目标单元回波，用本发明的方法分别在每一个脉冲重复频率和信杂比的组合下对仿真产生的1000个目标单元回波和1000个非目标单元回波进行特征检测，检测结果如图3所示。

从图3可以看出，当雷达脉冲重复频率为20khz时，本发明在信杂比大于等于-25db时能取得大于90％的检测概率,在信杂比大于等于-25db时能取得小于10％的虚警概率；当雷达脉冲重复频率为2khz时，本发明在信杂比大于等于-20db时能取得大于90％的检测概率,在信杂比大于等于-25db时能取得小于10％的虚警概率。

仿真实验结果证明，本发明在低重频、强杂波环境下能对悬停无人机进行有效的检测。