CN105353362B - 一种融合能量和变化率两种测量的船长估计方法 - Google Patents

Abstract

一种融合能量和变化率两种测量的船长估计方法，通过高分辨雷达获取目标舰船的高分辨距离像，并对取得的HRRP数据进行归一化预处理；将归一化预处理后的HRRP数据进行低通滤波处理；对低通滤波处理结果中每个点求微分，找出所有可能的端点位置；设定待测的起点和终点检测时能量比的双门限，并初始化待检测信号窗的大小；从HRRP波形图的两端向内逼近，并将能量定义为待检测信号窗内所有点的幅值的绝对值之和，在HRRP波形图中，船头(或船尾)处船身一侧的幅值比船身以外的幅值明显更大，因而能量更高，逐个比较待测点两侧能量的比值与门限的大小，直至找出船头和船尾在待检测信号窗中的相对位置；求船头和船尾的差值得到舰船长度。

Description

一种融合能量和变化率两种测量的船长估计方法

技术领域

本发明属于雷达目标识别领域，涉及一种融合能量和变化率两种测量的船长估计方法。

背景技术

高分辨距离像(HRRP)是用宽带雷达信号获取的目标散射点复子回波在雷达射线上投影的向量和，其特点是发出某一波长的高频信号，通过反射成像时间和位置，从而得出高分辨率距离像，它提供了目标散射点沿距离方向的分布信息,具有目标重要的结构特征，对目标识别与分类十分有价值。舰船的一维距离像反映出舰船主体结构分布，高峰一般出现在船楼、桅杆、烟囱等较高处，峰堆之间的波谷也一定程度上反映出舰船主体结构之间的相对位置。

对于HRRP数据幅度敏感性问题，通常采取对原始高分辨距离像归一化来解决。海上舰船作为一种移动目标，其回波是目标后向电磁散射和海杂波背景电磁散射的叠加，让预处理后的HRRP数据经过滤波处理，可以有效减少杂波的干扰。传统的船长估计算法总有一定的局限性，主要根据实际船长、投影船长、姿态角之间的几何关系，对HRRP数据进行单阈值处理，通过确定目标区域所占据的距离单元数来计算船长，但这种算法并没有验证其估测的精度。还有一种基于尾迹特征的SAR图像舰船长度估计方法，先提取波高谱，再根据波高谱和波幅函数之间的关系计算波幅函数，最后根据波幅函数曲线的周期特性估计舰船长度，这种算法忽略了实际环境中合成孔径雷达图像存在大量的噪声和杂波干扰，因此不能保证算法的鲁棒性，而且SAR与高分辨距离像的机理不同，导致其应用起来有困难。

数字滤波器根据冲激响应的不同，可分为有限冲激响应(FIR)滤波器和无限冲激响应(IIR)滤波器；从滤波特性上分类，则有低通、高通、带通、带阻几种类型。IIR滤波器采用递归性结构，引入反馈环节，通常由延时、乘以系数、相加等基本运算组合成直接型、正准型、级联型、并联型四种结构形式。切比雪夫滤波器是一种在通带或阻带上频率响应幅度等纹波动的滤波器。振幅特性在通带内是等波纹、在阻带内是单调的称为切比雪夫Ⅰ型滤波器；振幅特性在通带内是单调的、在阻带内是等波纹的称为切比雪夫Ⅱ型滤波器。

发明内容

本发明提供了一种可行性高、准确率高的融合能量和变化率两种测量的船长估计方法。

本发明采用的技术方案是：

一种融合能量和变化率两种测量的船长估计方法，其具体步骤如下：

(1)通过高分辨雷达获取目标舰船的高分辨距离像HRRP，并对取得的HRRP数据进行归一化预处理，设D_i(i＝1,2,...,m m是数据长度)表示原始HRRP中的每个数据，则归一化后的HRRP数据D'_i可用公式(1)计算：

式(1)中,max函数表示取{}内所有数值的最大者，min函数表示取{}内部数值的最小者；

(2)将归一化预处理后的HRRP数据D'_i进行低通滤波处理，输出结果为A_i(i＝1,2,...,m)；

