1.   引言

调频连续波(Frequency Modulated Continuous Wave, FMCW)毫米波雷达具有抗干扰能力强、集成度高、成本低等优势，涉及生命体征检测、水文监测和智能驾驶等多种场景^[1–4]。然而，天线孔径限制了角度分辨率，因此提高雷达系统前视图像分辨率至关重要。随着多输入多输出(Multiple-Input Multiple-Output, MIMO)技术发展，利用收发器件的空间位置信息来构造虚拟阵列，可以具备比真实阵列更高的孔径以实现高角度分辨率。MIMO雷达依赖于发射信号的可分离性，这需要正交波形生成技术支持，本文主要研究已成熟应用的时分复用(Time-Division Multiplexing, TDM)策略^[5]下的前视场景。

在传统车载MIMO毫米波雷达中，主要采用距离多普勒角度三维快速傅里叶变换(Fast Fourier Transform, FFT)方法完成目标估计^[6]，对每个脉冲做第1维距离FFT，接着对距离FFT结果做第2维多普勒FFT，然后对每个接收天线回波的二维FFT处理结果做第3维角度FFT，可以检测得到目标的距离、多普勒、角度和信号幅度等信息。然而随着车载毫米波雷达应用场景改变，急需实现高性能高分辨率前视成像，并增强目标探测和识别能力，因此传统的信号处理框架必须改变以满足现代汽车自动驾驶需求，本文从以下两个方面来阐述现有方法的不足。

第1个问题是相干处理间隔(Coherent Processing Interval, CPI)内距离多普勒走动问题。车载雷达已成为无人驾驶和高级驾驶员辅助系统的关键推动者，不同于传统的车载毫米波前视、角雷达等近距离跟踪、防撞等简单任务，其应用场景已经从密集的城市环境拓展到高速公路等复杂路况。然而高速运动目标限制了CPI的选取，从而限制了雷达的最大探测距离。此外，它也限制了长时间相干积分方法的实施。这里考虑相干处理的方法，距离走动补偿的主要方法为基于Keystone变换和基于Radon变换方法。Keystone变换原理主要是对每个快时间样本下的慢时间样本进行重新采样^[7]。由于其简单性，Keystone变换技术已成功应用于合成孔径雷达、动目标指示等^[8]。在这些应用中，雷达脉冲重复频率相对目标多普勒频移较大，目标多普勒走动变化缓慢。然而，该算法不适用于车载雷达，因为车载高速目标多普勒频率可能超过脉冲重复频率，发生多普勒模糊，使得Keystone变换无法使用。另一类基于Radon-变换的方法也被广泛应用，如Radon傅里叶变换^[9]、Radon傅里叶逆变换^[10]等。其沿着距离走动曲线相干积分，从而提高了处理增益。然而，在距离慢时间域中，每个不同距离速度的目标对应一个距离走动曲线。因此，基于Radon的方法必须在慢时间内完成所有候选的距离-速度对的相干积分，这会极大增加计算复杂性。

第2个问题是高分辨成像。车载毫米波雷达高分辨率环绕成像，远距离探测小型快速移动物体，检测坑洞、道路碎片和头顶障碍物，以及目标分类等各种任务都需要分辨率的前视图像作为输入。传统毫米波雷达主要采用数字波束形成(FFT, Capon或Bartlett滤波器)来实现角度估计，这有利于在嵌入式设备中通过单快拍实现。然而，这些方法都不具备高分辨率角度测量能力，无法突破系统参数决定的瑞利极限。虽然采用子空间的超分辨算法，如多重信号分类方法(Multiple Signal Classification, MUSIC)、迭代自适应方法(Iterative Adaptive Approach, IAA)^[11]等一定程度上突破了物理孔径的角度分辨率，但是在强噪声环境中依然难以获得令人满意的成像结果，甚至强点噪声会引起角度目标个数估计错误。因此为了得到更好的毫米波雷达成像效果，有学者将一维的角度分辨率算法扩展到了距离-角度等高维空间联合估计，其中多维FFT是最直接、最简单的方法，当然也是性能较差的方法。于是衍生出FFT和超分辨率算法结合的方法，例如DFT-ESPIRT、聚类ESPIRT^[12,13]等。然而，这些基于子空间的方法的结果依赖于多快拍估计的协方差矩阵的准确性，这在实际中很难得到满足。

根据雷达目标特性理论，可以通过多个聚类点来逼近雷达图像，以满足空间域中的稀疏特征。因此，压缩感知(Compressive Sensing, CS)有望成为最有前途的高分辨率成像技术之一，能提供超过瑞利极限的分辨率^[14,15]。文献[16]给出了不同L1-范数正则下雷达目标角度稀疏重建方法的比较，文献[17]借用了图像处理中常用的梯度约束的正则化方法，取得了较好的雷达目标轮廓成像效果。这些正则化方法无不例外都会出现噪声放大问题，导致出现大量虚警。文献[18]采用L1-范数正则和低秩双约束模型，Shu等人^[19]则是将截断奇异值分解和L1-范数正则结合，实现雷达前视超分辨成像的噪声抑制。虽然组合约束可以减弱全局噪声对单一正则的影响，但是对于分布稀疏的强点噪声依旧敏感，对强噪声附近目标存在虚警、漏警，同时增加正则参数会降低算法实用性。另一类通用方法是基于贝叶斯定理，利用场景中目标和噪声的统计特性，将稀疏重建问题转化为解空间估计问题，并构建相对应的目标函数。基于贝叶斯的方法可以根据实际情况，对噪声及目标散射信息参数化建模，可以有效利用噪声的先验统计信息，在低信噪比下信号重构性能好，与正则化方法相比，避免了手工设置正则参数的麻烦，具有比传统CS重构方法更小的估计误差，如稀疏贝叶斯学习(Bayesian Compressive Sensing, BCS)^[20,21]算法、快速贝叶斯匹配追踪算法(Fast Bayesian Matching Pursuit, FBMP)^[22,23]、梯度投影法(Gradient Projection for Sparse Reconstruction, GPSR)^[24,25]等。

