1.   引言

双基雷达，也称双站雷达，与单基雷达相比，由于收发系统分置而具有更强的生存能力和抗干扰性能，能够获取更加丰富的目标信息，并且有效降低系统成本^[1]。因此，双基雷达在现代电子战中具有独特优势与潜力，成为当前雷达技术的研究热点^[2,3]。其中，双基逆合成孔径雷达(Bistatic Inverse Synthetic Aperture Radar, Bistatic ISAR)成像在防空反导、国土防空、空间监测等领域更是具有广阔的应用前景。

文献[4]研究了空间高速目标的双基ISAR成像算法，通过构造补偿函数，应用Keystone变换解决越距离单元徙动问题，获得清晰的空间高速目标2维像。文献[5]基于最小熵准则利用分数阶傅里叶变换获得了高速目标的高分辨1维距离像，并提出一种基于频率平滑度的成像时间选择方法获取复杂运动目标的方位像。文献[6]研究了双基地角时变情况下的双基ISAR成像，利用1阶泰勒展开公式推导得到畸变项和散焦项表达式，提出一种畸变估计方法。文献[7,8]研究了双基ISAR成像的畸变校正以及距离向和方位向定标问题，实现对成像结果的距离向和方位向定标。文献[9,10]中研究了低信噪比条件下双基ISAR成像算法，分别提出基于贝叶斯框架的稀疏孔径ISAR成像方法和基于多量测向量模型的复数近似消息传递的双基ISAR 1维距离成像分辨率增强算法。文献[11]对数据缺失情况下的双基ISAR成像进行研究，提出稀疏恢复内插值的方法重构目标2维像。文献[12]对双基ISAR系统中的微多普勒效应进行研究，提出一种扩展Hough变换提取目标实际微动特征的方法。随着信号处理技术的不断发展，基于干涉处理的双基3维干涉成像技术也逐渐成为研究热点^[13-15]。这些研究成果对促进双基ISAR成像系统的实用化发挥了重要作用。

ISAR成像技术通过发射大带宽信号获得距离向的高分辨，利用目标相对雷达运动产生的多普勒频移，实现方位向上的高分辨成像。当目标在非合作情况下，或者凝视/近凝视的非理想观测几何条件下，ISAR成像效果急剧恶化^[16]。虽然现阶段可以通过运动补偿^[4]或者畸变校正^[6,7]的方法获取聚焦较好的成像结果，但这类方法过程仍然较为复杂。雷达关联成像(Radar Coincidence Imaging, RCI)通过阵列天线产生在时间和空间上随机分布的辐射场，依靠雷达信号本身的差异性实现对目标的高分辨成像^[17]。其方位向分辨率不再依赖于目标相对雷达运动产生的多普勒频移，而是与辐射场的差异性有关，即辐射场中相邻两点信号的差异性决定了这两点的分辨^[18,19]。与ISAR成像相比，方位向分辨不依赖于对多个回波信号做相干处理，因此，RCI不需要对回波信号进行相干积累，成像更快。

本文基于RCI成像原理，提出双基RCI成像方法，并针对运动目标成像问题展开研究，由于双基RCI系统与传统RCI系统存在结构差异，并且目标运动使得推演的双基RCI 2维随机辐射场与真实辐射场存在失配。因此，本文在构建双基RCI系统模型的基础上，首先建立双基RCI参数化稀疏表征模型，然后提出一种新的成像算法—基于稀疏贝叶斯学习的迭代关联成像(Iteration Coincidence Imaging based on Sparse Bayesian Learning, ICI-SBL)算法，利用该算法可以实现对运动目标进行成像，同时完成对目标运动参数的精确估计。

2.   回波信号模型

双基RCI系统发射机天线采用均匀线性阵列(Uniform Linear Array, ULA)，如图1所示，由${N_{\rm t}}$个阵元组成，每个阵元分别发送不相关的随机频率调制信号^[19]，每个随机频率调制信号由$M$个子脉冲构成，子脉冲时宽为${T_{\rm s}}$，第$m$个子脉冲发射时刻为${t_m} = \left( {m - 1} \right){T_{\rm s}}$, $m = 1,2, ·\!·\!· ,M$。以子脉冲发射时刻为起点的时间用${t_k}$表示。${t_k}$用来计量电波传播的时间，而${t_m}$是计量发射子脉冲的时刻，这两个时间与全时间的关系为${t_k} = t - {t_m}$。因此第$\rm nt$个阵元发送的信号可以写成

图 1  双基RCI系统示意图

下载: 全尺寸图片 幻灯片

$$ \begin{align} {s_{{\rm{nt}}}}\left( {{t_k},{t_m}} \right) =& {\rm{rect}}\left( {\frac{{{t_k}}}{{{T_{\rm s}}}}} \right)\\ {\rm{}}& \cdot{\rm{exp}}\left\{ {{\rm{j}}{2{\text{π}}}\left( {{f_{\rm c}} + C_{{\rm{nt}}}^m \cdot \Delta f} \right){t_k} + {\rm{j}}{\varphi _{0{\rm{nt}}}}} \right\},\\ {\rm{}}&{\rm{nt}} = 1,2, ·\!·\!· ,{N_{\rm t}} \end{align} $ $

(1)