(3)对A_i(i＝1,2,...,m)中每个点求微分，找出所有可能的端点位置，所有待测起点的位置用H_j(j＝1,2，…,n₁,n₁是H_j长度)表示，所有待测终点的位置用I_k(k＝1,2,…,n₂,n₂是I_k长度)表示；

(4)设定待测的起点和终点检测时能量比的双门限，分别用T_s和T_e表示；并初始化待检测信号窗的大小W，初始化低通滤波器延迟大小D；

(5)检测端点，从HRRP波形图的两端向内逼近，在H_j中从前往后搜索起点，在I_k中从后向前搜索终点，并将能量定义为A_i的一个窗内所有点的幅值的绝对值之和，在HRRP波形图中，船头(或船尾)处船身一侧的幅值比船身以外的幅值明显更大，因而能量更高，逐个比较待测点两侧能量的比值与门限的大小，直至找出船头和船尾在A_i中的相对位置，用P_s表示起点位置，P_e表示终点位置；

(6)求取舰船长度，设船长为Len，则Len＝P_e-P_s。

进一步，步骤(2)中的低通滤波器采用切比雪夫Ⅰ型的IIR低通滤波器，其差分方程可用式(2)表示：

式(2)中，x(n)是输入信号，y(n)是输出信号，a_k和b_k为滤波器系数；并设置低通滤波器的阶数、截止频率、通带纹波、采样率。

进一步，步骤(3)中端点位置的确定步骤包括：

1)由于相邻数据间的步长为1，求微分等同于差分化，即：

△A_i＝A_i+1-A_i 其中i＝1,2,…,m-1 (3)

2)再对△A_i归一化处理，使：

其中，i＝1,2,…,m-1 (4)

3)设定两个临界值g₁和g₂，H_j为△A'_i中所有大于g₁的点，即I_k为△A'_i中所有小于g₂的点，即

进一步，步骤(5)中起点P_s的检测步骤如下：

A、初始化起点，令P_s＝H_1+j-D，式中j＝0,1,2,...,l₁，检测前先令j＝0，令l₁＝W，即最多检测W次；

B、计算能量值，根据定义，待测起点右侧的能量用式(5)计算：

其中，i_min＝P_s，i_max＝min{P_s+W,m} (5)

式(5)中，m是A_i总长度，引入min函数是为了当右侧数据量少于一个窗时，取到最后一个点为止，下同；待测起点左侧的能量用式(6)计算：

其中，i_min＝max{P_s-W,v}，i_max＝P_s-1 (6)

式(6)中，v＝1，即A_i中第一个点所在位置，引入max函数是为了当左侧数据量少于一个窗时，从第一个点开始取，下同；

C、门限判决，将步骤B计算出的能量值之比E_s1/E_s0与设定的T_s比较，当E_s1/E_s0>T_s则退出循环，确定起点位置；否则，令j＝j+1并返回步骤A。

进一步，步骤(5)中终点P_e的检测步骤如下：

a、初始化终点，令式中k＝0,1,2,...,l₂，检测前先令k＝0，令l₂＝W；

b、计算能量值，采用式(7)计算待测终点右侧的能量值：

其中，i_min＝P_e+1，i_max＝min{P_e+W,m} (7)

采用式(8)计算待测终点左侧的能量值：

其中，i_min＝max{P_e-W,v},i_max＝P_e (8)

c、门限判决，将步骤b计算出的能量值之比E_e1/E_e0与设定T_e比较，当E_e1/E_e0>T_e则退出循环，确定终点位置；否则，令k＝k+1并返回步骤a。

本发明首先需要将HRRP数据进行归一化，以消除幅度敏感性的影响；将归一化后的HRRP数据进行低通滤波处理，可在一定程度上降低海杂波等噪声的干扰，从而提取到对舰船长度估计有用的信息；HRRP数据波形图中，端点(对应船头和船尾)处的曲线走势通常较“陡”，即斜率的绝对值较大，对低通滤波结果A_i中每个点求微分，找出所有可能的端点位置；然后通过对各确定端点位置的能量比与预设能量比的比较来确定起点位置和终点位置，最后求起点位置和终点位置差值得到舰船长度，得到的结果准确率高，具有较高的可行性。