在本文中，首先对车载平台与目标相对运动引起的瞬时多普勒信息和距离走动信息、通道相位跳变等问题进行描述和论证推导，同时提出一种有效简单的基于多普勒域的补偿方法。该方法首先通过对快时间多普勒图进行插值完成多普勒/距离走动校正，并且通过点对点补偿相位完成了多普勒距离解耦合和通道相位校正。更一步，针对现有算法在低信噪比下雷达成像性能较差的问题，研究了基于贝叶斯理论的复杂环境下雷达距离角度成像方法。上述CS方法中，FBMP在低信噪比下也具有更好的重构性能，但是其重复搜索策略导致时间开销较大，因此本文提出一种基于最大后验准则(Maximum A Posteriori, MAP)的改进的贝叶斯匹配追踪算法，解决TDM-MIMO毫米波雷达前视成像超分辨问题。

2.   TDM MIMO信号模型

常见的车载毫米波MIMO雷达前视成像场景，其具备单次快速成像能力，几何示意图如图1(a)所示。雷达平台采用收发一体${M_{\text{t}}}$发${M_{\text{r}}}$收阵列天线体制，一个具有两个发射天线的FMCW TDM MIMO雷达示例如图1(b)所示，因此对于载频为${f_0}$，调制斜率为$\gamma$，脉冲宽度为${T_{\text{p}}}$，脉冲重复周期为${T_{\text{c}}}$的FMCW信号表示为

图  1  系统配置示意图

Figure  1.  System configuration diagram

下载: 全尺寸图片 幻灯片

将第m发射阵元和第n接收阵元产生的回波信号与发射信号${y_{{\rm{RF}}}}(t)$进行混频处理，再经过低通滤波后即得到回波的差拍信号${y_{m,n}}(t)$为

不妨假设一运动目标q与雷达径向距离为${R_q}$，并以恒定的径向速度${v_q}$，从角度${\theta _q}$远离雷达，则回波时延为$\tau = {\left({2\left( {{R_q} + {v_q}{t_{\text{f}}} + {v_q}{t_{\text{s}}}} \right) + {\varDelta _q}_{,1} + {\varDelta _q}_{,2}}\right)/ {\rm{c}}}$，其中c为光速，${t_{\rm{s}}}$表示慢时间，${t_{\rm{f}}}$表示快时间，${\varDelta }_{q}{}_{,1}=2\left(m-1\right){v}_{q}{T}_{\text{p}}$表示由于TDM策略下发射天线切换导致的运动残留距离，${\varDelta }_{q}{}_{,2}=\left[\left(m-1\right){d}_{\text{t}}+ \left(n-1\right){d}_{\text{r}}\right]\text{sin}\theta$表示接收发射天线阵元的间距引起的波程差，其中${d_{\rm{t}}}$, ${d_{\rm{r}}}$分别表示发送和接收阵列的阵元间距，并且有${d_{\text{t}}} = {M_{\text{r}}}{d_{\text{r}}}$和${d_{\text{r}}} = {\lambda \mathord{\left/ {\vphantom {\lambda 2}} \right. } 2}$。

将时间延迟$\tau$代入式(2)，${\tau ^2}$项可以看作一个附加的固定相位项(由于${{\rm{c}}^{ - 2}}$非常小)，可得

(1) 第1个指数项表示由目标初始距离导致的差拍信号，被用来进行距离估计。第5个指数项表示由于相对运动带来的多普勒信息，被用来进行速度估计。第6个指数项表示由于MIMO雷达不同收发天线之间的波程差引起的相位差，一般用来进行目标角度估计。

(2) 第2个指数项$\exp \left\{ {{{\rm{j}}2\pi }}\left( \dfrac{{2\gamma {{\Delta }}{{{R}}_{q,m,n}}}}{{\rm{c}}} + \dfrac{{2{v_q}{f_0}}}{{\rm{c}}} \right){t_{\text{f}}} \right\}$，其中${{\Delta }}{{{R}}_{q,m,n}} = \dfrac{1}{2}\left( {{\varDelta _q}_{,1} + {\varDelta _q}_{,2}} \right)$。由车载雷达测量目标速度小，其最大不模糊为${V_{{\text{max}}}} = {\lambda \mathord{\left/ {\vphantom {\lambda {4{M_{\text{t}}}{T_{\text{p}}}}}} \right. } ({4{M_{\text{t}}}{T_{\text{p}}}})}$ ，代入$\Delta {R_{q,m,n}}$中，可知：

因此，${{\Delta }}{{{R}}_{q,m,n}}$所产生的距离偏移远小于车载雷达距离分辨率${\rho _{\text{r}}}$，其在此项中可以忽略不计，但${{2{v_q}{f_0}} \mathord{\left/ {\vphantom {{2{v_q}{f_0}} c}} \right. } {\rm{c}}}$使目标的距离估计有一个附加的偏移，影响估计结果。