式中，$\rm{j}$为虚数单位，${\rm{rect}}\left( {{{{t_k}}/{{T_{\rm s}}}}} \right)$为矩形窗函数，$\left| {{t_k}} \right| \le {T_{\rm s}}/2$时，${\rm{rect}}\left( {{{{t_k}}/{{T_{\rm s}}}}} \right) = 1$, $\left| {{t_k}} \right| > {T_s}/2$时，${\rm{rect}}\left( {{{{t_k}}/{{T_{\rm s}}}}} \right) = 0$; ${f_{\rm c}}$为信号载频，$\Delta f$为频率间隔，$\left\{ {C_{\rm{nt}}^m,m = 1,2, ·\!·\!· ,M} \right\}$为一服从均匀分布的随机序列，${\varphi _{0\rm{nt}}}$为第$\rm{nt}$个发射阵元发送信号的初始相位。

双基RCI系统接收机天线同样采用ULA，由${N_{\rm r}}$个阵元组成，建立雷达坐标系$OXY\;$和目标本地坐标系$O'xy$，双基雷达发射机位于雷达坐标系原点$O$处，接收机位于雷达坐标系$\left( {L,0} \right)$处，$L$为雷达基线长度。其中，目标本地坐标系$x$与$y$轴分别平行于雷达坐标系$X$轴和$Y\;$轴。

在随机频率调制信号回波信号建模过程中，采用了“stop-go”模型^[20]，忽略目标在子脉冲内的运动。因此，对一个有$P\;$个散射点的目标，可得其第$\rm{nr}$个接收阵元接收的回波信号${y_{\rm{nr}}}\left( {{t_k},{t_m}} \right)$, ${\rm{nr}} \!=\! 1,2, ·\!·\!· ,{N_{\rm R}}$。

$$ {y_{{\rm{nr}}}}\left( {{t_k},{t_m}} \right) = \sum\limits_{p = 1}^P {\sum\limits_{{\rm{nt}} = 1}^{{N_{\rm t}}} {{\sigma _p}{s_{{\rm{nt}}}}\left( {{t_k} - {{{r_{{\rm{tpr}}}}\left( {{t_m}} \right)}/{\rm c}},{t_m}} \right)} } $ $

(2)

式中，${\sigma _p}$为第$p$个散射点的前向散射系数，$p = 1,2, ·\!·\!· ,P\;$, $\rm c$表示光在真空中传播的速度，${r_{\rm{tpr}}}\left( {{t_m}} \right)$表示在${t_m}$时刻第$\rm{nt}$个发射阵元和第$\rm{nr}$个接收阵元到第$p$个散射点的距离和。

选择目标本地坐标系原点$O'$作为参考点，并已知其在雷达坐标系中的坐标，假设目标做平稳匀速运动，用极坐标形式表示的速度为${\text{v}} = {\left[ {v,\theta } \right]^{\rm{T}}}$, $v$为目标运动速度，$\theta$为目标运动方向，即速度方向与$x$轴的夹角，则

$$ {r_{{\rm{tpr}}}}\left( {{t_m}} \right) = \left| {{{\text{r}}_p} + {\text{v}}{t_m} + {{\text{R}}_{\rm t}}} \right| + \left| {{{\text{r}}_p} + {\text{v}}{t_m} + {{\text{R}}_{\rm r}}} \right| $ $

(3)

式中，${{\text{r}}_p}$为目标散射点在目标本地坐标系的位置矢量，${{\text{R}}_{\rm t}}$和${{\text{R}}_{\rm r}}$分别表示第$\rm{nt}$个发射阵元到参考点$O'$的位置矢量和第$\rm{nr}$个接收阵元到参考点$O'$的位置矢量。

由于$\left| {{{\text{r}}_p} + {{\text{R}}_{\rm t}}} \right| \gg v{t_m}$, $\left| {{{\text{r}}_p} + {{\text{R}}_{\rm r}}} \right| \gg v{t_m}$，因此可对${r_{\rm{tpr}}}\left( {{t_m}} \right)$求取近似值，得

$$ \begin{aligned} {r_{{\rm{tpr}}}}\left( {{t_m}} \right) =& \left| {{{\text{r}}_p} + {\text{v}}{t_m} + {{\text{R}}_{\rm t}}} \right| + \left| {{{\text{r}}_p} + {\text{v}}{t_m} + {{\text{R}}_{\rm r}}} \right|\\ \approx& \left| {{{\text{r}}_p} + {{\text{R}}_{\rm t}}} \right| + \left| {{{\text{r}}_p} + {{\text{R}}_{\rm r}}} \right| \\ {\rm{}}&+ v{t_m}\left[ {\cos \left( {\theta - {\varphi _{\rm t}}} \right) + \cos \left( {\theta - {\varphi _{\rm r}}} \right)} \right] \end{aligned} $ $

(4)

式中，${\varphi _{\rm t}}$和${\varphi _{\rm r}}$分别为双基雷达发射机和接收机天线雷达视线方向方位角，且${\varphi _{\rm t}},{\varphi _{\rm r}} \in \left( {0,{\text{π}}} \right)$。将式(4)代入式(2)，则得到双基RCI系统接收机接收到的目标回波信号。

3.   双基雷达关联成像参数化稀疏表征模型

