[发明专利]基于相关向量机线性阵列SAR三维成像快速贝叶斯压缩感知方法

权利要求书

1.一种基于相关向量机的线性阵列SAR三维成像快速贝叶斯压缩感知方法，其特征是它包括如下步骤：

步骤1、初始化SAR系统参数：

初始化SAR系统参数包括：平台速度矢量记为阵列天线各阵元初始位置矢量，记做其中l为天线各阵元序号，N_L为阵列天线的阵元总数；阵列天线长度，记做L；雷达发射信号载频为f_c；雷达发射信号的调频斜率为f_dr；脉冲重复时间记为PRI；雷达系统的脉冲重复频率为PRF；雷达发射信号带宽记做B_r；电磁波在空气中的传播速度记做C；距离向快时间记做n，n＝1,…,N_R，N_R为距离向快时刻总数，方位向慢时刻记做k，k＝1,…,K，K为方位向慢时刻总数；上述参数均为SAR系统标准参数，其中雷达信号载频f_c，雷达发射信号的调频斜率f_dr，脉冲重复时间PRI，雷达系统的脉冲重复频率PRF，雷达发射信号带宽B_r，阵列天线的阵元总数N_L，阵列天线长度L在线阵SAR系统设计过程中已经确定；平台速度矢量阵列天线各阵元初始位置矢量在SAR观测方案设计中已经确定；根据SAR成像系统方案和观测方案，SAR成像方法需要的初始化成像系统参数均为已知；

步骤2、划分SAR的成像场景空间：

以雷达波束照射场区域地平面和垂直于该地平面向上的单位向量所构成的空间直角坐标系作为线阵SAR的成像场景目标空间Ω，其中水平横向和水平纵向构成阵列维成像空间；初始化水平横向成像场空间长度为L_x，水平纵向成像场空间长度为L_y；将成像场景目标空间Ω均匀划分为大小相等的立体单元格，成像场景空间在水平横向、水平纵向单元格数分别为M_x，M_y；

采用公式计算得到水平横向、水平纵向的单元格大小，分别记为d_x和d_y；成像场景空间高度向单元格总数为N_G，高度向单元格的大小为线阵SAR成像系统距离向分辨率，记为d_z；

采用公式计算得到划分后的成像场景空间Ω中第n个等距离单元格的阵列维成像空间中第m_x个水平横向单元格第m_y个水平纵向单元格所对应的元素的位置，记为其中m_x＝1，…,M_x，m_y＝1，…,M_y；

采用公式得到成像场景空间散射系数矩阵，记做δ，其中m＝(m_y-1)M_x+m_x＝1，…,M，n＝1,…,N_R，N_R为距离向快时刻总数，为划分后的成像场景空间Ω中第n个等距离单元格的阵列维成像空间中第m个元素的散射系数，M＝M_x·M_y为第n个等距离单元格的阵列维成像空间中的等效单元格总数；

步骤3、建立线阵SAR的测量矩阵：

步骤3.1、在实际线阵SAR成像中，原始回波数据由数据接收机提供；线阵SAR中二维等效阵列中第l个天线阵元在第n个距离向快时刻的原始回波数据记做s(n,l,k)；采用标准合成孔径雷达距离向脉冲压缩方法对s(n,l,k)进行距离向脉冲压缩，得到距离向压缩后的线阵SAR数据，记做s_AC(n,l,k)；

采用公式S_n＝{S_n(i)＝s_AC(n,l,k)}，i＝k×N_L+l，k＝1,…,K，l＝1,…,N_L，计算得到第n个等距离单元格回波信号，记为S_n，n＝1,…,N_R，N_L为步骤1中初始化得到的阵列天线阵元总数，K为步骤1中初始化得到的方位向慢时刻总数，其中N_R为步骤1中初始化得到的距离向快时刻总数，N_A＝K×N_L为回波信号的行数，即线阵SAR中二维等效阵列中的元素的总数；

步骤3.2、采用公式计算得到第l个阵列天线在第k个方位向慢时刻的位置矢量，记为其中N_L为步骤1中初始化得到的阵列天线阵元总数，K为步骤1中初始化得到的方位向慢时刻总数，为步骤1中初始化得到的阵列天线各阵元初始位置，为步骤1中初始化得到的平台速度，PRF为步骤1中初始化得到的雷达系统的脉冲重复频率；

采用公式计算得到在第k个方位向慢时刻线阵SAR成像场景目标空间Ω中第n个等距离单元格中阵列维成像空间中第m个元素到第l个天线阵元的时间延时，记为

