CN103091674B9 - 基于hrrp序列的空间目标高分辨成像方法 - Google Patents

基于hrrp序列的空间目标高分辨成像方法 Download PDF

Info

Publication number
CN103091674B9
CN103091674B9 CN201210596470.XA CN201210596470A CN103091674B9 CN 103091674 B9 CN103091674 B9 CN 103091674B9 CN 201210596470 A CN201210596470 A CN 201210596470A CN 103091674 B9 CN103091674 B9 CN 103091674B9
Authority
CN
China
Prior art keywords
matrix
amplitude
scattering point
hrrp
time
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.)
Active
Application number
CN201210596470.XA
Other languages
English (en)
Other versions
CN103091674A (zh
CN103091674B (zh
Inventor
白雪茹
周峰
刘妍
保铮
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Xidian University
Original Assignee
Xidian University
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Xidian University filed Critical Xidian University
Priority to CN201210596470.XA priority Critical patent/CN103091674B9/zh
Publication of CN103091674A publication Critical patent/CN103091674A/zh
Application granted granted Critical
Publication of CN103091674B publication Critical patent/CN103091674B/zh
Publication of CN103091674B9 publication Critical patent/CN103091674B9/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明公开了一种基于HRRP序列的空间目标高分辨成像方法,其步骤包括:(1)雷达录取ISAR回波;(2)高速运动补偿;(3)获取高质量HRRP;(4)散射点航迹关联;(5)高分辨成像。本发明首先通过稀疏信号重构获得高质量HRRP序列,之后采用有效的航迹关联方法从HRRP序列中生成散射点航迹矩阵,进而通过带约束条件的矩阵分解实现高分辨成像。该方法克服了基于分段伪Keystone变换方法对PRF要求高、对目标运动模型要求严格、对高速旋转目标成像时精确平动补偿和初相校正难度大、运算量高等缺点,具有在低PRF甚至方位多普勒模糊条件下能够进行高分辨成像、对目标运动模型具有普适性、不需平动补偿、效率高、图像聚焦良好等优点。

Description

基于HRRP序列的空间目标高分辨成像方法
技术领域
本发明属于信号处理技术领域,更进一步涉及雷达成像领域中的基于高分辨一维 距离像(high-resolution range profiles,HRRP)序列的空间目标高分辨成像方法。本 发明可以在低重频条件下,有效地对空间目标进行高分辨成像,避免方位欠采样、平 动量和随机初相误差对图像质量的影响,为后续目标识别提供有力保障。
背景技术
为获得对空间目标散射点分布特征更为精细和准确的描述,目前对雷达图像分辨 率的要求在不断提高。其中距离向高分辨依靠发射大的时宽-带宽积信号获得;而方位 向高分辨则依靠成像时间内目标相对于雷达视线运动所引起的多普勒变化获得。为满 足距离向高分辨成像要求,可以通过信号处理方法实现多波段超宽带合成,或通过系 统硬件升级等方法使雷达系统发射超宽带信号。但是在获得距离向高分辨率的同时, 若仍采用传统的成像方法获得方位向高分辨率时会遇到新的问题:首先,距离分辨率 提高的同时会使散射点更易发生越距离单元走动,从而限制了距离-瞬时多普勒、时频 分析等成像方法的应用;其次,由于电磁波受大气扰动的影响,回波相位中会包含时 变的随机初相误差。当散射点存在越距离单元走动时,对该误差的校正及平动补偿会 变得更加困难;最后,由于地基雷达发射功率有限,当对相距较远的空间目标进行探 测时,需要降低雷达脉冲重复频率以保证雷达接收机对信噪比的要求。在这种情况下, 对于旋转频率较高的机动目标,其回波会出现方位欠采样,混迭的方位频谱将进一步 增加对空间目标高分辨成像的难度。因此,有必要研究在低重频条件下的新的高分辨 成像方法。
当采用宽带雷达对目标进行高分辨观测时,根据投影定理,目标的高分辨一维距 离像(HRRP:high-resolution range profiles)将反映出散射点分布沿雷达视线方向的 投影。当雷达视线连续变化时,不同观测时刻对应的HRRP将形成HRRP序列。当散 射点不存在越距离单元走动时,HRRP序列仅描述散射点沿距离向的分布;当采用超 宽带对机动目标进行观测时,由于发生越距离单元走动,实HRRP序列会反映出散射 点在观测时间内的运动轨迹,从而包含散射点的二维或三维分布信息。因此,通过合 适的算法设计能够实现基于HRRP序列的目标高分辨成像。同传统成像方法相比,这 类方法不会受到方位欠采样和随机初相误差的影响,同时降低了成像对雷达PRF的要 求。
Kai Huo等人在文献“A novel imaging method for fast rotating targets based on the segmental pseudo Keystone transform”(IEEE Trans.Geosci.Remote Sens.,vol.49,no.4, pp.1464–1472,Apr.2011)中提出基于分段伪Keystone变换的快速旋转目标二维成像 方法,该方法在平动可被准确补偿,并且目标的旋转频率较高且在观测时间内恒定前 提下可以精确估计出目标的位置。但是,该方法仍存在的不足是,当散射点存在严重 越距离单元走动而使传统相邻相关平动补偿方法失效时,没有考虑有效的平动补偿方 法;当目标存在变加速旋转、或雷达进行重频抖动(即脉冲重复频率随时间变化)以 实现ECCM(electronic counter-countermeasures)时,该方法假设的正弦模型已经不再 满足,从而导致图像模糊;至少需要接收到目标一个旋转周期内的回波才能得到低旁 瓣、聚焦良好的图像,当回波次数较少时,旁瓣很高;运算复杂度较高。
发明内容
本发明的目的在于克服上述现有技术的不足,提出一种基于HRRP序列的空间目 标高分辨成像方法。该方法弥补了基于分段伪Keystone变换方法具有对PRF要求高、 对目标运动模型要求严格、对高速旋转目标成像时精确平动补偿和初相校正难度大、 参数搜索运算量高的不足,不需要对回波进行精确的平动补偿和初相校正,在低重频 条件下仅利用散射点在距离-慢时间平面上的航迹特征进行非参数化成像即可得到目 标聚焦良好的高分辨图像。
实现本发明的基本思路是:首先,采用稀疏信号重构方法获得目标高质量的HRRP 序列;其次,根据散射点的幅度和运动特征,采用航迹关联方法从HRRP序列中提取 出目标强散射点对应的轨迹并生成散射点航迹矩阵;最后,从航迹矩阵中求解出散射 点坐标并得到高分辨图像。
本发明的具体步骤如下:
(1)雷达录取ISAR回波;
(2)高速运动补偿
2a)将参考距离为雷达到场景中心距离,频率、调频率与雷达发射信号相同的线 性调频信号作为参考信号;
2b)将参考信号取共轭后与雷达录取的ISAR回波相乘得到解线调频处理结果;
2c)按照下式构造相位补偿函数:
Figure GDA0003820581330000031
其中,
Figure GDA0003820581330000032
为相位补偿函数,exp(·)为指数函数运算符,j为虚数单位,c为光速, v为测定的目标速度,t为方位维的采样时间,R(t)为t时刻场景中心到雷达的参考 距离,f1为信号载频,γ为调频率,f2为距离频率;
2d)将相位补偿函数与步骤2b)得到的解线调频处理结果相乘得到高速运动补偿 后回波;
(3)获取高质量HRRP
3a)按照下式构造过冗余字典:
Figure GDA0003820581330000033
其中,Φ为包含N个基向量的过冗余字典,N取为距离频率点数的整数倍,exp(·) 为指数函数运算符,j为虚数单位,c为光速,f2为距离频率,
Figure GDA0003820581330000034
f1 为载频,γ为调频率,
Figure GDA0003820581330000035
为距离维采样时间,R为参考距离,RΔ(1),RΔ(2),…RΔ(N) 为构造的N个散射点斜距,
Figure GDA0003820581330000036
L为雷达观测场景长度,相邻 两个距离采样点的距离步进量满足关系
Figure GDA0003820581330000037
M为距离频率点数, B为信号带宽;
3b)令m=1,其中,m表示方位单元序号;
3c)将步骤2d)得到的高速运动补偿后回波的第m个方位单元的复回波作为残余 信号,并生成一个空矩阵Φ1
3d)用残余信号向步骤3a)得到的过冗余字典中的基向量作投影,并且计算内积;
3e)将内积最大的基向量记录在步骤3c)生成的矩阵Φ1的第m列中,按照投影 方法计算内积最大的基向量对应的投影系数,记录得到的投影系数,同时将内积最大的 基向量从过冗余字典中剔除;
3f)从残余信号中减去矩阵Φ1与记录的投影系数的乘积,得到更新的残余信号;
3g)重复步骤3d)至步骤3f)的操作,直至更新的残余信号能量低于噪声门限,将 步骤3e)中记录的投影系数取模值作为对应基向量所代表斜距处散射点的幅度,以基向 量对应的散射点斜距为横坐标,幅度为纵坐标绘制HRRP;
3h)将m递增1个序号,重复步骤3c)至步骤3g)的操作,直至绘制完所有方位 单元对应的HRRP;
(4)散射点航迹关联
4a)令k=1,其中,k表示方位单元序号,生成一个大小为N1×M的空的航迹矩 阵,其中,N1为第k个方位单元对应的HRRP中峰值点的个数,M表示散射点回波 次数,将第k个方位单元对应的HRRP中所有峰值点对应的位置坐标和幅度记录在航 迹矩阵第k列中,将k递增一个序号;
4b)在第k个方位单元对应的HRRP中寻找峰值点,并将所有峰值点对应的位置 坐标、幅度记录在第k时刻的位置幅度集合中,根据小转角内散射点幅度近似不变且 位置接近准则,将第k时刻位置幅度集合中的元素与航迹矩阵第k-1列中元素进行关 联,将关联得到的第k时刻散射点位置坐标和幅度记录在航迹矩阵第k列;
4c)将k递增一个序号,执行步骤4b);
4d)分别将各散射点在航迹矩阵中第k-1列记录的位置坐标与第k-2列记录的 位置坐标相减,并除以脉冲重复周期得到速度,将航迹矩阵中第k列记录的位置坐标 减去第k-1列记录的位置坐标的两倍,再加上第k-2列记录的位置坐标,将计算结 果除以脉冲重复周期的平方得到加速度;
4e)用速度乘以脉冲重复周期并加上加速度与脉冲重复周期平方乘积的0.5倍, 求得距离步进量;
4f)分别将各散射点在航迹矩阵中第k列记录的位置坐标、步骤4d)得到的速度 与加速度作为第k时刻散射点状态向量;
4g)根据Kalman预测准则,对第k+1、k+2及k+3时刻的状态向量进行预测, 并将航迹矩阵中记录的对应散射点幅度求均值作为预测幅度值;
4h)在第k+1个方位单元对应的HRRP中,分别以各散射点在航迹矩阵中第k列 记录的位置坐标为中心,以距离步进量的三倍为搜索宽度寻找峰值点,将所有峰值点 对应的位置坐标和幅度记录在第k+1时刻的位置幅度集合中;
4i)在第k+2个方位单元对应的HRRP中,分别以各散射点在第k+1时刻的位置 幅度集合中所有位置坐标的均值为中心,以距离步进量的三倍为搜索宽度寻找峰值 点,将所有峰值点对应的位置坐标和幅度记录在第k+2时刻的位置幅度集合中;
4j)在第k+3个方位单元对应的HRRP中,分别以各散射点在第k+2时刻的位置 幅度集合中所有位置坐标的均值为中心,以距离步进量的三倍为搜索宽度寻找峰值 点,将所有峰值点对应的位置坐标和幅度记录在第k+3时刻的位置幅度集合中;
4k)根据归一化的最近邻准则函数对散射点航迹进行关联,求得散射点在第k+1 时刻的真实位置坐标和幅度,并记录在航迹矩阵第k+1列中;
4l)将步骤4g)预测得到的k+1时刻状态向量中的速度和加速度按照步骤4e)进 行计算,得到距离步进量,将k递增一个序号,重复步骤4g)至步骤4k),直至完成 散射点在所有方位单元对应的HRRP中的位置关联,从航迹矩阵中抽取散射点的位置 坐标并且剔除对应幅度值,得到更新的航迹矩阵;
(5)高分辨成像
5a)将航迹矩阵中每个元素都减去该元素所在行所有元素的均值,得到平动校正 后的航迹矩阵;
5b)对平动校正后的航迹矩阵进行带约束的矩阵分解得到散射点分布矩阵和雷达 视线矩阵;
5c)用散射点分布矩阵重构目标图像。
与现有技术相比,本发明具有以下优点。
第一,本发明通过采用目标回波对应的实HRRP序列进行高分辨成像,克服了现 有技术中对雷达PRF要求高,在方位欠采样条件下无法得到聚焦良好的高分辨图像等 不足,使得本发明具有在低PRF甚至方位多普勒模糊条件下能够进行高分辨成像。
第二,本发明通过采用航迹关联方法从HRRP序列中提取出目标强散射点对应的 轨迹并生成散射点航迹矩阵,通过对航迹矩阵的简单处理实现对目标的平动补偿,克 服了现有技术中对高速旋转目标成像时精确平动补偿和初相校正难度大的不足,使得 本发明在回波存在剩余平动量和随机初相误差的条件下能够得到目标聚焦良好的高 分辨图像。
第三,本发明通过带约束条件的矩阵分解,采用非参数方法得到目标散射点坐标 和等效雷达视线的分布矩阵,克服了现有技术中对目标运动模型要求严格、高质量成 像时要求目标回波数较多,参数搜索时运算量较大的不足,使得本发明具有对目标运 动模型具有稳健性、成像要求观测回波次数少、运算效率高、图像聚焦良好的优点。
附图说明
图1为本发明的流程图;
图2为本发明的仿真图。
具体实施方式
下面结合附图1,对本发明具体实施方式作进一步的详细描述:
步骤1,雷达录取目标的逆合成孔径ISAR回波,雷达以脉冲重复频率发射并接 收脉冲,得到以距离为行向量以方位为列向量的ISAR回波。
步骤2,高速运动补偿
2a)将参考距离为雷达到场景中心距离,频率、调频率与雷达发射信号相同的线 性调频信号作为参考信号,其中参考信号如下:
Figure GDA0003820581330000061
其中,
Figure GDA0003820581330000062
为参考信号,
Figure GDA0003820581330000063
为距离维采样时间,t1为方位维采样时间,rect(·)为 矩形窗函数,
Figure GDA0003820581330000064
R′为参考距离,T为参考信号的脉宽,j为虚数, f为中心频率,t为距离维采样时间与方位维采样时间之和,c为光速,γ为调频率。
2b)将参考信号取共轭后与雷达录取的ISAR回波相乘得到解线调频处理结果;
2c)按照下式构造相位补偿函数:
Figure GDA0003820581330000065
其中,
Figure GDA0003820581330000066
为相位补偿函数,exp(·)为指数函数运算符,j为虚数单位,c为光速, v为测定的目标速度,t为方位维的采样时间,R(t)为t时刻场景中心到雷达的参考 距离,f1为信号载频,γ为调频率,f2为距离频率。
2d)将相位补偿函数与步骤2b)得到的解线调频处理结果相乘得到高速运动补偿 后回波。
步骤3,获取高质量HRRP
3a)按照下式构造过冗余字典:
Figure GDA0003820581330000071
其中,Φ为包含N个基向量的过冗余字典,N取为距离频率点数的整数倍,exp(·) 为指数函数运算符,j为虚数单位,c为光速,f2为距离频率,
Figure GDA0003820581330000072
f1 为载频,γ为调频率,
Figure GDA0003820581330000073
为距离维采样时间,R为参考距离,RΔ(1),RΔ(2),…RΔ(N) 为构造的N个散射点斜距,
Figure GDA0003820581330000074
L为雷达观测场景长度,相邻 两个距离采样点的距离步进量满足关系
Figure GDA0003820581330000075
M为距离频率点数, B为信号带宽。
3b)令m=1,其中,m表示方位单元的序号;
3c)将步骤2d)得到的高速运动补偿后回波的第m个方位单元的复回波作为残余 信号,并生成一个空矩阵Φ1
3d)用残余信号向步骤3a)得到的过冗余字典中的基向量作投影,并且计算内积;
3e)将内积最大的基向量记录在步骤3c)生成的矩阵Φ1的第m列中,按照投影 方法计算内积最大的基向量对应的投影系数,记录得到的投影系数,同时将内积最大的 基向量从过冗余字典中剔除。其中,投影方法可描述为如下表达式:
Figure GDA0003820581330000076
其中,α为投影系数,Φ1为记录最大投影系数对应基向量的矩阵,
Figure GDA0003820581330000077
为Φ1的共 轭转置,(·)-1为矩阵求逆运算,g为残余信号。
3f)从残余信号中减去矩阵Φ1与记录的投影系数的乘积,得到更新的残余信号;
3g)重复步骤3d)至步骤3f)的操作,直至更新的残余信号能量低于噪声门限,将 步骤3e)中记录的投影系数取模值作为对应基向量所代表斜距处散射点的幅度,以基向 量对应的散射点斜距为横坐标,幅度为纵坐标绘制HRRP。其中噪声门限可以设定为对 步骤3c)得到的残余信号进行逆傅里叶变换并取模值,将取模值后的图像中远离峰值 点的其余点对应的平均功率作为噪声门限。
3h)将m递增1个序号,重复步骤3c)至步骤3g)的操作,直至绘制完所有方位 单元对应的HRRP。
步骤4,散射点航迹关联
4a)令k=1,其中,k表示方位单元的序号,生成一个大小为N1×M的空的航迹 矩阵,其中,N1为第k个方位单元的HRRP中峰值点的个数,M表示散射点回波次 数,将第k个方位单元的HRRP中所有峰值点对应的位置坐标、幅度记录在航迹矩阵 的第k列中,将k递增一个序号;
4b)在第k个方位单元的HRRP中寻找峰值点,并将所有峰值点对应的位置坐标、 幅度记录在第k时刻的位置幅度集合中,根据小转角内散射点幅度近似不变且位置接 近准则,将第k时刻位置幅度集合中的元素与航迹矩阵的第k-1列中元素进行关联, 将关联得到的第k时刻散射点的位置坐标及幅度记录在航迹矩阵的第k列。其中小转 角内散射点幅度近似不变位置接近准则为当前时刻散射点幅度应该满足如下等式:
Ak=arg min||A′k-Ak-1||2
其中,Ak为当前散射点在第k时刻对应的最优幅度值,argmin(·)为求使目标函 数取最小值时对应变量值的运算,||·||2为求2范数,Ak′为第k时刻的位置幅度集合中 当前散射点对应的所有幅度值,Ak-1为航迹矩阵中第k-1列中当前散射点对应的幅 度值。
4c)将k递增一个序号,执行步骤4b);
4d)分别将各散射点在航迹矩阵中第k-1列记录的位置坐标与第k-2列记录的 位置坐标相减,并除以脉冲重复周期得到速度,将航迹矩阵中第k列记录的位置坐标 减去第k-1列记录的位置坐标的两倍,再加上第k-2列记录的位置坐标,将计算结 果除以脉冲重复周期的平方得到加速度;
4e)用速度乘以脉冲重复周期并加上加速度与脉冲重复周期平方乘积的0.5倍, 求得距离步进量;
4f)分别将各散射点在航迹矩阵中第k列记录的位置坐标、步骤4d)得到的速度 及加速度作为第k时刻散射点状态向量;
4g)根据Kalman预测准则,对第k+1、k+2及k+3时刻的状态向量进行预测, 并将航迹矩阵中记录的对应散射点幅度求均值作为预测幅度值。其中,Kalman预测 准则如下:
Figure GDA0003820581330000091
其中,
Figure GDA0003820581330000092
为由k时刻预测的k+1时刻的状态,Φ为状态转移矩阵,
Figure GDA0003820581330000093
T=1/PRF为脉冲重复间隔,PRF为雷达脉冲重复频率,
Figure GDA0003820581330000094
为 由k-1时刻预测的k时刻的状态,Gk为k时刻增益矩阵, ΦGk=ΦPk/k-1HT[HPk/k-1HT+Rk]-1,Pk/k-1为k时刻预测误差的协方差矩阵, Pk/k-1=ΦPk-1ΦT+Qk-1,Pk-1为k-1时刻平滑估计的协方差矩阵,ΦT为Φ的转置, Qk-1为k-1时刻动态噪声的协方差矩阵,
Figure GDA0003820581330000095
n(k-1)是均值为零 的高斯白噪声序列在k-1时刻的取值,H为观测矩阵,H=[100],HT为H的 转置,Rk为k时刻观测噪声的协方差矩阵,(·)-1为矩阵求逆运算,yk为k时刻的观 测值,
Figure GDA0003820581330000096
为k-1时刻状态的估计值。
4h)在第k+1个方位单元的HRRP中,分别以各散射点在航迹矩阵中第k列记录 的位置坐标为中心,以距离步进量的三倍为搜索宽度寻找峰值点,将所有峰值点对应 的位置坐标和幅度记录在第k+1时刻的位置幅度集合中;
4i)在第k+2个方位单元的HRRP中,分别以各散射点在第k+1时刻的位置幅度 集合中所有位置坐标的均值为中心,以距离步进量的三倍为搜索宽度寻找峰值点,将 所有峰值点对应的位置坐标和幅度记录在第k+2时刻的位置幅度集合中;
4j)在第k+3个方位单元的HRRP中,分别以各散射点在第k+2时刻的位置幅度 集合中所有位置坐标的均值为中心,以距离步进量的三倍为搜索宽度寻找峰值点,将 所有峰值点对应的位置坐标和幅度记录在第k+3时刻的位置幅度集合中;
4k)根据归一化的最近邻准则函数对散射点航迹进行关联,求散射点在第k+1时 刻的真实位置坐标和幅度,并记录在航迹矩阵的第k+1列中。其中,归一化的最近邻 准则函数如下:
Figure GDA0003820581330000101
其中,rk为当前散射点在第k时刻的最优位置坐标,argmin(·)为求使目标函数 取最小值时对应变量值的运算,Ak+1为第k+1时刻的位置幅度集合中当前散射点对 应的所有幅度值,
Figure GDA0003820581330000102
为当前散射点在航迹矩阵中前k列的幅度均值,rk+1为第k+1 时刻的位置幅度集合中当前散射点对应的所有位置坐标,
Figure GDA0003820581330000103
为由k时刻预测 k+1时刻的状态向量中的位置坐标,rk+2为第k+2时刻的位置幅度集合中当前散射 点对应的所有位置坐标,
Figure GDA0003820581330000104
为由k+1时刻预测k+2时刻的状态向量中的位 置坐标,rk+3为第k+3时刻的位置集合中当前散射点对应的所有位置坐标,
Figure GDA0003820581330000105
为由k+2时刻预测k+3时刻的状态向量中的位置坐标。
4l)将步骤4g)预测得到的k+1时刻状态向量中的速度和加速度按照步骤4f)进 行计算,得到距离步进量,将k递增一个序号,重复步骤4g)至步骤4k),直至完成 散射点在所有方位单元的HRRP中的位置关联,从航迹矩阵中抽取散射点的位置坐标 并且剔除对应幅度值,得到更新的航迹矩阵。
步骤5,高分辨成像
5a)将航迹矩阵中每个元素都减去该元素所在行所有元素的均值,得到平动校正 后的航迹矩阵;
5b)对平动校正后的航迹矩阵进行带约束的矩阵分解得到散射点分布矩阵和雷达 视线矩阵,其中带约束的矩阵分解步骤如下:
第一步,对平动校正后的航迹矩阵进行矩阵分解:
Figure GDA0003820581330000111
其中,W为平动校正后的航迹矩阵,
Figure GDA0003820581330000112
为对W矩阵分解后得到的两个酉矩阵;
第二步,对下列方程求最小二乘解,得到尺度变换矩阵:
Figure GDA0003820581330000113
其中,ln
Figure GDA0003820581330000114
的第n个行向量,
Figure GDA0003820581330000115
为第一步奇异值分解得到的酉矩阵,n∈[1,Na], Na为方位向采样点数,A为尺度变换矩阵,AT为A的转置,
Figure GDA0003820581330000116
为ln的转置,I为 单位矩阵;
第三步,对下列非线性最小二乘方程组求解,得到旋转矩阵:
Figure GDA0003820581330000117
其中,l为初始时刻雷达视线矢量,
Figure GDA0003820581330000118
为从
Figure GDA0003820581330000119
中得到的对初始时刻雷达视线矢量的 估计,
Figure GDA00038205813300001110
为第一步奇异值分解得到的酉矩阵,A1为旋转矩阵,
Figure GDA00038205813300001111
为A1的转置,I为 单位矩阵;
第四步,获取真实的等效雷达视线旋转矩阵与散射点分布矩阵:
Figure GDA00038205813300001112
其中,R为等效雷达视线旋转矩阵,S为散射点分布矩阵,
Figure GDA00038205813300001113
为第一步奇异 值分解得到的两个酉矩阵,A为第二步得到的尺度变换矩阵,A1为第三步得到的旋 转矩阵。
5c)用散射点分布矩阵重构目标图像。
下面结合附图2对本发明的效果做进一步说明。
附图2所示的仿真图在MATLAB7.0软件下进行的,仿真数据的参数如下:雷达带宽 为8GHz,载频为96GHz,PRF为150Hz,观测时间为0.1s,方位向回波为15次,接 收回波的信噪比为10dB,仿真数据包含5个强散射点,目标旋转频率在[3.5π,4.5π]rad/s 间随机变化,速度为4000m/s。按照仿真参数,方位无模糊采样所需的PRF应为9651Hz, 因此仿真的PRF远远小于传统成像方法所需的PRF。
图2(a)是采用基于HRRP序列的空间目标高分辨成像方法对散射点航迹进行关 联的结果,其中水平坐标表示回波序号,垂直坐标表示散射点航迹,分别用加号、乘 号、星号、圆圈及实点代表了不同的航迹,从图中可看出,共五条散射点航迹。
图2(b)是目标旋转频率随回波的变化示意图,其中水平坐标表示回波序号,垂 直坐标表示旋转角频率,单位为rad/s,由图中可看出旋转角频率随时间呈随机变化。
图2(c)是采用基于HRRP序列的空间目标高分辨成像方法对图2(a)进行高分 辨成像的结果图,其中,水平坐标表示散射点的水平位置,垂直坐标表示散射点的垂直 位置,单位均为米,实际的散射点位置用星号表示,估计出的散射点位置用圆圈表示, 由图中可看出估计的所有散射点位置与实际散射点位置均重合到一起,说明虽然目标存 在平动,旋转角频率随机变化并且脉冲重复频率较低,但本方法仍能够准确估计散射点 位置,从而可以获得高分辨的图像,证明了该方法的有效性。