雷达关联成像是一种不依赖相对运动的成像方法，可以对静止目标进行成像。然而，在实际应用中，目标总是运动的，因此，考虑目标做平稳匀速运动的情况^[21]，建立了运动目标的双基RCI参数化稀疏表征模型。首先在目标本地坐标系$O'xy$中，选取将整个目标包含在内的成像区域$I$，且成像区域$I$随目标一起以速度${\text{v}}$做平稳匀速运动。成像区域内任意一点的位置矢量为${\text{r}}$，利用ULA发送的随机频率调制信号，可以在成像区域内叠加得到探测信号${S_I}\left( {{\text{r}},{t_k},{t_m}} \right)$^[17], ${\text{r}} \in I$

$$ {S_I}\left( {{\text{r}},{t_k},{t_m}} \right) = \sum\limits_{{\rm{nt}} = 1}^{{N_{\rm t}}} {{s_{{\rm{nt}}}}\left( {{t_k} - {{{r_{\rm t}}\left( {{t_m}} \right)}/{\rm c}},{t_m}} \right)} $ $

(5)

式中，${r_{\rm t}}\left( {{t_m}} \right) = \left| {{\text{r}} + {\text{v}}{t_m} + {{\text{R}}_{\rm t}}} \right|$为${t_m}$时刻第$\rm{nt}$个发射阵元到成像区域任意一点的距离。探测信号被目标散射，再被接收机天线接收得到回波信号，可表示为成像区域内探测信号${S_I}\left( {{\text{r}},{t_k},{t_m}} \right)$的叠加

$$ \begin{align} {y_{{\rm{nr}}}}\left( {{t_k},{t_m}} \right) =& \int_I {\sigma \left( {\text{r}} \right){S_I}\left( {{\text{r}},{t_k} - {{{r_{\rm r}}\left( {{t_m}} \right)}/{\rm c}},{t_m}} \right){\rm{d}}{\text{r}}} \\ =& \int_I \sigma \left( {\text{r}} \right)\sum\limits_{{\rm{nt}} = 1}^{{N_{\rm t}}} {s_{{\rm{nt}}}}\left( {t_k} - \left( {r_{\rm t}}\left( {{t_m}} \right) \right.\right.\\ {\rm{}}& \left.\left.+ {r_{\rm r}}\left( {{t_m}} \right) \right)/{\rm c},{t_m} \right) {\rm{d}}{\text{r}} \end{align} $ $

(6)

式中，${r_{\rm r}}\left( {{t_m}} \right) = \left| {{\text{r}} + {\text{v}}{t_m} + {{\text{R}}_{\rm r}}} \right|$为${t_m}$时刻第$\rm{nr}$个接收阵元到成像区域任意一点的距离。根据文献[17]，可得2维随机辐射场为

$$ S\left( {{\text{r}},{t_k},{t_m}} \right) = \sum\limits_{{\rm{nt}} = 1}^{{N_{\rm t}}} {{s_{{\rm{nt}}}}\left( {{t_k} - {{{r_{{\rm{tr}}}}\left( {{t_m}} \right)}/{\rm c}},{t_m}} \right)} $ $

(7)

式中，

$$ \begin{align} {r_{{\rm{tr}}}}\left( {{t_m}} \right) =& {r_{\rm t}}\left( {{t_m}} \right) + {r_{\rm r}}\left( {{t_m}} \right) = \left| {{\text{r}} + {\text{v}}{t_m} + {{\text{R}}_{\rm t}}} \right| \\ {\rm{}}&+ \left| {{\text{r}} + {\text{v}}{t_m} + {{\text{R}}_{\rm r}}} \right|\\ \approx& \left| {{\text{r}} + {{\text{R}}_{\rm t}}} \right| + \left| {{\text{r}} + {{\text{R}}_{\rm r}}} \right| \\ {\rm{}} & + v{t_m}\left[ {\cos \left( {\theta - {\varphi _{\rm t}}} \right) + \cos \left( {\theta - {\varphi _{\rm r}}} \right)} \right] \end{align} $ $

(8)

然后，对成像区域离散化，将其划分为$Q$个成像单元，则成像区域可以表示为$V = \left\{ {{{\text{r}}_1},{{\text{r}}_2}, ·\!·\!· ,{{\text{r}}_Q}} \right\}$。其中，${{\text{r}}_q}$为第$q$个成像单元的位置矢量，$q = 1,2, ·\!·\!· ,Q$。所有成像单元的散射系数组合为散射强度矢量：${\text{σ}} = {\left[ {{\sigma _1},{\sigma _2}, ·\!·\!· ,{\sigma _Q}} \right]^{\rm{T}}}$，其中，${\sigma _q}$表示第$q$个成像单元的散射系数。

再对回波信号${y_{\rm{nr}}}\left( {{t_k},{t_m}} \right)$进行时域采样，对每个子脉冲回波采样一次，相应的采样时刻可以表示为${{\text{t}}_{\rm s}} = \left[ {t_{\rm s}^1,t_{\rm s}^2, ·\!·\!· ,t_{\rm s}^M} \right]$，其中，$t_{\rm s}^m$表示第$m$个采样时刻，$m = 1,2, ·\!·\!· ,M$。得到离散化回波信号为