本发明的有益效果：结果准确率高，具有较高的可行性。

附图说明

图1是本发明的流程示意图。

图2是本发明的低通滤波器的幅值与相位曲线图。

图3是某目标舰船的HRRP数据滤波前和滤波后的波形图。

具体实施方式

下面结合具体实施例来对本发明进行进一步说明，但并不将本发明局限于这些具体实施方式。本领域技术人员应该认识到，本发明涵盖了权利要求书范围内所可能包括的所有备选方案、改进方案和等效方案。

参照图1，一种融合能量和变化率两种测量的船长估计方法，其具体步骤如下：

(1)通过高分辨雷达获取目标舰船的高分辨距离像HRRP，并对取得的HRRP数据进行归一化预处理，设D_i(i＝1,2,...,m m是数据长度)表示原始HRRP中的每个数据，则归一化后的HRRP数据D'_i可用公式(1)计算：

式(1)中,max函数表示取{}内所有数值的最大者，min函数表示取{}内部数值的最小者；

(2)将归一化预处理后的HRRP数据D'_i进行低通滤波处理，输出结果为A_i(i＝1,2,...,m)；

(3)对A_i(i＝1,2,...,m)中每个点求微分，找出所有可能的端点位置，所有待测起点的位置用H_j(j＝1,2，…,n₁,n₁是H_j长度)表示，所有待测终点的位置用I_k(k＝1,2,…,n₂,n₂是I_k长度)表示；

(4)设定待测的起点和终点检测时能量比的双门限，分别用T_s和T_e表示，并令T_s＝T_e＝2；并初始化待检测信号窗的大小W＝50，对待检测信号进行定长截取，来计算这个窗内的能量；初始化低通滤波器延迟大小D＝0；

(5)检测端点，从HRRP波形图的两端向内逼近，在H_j中从前往后搜索起点，在I_k中从后向前搜索终点，并将能量定义为A_i的一个窗内所有点的幅值的绝对值之和，在HRRP波形图中，船头(或船尾)处船身一侧的幅值比船身以外的幅值明显更大，因而能量更高，逐个比较待测点两侧能量的比值与门限的大小，直至找出船头和船尾在A_i中的相对位置，用P_s表示起点位置，P_e表示终点位置；本发明假设起点在左，终点在右，与舰船的实际位置无关；

(6)求取舰船长度，设船长为Len，则Len＝P_e-P_s。

本实施例步骤(2)中的低通滤波器采用切比雪夫Ⅰ型的IIR低通滤波器，其差分方程可用式(2)表示：

式(2)中，x(n)是输入信号，y(n)是输出信号，a_k和b_k为滤波器系数；并将其阶数设为8，截止频率设为500HZ，通带纹波设为0.009db，采样率设为40050HZ。

本实施例步骤(3)中端点位置的确定步骤包括：

1)由于相邻数据间的步长为1，求微分等同于差分化，即：

△A_i＝A_i+1-A_i 其中i＝1,2,…,m-1 (3)

2)再对△A_i归一化处理，使：

其中，i＝1,2,…,m-1 (4)

3)设定两个临界值g₁＝0.03和g₂＝-0.03，H_j为△A'_i中所有大于g₁的点，即I_k为△A'_i中所有小于g₂的点，即

本实施例步骤(5)中起点P_s的检测步骤如下：

A、初始化起点，令P_s＝H_1+j-D，式中j＝0,1,2,...,l₁，检测前先令j＝0，令l₁＝50，即最多检测50次；

B、计算能量值，根据定义，待测起点右侧的能量用式(5)计算：

其中，i_min＝P_s，i_max＝min{P_s+W,m} (5)

式(5)中，m是A_i总长度，引入min函数是为了当右侧数据量少于一个窗时，取到最后一个点为止，下同；待测起点左侧的能量用式(6)计算：