(3) 第3个指数项$\exp \left\{ {{{{\rm{j}}2\pi }}\dfrac{{2\gamma {v_q}{t_{\text{f}}}{t_{\text{s}}}}}{{\rm{c}}}} \right\}$会导致目标在慢时间发生距离走动，使得目标在距离和多普勒频率偏移，会存在散焦和能量泄露问题。

(4) 第4个指数项$\exp \left\{ {{{{\rm{j}}2\pi }}\dfrac{{2\gamma {v_q}t_{\text{f}}^2}}{{\rm{c}}}} \right\}$表示目标在快时间的距离走动信息，在车载毫米波雷达中，由于${T_{\rm{p}}}$极短，可以认为目标回波是单频连续波信号，不予考虑。

(5) 第7个指数项$\exp \left\{ {{{{\rm{j}}2\pi }}\dfrac{{2{v_q}}}{\lambda }{T_{\text{p}}}\left( {m - 1} \right)} \right\}$会导致进行角度估计的时候发生不同通道间的相位跳变，影响角度估计精度。

传统的雷达信号处理将式(3)简化为3个维度的单频信号。忽略第2项和第3项，距离和速度估计问题可以分别在快时间和慢时间维度解耦，忽略第7项，角度估计可以和速度参数解耦。然而在高性能汽车雷达中，这些忽略项会使得目标回波发生距离走动、相位跳变等现象，导致估计结果性能急剧降低。

接下来，本文将采用一个更精确的信号描述模型。首先，将式(3)通过MIMO通道合成，得到联合角度距离多普勒三维信号模型，并由采样频率${f_{\text{s}}}$离散化得到：

其中，$\left[ \cdot \right]$表示向下取整， K是每个chirp的ADC采样总数，索引值为k；L是用于多普勒测量的TDM MIMO chirp周期数，索引值为l。

3.   基于多普勒域补偿的多维联合估计算法

3.1   多普勒相位补偿

由于传统毫米波方法缺乏对于距离和多普勒走动校正措施，会导致传统距离-多普勒谱发生频谱展宽，影响估计精度。为进一步分析距离/多普勒走动信息，首先对式(4)沿快时间进行FFT计算，得到距离慢时间图：

同样，也可以对式(4)沿慢时间进行FFT计算，得到快时间多普勒图：

其中，${\text{sinc}}\left( x \right) = \dfrac{{{{\sin\pi }}x}}{{{\pi}x}}$。

从式(5)可知，目标距离发生走动的条件为

其中，${\rho _{\text{r}}} = \dfrac{{\rm{c}}}{{2B}}$ 表示距离分辨率。

同样从式(6)可知目标多普勒发生走动的条件为

其中，${\rho _{\text{d}}} = \dfrac{1}{{L{T_{\text{c}}}}}$表示目标多普勒频率分辨率。

从式(8)可以得到，目标多普勒走动的条件与目标发生距离走动的条件一致。同样注意到式(5)中，距离慢时间图中的距离走动曲线由两个未知因素决定参数，即距离${R_q}$和速度${v_q}$。这意味着在给定一个距离和速度的情况下，只能得到一条对应的距离走动曲线。因此基于Radon和聚焦的方法，要想实现相干积分必须严格按照曲线进行，这会大大增加计算复杂度。相比之下，从式(6)中可以看到多普勒走动曲线仅由目标速度决定，而且运动导致的多普勒频偏、通道相位跳变也只与目标速度相关。因此，接下来采用一种多普勒补偿方法，完成三维信号的解耦合。

首先，采用最近邻域或线性插值算法完成多普勒走动校正。在本文中，通过在每个多普勒下，采用8点的sinc插值法完成，如下：

假设多普勒插值是精确的，式(6)可以重写为

接下来考虑补偿第2个和第7个指数项，从式(3)中直接进行相位补偿是不可能的。因此，利用其多普勒相关性，构建通道-快时间-多普勒相位补偿系数如下：

通过元素点乘完成相位校正，其计算成本可以忽略不计。因此，完全补偿后的相位为

3.2   多普勒模糊补偿

在本文所采用的TDM-MIMO策略下，所能估计的最大不模糊速度为

当发生模糊的速度${v_q}$所对应的真实速度为${v_{{\text{true}}}} = {v_q} \pm 2{V_{{\text{max}}}}$时，从式(6)可以得到：

比较式(6)和式(14)，可以看到多普勒${v_q}$和${v_{{\text{true}}}}$的走动曲线的起始速度(模糊后的)是一致的，但是两者具有不同的斜率。在不进行模糊补偿下，如果按照式(9)进行插值操作，会使得补偿后的走动曲线偏差更大，只有在正确的多普勒下走动曲线才会被校正。也正是根据这种特性，可以从最终的多普勒响应中，选取幅度更大的那个作为正确值。值得注意的一点，不正确的多普勒估计也会影响式(11)所进行的通道相位补偿，具体见第5节。

3.3   距离角度联合估计模型

从式(12)可以看到，回波信号在3个维度实现了解耦合，后续可按照常规方法完成距离FFT、二维恒虚警检测(Constant False Alarm Rate, CFAR)和角度FFT等操作，检测出目标的距离、多普勒、角度和信号幅度等信息。但是由于发射阵元/接收机阵元或发射带宽较窄等雷达参数限制了系统分辨率，导致传统FFT信号处理方法得到的图像质量较差，只能实现对强散射点的检测，远远满足不了前视成像的需求。因此基于式(12)，本文提出了一种距离角度联合估计算法。

