一种基于Hilbert-Huang变换的超宽带时域RCS测量信号处理方法与流程

本发明涉及针对雷达目标的rcs测量领域，具体涉及一种基于hilbert-huang变换的超宽带时域rcs测量信号处理方法。

背景技术

无论是研究物体的电磁散射特性还是研制具有突防能力的隐身武器系统，rcs测试具有非常重要的意义。通过rcs测试可以验证电磁散射计算的理论和算法，更重要的是，对复杂目标，电磁散射理论计算已非常困难，而通过测试可以直观的获得目标的电磁散射特性数据从而避开复杂的电磁仿真计算。研究目标雷达散射截面的方法分为理论探究和测试技术。虽然现如今有整套完整的电磁理论，以便分析典型的散射原理，但是对于复合材料的目标物体，通过我们现在的技术很难准确的计算其rcs，因此，测试技术是最准确、最有效的方法。通过测量各种目标物体，可以得到其基本的散射现象，通过对比检验理论分析，同时通过实际测量，获得各种目标物体的特征数据，可以建立目标特性的数据库，以便以后使用。测试目标物体rcs值的测试场地可以分为场外测试和室内测试，因为不同的气候都会影响场外的测试结果，并且不能做到有效果的保密性，所以一般都采用室内测试进行研制飞行器等保密性工作，这样的室内测试不仅使研究人员有一个相对舒适的工作环境，而且这也减少了三分之一以上的测试时间，提高了测试效率。

由于采用频谱的不同，可以把测试方法分为点频测试和扫频测试。点频测试(使用窄带接收机接收单极化窄带电磁波)是通过测试可以获得目标物体的具体位置信息和目标的时域波形，然后计算得到测试目标物体的各种信息，如目标物体方位、角速度、径向速度和各种粗略的特征信号等。扫频测试(使用宽带接收机接收宽带电磁波)不仅可以获得以上的信息，而且可以获得目标物体的径向方向上的一维距离像，以及一维结构特征。除此之外，使用全极化电磁波照射目标物体，可以从雷达的回波信号中获得目标的极化特征，然后再使用目标物体的极化特点识别目标的具体特征，进而对目标物体进行识别与分类。

现有的rcs测量方法多采用频域测试，针对所需了解的rcs频段范围，需要选择多个频点进行多次测试；同时为了降低外界干扰、噪声等影响，通常需要在微波暗室环境下进行测试，这样会导致在针对大型目标进行测试时，很难寻找到合适测试环境。因此频域测量方法存在测试步骤繁琐、测试成本较大、测试结果易波动等缺点。

因此需要一种新的rcs测量方法，来解决现有技术的缺陷。

技术实现要素：

本发明的目的是为了解决现有的rcs测量方法很难寻找到合适测试环境，因而导致频域测量方法存在测试步骤繁琐、测试成本较大、测试结果易波动等缺点，而提出一种基于hilbert-huang变换的超宽带时域rcs测量信号处理方法，包括：

步骤一、对原始雷达回波信号进行经验模态分解，得到所有imf分量以及残差分量；

步骤二、计算所述原始雷达回波信号的近似熵；

步骤三、将所述近似熵作为阈值，将步骤一中得到的imf分量进行后项叠加，以使原始雷达回波信号重构。

本发明的有益效果为：本发明通过进行经验模态分解将雷达回波信号分解为一组本征模态函数，接下来使用基于近似熵阈值的emd信号提取方法，对这一组本征模态函数进行提取，提取出的信号经过复原即可很好的降低信号不确定性，尽量消除其中的噪声对测量结果的影响。此方案对于雷达目标rcs的准确探测有重要意义。

附图说明

图1为本发明的基于hilbert-huang变换的超宽带时域rcs测量信号处理方法的流程图；

图2为后项叠加方法的示意图；stepn表示后项叠加过程中的第n次叠加；

图3为后项叠加滤波的示意图。

具体实施方式

具体实施方式一：本实施方式的基于hilbert-huang变换的超宽带时域rcs测量信号处理方法，如图1所示，包括：

步骤一、对原始雷达回波信号进行经验模态分解，得到所有imf分量以及残差分量。其中imf表示本征模态函数。

步骤二、计算所述原始雷达回波信号的近似熵。

步骤三、将所述近似熵作为阈值，将步骤一中得到的imf分量进行后项叠加，以使原始雷达回波信号重构。