$$ \begin{align} {\text{y}} \!=& \left[ {y_1}\left( 1 \right), ·\!·\!· ,{y_1}\left( M \right),{y_2}\left( 1 \right), ·\!·\!· ,{y_2}\left( M \right), ·\!·\!· ,\right.\\ {\rm{}}& \left.{y_{{N_{\rm r}}}}\left( 1 \right), ·\!·\!· ,{y_{{N_{\rm r}}}}\left( M \right) \right]^{\rm{T}} \end{align} $ $

(9)

通过对成像区域和时域离散化，可以得到离散化的2维随机辐射场

$$ {\text{S}}\left( {v,\theta } \right) = \left[ {\begin{array}{*{20}{c}} {{S_1}\left( {1,1} \right)}\!\!&\!\!{{S_1}\left( {1,2} \right)}\!\!&\!\! ·\!·\!· \!\!&\!\!{{S_1}\left( {1,Q} \right)} \\ \vdots \!\!&\!\! \vdots \!\!&\!\! \ddots \!\!&\!\! \vdots \\ {{S_1}\left( {M,1} \right)}\!\!&\!\!{{S_1}\left( {M,2} \right)}\!\!&\!\! ·\!·\!· \!\!&\!\!{{S_1}\left( {M,Q} \right)} \\ \vdots \!\!&\!\! \vdots \!\!&\!\! \ddots \!\!&\!\! \vdots \\ {{S_{{N_{\rm r}}}}\left( {M,1} \right)}\!\!&\!\!{{S_{{N_{\rm r}}}}\left( {M,2} \right)}\!\!&\!\! ·\!·\!· \!\!&\!\!{{S_{{N_{\rm r}}}}\left( {M,Q} \right)} \end{array}} \right] $ $

(10)

由于${\text{S}}\left( {v,\theta } \right)$中含有目标运动参数$v$和$\theta$，根据压缩感知(Compressive Sensing, CS)理论，考虑复高斯白噪声的影响，得到双基RCI参数化稀疏表征模型为

$${\text{y}} = {\text{S}}\left( {v,\theta } \right) \cdot {\text{σ}} + {\text{w}}$ $

(11)

式中，${\text{y}} \in {\mathbb{C}^{M \times 1}}$为离散化回波信号矢量，${\text{S}}\left( {v,\theta } \right)$$\in {\mathbb{C}^{M \times Q}}$为离散化的随机辐射场矩阵，${\text{σ}} \in {\mathbb{C}^{Q \times 1}}$为成像区域的散射强度矢量，${\text{w}} \in {\mathbb{C}^{M \times 1}}$为复高斯白噪声矢量。$M$为子脉冲数目，$Q$为成像单元个数。

4.   基于稀疏贝叶斯学习的迭代关联成像算法

为了求解双基RCI参数化稀疏表征模型，本文提出了ICI-SBL算法，可以获得目标成像结果，同时估计出目标运动参数。首先，利用跟踪雷达获取的目标速度，作为目标运动参数的先验信息。由于目标运动速度的方向初始值与真实值相差较小，故对2元函数${\text{S}}\left( {v,\theta } \right)$可利用泰勒(Taylor)公式在点$\left( {{v_0},{\theta _0}} \right)$展开并忽略高次项的影响，得到

$$ \begin{aligned} {\text{S}}\left( {v,\theta } \right) \approx& {\text{S}}\left( {{v_0},{\theta _0}} \right) + {\left. {\frac{{\partial {\text{S}}\left( {v,\theta } \right)}}{{\partial v}}} \right|_{v = {v_0}\atop \theta = {\theta _0}}} \cdot \left( {v - {v_0}} \right) \\ {\rm{}}& + {\left. {\frac{{\partial {\text{S}}\left( {v,\theta } \right)}}{{\partial \theta }}} \right|_{v = {v_0}\atop \theta = {\theta _0}}} \cdot \left( {\theta - {\theta _0}} \right)\\ =& {{\text{S}}_0} + {{\text{S}}_1} \cdot \Delta v + {{\text{S}}_2} \cdot \Delta \theta \end{aligned} $ $

(12)

式中，${v_0}$表示已知的目标运动速度初始值，${\theta _0}$表示已知的目标运动方向初始值。$\Delta v = v - {v_0}$表示目标运动速度真值与初始值之差，$\Delta \theta = \theta - {\theta _0}$表示目标运动方向真值与初始值之差。从而式(11)变为

$$ \begin{align} {\text{y}} &= \left( {{{\text{S}}_0} + {{\text{S}}_1} \cdot \Delta v + {{\text{S}}_2} \cdot \Delta \theta } \right) \cdot {\text{σ}} + {\text{w}} \\ {\rm{}}&= {\text{D}}\left( {v,\theta } \right) \cdot {\text{σ}} + {\text{w}} \end{align} $ $

(13)

基于贝叶斯框架，对各信号、参数建立合理的贝叶斯模型，通过贝叶斯推理获得对目标运动参数的估计，并反演得到稀疏重构信号。

4.1   贝叶斯模型

对于待重构的成像区域散射强度矢量${\text{σ}}$，建立一种2层贝叶斯稀疏先验模型。首先对于成像区域散射强度矢量${\text{σ}}$先验有${\text{σ}} \sim \rm{CN}\left( {{\text{σ}}\left| {0,{\text{Λ}}} \right.} \right)$，即散射强度矢量${\text{σ}}$先验服从均值为0，协方差为${\text{Λ}} =$$\rm{diag}\left( {\text{α}} \right) = \rm{diag}\left( {\left[ {{\alpha _1},{\alpha _2}, ·\!·\!· ,{\alpha _Q}} \right]} \right)$的复高斯分布，其中，${\text{α}}$为协方差矩阵对角线元素组成的矢量，${\alpha _q}$为${\text{α}}$第$q$个元素，其概率密度函数为