首先通过CFAR或者峰值检测方法提取出候选多普勒单元，然后对于特定的多普勒单元，进行距离角度联合估计。首先假设某一多普勒单元内，存在Q个可能目标，因此：

其中，${\varphi _q}$表示初始相位，w表示加性高斯白噪声。

使用信号模型(15)，将遍历${M_{\text{t}}}{M_{\text{r}}}$个虚拟通道的快时间序列的观测数据形成矩阵表示：

其中，${{\boldsymbol{Y}}_{\text{s}}}$, W分别表示回波信号和噪声信号的二维矩阵形式，${\boldsymbol{X}} \in {C^{{Q_1} \times {{{Q}}_2}}}$是Q个稀疏目标的距离角度空间散射矩阵，并且有${{\boldsymbol{d}}_q}{\text{,}}\;{{\boldsymbol{g}}_q}$分别表示角度为$\theta _q$和距离维$R_q$下的导向矢量：

利用Kronecker积的性质，可以得到

其中，${{\boldsymbol{Y}}} \in {C^{({M_{\rm{t}}}{M_{\rm{r}}}K) \times 1}}$是回波信号的矢量化，${{\boldsymbol{x}}} = {\text{vec}}\left( {{\boldsymbol{X}}} \right) = {\left( {{{{x}}_1}{\text{ }}{{{x}}_2}{\text{ }} \cdots {\text{ }}{{{x}}_Q}} \right)^{\text{H}}}$是目标散射系数矩阵矢量化，${{\boldsymbol{A}}} = {{{\boldsymbol{G}}}^{\text{T}}} \otimes {{\boldsymbol{D}}}$是一个${M_{\text{t}}}{M_{\text{r}}}K$行Q列的基矩阵，其中$\left( {Q\left( {u - 1} \right) + v} \right)$列由${\boldsymbol{d}}_u^{\rm{T}} \otimes {{\boldsymbol{g}}_v}$计算得到。

综上所述，式(17)可以写成${{\boldsymbol{Y}} = {\boldsymbol{Ax}}}$模型，因此可以引入CS来解决传统研究方法的不足。而过度依赖于预设参数的非贝叶斯方法的解决方案在非理想条件下会导致较差的性能，因此为了获得独立于被估计参数的性能度量，通常优选基于贝叶斯框架的方法，并可实现超出雷达参数(例如带宽和天线阵元的数量)设置的极限分辨率，因此本文提出基于MAP的贝叶斯匹配追踪算法来进行求解。

综上所述，针对车载毫米波雷达前视成像下的距离走动以及高性能成像问题，提出了基于多普勒补偿的多维联合估计算法，其算法流程图如图2所示。算法主要分为3个部分：(1)多普勒FFT以及用于速度解模糊的多普勒拓展；(2)采用插值完成多普勒走动校正以及用于解决距离、多普勒、角度解耦合的相位补偿；(3)采用基于MAP的贝叶斯匹配追踪算法来完成距离角度联合估计。

图  2  基于多普勒补偿的多维联合估计算法流程图

Figure  2.  The flowchart of joint estimation algorithm based on Doppler compensation

下载: 全尺寸图片 幻灯片

4.   基于MAP的贝叶斯匹配追踪算法

通过引入伯努利-高斯(Bernoulli-Gaussian, BG)模型，将二进制稀疏向量定义为${{\boldsymbol{s}}} = {\left[ {{s_1}{\text{ }}{s_2}{\text{ }} \cdots {\text{ }}{s_Q}} \right]^{\text{T}}} \in {\left\{ {0,1} \right\}^Q}$，其中${s_q}$指示原子${{\boldsymbol{a}}_q}$是激活或者失活状态(分别对应1和0)，并服从伯努利分布：$p\left( {{s_q}} \right) = {\rm{Ber}}\left( {{p_q}} \right)$。假设x服从如下模型：

其中，${x_q}|{s_q}{\text{～}}N\left( {0,\sigma _{{s_q}}^2} \right)$ ，为了简化模型令$\left\{ {{p_q}} \right\}_{q = 1}^Q = {p_1}$, $\sigma _0^2 = 0$ 。

通过式(13)可知，条件概率分布 $p\left( {{{\boldsymbol{x}}}|{{\boldsymbol{Y}}},{{\boldsymbol{s}}}} \right)$的均值如下：

其中，${{\boldsymbol{R}}}\left( {\boldsymbol{s}} \right)$为对角矩阵，满足${\left[ {{{\boldsymbol{R}}}\left( {\boldsymbol{s}} \right)} \right]_{q,q}} = \sigma _{{s_q}}^2$, ${{\boldsymbol{\varPhi}} }\left( {{\boldsymbol{s}}} \right) = {{\boldsymbol{AR}}}\left( {\boldsymbol{s}} \right){{{\boldsymbol{A}}}^{\text{H}}} + {\sigma ^2}{{{\boldsymbol{I}}}_Q}$。

为了得到联合矢量x的最大似然估计，使得

结合式(13)和式(14)，可知唯一的困难是${{{\boldsymbol{s}}}_{{\text{map}}}}$的求解。然而隐向量s有${2^Q}$种取值，贪婪遍历选择的方法是不太可能的。因此下文将使用基于树结构的快速搜索算法寻找最可能的一个较小支持集${{{\boldsymbol{s}}}_{{\text{map}}}}$，并同时计算该状态下的相关概率。