采用公式计算得到回波信号S_n与散射系数矩阵δ之间的测量矩阵，记做Θ_n；其中M为步骤2中得到的第n个等距离单元格的阵列维成像空间中的等效单元格总数，||·||₂表示向量L2范数，为步骤2中得到的划分后的成像场景空间Ω中第n个等距离单元格的阵列维成像空间中第m_x个水平横向单元格第m_y个水平纵向单元格所对应的元素的位置，C为步骤1中初始化得到的电磁波在空气中的传播速度，S_n为步骤3.1中得到的脉冲压缩后第n个等距离单元格回波信号，n＝1,…,N_R，N_R为距离向快时刻总数；

步骤4、对测量矩阵、回波信号、散射系数进行实数域表示：

采用公式和计算得到实数域回波信号和散射系数其中N_A为步骤3.1中得到的回波信号S_n的行数，M为步骤2中得到的第n个等距离单元格的阵列维成像空间中的等效单元格总数，Re(·)和Im(·)分别表示复数域数据的实部和虚部，为划分后的成像场景空间Ω中第n个等距离单元格的阵列维成像空间中第m个元素的散射系数；采用公式计算得到实数域测量矩阵，记为S_n为步骤3.1中得到的第n个等距离单元格回波信号，Θ_n为步骤3.2中回波信号S_n与散射系数矩阵δ之间的测量矩阵；

步骤5、初始化目标区域：

步骤5.1、初始化噪声方差β₀；

步骤5.2、采用公式计算得到的第m_c列向量，记为其中为步骤4中的实数域测量矩阵；采用公式计算得到第n个等距离单元格中第m_c个散射系数，记为其中为步骤4中第n个等距离单元格的实数域回波信号，||·||为范数，∑(·)为求和运算符号，(·)^T为矩阵转置的运算符号；

采用公式初始化第n个等距离单元格回波信号的散射系数，记为其中M为步骤2中得到的第n个等距离单元格的阵列维成像空间中的等效单元格总数；采用公式初始化超参数，记为其中为第m_c个散射系数对应的超参数；

步骤5.3、遍历步骤5.2中的散射系数得到最大元素对应的散射单元，记为Id₀；并根据公式得到初始化目标区域，记为

采用公式计算得到散射系数中的最大值对应的超参数，记为其中为的第Id₀列元素，为步骤4中实数域测量矩阵，为步骤4中第n个等距离单元格的实数域回波信号，为步骤5.2中第n个等距离单元格回波信号的散射系数的初始值，||·||为范数，β₀为步骤5.1初始化得到的噪声方差；

步骤6、计算最优超参数来提取目标区域：

步骤6.1、初始化迭代终止门限ε_e，初始化第i_v个超参数的迭代误差其中i_v＝1,…,2M，M为步骤2中得到的第n个等距离单元格的阵列维成像空间中的等效单元格总数；采用公式初始化超参数的迭代误差，记为E⁰；

在第t次迭代中，若t＝1，则为初始目标区域；若t1，为经过t-1迭代得到的目标区域；其中为步骤5.3中得到的初始目标区域；

步骤6.2、采用公式计算得到第t-1次迭代的超参数的对角矩阵，记为其中为第t-1次迭代中第i_v个散射系数对应的超参数，M为步骤2中得到的第n个等距离单元格的阵列维成像空间中的等效单元格总数；

采用公式计算协方差矩阵，记为C_t-1；其中Ι为单位矩阵，β为信号噪声的方差，为步骤4中得到的实数域测量矩阵，N_A为步骤3.1中得到的回波信号S_n的行数，(·)^T为步骤5.2中矩阵转置的运算符号；

采用公式更新第i_v个超参数对应的稀疏因子、质量因子，记为和其中为的第i列元素，为步骤4中第n个等距离单元格的实数域回波信号；

步骤6.3、利用步骤6.2中的和采用公式更新第t次迭代后的第i_v个超参数，记为

采用公式更新第t次迭代后的第i_v个超参数的估计误差，记为采用公式计算得到所有超参数的迭代误差，记为E^t；其中M为步骤2中得到的第n个等距离单元格的阵列维成像空间中的等效单元格总数；

步骤6.4、如果则迭代终止，采用公式得到最终目标区域，记为其中E⁰为步骤6.1中所有超参数的迭代误差的初始值，||·||为范数，ε_e为步骤6.1中的迭代终止门限；

否则，选取第t次迭代的中最大误差项对应的散射单元，记为Id_t；

采用公式更新步骤6.1中的目标区域记为其中为步骤6.2中第t次迭代中第Id_t个散射单元的为步骤6.2中第t次迭代中第Id_t个散射单元的为第t次迭代中第Id_t个散射单元的超参数，∪为添加元素运算符号，/为删除元素运算符号；