具体实施方式二：本实施方式与具体实施方式一不同的是：步骤一的具体过程为：

对于原始雷达回波信号x(t)，它的上、下包络线分别为u(t)和v(t)，则上、下包络的平均曲线为：

m(t)＝[u(t)+v(t)]/2(1-1)

用h1(t)代替x(t)，u1(t)和v1(t)分别表示了与h1(t)相适应的上下包络线，重复以上过程，即：

h2(t)＝x1(t)-m1(t)(1-3)

hk(t)＝hk-1(t)-mk-1(t)(1-5)

若hk(t)满足imf条件，则将hk(t)作为第一个imf分量c1(t)；将原始雷达波回波信号除去c1(t)后的剩余部分记为r1(t)，即：

c1(t)＝hk(t)(1-6)

r1(t)＝x(t)-c1(t)(1-7)

对信号的剩余部分r1(t)继续进行emd分解，直至分解所剩余部分为一单调信号或其值小于预先给定的值为止；最终分解将获得所有imf分量及残差分量：

x(t)＝c1(t)+c2(t)+…+cn(t)+rn(t)(1-8)

其它步骤及参数与具体实施方式一相同。

具体实施方式三：本实施方式与具体实施方式一或二不同的是：

步骤二具体为：

设原始雷达回波数据包含n点数据，即x＝[x(1),x(2),...,x(n)]；

(1)根据原始雷达回波数据，构建一组m维雷达回波子矢量：

x(i)＝[x(i),x(i+1),…,x(i+m-1)]，其中i＝1,…,(n-m+1)

(2)定义不同雷达回波子矢量x(i),x(j)间的距离d[x(i),x(j)]，其物理意义是两子矢量对应表示雷达回波能量元素中差值最大的一个，即：其中j＝1,…,(n-m+1)且j≠i；

(3)此时雷达回波子矢量x(i),x(j)中其他对应元素间差值小于d，并对每一个i值计算子矢量x(i)与其余子矢量x(j),j＝1,…,i-1,i+1…,(n-m+1)的距离d[x(i),x(j)]。

(4)在给定阈值r的条件下，计算每一个i值下不同雷达回波矢量间距离d[x(i),x(j)]小于r的数目q，再计算该值与与距离总数n-m的比值，记作即：

(5)对取对数，求其对所有i的平均值φ^m(r)：

(6)再把维数加1，变成m+1，重复前四个步骤，进一步得和

(7)得到该组雷达回波数据的近似熵为

当n为有限值时步骤(7)中近似熵的统计估计表示为:

fapen(m,r,n)＝φ^m(r)-φ^m+1(r)(3-7)

其它步骤及参数与具体实施方式一或二相同。

具体实施方式四：本实施方式与具体实施方式一至三之一不同的是：基于hilbert-huang变换的超宽带时域rcs测量信号处理方法还包括：

步骤四、对于步骤三中重构后的雷达回波信号进行以下计算得到目标物体的rcs；

其中eobject(f)是步骤三中重构后雷达回波信号对应的频谱；σobject(f)为所述目标物体的rcs，c(f)为校准函数，c(f)是用标准物体先测试一遍得到的；f是频率，

其它步骤及参数与具体实施方式一至三之一相同。

具体实施方式五：本实施方式是具体实施方式一至四的结合，下面对其中原理做具体描述。

1、概述

本发明所提出的针对时域rcs测量时对待测目标反射信号进行处理处理的算法以瞬时频率为基本研究对象，以信号固有模态信号为基本信号单元，与传统的时频分析方法相比有明显的区别，其中瞬时频率是该算法的核心内容。hht是一种自适应的信号处理方法，适用于非平稳信号的处理。通过进行经验模态分解将雷达回波信号分解为一组本征模态函数，接下来使用基于近似熵阈值的emd信号提取方法对这一组本征模态函数进行提取，提取出的信号经过复原即可很好的降低信号不确定性，尽量消除其中的噪声对测量结果的影响。

2、经验模态分解

经验模态分解(empiricalmodedecomposition,emd)是希尔伯特-黄变换(hilbert-huangtransformation,hht)中的重要部分，hht通过emd将一个时间序列分解成一组本征模态函数(intrinsicmodefunction,imf)，然后经hilbert变换后得到信号的瞬时频率特征。从与小波变换等常规方法相比，emd不需要预先设定基函数，可根据信号自身的特征进行分解，具有自适应性，所得的瞬时频率特征能够更好地突出了数据的局部特征，非常适用于非平稳性、非线性过程的信号处理中。