$$ \begin{aligned} p\left( {{\text{σ}}\left| {\text{α}} \right.} \right) =& \frac{1}{{{{\text{π}}^Q}\left| {\text{Λ}} \right|}}{\rm{exp}}\left( { - {{\text{σ}}^{\rm{H}}}{{\text{Λ}}^{ - 1}}{\text{σ}}} \right)\\ =& \prod\limits_{q = 1}^Q {{{\left( {2{\text{π}}{\alpha _q}} \right)}^{ - {1/2}}}{\rm{exp}}\left( { - \frac{{\sigma _q^2}}{{2{\alpha _q}}}} \right)} \end{aligned} $ $

(14)

式中，${{\text{σ}}^{\rm{H}}}$表示${\text{σ}}$的共轭转置矩阵，${\text{α}}$与高斯分布共轭，服从伽马分布

$$ p\left( {{\text{α}};\rho } \right) = \prod\limits_{i = 1}^Q {\varGamma \left( {{\alpha _i}\left| {1,\rho } \right.} \right)} $ $

(15)

式中，$\varGamma \left( {\xi \left| {{v_1},{v_2}} \right.} \right) = {\left[ {\varGamma \left( {{v_1}} \right)} \right]^{ - 1}}{v_2}^{{v_1}}{\xi ^{{v_1} - 1}}{e^{ - {v_2}\xi }}$, $\varGamma \left( \cdot \right)$表示伽马函数。

对于复高斯白噪声矢量${\text{w}}$，则${\text{w}} \sim \rm{CN}$$\left( {{\text{n}}\left| {0,\alpha _0^{ - 1}{{\text{I}}}} \right.} \right)$，即均值为0，方差为$\alpha _0^{ - 1}{{\text{I}}}$的复高斯分布，为了满足高斯分布的共轭特性，假设${\alpha _0}$服从伽马分布，$\alpha _0^{ - 1}$表示噪声功率，则

$$ p\left( {{\alpha _0};{\rm a,b}} \right) = \varGamma \left( {{\alpha _0}\left| {\rm a,b} \right.} \right) $ $

(16)

式中，${\rm a,b} \to 0$并且为一常数。

4.2   贝叶斯推理

利用期望最大化(EM)方法作为贝叶斯推理的准则，得到成像区域散射强度矢量${\text{σ}}$的后验分布是一个复高斯分布为$p\left( {{\text{σ}}\left| {{\text{y}},{\alpha _0},{\text{α}},v,\theta } \right.} \right) = {\rm{CN}}\left( {{\text{σ}}\left| {{{\text{μ}}_{{\text{σ}}}},{{\text{Σ}}_{{\text{σ}}}}} \right.} \right)$，即均值为${{\text{μ}}_{{\text{σ}}}}$，协方差为${{\text{Σ}}_{{\text{σ}}}}$的复高斯分布

$$ {{\text{Σ}}_{{\text{σ}}}} = {\left( {{\alpha _0}{{\text{D}}^{\rm{H}}}{\text{D}} + {{\text{Λ}}^{ - 1}}} \right)^{ - 1}} $ $

(17)

$$ {{\text{μ}}_{{\text{σ}}}} = {\alpha _0}{{\text{Σ}}_{{\text{σ}}}}{{\text{D}}^{\rm{H}}}{\text{y}} \hspace{42pt} $ $

(18)

式中，${\text{D}}$是随机辐射场${\text{S}}\left( {v,\theta } \right)$在点$\left( {{v_0},{\theta _0}} \right)$处泰勒展开1阶近似值${\text{D}}\left( {{v_0},{\theta _0}} \right)$的简写形式，利用${{\text{μ}}_{{\text{σ}}}}$可以对目标运动参数$v$和$\theta$进行更新，通过使$p\left( {{\text{y}}\left| {{\text{σ}},{\alpha _0},v,\theta } \right.} \right)$$p\left( v \right)p\left( \theta \right)$最大化，可以获得其估计值。也就是使

$$ \begin{align} {\rm{}}&\left\| {{\text{y}} - \left( {{{\text{S}}_0} + {{\text{S}}_1} \cdot \Delta v + {{\text{S}}_2} \cdot \Delta \theta } \right) \cdot {\text{σ}}} \right\|_2^2 \\ {\rm{}}&= \left\| {{\text{y}} - \left( {{{\text{S}}_0} + {{\text{S}}_1} \cdot \Delta v + {{\text{S}}_2} \cdot \Delta \theta } \right) \cdot {{\text{μ}}_{{\text{σ}}}}} \right\|_2^2 \end{align} $ $

(19)

最小，则

$$ v \!=\! \arg \min \left\{ {\left\| {{\text{y}} - \left( {{{\text{S}}_0} + {{\text{S}}_1} \cdot \Delta v + {{\text{S}}_2} \cdot \Delta \theta } \right) \cdot {{\text{μ}}_{{\text{σ}}}}} \right\|_2^2} \right\} $ $

(20)

