本发明公开了一种基于稀疏贝叶斯学习的干扰环境逆合成孔径雷达ISAR高分辨成像方法，本发明的实现步骤是：(1)录取目标有效回波；(2)对有效回波进行距离压缩并转置；(3)生成字典；(4)稀疏贝叶斯建模；(5)利用最大后验‑期望最大MAP‑EM算法求解稀疏贝叶斯模型，计算存在回波的每个距离单元的权向量；(6)重构权值矩阵；(7)转置并得到得到干扰环境下的二维高分辨逆合成孔径雷达ISAR成像结果。本发明实现了在干扰环境下聚焦良好的高分辨二维逆合成孔径雷达ISAR成像，可用于强干扰等复杂电磁环境下对空间与空中非合作目标的高分辨二维成像。
1.一种基于稀疏贝叶斯学习的干扰环境逆合成孔径雷达ISAR高分辨成像方法，其特征在于，包括如下步骤：(1)录取目标有效回波：
在干扰环境下，逆合成孔径雷达ISAR向运动目标发射线性调频信号，获取干扰环境下所发射线性调频信号的缺损回波矩阵Sr，缺损回波矩阵Sr的维数为Nr×Na，其中，Nr表示缺损回波矩阵Sr的距离向采样点数，Na表示缺损回波矩阵Sr的方位向采样点数，剔除缺损回波矩阵Sr中缺损的列向量，得到维数为Nr×N1的有效回波矩阵Se，其中N1表示方位向有效样本数；
(2a)对有效回波矩阵Se沿距离向进行脉冲压缩，得到距离向脉冲压缩后的矩阵S；
(2b)对距离向脉冲压缩后的矩阵S进行转置操作，得到转置矩阵S′；
(3a)以 为元素,构造维数为Na×Na的傅里叶字典Φ，其中e(·)表示以自然常数为底的指数操作，j表示虚数单位，π表示圆周率，m表示傅里叶字典Φ的行序号，n表示傅里叶字典Φ的列序号，行序号和列序号的取值范围均为[-Na/2,Na/2-1]；
(3b)剔除傅里叶字典Φ的第t1,…,tp行，得到维数为N1×Na的有效傅里叶字典Φe，其中t1,…,tp表示缺损回波矩阵Sr中缺损列向量对应的列序号；
(3c)利用公式 构造实傅里叶字典Φ′，其中Re(·)表示取实部操作，Im(·)表示取虚部操作；
(4a)按照式(1)，构建每个距离单元的稀疏信号模型：s′f＝Φeω′q+ε′ (1)
其中，s′f表示转置矩阵S′的第f个列向量，Φe表示有效傅里叶字典，ω′q表示第q个距离单元的复权向量，ε′表示复干扰向量，f与q的取值相同；
(4c)利用伽马Gamma-高斯Gaussian分层模型对每个距离单元的权向量建模；
(5)利用最大后验-期望最大MAP-EM算法求解稀疏贝叶斯模型，计算存在回波的每个距离单元的权向量：(5a)设置h1为距离脉冲压缩后的矩阵S中第一个存在回波的距离单元序号，h2为距离脉冲压缩后的矩阵S中最后一个存在回波的距离单元序号，距离单元的初始序号q＝h1；
(5b)设置初始迭代次数k＝1，初始化模型先验参数：在第q个距离单元中，令干扰向量的聚类数初始值为10，令类别矩阵初始值的维数为2N1×10，矩阵每行中仅有一个元素为1而其它元素为0，非零元素序号在整数[0,10]中随机产生，令精度向量所有元素的初始值为
0.01，令均值向量所有元素的初始值随机产生，其中精度向量和均值向量的元素个数均为
10，令精度矩阵的初始值为对角线元素全为100的对角阵，维数为2Na×2Na，令混合系数的超参数全为1，精度向量每个元素的超参数全为10-4，精度矩阵每个对角线元素的超参数全为
(5c)利用迭代公式，依次计算第k次迭代的各个参数的值；
(5d)以类别矩阵L(k)每一行中最大元素对应的列序号为元素构造标签向量(k)
其中 为正整数，表示类别矩阵L 第1,…,2N1行中最大元素对应的聚类序号；
(5e)删除类别矩阵L 中第l1,…,lr个列向量，l1,…,lr表示标签向量 中重复元素对应的序号，类别矩阵L(k)剩余列的个数为干扰向量ε的聚类数J(k)，对类别矩阵L(k)按行进行归一化操作：将类别矩阵L(k)的每个元素除以所在行元素之和，删除混合系数(k) (k) (k)π 、均值向量μ 和精度向量α 中第l1,…,lr个元素；
(5f)判断终止函数E是否小于终止阀值η＝10-3，若是，则停止迭代，执行步骤(5g)；否则，更新迭代次数k＝k+1，执行步骤(5c)；
(5g)更新距离单元序号q＝q+1，判断更新后的距离单元序号q是否大于h2，若是，则停止对距离单元的搜索，执行步骤(6)；否则，执行步骤(5b)；
生成维数为2Na×Nr的权值矩阵W，权值矩阵W的第y列为第q个距离单元的权向量ωq，距离单元序号q满足q＝h1,...,h2，y与q的取值相同，其余元素全部为零，再利用公式W′＝W1+jW2生成复权值矩阵W′，其中W1表示权值矩阵W第1,…,Na行元素构成的矩阵，W2表示权值矩阵W第Na+1,…,2Na行元素构成的矩阵；
对复权值矩阵W′进行转置操作，得到干扰环境下的二维高分辨逆合成孔径雷达ISAR成像结果。
2.根据权利要求1所述基于稀疏贝叶斯学习的干扰环境逆合成孔径雷达ISAR高分辨成像方法，其特征在于，步骤(2a)中所述的对有效回波矩阵Se沿距离向进行脉冲压缩是按照如下步骤实现的：第一步，将逆合成孔径雷至场景中心的距离作为参考距离，将与逆合成孔径雷达发射信号载频、调频率相同，距离为参考距离的线性调频信号作为参考信号Sref；
第二步，将参考信号Sref取共轭后与接收的有效回波矩阵Se相乘，得到解线频调后的信号Srd；
第三步，对解线频调后的信号Srd沿快时间维作傅里叶变换，得到距离脉冲压缩后的矩阵S。
3.根据权利要求1所述基于稀疏贝叶斯学习的干扰环境逆合成孔径雷达ISAR高分辨成像方法，其特征在于，步骤(4b)中所述利用高斯混合模型对干扰向量建模的具体步骤如下：第一步，将复干扰向量ε′的实部和虚部拼接为干扰向量ε，其中实部在上，虚部在下；
第二步，干扰向量ε中每个元素的概率分布按照式(2)给出：其中，p(εu|π,μ,α)表示干扰向量ε第u个元素的概率密度，J表示干扰向量ε的聚类数，Σ表示求和操作， 表示混合系数π的第j1个元素， 表示均值向量μ的第j2个元素， 表示精度向量α的第j3个元素，N(·)表示高斯分布的概率密度，j1、j2、j3的取值相同；
第三步，精度向量α中每个元素的先验分布按照式(3)给出：其中， 表示精度向量α第j3个元素的概率密度，Gam(·)表示伽马Gamma分布的概率密度，a和b表示精度向量α第j3个元素 的超参数；
p(π)＝Dir(π|γ1,...,γJ) (4)其中，p(π)表示混合系数π的概率密度，Dir(·)表示狄利克雷Dirichlet分布的概率密度，γ1,...,γJ表示混合系数π的超参数。
4.根据权利要求1所述基于稀疏贝叶斯学习的干扰环境逆合成孔径雷达ISAR高分辨成像方法，其特征在于，步骤(4c)中所述的利用伽马Gamma-高斯Gaussian分层模型对每个距离单元的权向量建模的具体步骤如下：第一步，将每个距离单元复权向量的实部和虚部拼接为对应的权向量，其中实部在上，虚部在下；
第二步，每个距离单元权向量的先验分布按照式(5)给出：p(ωq)＝N(0,Λ-1) (5)
其中，p(ωq)表示第q个距离单元权向量的概率密度，精度矩阵Λ为对角阵；
第三步，精度矩阵中每个对角线元素的先验分布按照式(6)给出：p(λd)＝Gam(λd|e,f) (6)
其中，p(λd)表示精度矩阵Λ第d个元素的概率密度，d＝1,…,2Na，e和f表示精度矩阵Λ第d个对角线元素λd的超参数。
5.根据权利要求1所述基于稀疏贝叶斯学习的干扰环境逆合成孔径雷达ISAR高分辨成像方法，其特征在于，步骤(5c)中所述利用迭代公式，依次计算第k次迭代的各个参数值的具体步骤如下：第一步，按照式(7)，计算干扰向量每个聚类的有效数目：其中， 表示干扰向量ε第c个聚类的有效数目，N1表示方位向有效样本数，Σ表示求和操作， 表示类别矩阵L(k-1)第i行的第j个元素，j＝1,…,J(k-1)，J(k-1)表示干扰向量ε的聚类数，c与j的取值相同；
其中， 表示混合系数π 的第j1个元素， 表示混合系数π的超参数，Na表示缺损回波矩阵Sr的方位向采样点数，j1、j4与j的取值相同；
其中， 表示精度向量α(k-1)的第j3个元素， 表示实傅里叶字典Φ′转置矩阵的第v个列向量，Λ(k)表示精度矩阵， 表示第f个距离单元回波向量sf的第z个元素， 表示均值向量μ(k-1)的第j2个元素，[·]-1表示矩阵求逆操作，[·]T表示转置操作，v、z与i的取值相同，f与q的取值相同，j2、j3与j的取值相同；
其中，a和b表示精度向量α第j3个元素 的超参数；
第七步，按照式(13)，计算精度矩阵的每个对角线元素：其中， 表示精度矩阵Λ(k)的第d个对角线元素，e和f表示精度矩阵Λ第d个对角线元素λd的超参数， 表示第q个距离单元权向量 的第x个元素，d＝1,…,2Na，x与d的取值相同；
6.根据权利要求1所述基于稀疏贝叶斯学习的干扰环境逆合成孔径雷达ISAR高分辨成像方法，其特征在于，步骤(5f)中所述的终止函数E如式(14)所示：其中，下界函数 的表达式如式(15)所示：
其中，N1表示方位向有效样本数，J(k)表示干扰向量ε的聚类数， 表示类别矩阵L(k) (k) (k)第i行的第j个元素， 表示混合系数π 的第j1个元素， 表示精度向量α 的第j3个元素， 表示第f个距离单元回波向量sf的第z个元素， 表示实傅里叶字典Φ′转置矩阵的第v个列向量， 表示第q个距离单元的权向量， 表示均值向量μ(k)的第j2个元素，Λ(k)表示精度矩阵，a和b表示精度向量α第j3个元素 的超参数， 表示混合系数π的超参数，In表示以自然常数为底的对数操作，v、z与i的取值相同，j1、j2、j3、j4与j的取值相同，f与q的取值相同。
基于稀疏贝叶斯学习的干扰环境ISAR高分辨成像方法
技术领域
[0001] 本发明属于信息技术领域，更进一步涉及雷达信号处理技术领域中的一种基于稀疏贝叶斯学习的干扰环境逆合成孔径雷达ISAR(Inverse Synthetic Aperture Radar)高分辨成像方法。本发明可用于强干扰等复杂电磁环境下对空间与空中非合作目标的高分辨二维成像。
背景技术
[0002] 由于具有全天时、全天候、高分辨率和远距离等特点，高分辨率逆合成孔径雷达ISAR在航空与航天目标观测中发挥着重要作用。目前，目标所处的电磁环境愈加复杂，而干扰种类与强度也日益增加，甚至会完全遮盖目标而无法成像。此外，当受雷达工作模式限制而无法对目标进行连续观测时，会产生方位缺损回波。因此，回波缺损、强干扰情况下的空间目标高分辨成像是提高现有雷达对空间、空中目标探测与监视能力的关键技术。
[0003] 哈尔滨工业大学在其申请的发明专利文献“基于参数估计的压缩感知成像方法”(公开号：CN106772375A，申请号：201611226167.5)中公开了一种基于参数估计的压缩感知成像方法。该方法的具体步骤为：对回波数据进行距离压缩和运动补偿，运用压缩感知原理进行图像对比度检索并进行信号参数估计，依据所得估计参数进行目标高分辨成像。该方法虽然能够在回波缺损情况下实现高分辨成像，但该方法存在的不足之处是，针对L1范数的数值优化方法运算量较大，在干扰环境下无法高效地获得聚焦良好的ISAR像。
[0004] 成萍，司锡才，姜义成，许荣庆在其发表的论文“基于稀疏贝叶斯学习的稀疏信号表示ISAR成像方法”(《电子学报》2008，36(3)：547-550)中提出一种基于稀疏贝叶斯学习的稀疏信号表示ISAR成像方法。该方法基于稀疏信号表示理论，将ISAR高分辨成像问题转化为稀疏信号表示问题，假设噪声服从高斯分布，利用EM算法求解模型参数，最终实现高分辨二维ISAR成像。该方法虽然能够在干扰环境下实现二维ISAR成像，但该方法存在的不足之处是，高斯模型不适用于某些干扰信号，由于建模不准确而无法充分利用环境先验信息，在强干扰环境下无法获得聚焦良好的ISAR像。
发明内容
[0005] 本发明针对上述现有技术中对L1范数的数值优化方法运算量较大问题，对干扰信号建模不准确，无法充分利用环境先验信息的局限性，提出一种基于稀疏贝叶斯学习的干扰环境逆合成孔径雷达ISAR高分辨成像方法。该方法通过概率建模引入目标与环境的先验信息，进而对模型参数进行有效求解，最终实现强干扰等复杂电磁环境下的空间、空中目标高分辨逆合成孔径雷达ISAR成像。
[0006] 为实现本发明的目的，本发明的技术思路是：基于稀疏信号重构理论，将ISAR高分辨成像问题转化为稀疏线性回归问题，利用高斯混合模型对未知干扰信号进行建模，进而采用最大后验-期望最大MAP-EM算法求解权值向量，实现干扰环境下的高分辨二维成像。
[0007] 为实现上述目的，本发明的主要步骤如下：
[0008] (1)录取目标有效回波：
[0009] 在干扰环境下，逆合成孔径雷达ISAR向运动目标发射线性调频信号，获取干扰环境下所发射线性调频信号的缺损回波矩阵Sr，缺损回波矩阵Sr的维数为Nr×Na，其中，Nr表示缺损回波矩阵Sr的距离向采样点数，Na表示缺损回波矩阵Sr的方位向采样点数，剔除缺损回波矩阵Sr中缺损的列向量，得到维数为Nr×N1的有效回波矩阵Se，其中N1表示方位向有效样本数；
[0010] (2)距离压缩并转置：
[0011] (2a)对有效回波矩阵Se沿距离向进行脉冲压缩，得到距离向脉冲压缩后的矩阵S；
[0012] (2b)对距离向脉冲压缩后的矩阵S进行转置操作，得到转置矩阵S′；
[0013] (3)生成字典：
[0014] (3a)以 为元素,构造维数为Na×Na的傅里叶字典Φ，其中e(·)表示以自然常数为底的指数操作，j表示虚数单位，π表示圆周率，m表示傅里叶字典Φ的行序号，n表示傅里叶字典Φ的列序号，行序号和列序号的取值范围均为[-Na/2,Na/2-1]；
[0015] (3b)剔除傅里叶字典Φ的第t1,…,tp行，得到维数为N1×Na的有效傅里叶字典Φe，其中t1,…,tp表示缺损回波矩阵Sr中缺损列向量对应的列序号；
[0016] (3c)利用公式 构造实傅里叶字典Φ′，其中Re(·)表示取实部操作，Im(·)表示取虚部操作；
[0017] (4)稀疏贝叶斯建模：
[0018] (4a)按照式(1)，构建每个距离单元的稀疏信号模型：
[0019] s′f＝Φeω′q+ε′ (1)
[0020] 其中，s′f表示转置矩阵S′的第f个列向量，Φe表示有效傅里叶字典，ω′q表示第q个距离单元的复权向量，ε′表示复干扰向量，f与q的取值相同；
[0021] (4b)利用高斯混合模型对干扰向量建模；
[0022] (4c)利用伽马Gamma-高斯Gaussian分层模型对每个距离单元的权向量建模；
[0023] (5)利用最大后验-期望最大MAP-EM算法求解稀疏贝叶斯模型，计算存在回波的每个距离单元的权向量：
[0024] (5a)设置h1为距离脉冲压缩后的矩阵S中第一个存在回波的距离单元序号，h2为距离脉冲压缩后的矩阵S中最后一个存在回波的距离单元序号，距离单元的初始序号q＝h1；
[0025] (5b)设置初始迭代次数k＝1，初始化模型先验参数：在第q个距离单元中，令干扰向量的聚类数初始值为10，令类别矩阵初始值的维数为2N1×10，矩阵每行中仅有一个元素为1而其它元素为0，非零元素序号在整数[0,10]中随机产生，令精度向量所有元素的初始值为0.01，令均值向量所有元素的初始值随机产生，其中精度向量和均值向量的元素个数均为10，令精度矩阵的初始值为对角线元素全为100的对角阵，维数为2Na×2Na，令混合系数的超参数全为1，精度向量每个元素的超参数全为10-4，精度矩阵每个对角线元素的超参数全为10-4；
[0026] (5c)利用迭代公式，依次计算第k次迭代的各个参数的值；
[0027] (5d)以类别矩阵L(k)每一行中最大元素对应的列序号为元素构造标签向量其中 为正整数，表示类别矩阵L(k)第1,…,2N1行中最大元素对应的聚类序号；
[0028] (5e)删除类别矩阵L(k)中第l1,…,lr个列向量，l1,…,lr表示标签向量(k)中重复元素对应的序号，类别矩阵L 剩余列的个数为干扰向量ε的聚类数J
(k)，对类别矩阵L(k)按行进行归一化操作：将类别矩阵L(k)的每个元素除以所在行元素之和，删除混合系数π(k)、均值向量μ(k)和精度向量α(k)中第l1,…,lr个元素；
[0029] (5f)判断终止函数E是否小于终止阀值η＝10-3，若是，则停止迭代，执行步骤(5g)；否则，更新迭代次数k＝k+1，执行步骤(5c)；
[0030] (5g)更新距离单元序号q＝q+1，判断更新后的距离单元序号q是否大于h2，若是，则停止对距离单元的搜索，执行步骤(6)；否则，执行步骤(5b)；
[0031] (6)重构权值矩阵：
[0032] 生成维数为2Na×Nr的权值矩阵W，权值矩阵W的第y列为第q个距离单元的权向量ωq，距离单元序号q满足q＝h1,...,h2，y与q的取值相同，其余元素全部为零，再利用公式W′＝W1+jW2生成复权值矩阵W′，其中W1表示权值矩阵W第1,…,Na行元素构成的矩阵，W2表示权值矩阵W第Na+1,…,2Na行元素构成的矩阵；
[0033] (7)转置并二维高分辨成像：
[0034] 对复权值矩阵W′进行转置操作，得到干扰环境下的二维高分辨逆合成孔径雷达ISAR成像结果。
[0035] 本发明与现有的技术相比，具有以下优点：
[0036] 第一，本发明基于稀疏贝叶斯模型，利用高斯混合模型对未知干扰信号进行建模，充分利用了环境的先验信息，克服了现有技术对干扰信号建模不准确的局限性，使得本发明提高了逆合成孔径雷达ISAR的成像质量。
[0037] 第二，本发明采用最大后验-期望最大MAP-EM算法，求解稀疏线性回归问题的权向量，能够有效降低运算量，提高运算效率，克服了现有技术中针对L1范数数值优化方法运算量较大等问题，使得本发明在干扰环境下，可以高效地获得聚焦良好的高分辨二维逆合成孔径雷达ISAR图像。
附图说明
[0038] 图1为本发明的流程图；
[0039] 图2(a)为缺损回波脉压后的距离慢时间图；
[0040] 图2(b)为利用正交匹配追踪OMP方法对缺损回波重构后的距离多普勒图；
[0041] 图2(c)为利用相关向量机RVM对缺损回波重构后的距离多普勒图；
[0042] 图2(d)为利用本发明对缺损回波重构后的距离多普勒图。
具体实施方式
[0043] 下面结合附图对本发明做进一步的描述。
[0044] 参照附图1，本发明的具体实施步骤如下：
[0045] 步骤1，录取目标有效回波。
[0046] 逆合成孔径雷达发射的电磁波在传播过程中遇到目标后，所反射的回波在传播过程中受到干扰的影响，最终被雷达接收机接收，在雷达显示器上显示出目标的缺损回波矩阵Sr，剔除缺损回波矩阵Sr中缺损的列向量，得到有效回波矩阵Se。
[0047] 步骤2，距离压缩并转置。
[0048] 将逆合成孔径雷达到场景中心的距离作为参考距离，将与逆合成孔径雷达发射信号载频、调频率相同，距离为参考距离的线性调频信号作为参考信号Sref，将参考信号Sref取共轭后与有效回波矩阵Se相乘，得到解线频调后的信号Srd：
[0049]
[0050] 其中，为距离快时间，tm为方位慢时间，Sref(·)为参考信号，Srd为解线频调后的信号，*表示共轭操作。
[0051] 对解线频调后的信号Srd沿快时间维作傅里叶变换，得到距离脉冲压缩后的矩阵S。
[0052] 对距离脉冲压缩后的矩阵S进行转置操作，得到转置矩阵S′。
[0053] 步骤3，生成字典。
[0054] 第1步，以 为元素,构造维数为Na×Na的傅里叶字典Φ，其中e(·)表示以自然常数为底的指数操作，j表示虚数单位，π表示圆周率，m表示傅里叶字典Φ的行序号，n表示傅里叶字典Φ的列序号，行序号和列序号的取值范围均为[-Na/2,Na/2-1]。
[0055] 第2步，剔除傅里叶字典Φ的第t1,…,tp行，得到维数为N1×Na的有效傅里叶字典Φe，其中t1,…,tp表示缺损回波矩阵Sr中缺损列向量对应的列序号。
[0056] 第3步，利用公式 构造实傅里叶字典Φ′，其中Re(·)表示取实部操作，Im(·)表示取虚部操作。
[0057] 步骤4，稀疏贝叶斯建模。
[0058] A.按照式(1)，构建每个距离单元的稀疏信号模型：
[0059] s′f＝Φeω′q+ε′ (1)
[0060] 其中，s′f表示转置矩阵S′的第f个列向量，Φe表示有效傅里叶字典，ω′q表示第q个距离单元的复权向量，ε′表示复干扰向量，f与q的取值相同；
[0061] B.利用高斯混合模型对干扰向量建模：
[0062] 第1步，将复干扰向量ε′的实部和虚部拼接为干扰向量ε，其中实部在上，虚部在下；
[0063] 第2步，干扰向量ε中每个元素的概率分布按照式(2)给出：
[0064]
[0065] 其中，p(εu|π,μ,α)表示干扰向量ε第u个元素的概率密度，J表示干扰向量ε的聚类数，Σ表示求和操作， 表示混合系数π的第j1个元素， 表示均值向量μ的第j2个元素，表示精度向量α的第j3个元素，N(·)表示高斯分布的概率密度，j1、j2、j3的取值相同；
[0066] 第3步，精度向量α中每个元素的先验分布按照式(3)给出：
[0067]
[0068] 其中， 表示精度向量α第j3个元素的概率密度，Gam(·)表示伽马Gamma分布的概率密度，a和b表示精度向量α第j3个元素 的超参数；
[0069] 第4步，混合系数的先验分布按照式(4)给出：
[0070] p(π)＝Dir(π|γ1,...,γJ) (4)
[0071] 其中，p(π)表示混合系数π的概率密度，Dir(·)表示狄利克雷Dirichlet分布的概率密度，γ1,...,γJ表示混合系数π的超参数。
[0072] C.利用伽马Gamma-高斯Gaussian分层模型对每个距离单元的权向量建模：
[0073] 第1步，将每个距离单元复权向量的实部和虚部拼接为对应的权向量，其中实部在上，虚部在下；
[0074] 第2步，每个距离单元权向量的先验分布按照式(5)给出：
[0075] p(ωq)＝N(0,Λ-1) (5)
[0076] 其中，p(ωq)表示第q个距离单元权向量的概率密度，精度矩阵Λ为对角阵；
[0077] 第3步，精度矩阵中每个对角线元素的先验分布按照式(6)给出：
[0078] p(λd)＝Gam(λd|e,f) (6)
[0079] 其中，p(λd)表示精度矩阵Λ第d个元素的概率密度，d＝1,…,2Na，e和f表示精度矩阵Λ第d个对角线元素λd的超参数。
[0080] 步骤5，利用最大后验-期望最大MAP-EM算法求解稀疏贝叶斯模型，计算存在回波的每个距离单元的权向量。
[0081] A.设置h1为第一个存在回波数据的距离单元序号，h2为最后一个存在回波数据的距离单元序号，距离单元的初始序号q＝h1。
[0082] B.设置初始迭代次数k＝1，初始化模型先验参数：在第q个距离单元中，令干扰的聚类数初始值为10，令类别矩阵初始值的维数为2N1×10，矩阵每行中仅有一个元素为1而其它元素为0，非零元素序号在整数[0,10]中随机产生，令精度向量所有元素的初始值为0.01，令均值向量所有元素的初始值随机产生，其中精度向量和均值向量的元素个数均为
10，令精度矩阵的初始值为对角线元素全为100的对角阵，维数为2Na×2Na，令混合系数的超参数全为1，精度向量每个元素的超参数全为10-4，精度矩阵每个对角线元素的超参数全为
10-4。
[0083] C.利用迭代公式，依次计算第k次迭代的各个参数的值：
[0084] 第1步，按照式(7)，计算干扰向量每个聚类的有效数目：
[0085]
[0086] 其中， 表示干扰向量ε第c个聚类的有效数目，N1表示方位向有效样本数，Σ表(k-1) (k-1) (k-1)示求和操作， 表示类别矩阵L 第i行的第j个元素，j＝1,…,J ，J 表示干扰向量ε的聚类数，c与j的取值相同；
[0087] 第2步，按照式(8)，计算混合系数的每个元素：
[0088]
[0089] 其中， 表示混合系数π(k)的第j1个元素， 表示混合系数π的超参数，Na表示缺损回波矩阵Sr的方位向采样点数，c、j1、j4与j的取值相同；
[0090] 第3步，按照式(9)，计算第q个距离单元的权向量
[0091]
[0092] 其中， 表示精度向量α(k-1)的第j3个元素， 表示实傅里叶字典Φ′转置矩阵的第v个列向量，Λ(k)表示精度矩阵， 表示第f个距离单元回波向量sf的第z个元素，表示均值向量μ(k-1)的第j2个元素，[·]-1表示矩阵求逆操作，[·]T表示转置操作，v、z与i的取值相同，f与q的取值相同，j2、j3与j的取值相同；
[0093] 第4步，按照式(10)，计算精度向量的每个元素：
[0094]
[0095] 其中，a和b表示精度向量α第j3个元素 的超参数，v、z与i的取值相同，c、j2、j3与j的取值相同，f与q的取值相同；
[0096] 第5步，按照式(11)，计算均值向量的每个元素；
[0097]
[0098] 其中，v、z与i的取值相同，f与q的取值相同，c、j2、j3与j的取值相同；
[0099] 第6步，按照式(12)，计算类别矩阵的每个元素：
[0100]
[0101] 其中，i＝1,...,2N1，v与i的取值相同，j1、j2、j3与j的取值相同；
[0102] 第7步，按照式(13)，计算精度矩阵的每个对角线元素：
[0103]
[0104] 其中， 表示精度矩阵Λ(k)的第d个对角线元素，e和f表示精度矩阵Λ第d个对角线元素λd的超参数， 表示第q个距离单元权向量 的第x个元素，d＝1,…,2Na，x与d的取值相同；
[0105] 第8步，更新精度矩阵Λ(k)。
[0106] D.以类别矩阵L(k)每一行中最大元素对应的列序号为元素构造标签向量其中 为正整数，表示类别矩阵L(k)第1,…,2N1行中最大元素对应的聚类序号。
[0107] E.删除类别矩阵L(k)中第l1,…,lr个列向量，l1,…,lr表示标签向量(k) (k)中重复元素对应的序号，类别矩阵L 剩余列的个数为干扰向量ε的聚类数J ，对类别矩阵L(k)按行进行归一化操作：将类别矩阵L(k)的每个元素除以所在行元素之和，删除混合系数π(k)、均值向量μ(k)和精度向量α(k)中第l1,…,lr个元素。
[0108] F.按照式(14)计算终止函数E：
[0109]
[0110] 其中，下界函数 的表达式如式(15)所示：
[0111]
[0112] 其中，N1表示方位向有效样本数，J(k)表示干扰向量ε的聚类数， 表示类别矩阵L(k)第i行的第j个元素， 表示混合系数π(k)的第j1个元素， 表示精度向量α(k)的第j3个元素， 表示第f个距离单元回波向量sf的第z个元素， 表示实傅里叶字典Φ′转置矩阵的第v个列向量， 表示第q个距离单元的权向量， 表示均值向量μ(k)的第j2个元素，Λ(k)表示精度矩阵，a和b表示精度向量α第j3个元素 的超参数， 表示混合系数π的超参数，In表示以自然常数为底的对数操作，v、z与i的取值相同，j1、j2、j3、j4与j的取值相同，f与q的取值相同。
[0113] 判断终止函数E是否小于终止阀值η＝10-3，若是，则停止迭代，执行本步骤的第G步；否则，更新迭代次数k＝k+1，执行本步骤的第C步。
[0114] G.更新距离单元序号q＝q+1，判断更新后的距离单元序号q是否大于h2，若是，则停止对距离单元的搜索，执行步骤6；否则，执行本步骤的第B步。
[0115] 步骤6，重构权值矩阵。
[0116] 生成维数为2Na×Nr的权值矩阵W，权值矩阵W的第y列为第q个距离单元的权向量ωq，距离单元序号q满足q＝h1,...,h2，y与q的取值相同，其余元素全部为零，再利用公式W′＝W1+jW2生成复权值矩阵W′，其中W1表示权值矩阵W第1,…,Na行元素构成的矩阵，W2表示权值矩阵W第Na+1,…,2Na行元素构成的矩阵。
[0117] 步骤7，转置并二维高分辨成像。
[0118] 对复权值矩阵W′进行转置操作，得到干扰环境下的二维高分辨逆合成孔径雷达ISAR成像结果。
[0119] 通过以下仿真对本发明的效果做进一步的说明：
[0120] 1.仿真条件及其参数。
[0121] 本发明的仿真实验采用工作在C波段的雷达，对应载频为10GHZ，带宽为0.4GHZ。目标长度为8.75米，翼展宽度为6米，回波数据的缺损率为50％。
[0122] 2.仿真内容。
[0123] 本发明的仿真实验有四个，其中：
[0124] 仿真实验1：对距离脉压后回波数据的第1～80列、207～320列以及451～512列产生缺损，并添加随机干扰，绘制其距离慢时间图，结果如图2(a)。
[0125] 仿真实验2：利用现有技术中的正交匹配追踪OMP方法，对缺损回波进行重构，绘制其距离多普勒图像，结果如图2(b)。
[0126] 仿真实验3：利用现有技术中的相关向量机RVM对缺损回波进行重构，绘制其距离多普勒图像，结果如图2(c)。
[0127] 仿真实验4：利用本发明对缺损回波进行重构，绘制其距离多普勒图像，结果如图2(d)。
[0128] 图2(a)给出了缺损率为50％的脉压后回波数据的距离慢时间图，图2(a)中的横坐标方向表示慢时间，纵坐标方向表示距离。
[0129] 图2(b)给出了利用现有技术中的正交匹配追踪OMP方法对缺损回波数据进行重构的距离多普勒图，图2(b)中的横坐标方向表示多普勒频率，纵坐标方向表示距离，可以看出利用现有的正交匹配追踪OMP方法得到的距离多普勒图在第240到280列附近具有强干扰，无法清楚地呈现目标的几何结构。
[0130] 图2(c)给出了利用现有技术中的相关向量机RVM对缺损回波数据进行重构的距离多普勒图，图2(c)中的横坐标方向表示多普勒频率，纵坐标方向表示距离，可以看出利用现有的相关向量机RVM得到的距离多普勒图受到残留干扰信号的影响，同样无法清楚地呈现目标的几何结构。
[0131] 图2(d)给出了利用本发明对缺损回波数据进行重构的距离多普勒成像图，图2(d)中的横坐标方向表示多普勒频率，纵坐标方向表示距离，与图2(b)和图2(c)对比可得，相比于现有技术中的正交匹配追踪OMP方法和相关向量机RVM重构后的距离多普勒图，利用本发明得到重构后的距离多普勒图能够清楚地呈现目标的几何结构，虚假点少，分辨力显著提高且图像聚焦良好。
[0132] 仿真结果表明，本发明基于稀疏信号重构理论将逆合成孔径雷达ISAR高分辨成像问题转化为稀疏线性回归问题，利用高斯混合模型对未知的干扰信号进行建模，采用最大后验-期望最大MAP-EM算法求解权值向量，充分利用了雷达回波数据的先验信息，在干扰环境下得到高分辨且聚焦良好的距离多普勒图像。
