CN108415017B - 复杂目标雷达散射特性稀疏表征的一维增广状态空间方法 - Google Patents

复杂目标雷达散射特性稀疏表征的一维增广状态空间方法 Download PDF

Info

Publication number
CN108415017B
CN108415017B CN201810441778.4A CN201810441778A CN108415017B CN 108415017 B CN108415017 B CN 108415017B CN 201810441778 A CN201810441778 A CN 201810441778A CN 108415017 B CN108415017 B CN 108415017B
Authority
CN
China
Prior art keywords
matrix
pole
dimensional
solving
state space
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
CN201810441778.4A
Other languages
English (en)
Other versions
CN108415017A (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.)
Beihang University
Original Assignee
Beihang 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 Beihang University filed Critical Beihang University
Priority to CN201810441778.4A priority Critical patent/CN108415017B/zh
Publication of CN108415017A publication Critical patent/CN108415017A/zh
Application granted granted Critical
Publication of CN108415017B publication Critical patent/CN108415017B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9021—SAR image post-processing techniques
    • G01S13/9027—Pattern recognition for feature extraction
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9004—SAR image acquisition techniques
    • G01S13/9011—SAR image acquisition techniques with frequency domain processing of the SAR signals in azimuth
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Physics & Mathematics (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • General Physics & Mathematics (AREA)
  • Electromagnetism (AREA)
  • Signal Processing (AREA)
  • Artificial Intelligence (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明公开了一种用于复杂目标雷达散射特性稀疏表征的一维增广状态空间方法，包括步骤：横向极点矩阵的求解——构造一维增广汉克尔矩阵，采用状态空间法求解横向极点矩阵；二维极点自适应配对——最小二乘求解纵向配对矩阵，对配对矩阵逐列进行一维状态空间求解纵向极点参数。同已有方法相比，本发明所提出方法具有更高的特征重构精度和算法运行效率。

Description

复杂目标雷达散射特性稀疏表征的一维增广状态空间方法
技术领域
本发明涉及目标雷达散射特性表征的技术领域，具体涉及复杂目标雷达散射特性稀疏表征的一维增广状态空间方法。
背景技术
电磁散射中心理论表明，目标的雷达散射特性可由少量孤立的散射中心来表征，即SAR成像获取的感兴趣散射中心的数目往往远小于反演这些散射中心所需的采样数据，因此，对复杂目标的雷达散射特性进行稀疏表征，一方面可以将高采样率、大数据量的雷达回波进行降噪压缩重构，缓解传统高分辨率雷达数据传输和储存的压力。另一方面稀疏表征的参数包含目标主散射点的位置、强度信息，这些二维结构特征可作为后续目标识别分类的重要依据。此外，由于稀疏重构的信号模型是对现有回波信号间变化规律的表征，将数据间的递推规律进一步延拓，可用于实现目标的超分辨率成像。
现有的雷达散射特性稀疏表征方法从理论上主要可分为基于现代谱估计理论的方法和基于压缩感知理论的方法，基于谱估计理论方法通常假设雷达回波符合某一特定的散射特征模型，将基于该模型求解出的参数作为稀疏重构的特征量；基于压缩感知理论方法的主要思路在于构造过完备的字典，将原始回波投影在字典空间的测量矢量和过完备字典作为稀疏重构的特征量。研究表明，这两类算法均能较为准确地重构出目标的强散射特征，但仍存在对复杂目标重构效率不高以及重构精度不高等问题，基于此，本发明提出一种适用于二维回波信号的一维增广状态空间法，旨在为提高稀疏表征参数求解效率和模型重构精度上提供一种新的思路。
与本发明相关的现有技术介绍如下：
1现有技术一的技术方案
基于旋转不变技术的二维参数估计法(Two Dimensional Estimation of SignalParameters via Rotation Invariant Technique，2D-ESPRIT)是目前常用于估计二维散射中心参数的方法(参见文献[1]Rouquette S.,Najim M.Estimation of frequenciesand damping factors by two-dimensional ESPRIT type methods[J].IEEETransactions on Signal Processing,2001,49(1):237-245.)。2D-ESPRIT的主要思路在于将基于极点模型的二维散射中心提取问题转化为两个一维问题分别进行估计，利用转换矩阵同时对横向和纵向极点矩阵进行特征值提取实现各向极点配对。该方法的优势在于不但能够对原始雷达回波进行稀疏表征，而且提取的极点对参数能够较为准确地估计出目标强散射中心的位置强度信息，实现高分辨率成像。
现有技术一的缺点
(1)由于2D-ESPRIT在极点参数估计时需构造两个方向的二次增强汉克尔矩阵(约为原始回波数据的MN/4倍)，过大的二次增强汉克尔矩阵导致奇异值分解计算量显著提高，故而整体运行效率不高，
(2)由于该方法的极点幅值估计精度会在横向纵向位置参数重复或相近时有明显下降，导致对散射中心密集分布的目标重构精度不高，且对于复杂目标易出现极点对数目增多时重构结果发散的问题。
2现有技术二的技术方案
二维状态空间法(Two Dimensional State Space Approach，2D-SSA)是另一种用于目标散射特性稀疏表征的方法。该方法基于对二维极点模型的参数进行估计的思路提取散射中心特征。与2D-ESPRIT算法不同的是，2D-SSA(参见文献[2]Piou J.E.,SystemRealization Using 2-D Output Measurements,American Control Conference[C],Jun.2004.2840-2844.)采用二维状态空间的方式分别对纵向和横向的极点进行求解，将两个极点矩阵(纵向、横向)投影至相同子空间的相位作为极点配对的依据。该方法的流程图如图2所示，该方法的优势在于不仅可以用极点对参数对目标雷达散射特性进行稀疏表征，还可以通过系统矩阵递推的方式实现对原始雷达回波的重构，进而实现超分辨成像。
现有技术二的缺点
(1)运行效率不高，同2D-ESPRIT方法一样，2D-SSA在极点参数估计时需构造二次增强汉克尔矩阵，过大的增强汉克尔矩阵导致运行缓慢。
(2)极点配对不准确，由于2D-SSA在极点配对时仅考虑相角因素，易出现相角相近的极点错配的现象，进而导致重构精度下降。
发明内容
本发明所需要解决的技术问题：针对现有复杂目标稀疏表征方法运行效率不高以及极点对配对不准确等问题，本发明提出一种一维增广状态空间方法，包括：(1)构造维数较小的一维增广汉克尔矩阵替代二次增强汉克尔矩阵，使参数求解效率相较于2D-ESPRIT、2D-SSA等算法得到显著提高；(2)采用逐方位自适应搜索的配对策略替代各向极点同时配对，在一定程度上解决横向纵向位置参数重复或相近时极点配对不准确的问题，提高了模型重构精度。
本发明解决上述技术问题采用的技术方案为：一种用于复杂目标雷达散射特性稀疏表征的一维增广状态空间方法，具体的技术方案如下：
(1)横向极点矩阵的求解——构造一维增广汉克尔矩阵，采用状态空间法求解横向极点矩阵；
(2)二维极点自适应配对——最小二乘求解配对矩阵，配对矩阵逐列进行一维状态空间求解纵向极点参数。
其中，所述的状态空间法求解横向极点矩阵，具体步骤如下：
步骤1-1：获取二维雷达回波数据的M×N雷达回波数据Y：
Figure BDA0001656169930000031
其中M代表方位点数，N代表频率点数，y(l,r)表示第l个方位下第r个频点的回波信号样本。
如果目标的散射特性可由其多散射中心来表征，则来自于目标原始二维信号可表示为：
Figure BDA0001656169930000032
式中，y(l,r)表示矩阵Y第l行第r列元素，g表示雷达回波数据中散射中心的数目，极点对(ai,si,pi)表示第i个散射中心的散射幅度与空间位置信息，其中ai表示第i个散射中心的复数幅度信息，si表示第i个散射中心的纵向极点，pi表示第i个散射中心的横向极点，w(l,r)表示信号噪声。
步骤1-2：将雷达回波数据Y的任意第k列的元素定义为Y(k)(k＝1,...,N)，构造一维增广汉克尔矩阵H：
Figure BDA0001656169930000033
Y(k)＝[y(1,k) y(2,k) ... y(M,k)]T (4)
L的取值范围：
Figure BDA0001656169930000041
式中[]表示向下取整。
步骤1-3：对H矩阵进行奇异值分解，得到矩阵中的信号空间分量：
Figure BDA0001656169930000042
式中上标*表示共轭转置，Usn、
Figure BDA0001656169930000043
表示信号空间分量，Uno、
Figure BDA0001656169930000044
表示噪声空间分量，式中信号空间分量Usn、
Figure BDA0001656169930000045
的模型阶MDL(x)选用最小描述长度准则进行判别：
Figure BDA0001656169930000046
步骤1-4：根据信号空间分量计算观察矩阵
Figure BDA0001656169930000047
Figure BDA0001656169930000048
步骤1-5：根据式(8)求解系统矩阵P：
Figure BDA0001656169930000049
其中：
Figure BDA00016561699300000410
Figure BDA00016561699300000411
步骤1-6：对系统矩阵P做特征值分解：
P＝MΛM-1 (12)
其中，M为特征向量矩阵，M-1为对应逆矩阵，Λ为特征值矩阵。
步骤1-7：获取特征值矩阵Λ的主对角元素[p1 p2 ... pn](n与式(7)中信号分量阶数一致)，求解横向极点矩阵：
Figure BDA00016561699300000412
其中，所述的配对矩阵逐列进行一维状态空间求解纵向极点参数，具体步骤如下：
步骤2-1：根据(13)所得横向极点矩阵Q和原始雷达回波矩阵Y，最小二乘求解配对矩阵K：
K＝YQ*(QQ*)-1 (14)
步骤2-2：对矩阵K逐列进行一维状态空间求解，以j＝1列为例，令：
T＝K(:,j) (15)
步骤2-3：构造矩阵T的汉克尔矩阵Ht：
Figure BDA0001656169930000051
式中，[]表示向下取整。
步骤2-4：同步骤1-3，对矩阵Ht进行奇异值分解，获取信号空间分量：
Ht＝URV* (17)
Ut＝U(:,1:e) (18)
Rt＝R(1:e,1:e) (19)
Vt＝V(:,1:e) (20)
式中上标*表示共轭转置，U、R、V为奇异值分解矩阵，Ut、Rt、Vt为信号空间分量，e为式(7)准则判别后的模型阶数。
步骤2-5：根据式(18-20)计算观察矩阵
Figure BDA0001656169930000052
和控制矩阵
Figure BDA0001656169930000053
步骤2-6：根据观察矩阵计算得到矩阵A和C：
Figure BDA0001656169930000054
Figure BDA0001656169930000055
式中
Figure BDA0001656169930000056
表示删除
Figure BDA0001656169930000057
第e行后的矩阵，
Figure BDA0001656169930000058
表示删除
Figure BDA00016561699300000513
第一行后的矩阵。
步骤2-7：构造矩阵
Figure BDA00016561699300000510
Figure BDA00016561699300000511
解方程：
Figure BDA00016561699300000512
最小二乘求解式(24)可以得到：
Figure BDA0001656169930000061
步骤2-8：对矩阵A进行特征值分解：
Figure BDA0001656169930000062
步骤2-9：根据式(22)、(25)、(26)中得到的矩阵C、B、M1和Λ1，计算系数向量A1和纵向极点向量S1：
Figure BDA0001656169930000063
S1＝[Λ1(1,1) Λ1(2,2) ... Λ1(e,e)] (28)
式中·表示点乘运算。
根据式(13)、(27)、(28)，完成横向极点pj的纵向极点配对：
[(A1(1),S1(1),pj) (A1(2),S1(2),pj) ... (A1(e),S1(e),pj)] (29)
步骤2-10：依次令j＝2,...,n，重复步骤202-209，直至求解出所有极点对(ai,si,pi)(式2所示，i＝1,2...,g)。对于每一极点对，所表征的散射中心位置为：
Figure BDA0001656169930000064
其中，Xi,Yi分别表示散射中心的横向和纵向距离，angle()表示幅角，Δθ表示角度间隔，Δf表示步进频率间隔，c＝3×108m/s表示光速。
根据极点对(ai,si,pi)重构的频域数据为：
Figure BDA0001656169930000065
式中
Figure BDA0001656169930000066
表示第l个方位角下第r个频点的目标回波信号样本。
步骤1-3判别信号分量模型阶数的最小长度准则可采用奇异值分析法、信息量准则法等方法和准则替换。
本发明技术方案带来的有益效果为：
(1)更高的运行效率。由于将一维增广汉克尔矩阵替代二次增强汉克尔矩阵，并在极点配对环节将求解高维数增强矩阵的奇异值分解问题转化为若干个维数较小矩阵依次求解的问题，所提方法的运行效率约为2D-ESPRIT、2D-SSA等方法的10倍以上。
(2)更高的极点配对正确率，由于在极点配对环节采用最小二乘自适应配对的方法，较好解决了横向纵向位置参数重复或相近时极点配对正确率下降的问题。
(3)更强的复杂目标稀疏表征鲁棒性，由于求解的极点对均为逐方位自适应搜索的结果，本发明所提方法解决了复杂目标随着极点对阶数升高稀疏表征结果不稳定的问题。
附图说明
图1为本发明提出的一维增广状态空间方法流程图；
图2为2D-SSA方法的流程图。
图3为典型复杂目标飞机ISAR数据进行2DFFT成像算法得到的二维距离像(动态方位50dB)；
图4为三种方法估计的散射中心位置示意图(极点对数目均取90)，其中，图4(a)为2D-SSA估计的散射中心位置示意图；图4(b)为2D-ESPRIT估计的散射中心位置示意图；图4(c)为一维增广状态空间法估计的散射中心位置示意图；
图5为极点对数目均增至307阶时三种方法估计的散射中心位置示意图：其中，图5(a)为原始图像；图5(b)为一维增广状态空间法估计的散射中心位置示意图；图5(c)为2D-SSA估计的散射中心位置示意图，图5(d)为2D-ESPRIT估计的散射中心位置示意图；
图6为三种方法重构的一维距离像的比较结果(极点对数目均取90)。
具体实施方式
下面结合附图和具体实例对本发明做进一步说明。
一种用于复杂目标雷达散射特性稀疏表征的一维增广状态空间方法，以实测典型复杂目标飞机的逆合成孔径雷达(ISAR)回波数据为例，其中数据频段选用2.25GHz～3.75GHz，步进频率为6MHz，方位角范围-5°-5°(鼻锥向)，方位步长为0.1°。实例选用2D-SSA、2D-ESPRIT两种方法作对比验证，验证范围包括极点对位置，重构精度以及运行效率。具体的技术方案如下：
(1)横向极点矩阵的求解——构造一维增广汉克尔矩阵，采用状态空间法求解横向极点矩阵；
步骤1-1：输入典型复杂目标的二维雷达回波101(方位)×251(频率)；
步骤1-2：根据式(3)(4)，构造一维增广汉克尔矩阵H(L取[N/2]＝125)，H矩阵大小为12827×125；
步骤1-3：对H矩阵进行奇异值分解，最小描述长度准则判别信号分量阶数为63；
步骤1-4：根据式(8)计算观察矩阵
Figure BDA0001656169930000071
步骤1-5：根据式(9)～(11)求解系统矩阵P；
步骤1-6：根据式(12)对系统矩阵P做特征值分解；
步骤1-7：根据式(13)求解横向极点矩阵Q；
(2)二维极点自适应配对——最小二乘求解配对矩阵，配对矩阵逐列进行一维状态空间求解纵向极点参数。
步骤2-1：根据式(14)求解配对矩阵K；
步骤2-2：对矩阵K逐列进行一维状态空间求解，式(15)得到T向量；
步骤2-3：根据式(16)构造矩阵T的汉克尔矩阵Ht(矩阵大小为50×50)；
步骤2-4：根据式(17)～(20)对矩阵Ht进行奇异值分解，获取信号空间分量，最小描述长度准则判别信号分量阶数为2；
步骤2-5：根据信号分量计算观察矩阵
Figure BDA0001656169930000081
和控制矩阵
Figure BDA0001656169930000082
步骤2-6：根据式(21)～(22)计算矩阵A和C；
步骤2-7：根据式(23)～(25)，最小二乘求解矩阵B；
步骤2-8：对矩阵A进行特征值分解，获取对应的特征值和特征值向量；
步骤2-9：根据式(27)～(28)计算系数向量A1和纵向极点向量S1，根据式(29)完成横向极点pj的纵向极点配对；
步骤2-10：依次令j＝2,...,n，重复步骤202～209，直至求解出所有横纵向极点对(ai,si,pi)。根据式(30)求解所有极点对对应的散射中心位置(Xi,Yi)，根据式(31)得到极点对重构的频域数据
Figure BDA0001656169930000083
实施例结果如下：
(1)极点对位置对比分析：图4(a)、图4(b)、图4(c)分别表示三种方法重构极点在目标ISAR像的几何位置关系(图中黑点为估计的散射中心位置)，结果表明：由于实测目标包含较为密集的散射中心分布，2D-SSA和2D-ESPRIT均能够较为准确地配对目标主要散射中心的位置(图4中等高线颜色最深的位置)，而对于密集分布且相对较弱的散射中心，由于强度接近且相位相近，这两种算法均出现明显的极点错配和漏配现象。对比这两种方法，由于采用逐列自适应搜索的配对策略，本发明所提方法极点配对的正确率显著高于前两种。图5为极点对数目增至307阶的结果，图中可以看出随着极点对数的进一步增多，本发明提出的方法会优先寻找较弱的散射中心，所配对极点的结果仍然正确，而2D-SSA错配的极点数目显著增多，2D-ESPRIT相对于2D-SSA配对正确率更高，但也会出现部分错配和漏配的现象。
(2)重构精度分析：表1为三种方法重构精度的重构误差均值REM和图像相似度β评价结果(计算公式如式(32)(33)所示)，对比表中的数据，由于错配和漏配极点带来的误差积累，2D-SSA和2D-ESPRIT重构误差均高于一维增广状态空间法，进而导致重构二维距离像的图像相似度显著低于本发明所提方法，图6为三种算法重构方位角为0°的一维距离像对比，可以看出，错配或漏配极点导致2D-SSA和2D-ESPRIT重构主散射中心出现不同程度相位偏移和幅值下降，相较于这两种方法，本发明重构的RCS主散射中心和原始数据的拟合程度更高，但也会存在较弱散射中心上(低于40dB)出现差异的情况。
Figure BDA0001656169930000091
Figure BDA0001656169930000092
式中REM为重构误差均值，y(i,j)表示原始频域数据，y(i,j)为不同方法根据估计极点对重构频域数据。β为图像相似系数，x(i,j)为y(i,j)经二维傅里叶变换得到的原始时域数据，x(i,j)为y(i,j)经二维傅里叶变换得到的估计值。
表1 三种方法的重构误差和图像相似度比较
Figure BDA0001656169930000093
(3)运行效率：表2所示为三种方法的运行效率的对比，表中本文算法的运行效率显著高于2D-SSA和2D-ESPRIT两种方法，原因在于，2D-SSA和2D-ESPRIT均需对原始数据矩阵进行二次增强，二次增强后的汉克尔矩阵约为本发明采用的一维增广汉克尔矩阵的M/4倍，在极点配对环节，本发明采用逐列最小二乘自适应配对的方式，将2D-SSA和2D-ESPRIT需求解的纵向二次增强矩阵转化为若干个维数较小矩阵依次求解。
表2 三种方法的运行效率对比
Figure BDA0001656169930000094

Claims (2)

1.一种复杂目标雷达散射特性稀疏表征的一维增广状态空间方法，其特征在于，该方法的具体如下：
(1)横向极点矩阵的求解——构造一维增广汉克尔矩阵，采用状态空间法求解横向极点矩阵；
(2)二维极点自适应配对——最小二乘求解配对矩阵，配对矩阵逐列进行一维状态空间求解纵向极点参数；
其中，所述的状态空间法求解横向极点矩阵，具体步骤如下：
步骤1-1：获取二维雷达回波数据的M×N雷达回波数据Y：
Figure FDA0003111786680000011
其中M代表方位点数，N代表频率点数；y(l,r)表示第l个方位下第r个频点的回波信号样本；
如果目标的散射特性可由其多散射中心来表征，则来自于目标原始二维信号可表示为：
Figure FDA0003111786680000012
式中，y(l,r)表示矩阵Y第l行第r列元素，g表示雷达回波数据中散射中心的个数，极点对(ai,si,pi)表示第i个散射中心的散射幅度与空间位置信息，其中ai表示第i个散射中心的复数幅度信息，si表示第i个散射中心的纵向极点，pi表示第i个散射中心的横向极点，w(l,r)表示信号噪声；
步骤1-2：将雷达回波数据Y的任意第k列的元素定义为Y(k)，k＝1,...,N，构造一维增广汉克尔矩阵H：
Figure FDA0003111786680000013
Y(k)＝[y(1,k) y(2,k) ... y(M,k)]T (4)
L的取值范围：
Figure FDA0003111786680000021
式中[]表示向下取整；
步骤1-3：对H矩阵进行奇异值分解，得到矩阵中的信号空间分量：
Figure FDA0003111786680000022
式中上标*表示共轭转置，Usn、
Figure FDA0003111786680000023
表示信号空间分量，Uno、
Figure FDA0003111786680000024
表示噪声空间分量，式中信号空间分量Usn、
Figure FDA0003111786680000025
的模型阶数MDL(x)选用最小描述长度准则进行判别：
Figure FDA0003111786680000026
式中x＝0,1,...,N-1，lm为H矩阵的奇异值，m＝1,2,...N，Π和∑分别为求积和求和运算；
步骤1-4：根据信号空间分量计算观察矩阵
Figure FDA0003111786680000027
Figure FDA0003111786680000028
步骤1-5：根据式(8)求解系统矩阵P：
Figure FDA0003111786680000029
其中，
Figure FDA00031117866800000210
Figure FDA00031117866800000211
步骤1-6：对系统矩阵P做特征值分解：
P＝MΛM-1 (12)
其中，M为特征向量矩阵，M-1为对应逆矩阵，Λ为特征值矩阵；
步骤1-7：获取特征值矩阵Λ的主对角元素[p1 p2 … pn]，n与式(7)中信号分量阶数一致，求解横向极点矩阵：
Figure FDA00031117866800000212
2.根据权利要求1所述的一种复杂目标雷达散射特性稀疏表征的一维增广状态空间方法，其特征在于：所述二维极点自适应配对，即通过配对矩阵逐列进行一维状态空间求解纵向极点参数，具体步骤如下：
步骤2-1：根据(13)所得横向极点矩阵Q和原始雷达回波矩阵Y，最小二乘求解配对矩阵K：
K＝YQ*(QQ*)-1 (14)
步骤2-2：对矩阵K逐列进行一维状态空间求解，以j＝1列为例，令：
T＝K(:,j) (15)
步骤2-3：构造矩阵T的汉克尔矩阵Ht：
Figure FDA0003111786680000031
式中，[]表示向下取整；
步骤2-4：同步骤1-3，对矩阵Ht进行奇异值分解，获取信号空间分量：
Ht＝URV* (17)
Ut＝U(:,1:e) (18)
Rt＝R(1:e,1:e) (19)
Vt＝V(:,1:e) (20)
式中上标*表示共轭转置，U、R、V为奇异值分解矩阵，Ut、Rt、Vt为信号空间分量，e为式(7)准则判别后的模型阶数；
步骤2-5：根据式(18)至(20)计算观察矩阵
Figure FDA0003111786680000032
和控制矩阵
Figure FDA0003111786680000033
步骤2-6：根据观察矩阵计算得到矩阵A和C：
Figure FDA0003111786680000034
Figure FDA0003111786680000035
式中
Figure FDA0003111786680000036
表示删除
Figure FDA0003111786680000037
第e行后的矩阵，
Figure FDA0003111786680000038
表示删除
Figure FDA0003111786680000039
第一行后的矩阵；
步骤2-7：构造矩阵
Figure FDA00031117866800000310
Figure FDA00031117866800000311
解方程：
Figure FDA0003111786680000041
最小二乘求解式(24)可以得到：
Figure FDA0003111786680000042
步骤2-8：对矩阵A进行特征值分解：
Figure FDA0003111786680000043
步骤2-9：根据式(22)、(25)、(26)中得到的矩阵C、B、M1和Λ1，计算系数向量A1和纵向极点向量S1：
Figure FDA0003111786680000044
S1＝[Λ1(1,1) Λ1(2,2) ... Λ1(e,e)] (28)
式中.表示点乘运算；
根据(13)、(27)、(28)，完成横向极点pj的纵向极点配对：
[(A1(1),S1(1),pj) (A1(2),S1(2),pj) ... (A1(e),S1(e),pj)] (29)
步骤2-10：依次令j＝2,...,n，重复步骤2-2至步骤2-9，直至求解出所有极点对(ai,si,pi)，i＝1,2...,g，对于每一极点对，所表征的散射中心位置为：
Figure FDA0003111786680000045
其中，Xi,Yi分别表示散射中心的横向和纵向距离，angle()表示幅角，Δθ表示角度间隔，Δf表示步进频率间隔，c＝3×108m/s表示光速；
根据极点对(ai,si,pi)重构的频域数据为：
Figure FDA0003111786680000046
式中
Figure FDA0003111786680000047
表示第l个方位角下第r个频点的目标回波信号样本。
CN201810441778.4A 2018-05-10 2018-05-10 复杂目标雷达散射特性稀疏表征的一维增广状态空间方法 Active CN108415017B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810441778.4A CN108415017B (zh) 2018-05-10 2018-05-10 复杂目标雷达散射特性稀疏表征的一维增广状态空间方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810441778.4A CN108415017B (zh) 2018-05-10 2018-05-10 复杂目标雷达散射特性稀疏表征的一维增广状态空间方法

Publications (2)

Publication Number Publication Date
CN108415017A CN108415017A (zh) 2018-08-17
CN108415017B true CN108415017B (zh) 2021-07-30

Family

ID=63138388

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810441778.4A Active CN108415017B (zh) 2018-05-10 2018-05-10 复杂目标雷达散射特性稀疏表征的一维增广状态空间方法

Country Status (1)

Country Link
CN (1) CN108415017B (zh)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109444842B (zh) * 2019-01-04 2020-09-25 北京环境特性研究所 一种目标电磁散射特性数据重构方法和装置
CN110098488B (zh) * 2019-05-17 2020-11-24 西安电子科技大学 一种基于特征模理论的低rcs超表面的模式转化方法
CN111308438B (zh) * 2020-03-06 2021-11-02 中国人民解放军海军航空大学 一种回波数据散射特征提取方法及系统

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6529794B1 (en) * 1997-08-22 2003-03-04 Siemens Aktiengesellschaft Method and device for measuring distance and speed
CN101776748A (zh) * 2009-01-12 2010-07-14 哈尔滨工业大学 复杂目标极点特征提取方法和装置
CN103412306A (zh) * 2013-07-10 2013-11-27 北京航空航天大学 基于多站多波段带宽融合的超宽带相参处理方法
CN106443629A (zh) * 2016-09-30 2017-02-22 西安空间无线电技术研究所 一种基于傅里叶基函数的雷达目标识别方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6529794B1 (en) * 1997-08-22 2003-03-04 Siemens Aktiengesellschaft Method and device for measuring distance and speed
CN101776748A (zh) * 2009-01-12 2010-07-14 哈尔滨工业大学 复杂目标极点特征提取方法和装置
CN103412306A (zh) * 2013-07-10 2013-11-27 北京航空航天大学 基于多站多波段带宽融合的超宽带相参处理方法
CN106443629A (zh) * 2016-09-30 2017-02-22 西安空间无线电技术研究所 一种基于傅里叶基函数的雷达目标识别方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
A Robust State Space Model for the Characterization of Extended Returns in Radar Target Signatures;Krishna Naishadham 等;《IEEE TRANSACTIONS ON ANTENNAS AND PROPAGATION》;20080630;第56卷(第6期);1742-1751 *
基于状态空间模型(SSM)的雷达目标二维散射中心参数提取;袁仕继 等;《舰船电子工程》;20100831;第30卷(第8期);54-56 *

Also Published As

Publication number Publication date
CN108415017A (zh) 2018-08-17

Similar Documents

Publication Publication Date Title
CN108415017B (zh) 复杂目标雷达散射特性稀疏表征的一维增广状态空间方法
CN110780298B (zh) 基于变分贝叶斯学习的多基isar融合成像方法
CN104076360B (zh) 基于压缩感知的二维sar稀疏目标成像方法
CN107085835B (zh) 基于四元数加权核范数最小的彩色图像去噪方法
CN107493106B (zh) 一种基于压缩感知的频率和角度联合估计的方法
CN110728706A (zh) 基于深度学习的sar图像精细配准方法
CN116908853B (zh) 高相干点选取方法、装置和设备
Hu et al. Inverse synthetic aperture radar imaging exploiting dictionary learning
CN110596706B (zh) 一种基于三维图像域投射变换的雷达散射截面积外推方法
CN109543717B (zh) 基于自适应邻域及字典的联合协作表达高光谱分类方法
CN112230221A (zh) 一种基于三维稀疏成像的rcs测量方法
CN108931770B (zh) 基于多维贝塔过程线性回归的isar成像方法
Rößing et al. Multi-antenna SAR tomography using superresolution techniques
CN107861125B (zh) 基于广义帕累托分布的贝叶斯压缩感知成像方法
CN107561536B (zh) 压缩感知逆合成孔径雷达偏离栅格的快速成像方法
CN112684445B (zh) 基于md-admm的mimo-isar三维成像方法
CN113240081B (zh) 针对雷达载频变换的高分辨距离像目标稳健识别方法
CN113093164B (zh) 一种平移不变和噪声稳健的雷达图像目标识别方法
CN111538003B (zh) 一种单比特压缩采样合成孔径雷达成像方法
CN112965066A (zh) 一种单比特带通采样合成孔径雷达成像方法
CN110263777B (zh) 基于空谱结合的局部保持投影算法的目标检测方法及系统
CN109633635B (zh) 基于结构化递归最小二乘的米波雷达测高方法
CN107704724B (zh) 基于Meridian分布的贝叶斯压缩感知的参数选取方法
CN112330578A (zh) 一种压缩感知与洁化算法结合的图像重建方法及系统
Zeng et al. Unsupervised 3D array-SAR imaging based on generative model for scattering diagnosis

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
GR01 Patent grant
GR01 Patent grant