$$ \theta \!=\! \arg \min \left\{ {\left\| {{\text{y}} - \left( {{{\text{S}}_0} + {{\text{S}}_1} \cdot \Delta v + {{\text{S}}_2} \cdot \Delta \theta } \right) \cdot {{\text{μ}}_{{\text{σ}}}}} \right\|_2^2} \right\} $ $

(21)

基于最小二乘准则，得

$${v^{{\rm{new}}}} = {\left[ {{{\text{A}}^{\rm{H}}}{\text{A}}} \right]^{ - 1}} \cdot {{\text{A}}^{\rm{H}}}{{\text{Y}}_{\rm{A}}} + {v_0}$ $

(22)

$${\theta ^{{\rm{new}}}} = {\left[ {{{\text{B}}^{\rm{H}}}{\text{B}}} \right]^{ - 1}} \cdot {{\text{B}}^{\rm{H}}}{{\text{Y}}_{\rm{B}}} + {\theta _0}$ $

(23)

式中，${v^{\rm{new}}}$, ${\theta ^{\rm{new}}}$分别表示对目标运动参数$v$和$\theta$的更新值，${\text{A}} \!=\! {{\text{S}}_1} \cdot {{\text{μ}}_{{\text{σ}}}}$, ${{\text{Y}}_{\rm{A}}} \!=\! {\text{y}} \!-\! \left( {{{\text{S}}_0} \!+\! {{\text{S}}_2} \cdot \Delta \theta } \right) \cdot {{\text{μ}}_{{\text{σ}}}}$, ${\text{B}} = {{\text{S}}_2} \cdot {{\text{μ}}_{{\text{σ}}}}$, ${{\text{Y}}_{\rm{B}}} = {\text{y}} - \left( {{{\text{S}}_0} + {{\text{S}}_1} \cdot \Delta v} \right) \cdot {{\text{μ}}_{{\text{σ}}}}$。

当成像区域散射强度矢量${{\text{σ}}}$后验分布的均值${{\text{μ}}_{{\text{σ}}}}$、目标运动参数$v$和$\theta$的估计值满足精度要求时，就完成了对目标的成像和运动参数的估计，若不满足，则需要对环境参数即噪声功率$\alpha _0^{ - 1}$，以及成像区域散射强度矢量${{\text{σ}}}$、先验协方差矩阵${{\text{Λ}}}$进行更新，对于噪声功率的倒数${\alpha _0}$和协方差矩阵对角线元素矢量${{\text{α}}}$，可以对其进行以下更新^[22]

$$ \alpha _0^{{\rm{new}}} = \frac{{M + {\rm{a}} - 1}}{{\left\| {{\text{y}} - {\text{D}} \cdot {\text{σ}}} \right\|_2^2 + {\rm{b}}}} \hspace{67pt} $ $

(24)

$$ \alpha _q^{{\rm{new}}} = \frac{{\sqrt {1 + 4\rm{\rho }{{\left| {{\sigma _q}} \right|}^2}} - 1}}{{2\rm{\text{ρ} }}},q = 1,2, ·\!·\!· ,Q $ $

(25)

式(24)中满足如下关系，$\left\| {{\text{y}} - {\text{D}} \cdot {\text{σ}}} \right\|_2^2 = \left\| {{\text{y}} - }\right.$$\alpha _0^{ - 1}\displaystyle\sum\nolimits_{q = 1}^Q{\left( {1 - \alpha _q^{ - 1}{\text{Σ}_{{\text{σ}}qq}}} \right)}$, ${\left| {{\sigma _q}} \right|^2} = {\left| {{\mu _{{\text{σ}}q}}} \right|^2} + {{\text{Σ}}_{{\text{σ}}qq}}$。$\rm{a}$, $\rm{b}$和$\text{ρ}$ 均为常数，${\left\| \cdot \right\|_2}$ 表示求2-范数。

再重新对目标进行成像和运动参数估计。

4.3   算法流程

综上所述，利用ICI-SBL算法对运动目标进行成像和运动参数估计时，首先应进行参数初始化^[22]，在迭代次数$k = 0$时，获得${\alpha _0}$, ${\text{α}}$的初始估计$\alpha _0^0$, ${{\text{α}}^0}$: $\alpha _0^0 = {{100}/{\rm{var} \left( {\text{y}} \right)}}$, ${{\text{α}}^0} = {{\left| {{\text{S}}_0^{\rm{H}}{\text{y}}} \right|}/{M{N_{\rm t}}Q}}$。目标运动参数$v$和$\theta$通过跟踪雷达获取其先验信息，${\text{ρ} } = 0.01$, $\rm{a} = \rm{b} = {10^{ - 4}}$。然后，令迭代次数$k = k + 1$，通过更新各参数，实现对未知参数的迭代求解。其终止条件是：循环次数达到预先设定值，或相邻2次成像区域散射强度矢量重构结果、运动参数估计达到预计精度。因此，本文所提ICI-SBL算法进行运动目标成像的具体步骤如表1。

表 1  ICI-SBL算法步骤