首先通过选择具有最大概率$p\left( {{{\boldsymbol{s}}}|{{\boldsymbol{Y}}}} \right)$的s∈S来定义${{{\boldsymbol{s}}}_{{\text{map}}}}$：

根据贝叶斯定理，其后验概率可以写成

所以，为了简洁，采用对数后验来定义一个搜索过程的度量，如下所示：

因此，确定s的搜索策略如下：

步骤1　初始化参数向量s，即从${{\boldsymbol{s}}} = {\mathbf{0}}$。

步骤2　依次从0到1激活s向量中的一个元素，并得到最大度量$\partial \left( {{\boldsymbol{s}}} \right)$的下标和相对应的矢量${{{\boldsymbol{s}}}_1} = \left[ {i_1^*} \right]$。

步骤3　选择${{\boldsymbol{s}}_1}$中的剩余0元素，产生$Q - 1$个大小为2的可选项，产生导致最大$\partial \left( {{\boldsymbol{s}}} \right)$的${{{\boldsymbol{s}}}_2} = \left[ {i_1^*,i_2^*} \right]$。

步骤4　重复进行，直到产生${{\boldsymbol{s}}_D} = \left[ {i_1^*,i_2^*, \cdots ,i_D^*} \right]$，其中D可以选择以稍大于基数${\left| {\boldsymbol{s}} \right|_0}$ 。

最终，${{{\boldsymbol{s}}}_{{\text{map}}}} = {{{\boldsymbol{s}}}_D}$，基于树结构搜索策略如图3所示，其中Q=5和D＝3。

图  3  基于树结构的搜索策略

Figure  3.  Tree-search-based search strategy

下载: 全尺寸图片 幻灯片