其中，i_min＝max{P_s-W,v}，i_max＝P_s-1 (6)

式(6)中，v＝1，即A_i中第一个点所在位置，引入max函数是为了当左侧数据量少于一个窗时，从第一个点开始取，下同；

C、门限判决，将步骤B计算出的能量值之比E_s1/E_s0与设定的T_s比较，当E_s1/E_s0>T_s则退出循环，确定起点位置；否则，令j＝j+1并返回步骤A。

本实施例步骤(5)中终点P_e的检测步骤如下：

a、初始化终点，令式中k＝0,1,2,...,l₂，检测前先令k＝0，令l₂＝50；

b、计算能量值，采用式(7)计算待测终点右侧的能量值：

其中，i_min＝P_e+1，i_max＝min{P_e+W,m} (7)

采用式(8)计算待测终点左侧的能量值：

其中，i_min＝max{P_e-W,v},i_max＝P_e (8)

c、门限判决，将步骤b计算出的能量值之比E_e1/E_e0与设定T_e比较，当E_e1/E_e0>T_e则退出循环，确定终点位置；否则，令k＝k+1并返回步骤a。

本发明首先需要将HRRP数据进行归一化，以消除幅度敏感性的影响；将归一化后的HRRP数据进行低通滤波处理，可在一定程度上降低海杂波等噪声的干扰，从而提取到对舰船长度估计有用的信息；HRRP数据波形图中，端点(对应船头和船尾)处的曲线走势通常较“陡”，即斜率的绝对值较大，对低通滤波结果A_i中每个点求微分，找出所有可能的端点位置；然后通过对各确定端点位置的能量比与预设能量比的比较来确定起点位置和终点位置，最后求起点位置和终点位置差值得到舰船长度，得到的结果准确率高，具有较高的可行性。

本发明实验随机选取了20个待测目标，它们是一批不同姿态角下的高分辨率雷达一维距离像实测数据，对应海面上行驶的军船或者民船。首先对每个目标的HRRP数据归一化预处理，然后利用设计好的IIR滤波器对其滤波，图2为滤波器的幅值和相位曲线图。图3为某目标的HRRP数据经过预处理和滤波后的波形图。实验对每组数据的两个端点分别循环检测，最后得出船长，并计算出了估测结果的相对误差和总的平均相对误差，如表1所示。其中，相对误差＝|估测船长-实际船长|/实际船长，平均相对误差是取20个目标相对误差的均值。

表1 20个舰船目标的船长估计结果

从表1的结果不难看出，使用本发明提到的算法估计船长具有很高的准确性，大部分舰船目标长度估计的相对误差都在10％以下，平均相对误差也只有6.28％。另外，实验使用的都是高分辨雷达实测到的回波数据，鲁棒性良好，具有足够的说明力。

Claims (5)

1.一种融合能量和变化率两种测量的船长估计方法，其具体步骤如下：

(1)通过高分辨雷达获取目标舰船的高分辨距离像HRRP，并对取得的HRRP数据进行归一化预处理，设D_i表示原始HRRP中的每个数据，其中i＝1,2,...,m，m是数据长度，则归一化后的HRRP数据D_i'可用公式(1)计算：

其中i＝1,2,...,m (1)

式(1)中,max函数表示取{}内所有数值的最大者，min函数表示取{}内部数值的最小者；

(2)将归一化预处理后的HRRP数据D′_i进行低通滤波处理，输出结果为A_i，其中i＝1,2,...,m；

(3)对A_i中每个点求微分，找出所有可能的端点位置，所有待测起点的位置用H_j表示，j＝1,2，…,n₁,n₁是H_j长度，所有待测终点的位置用I_k表示，其中k＝1,2,…,n₂,n₂是I_k长度；

(4)设定待测的起点和终点检测时能量比的双门限，分别用T_s和T_e表示；并初始化待检测信号窗的大小W，初始化低通滤波器延迟大小D；