Claims (7)

1.基于HRRP序列的空间目标高分辨成像方法,包括如下步骤:
(1)雷达录取ISAR回波;
(2)高速运动补偿
2a)将参考距离为雷达到场景中心距离,频率、调频率与雷达发射信号相同的线 性调频信号作为参考信号;
2b)将参考信号取共轭后与雷达录取的ISAR回波相乘得到解线调频处理结果;
2c)按照下式构造相位补偿函数:
Figure FDA0003820581320000011
其中,
Figure FDA0003820581320000012
为相位补偿函数,exp(·)为指数函数运算符,j为虚数单位,c为光速, v为测定的目标速度,t为方位维的采样时间,R(t)为t时刻场景中心到雷达的参考 距离,f1为信号载频,γ为调频率,f2为距离频率;
2d)将相位补偿函数与步骤2b)得到的解线调频处理结果相乘得到高速运动补偿 后回波;
(3)获取高质量HRRP
3a)按照下式构造过冗余字典:
Figure FDA0003820581320000013
其中,Ψ为包含N个基向量的过冗余字典,N取为距离频率点数的整数倍,exp(·) 为指数函数运算符,j为虚数单位,c为光速,f2为距离频率,
Figure FDA0003820581320000014
f1 为载频,γ为调频率,
Figure FDA0003820581320000015
为距离维采样时间,R为参考距离,RΔ(1),RΔ(2),…RΔ(N) 为构造的N个散射点斜距,
Figure FDA0003820581320000016
L为雷达观测场景长度,相邻 两个距离采样点的距离步进量满足关系
Figure FDA0003820581320000017
M为距离频率点数, B为信号带宽;
3b)令m=1,其中,m表示方位单元序号;
3c)将步骤2d)得到的高速运动补偿后回波的第m个方位单元的复回波作为残余 信号,并生成一个空矩阵Φ1
3d)用残余信号向步骤3a)得到的过冗余字典中的基向量作投影,并且计算内积;
3e)将内积最大的基向量记录在步骤3c)生成的矩阵Φ1的第m列中,按照投影 方法计算内积最大的基向量对应的投影系数,记录得到的投影系数,同时将内积最大的 基向量从过冗余字典中剔除;
3f)从残余信号中减去矩阵Φ1与记录的投影系数的乘积,得到更新的残余信号;
3g)重复步骤3d)至步骤3f)的操作,直至更新的残余信号能量低于噪声门限,将 步骤3e)中记录的投影系数取模值作为对应基向量所代表斜距处散射点的幅度,以基向 量对应的散射点斜距为横坐标,幅度为纵坐标绘制HRRP;
3h)将m递增1个序号,重复步骤3c)至步骤3g)的操作,直至绘制完所有方位 单元对应的HRRP;
(4)散射点航迹关联
4a)令k=1,其中,k表示方位单元序号,生成一个大小为N1×M1的空的航迹矩 阵,其中,N1为第k个方位单元对应的HRRP中峰值点的个数,M1表示散射点回波 次数,将第k个方位单元对应的HRRP中所有峰值点对应的位置坐标和幅度记录在航 迹矩阵第k列中,将k递增一个序号;
4b)在第k个方位单元对应的HRRP中寻找峰值点,并将所有峰值点对应的位置 坐标、幅度记录在第k时刻的位置幅度集合中,根据小转角内散射点幅度近似不变且 位置接近准则,将第k时刻位置幅度集合中的元素与航迹矩阵第k-1列中元素进行关 联,将关联得到的第k时刻散射点位置坐标和幅度记录在航迹矩阵第k列;
4c)将k递增一个序号,执行步骤4b),直至获得所有方位时刻的散射点位置坐 标和幅度;
4d)分别将各散射点在航迹矩阵中第k-1列记录的位置坐标与第k-2列记录的 位置坐标相减,并除以脉冲重复周期得到速度,将航迹矩阵中第k列记录的位置坐标 减去第k-1列记录的位置坐标的两倍,再加上第k-2列记录的位置坐标,将计算结 果除以脉冲重复周期的平方得到加速度;
4e)用速度乘以脉冲重复周期并加上加速度与脉冲重复周期平方乘积的0.5倍, 求得距离步进量;
4f)分别将各散射点在航迹矩阵中第k列记录的位置坐标、步骤4d)得到的速度 与加速度作为第k时刻散射点状态向量;
4g)根据Kalman预测准则,对第k+1、k+2及k+3时刻的状态向量进行预测, 并将航迹矩阵中记录的对应散射点幅度求均值作为预测幅度值;
4h)在第k+1个方位单元对应的HRRP中,分别以各散射点在航迹矩阵中第k列 记录的位置坐标为中心,以距离步进量的三倍为搜索宽度寻找峰值点,将所有峰值点 对应的位置坐标和幅度记录在第k+1时刻的位置幅度集合中;
4i)在第k+2个方位单元对应的HRRP中,分别以各散射点在第k+1时刻的位置 幅度集合中所有位置坐标的均值为中心,以距离步进量的三倍为搜索宽度寻找峰值 点,将所有峰值点对应的位置坐标和幅度记录在第k+2时刻的位置幅度集合中;
4j)在第k+3个方位单元对应的HRRP中,分别以各散射点在第k+2时刻的位置 幅度集合中所有位置坐标的均值为中心,以距离步进量的三倍为搜索宽度寻找峰值 点,将所有峰值点对应的位置坐标和幅度记录在第k+3时刻的位置幅度集合中;
4k)根据归一化的最近邻准则函数对散射点航迹进行关联,求得散射点在第k+1 时刻的真实位置坐标和幅度,并记录在航迹矩阵第k+1列中;
4l)将步骤4g)预测得到的k+1时刻状态向量中的速度和加速度按照步骤4e)进 行计算,得到距离步进量,将k递增一个序号,重复步骤4g)至步骤4k),直至完成 散射点在所有方位单元对应的HRRP中的位置关联,从航迹矩阵中抽取散射点的位置 坐标并且剔除对应幅度值,得到更新的航迹矩阵;
(5)高分辨成像
5a)将航迹矩阵中每个元素都减去该元素所在行所有元素的均值,得到平动校正 后的航迹矩阵;
5b)对平动校正后的航迹矩阵进行带约束的矩阵分解得到散射点分布矩阵和雷达 视线矩阵;
5c)用散射点分布矩阵重构目标图像。
2.根据权利要求1所述的基于HRRP序列的空间目标高分辨成像方法,其特征 在于,步骤3e)中所述的投影方法如下:
Figure FDA0003820581320000041
其中,α为投影系数,Φ1为记录最大投影系数对应基向量的矩阵,
Figure FDA0003820581320000042
为Φ1的共 轭转置,(·)-1为矩阵求逆运算,g为残余信号。
3.根据权利要求1所述的基于HRRP序列的空间目标高分辨成像方法,其特征 在于,步骤3g)中所述的噪声门限为,对步骤3c)得到的残余信号进行逆傅里叶变 换并取模值,将取模值后的图像中远离峰值点的其余点对应的平均功率作为噪声门 限。
4.根据权利要求1所述的基于HRRP序列的空间目标高分辨成像方法,其特征 在于,步骤4b)中所述的小转角内散射点幅度近似不变且位置接近准则为,当前时刻 散射点幅度应该满足如下等式:
Ak=arg min||A′k-Ak-1||2
其中,Ak为当前散射点在第k时刻对应的最优幅度值,arg min(·)为求使目标函 数取最小值时对应变量值的运算,||·||2为求2范数,A′k为第k时刻的位置幅度集合中 当前散射点对应的所有幅度值,Ak-1为航迹矩阵中第k-1列中当前散射点对应的幅 度值。
5.根据权利要求1所述的基于HRRP序列的空间目标高分辨成像方法,其特征 在于,步骤4g)中所述的Kalman预测准则如下:
Figure FDA0003820581320000043
其中,
Figure FDA0003820581320000044
为由k时刻预测的k+1时刻的状态,Φ为状态转移矩阵,
Figure FDA0003820581320000045
T=1/PRF为脉冲重复间隔,PRF为雷达脉冲重复频率,
Figure FDA0003820581320000046
为 由k-1时刻预测的k时刻的状态,Gk为k时刻增益矩阵, ΦGk=ΦPk/k-1HT[HPk/k-1HT+Rk]-1,Pk/k-1为k时刻预测误差的协方差矩阵, Pk/k-1=ΦPk-1ΦT+Qk-1,Pk-1为k-1时刻平滑估计的协方差矩阵,ΦT为Φ的转置, Qk-1为k-1时刻动态噪声的协方差矩阵,
Figure FDA0003820581320000051
n(k-1)是均值为零 的高斯白噪声序列在k-1时刻的取值,H为观测矩阵,H=[1 0 0],HT为H的 转置,Rk为k时刻观测噪声的协方差矩阵,(·)-1为矩阵求逆运算,yk为k时刻的观 测值,
Figure FDA0003820581320000052
为k-1时刻状态的估计值。
6.根据权利要求1所述的基于HRRP序列的空间目标高分辨成像方法,其特征 在于,步骤4k)中所述的归一化的最近邻准则函数如下:
Figure FDA0003820581320000053
其中,rk为当前散射点在第k时刻的最优位置坐标,arg min(·)为求使目标函数 取最小值时对应变量值的运算,A′k+1为第k+1时刻的位置幅度集合中当前散射点对 应的所有幅度值,
Figure FDA0003820581320000054
为当前散射点在航迹矩阵中前k列的幅度均值,r′k+1为第k+1 时刻的位置幅度集合中当前散射点对应的所有位置坐标,
Figure FDA0003820581320000055
为由k时刻预测 k+1时刻的状态向量中的位置坐标,r′k+2为第k+2时刻的位置幅度集合中当前散射 点对应的所有位置坐标,
Figure FDA0003820581320000056
为由k+1时刻预测k+2时刻的状态向量中的位 置坐标,r′k+3为第k+3时刻的位置集合中当前散射点对应的所有位置坐标,
Figure FDA0003820581320000057
为由k+2时刻预测k+3时刻的状态向量中的位置坐标。
7.根据权利要求1所述的基于HRRP序列的空间目标高分辨成像方法,其特征 在于,步骤5b)中所述的带约束的矩阵分解步骤如下:
第一步,对平动校正后的航迹矩阵进行分解:
Figure FDA0003820581320000058
其中,W为平动校正后的航迹矩阵,
Figure FDA0003820581320000059
为对W矩阵分解后得到的两个酉矩 阵;
第二步,对下列方程求最小二乘解,得到尺度变换矩阵:
Figure FDA0003820581320000061
其中,ln
Figure FDA0003820581320000062
的第n个行向量,
Figure FDA0003820581320000063
为第一步奇异值分解得到的酉矩阵,n∈[1,Na], Na为方位向采样点数,A为尺度变换矩阵,AT为A的转置,
Figure FDA0003820581320000064
为ln的转置,I为 单位矩阵;
第三步,对下列非线性最小二乘方程组求解,得到旋转矩阵:
Figure FDA0003820581320000065
其中,l为初始时刻雷达视线矢量,
Figure FDA0003820581320000066
为从
Figure FDA0003820581320000067
中得到的对初始时刻雷达视线矢量的 估计,
Figure FDA0003820581320000068
为第一步奇异值分解得到的酉矩阵,A1为旋转矩阵,
Figure FDA0003820581320000069
为A1的转置,I为 单位矩阵;
第四步,获取真实的等效雷达视线旋转矩阵与散射点分布矩阵:
Figure FDA00038205813200000610
其中,
Figure FDA00038205813200000611
为等效雷达视线旋转矩阵,S为散射点分布矩阵,
Figure FDA00038205813200000612
为第一步奇异 值分解得到的两个酉矩阵,A为第二步得到的尺度变换矩阵,A1为第三步得到的旋 转矩阵。
CN201210596470.XA 2012-12-14 2012-12-14 基于hrrp序列的空间目标高分辨成像方法 Active CN103091674B9 (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201210596470.XA CN103091674B9 (zh) 2012-12-14 2012-12-14 基于hrrp序列的空间目标高分辨成像方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201210596470.XA CN103091674B9 (zh) 2012-12-14 2012-12-14 基于hrrp序列的空间目标高分辨成像方法

Publications (3)

Publication Number Publication Date
CN103091674A CN103091674A (zh) 2013-05-08
CN103091674B CN103091674B (zh) 2014-09-17
CN103091674B9 true CN103091674B9 (zh) 2023-01-31

Family

ID=48204476

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201210596470.XA Active CN103091674B9 (zh) 2012-12-14 2012-12-14 基于hrrp序列的空间目标高分辨成像方法

Country Status (1)

Country Link
CN (1) CN103091674B9 (zh)

Families Citing this family (16)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103576130B (zh) * 2013-11-05 2015-09-09 西安电子科技大学 一种进动锥体的三维成像方法
CN103760558B (zh) * 2014-01-23 2017-02-08 电子科技大学 一种太赫兹雷达isar成像方法
CN104459668B (zh) * 2014-12-03 2017-03-29 西安电子科技大学 基于深度学习网络的雷达目标识别方法
CN104535970B (zh) * 2014-12-25 2017-01-11 西安电子工程研究所 基于最大值的频率步进雷达信号目标抽取方法
CN104865568B (zh) * 2015-06-02 2017-05-24 西安电子科技大学 基于稀疏重构的宽带雷达高速群目标分辨方法
CN105022060A (zh) * 2015-07-20 2015-11-04 合肥工业大学 针对快速空天目标的步进isar成像方法
CN105894117B (zh) * 2016-03-31 2020-02-18 北京航空航天大学 航迹预测方法及装置
CN108509989B (zh) * 2018-03-26 2020-04-21 西安电子科技大学 基于高斯选控玻尔兹曼机的hrrp识别方法
CN109709549A (zh) * 2019-02-28 2019-05-03 电子科技大学 一种前视雷达超分辨率成像方法
CN110082765B (zh) * 2019-04-29 2023-04-07 西安电子科技大学 基于三维重构的空间目标姿态外推方法
CN111830501B (zh) * 2020-06-28 2023-04-28 中国人民解放军战略支援部队信息工程大学 Hrrp历史特征辅助的信号模糊数据关联方法及系统
CN112069651B (zh) * 2020-07-23 2024-04-09 西安空间无线电技术研究所 一种基于isar成像的自旋稳定目标旋转轴估计方法
CN112014817B (zh) * 2020-08-24 2023-06-02 中国电子科技集团公司第三十八研究所 一种空间自旋目标三维重构方法
CN113640791B (zh) * 2021-06-09 2023-12-26 西安电子科技大学 一种基于距离和瞬时速度的空间目标三维姿态重构方法
CN113687325B (zh) * 2021-07-08 2024-02-06 西安电子科技大学 基于lp和hrrp模型的遮蔽小目标检测方法
CN114509736B (zh) * 2022-01-19 2023-08-15 电子科技大学 一种基于超宽带电磁散射特征的雷达目标识别方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101509972B (zh) * 2009-03-30 2011-06-29 西安电子科技大学 基于高分辨目标距离像修正相关矩阵的宽带雷达检测方法

Also Published As

Publication number Publication date
CN103091674A (zh) 2013-05-08
CN103091674B (zh) 2014-09-17

Similar Documents

Publication Publication Date Title
CN103091674B9 (zh) 基于hrrp序列的空间目标高分辨成像方法
CN107229048B (zh) 一种高分宽幅sar动目标速度估计与成像方法
Zhang et al. Translational motion compensation for ISAR imaging under low SNR by minimum entropy
CN111157985B (zh) 基于多站一维距离像序列的空间刚体目标三维重构方法
Yang et al. Compressed sensing radar imaging with compensation of observation position error
CN110082764B (zh) 基于稳健正则化层析方法的sar图像成像方法
Liu et al. A novel channel phase bias estimation method for spaceborne along-track multi-channel HRWS SAR in time-domain
CN110187347B (zh) 一种地球同步轨道星机双基合成孔径雷达大幅宽成像方法
Yeh et al. Rotational motion estimation for ISAR via triangle pose difference on two range-Doppler images
JP2010014700A (ja) 画像レーダ装置
CN104122551B (zh) 基于二维酉esprit的isar横向定标方法
Berizzi et al. Performance analysis of a contrast-based ISAR autofocusing algorithm
CN109188436B (zh) 适用于任意平台轨迹的高效双基sar回波生成方法
CN110988873A (zh) 基于能量中心提取的单通道sar舰船速度估计方法及系统
Wei et al. Sparse autofocus via Bayesian learning iterative maximum and applied for LASAR 3-D imaging
CN102012510A (zh) 基于时间——相位导数分布的逆合成孔径雷达成像方法
Wei et al. LASAR autofocus imaging using maximum sharpness back projection via semidefinite programming
CN104181514A (zh) 一种合成孔径雷达高精度运动补偿方法
Bucciarelli et al. Multi-grazing ISAR for side-view imaging with improved cross-range resolution
Li et al. Doppler keystone transform for SAR imaging of moving targets
CN112684446B (zh) 基于最小熵准则的Bi-ISAR横向定标与畸变校正方法
Thammakhoune et al. Moving target imaging for synthetic aperture radar via RPCA
CN105974413B (zh) 多基地外辐射源雷达成像系统的自聚焦方法
CN113484862A (zh) 一种自适应的高分宽幅sar清晰重构成像方法
CN112859074A (zh) 多频带多视角isar融合成像方法

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
CI03 Correction of invention patent
CI03 Correction of invention patent

Correction item: Claims|Description

Correct: correct

False: error

Number: 38

Page: ??

Volume: 30