根据hilbert-huang变换理论，它主要是用于给出信号的模态，稳定的单频率正弦信号可以理解为一种模态。然而绝大部分的自然界信号并不具备固有的模态函数，而是在任意时刻的观测数据都隐含了多个振荡的模式。而模态在物理意义上包含了与频率对应的关系，然而区别于傅里叶变换理论的是hilbert-huang变换理论中的模态并不仅限于傅里叶级数展开时的表征不同频率分量的正弦信号，hilbert-huang变换理论中的模态具有更宽松的定义，只需要在函数对称的对称性与局部零均值分布上进行考察。据此，huang提出了固有模态函数的定义。一个固有模态函数需要满足以下两个条件：

1.该信号极值点数与过零点数差值的绝对值小于等于1；

2.相同时刻，该信号的局部极大值和局部极小值的包络平均值为0。

设雷达回波信号为x(t)，它的上、下包络线分别为u(t)和v(t)，则上、下包络的平均曲线为：

m(t)＝[u(t)+v(t)]/2(1-1)

用h1(t)代替x(t)，u1(t)和v1(t)分别表示了与h1(t)相适应的上下包络线，重复以上过程，即：

h2(t)＝x1(t)-m1(t)(1-3)

hk(t)＝hk-1(t)-mk-1(t)(1-5)

在hk(t)满足imf条件下，这样就分解得到第一个imf。此时c1(t)和信号的剩余部分为r1(t)，即：

c1(t)＝hk(t)(1-6)

r1(t)＝x(t)-c1(t)(1-7)

对信号的剩余部分r1(t)继续进行emd分解，分解结束的条件是分解所剩余部分为一单调信号或其值小于预先给定的值。最终分解将获得所有imf分量及残差分量：

x(t)＝c1(t)+c2(t)+…+cn(t)+rn(t)(1-8)

3、近似熵阈值设定

用emd进行雷达回波信号噪声滤除的基本思路在于将雷达回波信号分别通过emd分解后得到各层本征函数imf；雷达回波信号低频噪声存在于小尺度的imf中，高频噪声则存在于大尺度的imf中,对应emd分解过程中先筛选出的分量，则认定在某一个imf阈值处，其筛选的分量中含有较多的噪声能量，将分解后的imf分量由残差，大尺度分量到小尺度分量逐一相加至该阈值对应imf分量，即可重建原始函数，得到处理后的雷达回波信号(信号中的有用成分)。此算法的目的是将估计的雷达回波的误差消除掉，将剩余部分作为雷达回波导航或相关应用的输入。下面给出理论说明。

设信号为x(t)，而x(t)是所有imf及残余量之和：

对分解后的imf函数叠加可以得到原始的雷达回波信号，这里叠加的方法可以按照后项叠加进行。后项叠加主要是从残差与后项本征模态函数累加开始，不断依次加入各阶本征模态函数，寻找合适的本征模态函数的集合重构出原始的雷达回波信号。从信号处理的角度可以分析得出，对emd后的雷达回波数据进行后项叠加还原也是一种低通滤波行为。利用emd进行滤波降噪的前提就是对很多信号而言，其大部分信息主要集中在低频段，相对来说高频信息较少，于是可以利用低频段的若干imf对雷达回波信号进行部分重建。如图2所示。

对于上述理论可以通过式(3-2)进行描述：

其中x(t)为信号观测值，s(t)为信号值，n(t)为噪声，为信号估计值。

根据上述方法，对上文中待测目标的数据进行emd处理后，得到了若干阶imf和一个残余分量。然后从高阶到低阶依次将信号叠加最后还原成原始信号，同样从低阶到高阶依次相加也能得到原始信号，没有信息丢失。残余分量则代表了信号的大致趋势。通过采用后项叠加的方法进行原始信号的复原。根据本征特征函数分量的特性，能够表现出雷达回波信号分量近似于原始函数的时间渐变特性

方法进行信号滤波降噪的关键就是寻找合适的imf阈值k，使得以该分界值为起点的imf进行信号重构获得的信噪比最大或均方误差最小。

关于imf阈值分解点的确定，采用了连续均方误差的准则，首先计算由下式定义的相邻两个估计的均方误差：