ICI-SBL算法步骤
输入：输入${\text{y}}$, ${{\text{S}}_0}$, ${{\text{S}}_1}$和${{\text{S}}_2}$；
(1) 迭代次数$k = 0$，初始化相关参数$\alpha _0^0$和${{\text{α}}^0}$，令$\rm{\rho } = 0.01$, $\rm{a} = \rm{b} = {10^{ - 4}}$, $v = {v_0}$, $\theta = {\theta _0}$，并设置终止迭代次数$K = 2000$，成像区域散射 强度矢量重构精度${\varepsilon _{{\text{σ}}}} = {10^{ - 2}}$，目标运动速度估计精度${\varepsilon _v} = {10^{ - 4}}$，目标运动方向估计精度${\varepsilon _\theta } = {10^{ - 6}}$；
(2) 令$k = k + 1$；
(3) 根据式(17)和式(18)，计算并更新第$k$次的成像区域散射强度矢量${{\text{σ}}^k}$后验分布的协方差${\text{Σ}}_{{\text{σ}}}^k$和均值${\text{μ}}_{{\text{σ}}}^k$；
(4) 目标运动参数更新过程：根据式(22)和式(23)，计算并更新第$k$次的目标运动参数${v^k}$和${\theta ^k}$；
(5) 判断是否满足终止条件：若$k > K$，或者$\left| {{\text{μ}}_{{\text{σ}}}^k - {\text{μ}}_{{\text{σ}}}^{k - 1}} \right| < {\varepsilon _{{\text{σ}}}}$，或者$\left| {{v^k} - {v^{k - 1}}} \right| < {\varepsilon _v}$且$\left| {{\theta ^k} - {\theta ^{k - 1}}} \right| < {\varepsilon _\theta }$，输出结果。否则，继续步骤(6)；
(6) 环境参数即噪声功率更新过程：根据式(24)，计算并更新噪声功率的倒数$\alpha _0^k$；
(7) 成像区域散射强度矢量${{{\text{σ}}}^k}$先验协方差矩阵更新过程：根据式(25)，计算并更新参数${{{\text{α}}}^k}$，跳转步骤(2)。
输出：输出成像区域散射强度矢量重构结果$\hat {{\text{σ}}} = {\text{μ}}_{{\text{σ}}}^k$，目标运动速度$\hat v = {v^k}$，目标运动方向$\hat \theta = {\theta ^k}$。

下载: 导出CSV 

| 显示表格

5.   仿真结果与分析

为了证明双基RCI参数化稀疏表征模型和ICI-SBL算法的有效性，通过计算机仿真进行验证，并对仿真结果进行了分析。雷达系统的具体仿真参数设置如下：发射机和接收机天线数目分别设置为11和4，基线长度为20 km，发射子脉冲时宽和数目分别为0.390625 μs和128，载频为10 GHz，带宽为3 GHz，频率跳变间隔为23.4375 MHz。假设一个以速度${{\text{v}}_t} = {\left[ {3\;{\rm{km/s}},{{\text{π}}/2}\;{\rm{rad}}} \right]^{\rm{T}}}$做平稳匀速运动的弹道目标有3个散射点，其在目标本地坐标系中的位置分别为${P_1}\left( {0\;{\rm{m}},2\;{\rm{m}}} \right)$, ${P_2}\left( {0.5\;{\rm{m}}, - 1.0\;{\rm{m}}} \right)$和${P_3}\left( { - 0.5\;{\rm{m}}, - 1.0\;{\rm{m}}} \right)$。参考点$O'$在雷达坐标系中的坐标为$\left( {23\;{\rm{km}},5\;{\rm{km}}} \right)$，成像区域大小为$6\;{\rm{m}} \times 6\;{\rm{m}}$，以$0.25\;{\rm{m}} \times 0.25\;{\rm{m}}$分辨率对成像区域进行划分，得到$25 \times 25$个成像单元，并假设所有散射点位于成像单元中心。

仿真1　收发天线阵元数目对随机辐射场性能的影响

为了能够定量表示随机辐射场的相关性，引用文献[23]中“条件数”的概念衡量随机辐射场相关性大小。根据条件数的定义，随机辐射场矩阵${{\text{S}}_0}$的条件数可表示为

$$\rm{cond}\left( {{{\text{S}}_0}} \right) = {\left\| {{{\text{S}}_0}} \right\|_2} \cdot {\left\| {{{\text{S}}_0}^{ - 1}} \right\|_2}$ $

(26)

当$\rm{cond}\left( {{{\text{S}}_0}} \right)$越小时，${{\text{S}}_0}$不相关性越显著，相反地，${{\text{S}}_0}$的相关性越强。在目标运动参数为初始值，对回波信号采样数目一定的条件下，通过对随机辐射场矩阵${{\text{S}}_0}$的条件数进行计算和比较，分析不同收发天线阵元数目对随机辐射场相关性的影响。

首先，发射机天线阵元以1 m的距离等间隔排布，数目分别为3, 5, 7, 9, 11，其它参数同具体仿真参数设置，计算条件数的结果如表2所示。

表 2  不同发射机天线阵元数目的条件数结果

发射机天线阵元数目${N_{\rm t}}$	3	5	7	9	11
条件数($\times {10^8}$)	2.1451	1.3506	0.8652	0.7393	0.5546

下载: 导出CSV 

| 显示表格

然后，接收机天线阵元以1 m的距离等间隔排布，数目分别为1, 2, 4, 8, 16，在回波信号采样数目一致的条件下，对回波信号时域采样分别为512, 256, 128, 64, 32，其它参数同具体仿真参数设置，计算条件数的结果如表3所示。