从上述迭代过程观察到，每次更新过程只涉及s中一个元素从0到1的变化。因此，只需要更新增量而不是整个度量来提高计算速度。首先定义增量为：${\varDelta _{{q}}}\left( {\boldsymbol{s}} \right) = \partial \left( {{{\boldsymbol{s}}'}} \right) - \partial \left( {{\boldsymbol{s}}} \right)$，其中 ${{\boldsymbol{s}}'}$和s除了q坐标元素外其余元素都一样，并且满足${\left[ {{{\boldsymbol{s}}'}} \right]_q} = 1, {\left[ {{\boldsymbol{s}}} \right]_q} = 0$。

令$\sigma _\Delta ^2 = \sigma _1^2 - \sigma _0^2$，有

然后根据矩阵逆定理，可得

其中

所以结合式(18)和式(19)，可得

通过式(19)可知，${{\boldsymbol{c}}}_q^{}$的每一次计算复杂度是${{O}}\left( {{Q^2}} \right)$，不利于实时计算。因此，提出基于迭代计算的一种复杂度${{O}}\left( Q \right)$的更新方法，如图4所示：

图  4  ${{\boldsymbol{c}}}_q$的迭代更新示意图

Figure  4.  The iterative process of calculating ${{\boldsymbol{c}}}_q$

下载: 全尺寸图片 幻灯片

在图4中，假设${{{\boldsymbol{s}}}^{{\text{pre}}}}$是s的前一状态，只有第${q^{{\text{pre}}}}$的元素修改。显然，可以计算并存储由状态${{{\boldsymbol{s}}}^{{\rm{pre}}}}$变化所产生的相应$\left\{ {\beta _q^{{\text{pre}}}} \right\}_{q = 1}^Q$ 和$\left\{ {{{\boldsymbol{c}}}_q^{{\text{pre}}}} \right\}_{q = 1}^Q$。所以${{\boldsymbol{c}}}_q$的更新方程如下：

整个算法开始于${{\boldsymbol{s}}} = {\mathbf{0}}$，初始化如下：

该方法与FBMP相比，减少了贪婪寻找次数，计算复杂度下降D倍(D为重复搜索次数)。因此，改进后的贝叶斯匹配追踪方法(IBMP)放弃寻找多个支撑域下的基于MMSE的估计量，反而采用基于MAP的单一估计，在保持重构精度的同时可以提高计算速度，具体的算法处理步骤如算法1所示。

1  基于MAP的贝叶斯匹配追踪方法

1.  The Bayesian Matching Pursuit (BMP) based on MAP

初始化： 　$\partial \left( {\mathbf{0} } \right) = \left( { - M{\text{In} }\;\pi - M{\text{In} } \;{\sigma ^2} - \dfrac{1}{ { {\sigma ^2} } }\left\| { {\boldsymbol{Y} } } \right\|_2^2} \right) + Q{\text{In} }\left( {1 - {p_1} } \right)$
${ { {\boldsymbol{c} } }_q} = \dfrac{ { { { {\boldsymbol{a} } }_q} } }{ { {\sigma ^2} } },{\beta _q} = \sigma _1^2{\left( {1 + \sigma _1^2{ {\boldsymbol{a} } }_q^{\text{H} }{ { {\boldsymbol{c} } }_q} } \right)^{ - 1} }$
${\partial _q} = \partial \left( {\mathbf{0} } \right) + {\text{In} }\dfrac{ { {\beta _q} } }{ { \sigma _1^2} } + {\beta _q}\left\| {\left( { { { {\boldsymbol{Y} } }^{\text{H} } }{ { {\boldsymbol{c} } }_q} } \right)} \right\|_2^2 + {\text{In} }\left( {\dfrac{ { {p_1} } }{ {1 - {p_1} } } } \right)$
$\varOmega = [\;],\;{ { {\boldsymbol{s} } }^{\left( 0 \right)} } = {\mathbf{0} }$
更新迭代过程：
${q^*} = {{{\rm{max}}\_{\rm{index}}} }\left( { {\partial _q} } \right)$
$\varOmega = \varOmega \cup {q^*},{ { {\boldsymbol{s} } }^{\left( { {d} } \right)} } = { { {\boldsymbol{s} } }^{\left( { { {d} } - 1} \right)} } + \delta \left[ { {q^*} } \right],{\partial ^{\left( d \right)} } = {\partial _{ {q^*} } }$
${ { {\boldsymbol{c} } }_q} = { { {\boldsymbol{c} } }_q} - {\beta _{ {q^*} } }{ { {\boldsymbol{c} } }_{ {q^*} } }({ {\boldsymbol{c} } }_{ {q} })^{\text{H} }{ { {\boldsymbol{a} } }_q}$
${\beta _q} = \sigma _1^2{\left( {1 + \sigma _1^2({ {\boldsymbol{a} } }_q)^{ {\rm{H} } }{ { {\boldsymbol{c} } }_q} } \right)^{ - 1} }$
${\partial _q} = {\partial ^{\left( d \right)} } + {\rm{In} }\dfrac{ { {\beta _q} } }{ {\sigma _\Delta ^2} } + {\beta _q}\left\| { { { {\boldsymbol{Y} } }^{\text{H} } }{ { {\boldsymbol{c} } }_q} } \right\|_2^2 + {\text{In} }\left( {\dfrac{ { {p_1} } }{ {1 - {p_1} } } } \right)$
输出：
${\hat { {\boldsymbol{x} } }_{ {\text{map} } } } = \sigma _1^2 \displaystyle\sum\limits_{i \in \varOmega } { {\delta _i}{ {\left( { {\boldsymbol{c} }_i^{} } \right)}^{\text{H} } }{ {\boldsymbol{Y} } } }$

下载: 导出CSV 

| 显示表格

5.   实验与分析

本节通过简单的仿真信号，以验证所提出的方法的有效性，重点讨论针对不同SNR的估计性能比较，最后通过实测数据进一步验证该方法的高性能前视成像能力。所有实验结果均与常规2D FFT方法、GPSR方法、BCS方法以及FBMP方法进行对比。GPSR算法中的L1-惩罚选择为$\tau = 0.1{\left\| {{{{\boldsymbol{A}}}^{\text{H}}}{{\boldsymbol{x}}}} \right\|_\infty }$，FBMP的重复搜索次数D=5。

5.1   复杂度分析

本文是通过线性插值以及元素点乘完成相位多普勒相位补偿，其计算复杂度为$O\left( { {M_{\rm{t}}}{M_{\rm{r}}}KL} \right)$，因此基于FFT的多普勒处理计算复杂度为$O\left( { {M_{\text{t}}}{M_{\text{r}}}KL\log \left( L \right) + {M_{\text{t}}}{M_{\text{r}}}KL} \right)$。借助于贝叶斯匹配追踪算法完成距离角度联合估计，其计算复杂度与后向投影(Back Projection, BP)算法相当，因此整个算法总体计算复杂度是$O\left( {M_{\text{t}}}{M_{\text{r}}}KL\left( \log \left( L \right) + {M_{\text{t}}}{M_{\text{r}}}K \right) + {M_{\text{t}}}{M_{\text{r}}}KL \right)$。 相比较于常规3D FFT的计算复杂度为$O\left( { {M_{\text{t}}}{M_{\text{r}}}KL\log \left( {{M_{\text{t}}}{M_{\text{r}}}KL} \right)} \right)$，本文算法的计算复杂度较大，即以牺牲算力为代价来谋取较好的算法性能。

5.2   仿真分析

本节实验中，利用仿真点目标来验证所提方法的性能。雷达信号的参数配置见表1，其中信噪比(Signal to-Noise Ratio, SNR)设置为5 dB。此处采用2发4收MIMO阵列，可以计算出系统角度分辨率15.7º，距离分辨率为0.1 m。

表  1  仿真实验雷达参数

Table  1.  Radar parameters for simulation experiment

雷达参数	数值
载频(GHz)	77
带宽(MHz)	1500
MIMO天线配置	2T × 4R
脉冲宽度(μs)	55
脉冲采样点数K	256
脉冲个数L	256
角度分辨率(°)	15.7
距离分辨率(m)	0.1
最大不模糊速度(m/s)	8.85

下载: 导出CSV 

| 显示表格

为了初步验证多普勒补偿方法的有效性，首先分别选取了目标在距离[10, 10, 10, 15] m，速度为[16, 16, 7, –15] m/s，角度为[–10°, 10°, 20°, 20°]的4个目标。由系统参数可知，所有目标都会发生距离走动现象，可以观察到3条快时间多普勒曲线都发生了倾斜，从曲线倾斜斜率和初始多普勒判断，两条多普勒曲线发生了多普勒模糊。首先对所有候选多普勒单元通过多普勒模糊拓展，使得原始的多普勒测量范围扩大一倍。在经过本文提出的多普勒走动补偿后，可以观察到在只有真实目标的多普勒曲线才能被正确校正，而多普勒曲线和初始多普勒不匹配时，多普勒走动现象不会被补偿甚至会加剧，如图5(b)所示。与此同时，原始距离多普勒算法得到的结果目标在距离和多普勒维度都发生了能量泄露，二维频谱发生明显展宽，如图5(c)所示。在经过多普勒补偿后，目标在正确的多普勒和距离上能量更为集中，频谱走动得到较好的抑制，图5(d)中红色圆圈标注正确目标位置，白色圆圈标注多普勒模糊假目标。因此，可以很简单地通过多普勒幅度进行比较，选取的幅度较大的多普勒目标作为无模糊的目标。值得说明的一点是图5(d)中的目标2是包含预设信号的前两组信号，他们有相同的速度和距离，但是具有不同的角度。

图  5  点目标仿真

Figure  5.  Point target simulation

下载: 全尺寸图片 幻灯片

以目标1为例，图6(a)显示了目标1的速度剖面，即图5(d)在15 m处的水平切片。可以看到一个强烈的虚假峰值在2.56 m/s处，这是真实目标–15 m/s的模糊目标。同理，对于目标2和目标3，图6(b)显示了图5(d)在10 m处的水平切片，同样可以观察两个目标的虚假模糊位置。可以计算出虚假峰值比真实峰值低约7 dB。因此，通过在两个可能性中选择更强的峰值能有效地确定真实的目标速度。

图  6  仿真目标速度维剖面

Figure  6.  The velocity profiles of simulation target

下载: 全尺寸图片 幻灯片

上述结果展示了多普勒走动校正对距离多普勒图成像的影响，接下来可以看到多普勒补偿对于通道相位跳变的影响，如图7所示。在上文多普勒距离正确目标位置图中，首先选取距离为10 m、速度为7 m/s的无多普勒模糊目标。图7(a)展示了原始8个阵列相位以及FFT角度估计结果，可以看到原始相位由于TDM策略的实施会导致虚拟接收通道之间发现相位跳变，再通过式(12)进行点对点的补偿后，得到的FFT的频谱估计结果更为准确、旁瓣更低。观察到上述前两个仿真目标，位于同一个距离多普勒单元，但是此时目标的速度为16 m/s，发生了多普勒模糊，同样看到进行模糊校正和无模糊校正下的相位补偿结果对比，如图7(b)所示，与单个角度下的结果一致， 补偿后的FFT估计结果能够正确分辨两个目标。从图7(b)可以看到，如果不经过多普勒解模糊，错误地进行多个通道之间补偿，会导致更为剧烈的相位跳变，加剧空间谱估计误差，影响角度分辨能力。

图  7  多普勒通道补偿影响

Figure  7.  Doppler channel compensation effect

下载: 全尺寸图片 幻灯片

需要说明的是，本文并不是利用FFT进行空间谱估计，这里采用FFF进行分析是为了更清晰地展现本文提出的多普勒补偿的方法对于空间相位的影响，从而为后续进行高分辨率距离角度成像提供输入。

最后，为了初步验证IBMP算法的距离角度联合估计性能，不妨假设所有的目标位于同一多普勒单元内，仿真目标在距离角度空间中呈现X形分布，如图8(a)所示。在相同距离内，最小间隔角为12º，低于系统角度分辨率15.7º，因此利用传统2DFFT结合CFAR的方法，无法将其分开，如图8(b)所示。图8分别给出了GPSR方法、BCS方法、FBMP方法以及本文提出的IBMP方法的成像结果对比。可以明显看出，强噪声环境使得GPSR方法、BCS方法无法清晰地对场景目标进行重构，在成像的散射目标附近，存在虚警、漏警的情况，导致成像效果较差，FBMP方法虽然能恢复大部分点目标，但对低于角度分辨率的两个目标也无法完全重建。相比而言，IBMP方法得到的点目标清晰可分辨，且目标能量基本得到聚集，而且相较于FBMP方法具有更低的计算复杂度。

图  8  距离角度仿真点目标分布以及成像结果对比

Figure  8.  Range-angle simulation point target distribution and comparison of imaging results

下载: 全尺寸图片 幻灯片

5.3   性能分析

为说明所提方法的有效性，使用归一化均方误差(Normalized Mean Square Error, NMSE)量化真实目标场景信号${{\boldsymbol{x}}_{{t}}}$与稀疏恢复后信号${{{\hat {\boldsymbol{x}}}}_{{t}}}$的误差：

其中，T为重复实验次数，在本文中取T=200。设置SNR的变化范围从–5～15 dB，更小的NMSE意味着更好的重构性能。

针对距离角度模型，不同方法的成像结果的NMSE如图9所示。可以看出在较高SNR条件下，各种方法的RMSE均较小，但是随着SNR的下降，本文所提出的IBMP和FBMP方法依旧具有较低的RMSE，其他方法结果的NMSE升高明显较快，表明了提出的方法具有稳健的高性能成像性能。在SNR为15 dB时，IBMP和FBMP算法的NMSE约比其他算法低9 dB；在SNR为5 dB时，比2DFFT+CFAR方法约低8 dB，比其他所有算法低6 dB；在SNR为–5 dB时，约比2DFFT+CFAR方法算法低11 dB，并且比其他算法低6 dB。也可以看出FBMP由于其重复搜索在较低信噪比下比IBMP算法更好，但是时间开销也增加了D倍，两者的重构误差精度大致相同。

图  9  距离角度估计的NMSE对比

Figure  9.  NMSE comparison of range-angle estimation

下载: 全尺寸图片 幻灯片

5.4   实测分析

为验证所提方法对实测数据的成像性能，在公路场景利用TI AWR1642雷达进行前视成像，系统关键参数配置见表2，实验场景如图10(a)所示。

表  2  AWR1642雷达关键参数

Table  2.  Radar parameters of AWR1642

雷达参数	值
载频(GHz)	77
带宽(MHz)	1500
MIMO天线配置	2T × 4R
脉冲宽度(μs)	160
脉冲采样点数K	128
脉冲个数L	15.7
角度分辨率(°)	0.1
距离分辨率(m)	0.095

下载: 导出CSV 

| 显示表格

图  10  多普勒走动、多普勒模糊实测目标成像结果对比

Figure  10.  Comparison of imaging results of measured targets under Doppler walking and Doppler ambiguity

下载: 全尺寸图片 幻灯片

从常规的距离多普勒图可以看到，1～6号目标均存在距离走动现象，同时1,3,4,6号目标都发生了多普勒模糊，如图10(b)所示。在经过本文提出的基于多普勒补偿的校正处理后，从得到的多普勒距离图中可以就看到，相比较于原始估计结果，本文提出的方法得到真实目标频谱泄露更小，而且目标能量更为集中，信噪比得到了明显提升，如图10(c)所示(其中红色圆圈标注正确目标位置)。由于真实目标的多普勒距离走动得到了正确的补偿使得FFT处理增益最大化，而发生多普勒模糊的假目标由于不正确补偿，导致其走动现象加剧，散焦现象严重。从距离多普勒频谱来看，真实目标的能量要比其发生多普勒模糊的假目标的能量大，这也是本文提出的基于多普勒补偿方法能解决多普勒模糊的原因。

图11显示了目标1和目标2的速度剖面，即图10(c)分别在22.26 m和15.62 m处的水平切片。可以计算出虚假峰值比真实峰值低约7 dB，和仿真结果一致，顺利完成实测信号多普勒解模糊。

图  11  实测目标速度维剖面

Figure  11.  The velocity profiles of measured target

下载: 全尺寸图片 幻灯片

本文在图10(a)所示公路场景采集了400帧数据，搭载雷达设备的运动车辆速度为10 m/s，传统的距离多普勒角度三维FFT处理如图12(a)所示，可以看到车辆左侧的树木等静态目标由于相对运动的关系都产生了多普勒速度，根据图10(c)得知其存在多普勒走动和多普勒模糊。因此，常规方法得到的最终估计速度存在模糊。在经过本文提出的多普勒走动校正和多普勒解模糊措施后，正确地解析出目标的真实速度，通过这种速度对应关系，很容易通过后续滤波实现静态目标的去除。同时由于正确的通道补偿系数，使得目标角度的估计结果产生的偏差更小，虚假散射点更少，但由于较低的角度分辨率，使得目标处点云松散，能量不聚焦，如图12(b)所示。在利用本文提出的方法距离角度超分辨率算法处理后，得到的目标点云轨迹更精确，而且在目标散射点处更为集中，如图12(c)所示。

图  12  400帧实测数据点云成像结果对比(颜色表示速度)

Figure  12.  Comparison of point cloud imaging results from 400 frames of measured data (the color indicates the velocity of the target)

下载: 全尺寸图片 幻灯片

在上述相同系统参数下，实验场景如图13(a)所示，两人并列间隔1 m，以相同速度从两辆汽车中间同步径直走向雷达，使得目标在距离多普勒域无法分开，进一步验证所提方法对实测数据的前视距离角度成像性能。

图  13  实测成像结果对比(颜色表示速度)

Figure  13.  Measured imaging results comparison (the color indicates the velocity of the target)

下载: 全尺寸图片 幻灯片

实测数据实验中不同方法成像结果如图13所示。可以看到由于系统分辨率的限制，在较远处，距离为10 m 左右，相邻目标间隔8°，低于系统角度分辨率，导致传统方法在距离无法有效分离，因此两人轨迹无法分离，如图13(b)所示。经过本文提出的距离角度联合求解模型，4种方法都能有效改善的成像质量，如图13(c)—图13(f)所示。但是GPSR方法、BCS方法、FBMP方法得到的点云数据分散且轨迹不连续，本文提出的IBMP方法能准确地对目标的运动点云数据进行分辨，更有利于后续目标跟踪及更高级的信号处理。

6.   结语

本文提出了一种适用于高性能前视雷达图像的重建和成像方法，克服传统方法及相关CS方法在低信噪比情况下的大量虚警、漏警情况。本文提出的多普勒域补偿和点对点回波校正方法，可以实现距离走动和多普勒走动校正，并且同步完成多普勒解模糊，完成多域信号解耦合。为解决传统单维度估计方法分辨率不足的缺点，本文采用贝叶斯框架下基于最大后验准则的快速迭代更新方法，利用毫米波雷达TDM-MIMO前视成像下的多域联合回波估计模型，能够获得噪声鲁棒性较好的前视稀疏图像，能够准确重建出动静目标的空间分布特性。为了评估所提出方法的性能，进行了一系列仿真和实测结果，证明该方法结果准确性以及超分辨能力，并表明其在毫米波雷达前视成像估计性能上有明显的优势。