步骤6.5、采用公式计算得到第t次迭代后的目标区域对应的测量矩阵，记为Φ_t；其中为步骤4中实数域测量矩阵的第列元素，i_m＝1,…,t；

步骤7、目标区域转换，得到复数域目标区域；

采用公式计算得到复数域目标区域，记为G_n；其中N_n为目标区域中的元素个数，ω_r为目标区域G_n中第r个元素在成像场景中的序号，为实值目标区域中第r个元素，round(·)为四舍五入取整运算符号；

步骤8、利用目标区域简化测量矩阵和散射系数：

采用公式计算得到Θ_n的第ω_r列向量，记为其中Θ_n为步骤4中的复数域测量矩阵，N_A是步骤3.1中得到的回波信号S_n的行数，ω_r为步骤7中的目标区域G_n中第r个元素在成像场景中的序号；

采用公式计算得到目标区域的简化散射系数向量，记为α_ω；其中为序号ω_r对应的散射系数；采用公式计算得到目标区域的简化测量矩阵，记为Θ_ω；其中N_n为步骤7中目标区域G_n中的元素个数；采用公式初始化目标区域的散射系数估计值记为其中(·)^H为Hermite矩阵运算符号，其中S_n为步骤3.1中得到的第n个等距离单元格回波信号；

步骤9、估计目标区域的散射系数：

步骤9.1、初始化正则化参数λ、平滑因子η、范数项p、迭代次数t、最大迭代次数N_max；误差终止迭代门限ε₀；

步骤9.2、采用公式计算得到第t次迭代的对角矩阵，记为ω_r为步骤7中目标区域G_n中第r个元素在成像场景中的序号，为步骤8中Θ_n的第ω_r列向量，为第t-1次迭代后信号噪声方差的估计值，为第t-1次迭代后估计值的第r个元素，N_A是步骤3.1中得到的回波信号S_n的行数，N_n为步骤7中目标区域G_n中的元素个数，diag{·}为对角矩阵符号，(·)^H为步骤8中Hermite矩阵运算符号，λ为步骤9.1中初始化得到的正则化参数，p为步骤9.1中初始化得到范数项；η为步骤9.1中初始化得到的平滑因子；

步骤9.3、采用公式计算得到的秩，记为K_c；rank(·)为矩阵的秩运算符号，为步骤9.2中第t次迭代的对角矩阵；

如果满足K_c＝N_n；采用公式计算得到第t次迭代的散射系数和信号噪声方差的估计值，分别记为和Θ_ω为步骤8中的简化测量矩阵，N_A是步骤3.1中得到的回波信号S_n的行数，S_n为步骤3.1中得到的第n个等距离单元格回波信号，(·)^H为步骤8中Hermite矩阵运算符号，||·||为范数，其中K_c为步骤9.2中得到对角矩阵的秩；

否则，采用公式对对角矩阵进行奇异值分解，其中是的第i_s1个奇异值，i_s1＝1,…,N_n，并且满足和为对应于的奇异向量，N_n为步骤7中目标区域G_n中的元素个数，Σ(·)为步骤5.2中求和运算符号；采用公式计算得到广义逆矩阵，记为其中K_c为步骤9.2中的秩；采用公式计算得到第t次迭代的散射系数的估计值，记为其中Θ_ω为步骤8中的简化测量矩阵，(·)^H为步骤8中Hermite矩阵运算符号；采用降序排列对按照幅值以降序排列，并得到其中第K_m个元素，记为采用公式重新提取目标区域，记为G_n；其中为第t次迭代中序号为ω_r的散射系数的估计值，ω_r为步骤7中目标区域G_n中第r个元素在成像场景中的序号；

步骤10、判断是否终止迭代：

如果且tN_max，则继续执行步骤9，且t＝t+1，其中为步骤9.3中得到的第t次迭代的散射系数的估计值，N_max为步骤9.1中初始化得到的高分辨率成像的最大迭代次数，ε₀为步骤9.1中初始化得到的误差迭代终止门限；

若不满足和t≥N_max任一条件，算法迭代终止，则输出得到的高分辨率成像中第n个等距离单元格的散射系数向量估计值

步骤11、全场景三维成像：

采用公式将各等距离单元格散射系数向量排列为三维矩阵形式，得到线阵SAR成像场景空间的三维高分辨率的成像结果，记为其中N_R为步骤1中初始化得到的距离向快时刻总数，其中为步骤10中得到的高分辨率成像中第n个等距离单元格散射系数估计结果；

至此，整个方法结束。