可以将imf能量首次发生转折的位置作为噪声起主导作用imf与信号起主导作用imf的分界，同时考虑到不同的信号包含不同的频率成分，在信号分量起主导作用的imf中，某些频率段的imf能量可能低于第一个能量转折处imf的能量而成为全局最小值，因此进一步修正了该结果为：

(1)若相邻两个估计的均方误差在全局极小值之前存在局部极小值，则k应该取第一个局部极小值所对应的位置加1；

(2)如果不存在局部极小值，则k取全局最小值所对应的位置加1。

但是针对雷达回波的信号而言，其主要的信号能量集中在其波峰所在的极短的时间内，按照平均能量的方法很难反映出其噪声与信号间的关系，信号的细节会在平均化这种粗粒度处理过程中消失殆尽。在利用imf分量进行后项叠加对原始信号进行复原的过程中，满足复原后信噪比最大或者均方误差最小也可以理解为估计信号与原始信号具有更强的契合度，包含阈值imf分量在内的后项imf叠加取得的信号与原始信号应该具有很强的一致性，而在无噪声的条件下，一致性最大，随着噪声的加强则相似性相应降低，可见噪声是影响一致性的充分条件。噪声给信号带来的影响还体现在信号的不确定性上，随着后项imf不断叠加，当叠加的分界线来到阈值时，认定后项imf受到噪声干扰较小，而待叠加的imf分量存在较大的噪声干扰。原始的雷达回波信号多为时域上的窄带脉冲，信号含义清晰，具有很强的确定性，而噪声表现为信号的不确定性，因此一旦叠加的imf分量越过阈值点，则应该在信号的不确定性上得到明显的变化，因此选取表示信号不确定性的熵作为阈值判定标准。本段的后项滤波方法示意图如图3所示。

熵值是衡量信号信息丰富程度的一个重要指标，熵值的大小表示信号所包含的平均信息量的多少。对于一组雷达回波信号而言，可以认为其时域上表现出来的能量的急剧变化是一种信息的表现，显然该信息样本与噪声样本代表的不确定信息是有本质区别的，是相互独立的。叠加前后的估计信号其信息量必然会发生明显变化，计算信息熵可以客观地评价雷达回波信号在临近阈值imf时叠加前后信息量的变化。

一般地，信息熵的含义表示为：

式中xi等于该信号点出现的概率。

在实际应用过程中选择近似熵作为简化模型。近似熵的定义如下

设原始雷达回波数据x＝[x(1),x(2),...,x(n)]，共n点数据。

(1)按序号连续顺序组成一组m维雷达回波子矢量：

x(i)＝[x(i),x(i+1),...,x(i+m-1)]，其中i＝1,…,(n-m+1)

(2)定义不同雷达回波子矢量x(i),x(j)间的距离d[x(i),x(j)]，其物理意义是两子矢量对应表示雷达回波能量元素中差值最大的一个，即：

(3)此时雷达回波子矢量x(i),x(j)中其他对应元素间差值小于d，并对每一个i值计算子矢量x(i)与其余子矢量x(j),j＝1,…,i-1,i+1…,(n-m+1)的距离d[x(i),x(j)]。

(4)在给定阈值r的条件下，计算每一个i值下不同雷达回波矢量间距离d[x(i),x(j)]小于r的数目，再计算该值与与距离总数n-m的比值，记作即：

(5)对取对数，求其对所有i的平均值φ^m(r)：

(6)再把维数加1，变成m+1，重复前四个步骤，进一步得和φⁿ⁺¹(r)。

(7)理论上该组雷达回波数据的近似熵为

步骤7给出的极限值趋近于1。然后这种极限关系的前提是n趋近于∞。当n为有限值时步骤7的统计估计表示为:

fapen(m,r,n)＝φ^m(r)-φ^m+1(r)(3-7)

近似熵的值取决于m、r值的大小，指出可行的取值范围为：m＝2，r＝0.1～0.2·std[x(i)],i＝1,...,n，其中std[x(i)]表示为x(i)的标准差。

实际上，这里以近似熵作为阈值主要是利用imf分量叠加时，信号不确定性增加的大小。近似熵愈大，说明不确定性增大愈大，因此噪声的影响也越大。

本发明还可有其它多种实施例，在不背离本发明精神及其实质的情况下，本领域技术人员当可根据本发明作出各种相应的改变和变形，但这些相应的改变和变形都应属于本发明所附的权利要求的保护范围。