表 3  不同接收机天线阵元数目的条件数结果

接收机天线阵元数目${N_{\rm r}}$	1	2	4	8	16
条件数($\times {10^{10}}$)	1.3644	0.0695	0.0057	0.0027	2.9464

下载: 导出CSV 

| 显示表格

分析表2和表3条件数计算结果可以得出，随着发射机天线阵元数目的增加，随机辐射场条件数减小，表明随机辐射场的不相关性愈加显著；而随着接收机天线阵元数目的增加，随机辐射场的条件数先减小后增大。这是由于在一定范围内，可以通过增加接收机天线阵元数目降低回波信号时域采样数，并提升随机辐射场的不相关性能。而随着接收机天线阵元数目的进一步增加，随机辐射场中相同的发射信号带来的相关性日益显著，从而造成条件数的增大，使随机辐射场的相关性增强。因此，对于接收机天线阵元数目和回波信号时域采样数的设置应予以综合考虑。

仿真2　成像与目标运动参数估计结果

为了验证模型和算法的有效性，在没有噪声的情况下，利用本文提出的ICI-SBL算法对成像区域散射强度矢量进行重构和目标运动参数进行估计，并且选择OGSBI^[22]和SBL^[23]算法进行比较。分别进行100次蒙特卡洛实验，其目标结果如图2所示。

图 2  3种算法目标成像结果对比

下载: 全尺寸图片 幻灯片

定义均方误差MSE来衡量成像质量：$\rm{MSE} =$$10\lg \left( {{{\left\| {\tilde {{\text{σ}}} - \hat {{\text{σ}}}} \right\|_2^2}/{\left\| {\tilde {{\text{σ}}}} \right\|_2^2}}} \right)$，式中，$\tilde {{\text{σ}}}$为成像区域散射强度矢量真实值，$\hat {{\text{σ}}}$为成像区域散射强度矢量估计值。分别计算OGSBI, SBL和ICI-SBL 3种算法目标成像结果的MSE，其结果分别为–8.03 dB, 1.44 dB和–19.29 dB。通过比较图2中采用不同算法的成像结果，结合计算得到的MSE，可以得出，在已知目标运动参数先验信息的条件下，由于SBL算法对目标运动参数初始值误差敏感，因此利用SBL算法无法对成像区域散射强度矢量进行重构；而OGSBI算法和本文提出的ICI-SBL算法可以对成像区域散射强度矢量进行重构，但本文提出的ICI-SBL算法获得的目标成像结果更为理想。并且，可以对目标运动参数进行更加精确的估计，图3为ICI-SBL算法对目标运动速度和方向的估计结果随迭代次数变化的曲线。

图 3  ICI-SBL算法目标运动参数估计结果

下载: 全尺寸图片 幻灯片

分别定义目标运动速度和方向相对误差${v_{\rm{err}}}$和${\theta _{\rm{err}}}$来衡量目标运动参数估计性能：${v_{\rm{err}}} = {{\left| {\tilde v - \hat v} \right|}/{\tilde v}} \!\times\!$100%和${\theta _{\rm{err}}} = {{\left| {\tilde \theta - \hat \theta } \right|}/{\tilde \theta }} \times 100\%$，式中，$\tilde v$和$\tilde \theta$为目标运动速度和方向真实值，$\hat v$和$\hat \theta$为目标运动速度和方向估计值。计算采用本文提出的ICI-SBL算法对目标运动参数估计结果的${v_{\rm{err}}}$和${\theta _{\rm{err}}}$，分别为0.06%和0.05%，与目标运动参数初始值的0.84%和0.40%相比，目标运动参数估计相对误差明显减小，运动参数更加精确。

仿真3　噪声对成像结果的影响

考虑环境噪声对所提算法成像性能产生的影响，假设在环境噪声服从复高斯分布前提下，信噪比SNR从0 dB到25 dB，以5 dB等间隔变化，分别采用OGSBL, SBL算法和本文提出的ICI-SBL算法对目标成像，并进行了100次蒙特卡洛实验，得到目标成像结果的MSE随SNR的变化曲线如图4所示。

图 4  3种算法目标成像结果MSE随信噪比变化曲线

下载: 全尺寸图片 幻灯片

图4中虚线表示采用OGSBL算法目标成像结果MSE随SNR变化曲线，随着SNR的提高，MSE下降不明显，且略有起伏；实线表示采用SBL算法目标成像结果MSE随SNR变化曲线，随着SNR的提高，MSE有所上升，这是因为SBL算法对目标运动参数初始值误差敏感，即使SNR提升，该算法仍然受运动参数误差影响而无法对目标进行有效成像；点划线表示本文所提的ICI-SBL算法对目标成像结果的MSE随SNR变化曲线，随着SNR的提升，MSE下降趋势明显，即成像性能得到提升。

6.   结束语

本文研究了运动目标的双基雷达关联成像问题，在构建双基雷达关联成像参数化稀疏表征模型的基础上，提出基于稀疏贝叶斯学习的迭代关联成像算法(ICI-SBL)，实现对运动目标的成像，同时完成对目标运动参数的精确估计。本文提供了一种新的双基雷达成像方法，对研究双基雷达成像系统具有一定的理论和实际意义。然而，弹道目标具有自旋、进动等复杂运动形式，会对参数化稀疏表征模型产生影响，因此仍然需要进一步研究。