(5)检测端点，从HRRP波形图的两端向内逼近，在H_j中从前往后搜索起点，在I_k中从后向前搜索终点，并将能量定义为A_i的一个窗内所有点的幅值的绝对值之和，在HRRP波形图中，船头或船尾处船身一侧的幅值比船身以外的幅值明显更大，因而能量更高，逐个比较待测点两侧能量的比值与门限的大小，直至找出船头和船尾在A_i中的相对位置，用P_s表示起点位置，P_e表示终点位置；

(6)求取舰船长度，设船长为Len，则Len＝P_e-P_s。

2.如权利要求1所述的一种融合能量和变化率两种测量的船长估计方法，其特征在于：步骤(2)中的低通滤波器采用切比雪夫Ⅰ型的IIR低通滤波器，其差分方程可用式(2)表示：

<mrow> <mi>y</mi> <mrow> <mo>(</mo> <mi>n</mi> <mo>)</mo> </mrow> <mo>=</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>K</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>N</mi> </munderover> <msub> <mi>a</mi> <mi>K</mi> </msub> <mi>y</mi> <mrow> <mo>(</mo> <mi>n</mi> <mo>-</mo> <mi>K</mi> <mo>)</mo> </mrow> <mo>+</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>K</mi> <mo>=</mo> <mn>0</mn> </mrow> <mi>M</mi> </munderover> <msub> <mi>b</mi> <mi>K</mi> </msub> <mi>x</mi> <mrow> <mo>(</mo> <mi>n</mi> <mo>-</mo> <mi>K</mi> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> </mrow>

式(2)中，x(n)是输入信号，y(n)是输出信号，a_K和b_K为滤波器系数；并设置低通滤波器的阶数、截止频率、通带纹波、采样率。

3.如权利要求1或2所述的一种融合能量和变化率两种测量的船长估计方法，其特征在于：步骤(3)中端点位置的确定步骤包括：

1)由于相邻数据间的步长为1，求微分等同于差分化，即：

△A_i＝A_i+1-A_i 其中i＝1,2,…,m-1 (3)

2)再对△A_i归一化处理，使：

其中，i＝1,2,…,m-1 (4)

3)设定两个临界值g₁和g₂，H_j为△A′_i中所有大于g₁的点，即I_k为△A′_i中所有小于g₂的点，即

4.如权利要求3所述的一种融合能量和变化率两种测量的船长估计方法，其特征在于：步骤(5)中起点P_s的检测步骤如下：

A、初始化起点，令P_s＝H_1+j-D，式中j＝0,1,2,...,l₁，检测前先令j＝0，令l₁＝W，即最多检测W次；

B、计算能量值，根据定义，待测起点右侧的能量用式(5)计算：

其中，i_min＝P_s，i_max＝min{P_s+W,m} (5)

式(5)中，m是A_i总长度，引入min函数是为了当右侧数据量少于一个窗时，取到最后一个点为止，下同；待测起点左侧的能量用式(6)计算：

其中，i_min＝max{P_s-W,v}，i_max＝P_s-1 (6)

式(6)中，v＝1，即A_i中第一个点所在位置，引入max函数是为了当左侧数据量少于一个窗时，从第一个点开始取，下同；

C、门限判决，将步骤B计算出的能量值之比E_s1/E_s0与设定的T_s比较，当E_s1/E_s0>T_s则退出循环，确定起点位置；否则，令j＝j+1并返回步骤A。

5.如权利要求4所述的一种融合能量和变化率两种测量的船长估计方法，其特征在于：步骤(5)中终点P_e的检测步骤如下：

a、初始化终点，令式中k＝0,1,2,...,l₂，检测前先令k＝0，令l₂＝W；

b、计算能量值，采用式(7)计算待测终点右侧的能量值：

其中，i_min＝P_e+1，i_max＝min{P_e+W,m} (7)

采用式(8)计算待测终点左侧的能量值：

其中，i_min＝max{P_e-W,v},i_max＝P_e (8)

c、门限判决，将步骤b计算出的能量值之比E_e1/E_e0与设定T_e比较，当E_e1/E_e0>T_e则退出循环，确定终点位置；否则，令k＝k+1并返回步骤a。