CN109917361A - 一种基于双基雷达的三维未知场景成像方法 - Google Patents

Abstract

该发明公开了一种基于双基雷达的三维未知场景成像方法，涉及雷达成像技术，特别涉及基于双基雷达的三维未知场景成像技术。首先基于双基雷达对未知区域进行多位置、多角度的扫描，使得测量的数据能够包含场景的全部信息。该扫描方案仅记录接收信号的衰减值，因此也适用于窄带雷达，能够有效降低***成本。然后通过本专利提出的TV‑MF‑ART稀疏重建算法，对测量向量分别执行代数迭代重建、正约束、总变分最小化约束和中值滤波操作后，最终得到了高精度的三维场景图像。因此，本发明具有扫描速率快、计算量小和重建精度高的优点。

Description

一种基于双基雷达的三维未知场景成像方法

技术领域

本发明涉及雷达成像技术，特别涉及基于双基雷达的三维未知场景成像技术。

背景技术

未知场景成像技术是通过发射天线发射特定频段的电磁波信号，接收天线接收场景的回波信号或者透射信号，并令收发天线移动完成多位置、多视角的探测，从而获取全部场景的回波并形成包含全部场景图像的未知场景全景图像的技术。该技术可以用于获取未知场景的完整图像，为目标精确定位、多径抑制提供了先验信息，在城市感知、灾难救援等领域发挥着重要作用。

传统的未知场景成像技术一般是基于MIMO雷达或SAR获取回波数据。然而为了得到场景物体准确的位置信息，上述的探测手段需要精确的相位信息和大的信号带宽进行保障，这会导致高的***复杂度和硬件成本。此外，该种方式利用的是回波信号，当场景环境复杂时，内部的多径信号会占主导地位，从而影响成像质量。因此，寻找一种新的成像方案迫在眉睫。

近年来，许多国内外研究机构尝试将医学基于X射线的计算机断层扫描成像理论(CT)应用于微波频段，也开展了相关理论与技术研究工作。Karanam等人根据电磁波在不同介质、不同厚度中衰减不同的这一基本理论，基于Wi-Fi的接收信号强度值(RSSI)和TVAL3算法实现了三维未知场景成像(C.R.Karanam and Y.Mostofi,“3D Through-Wall Imagingwith Unmanned Aerial Vehicles Using Wi-Fi,”in 2017 16th ACM/IEEEInternational Conference on Information Processing in Sensor Networks(IPSN),April 2017,pp.131–142)。Fhager A等人通过FDTD方法仿真验证了基于线性调频信号的微波时延谱测量，并采用收发分置同步扫描方式实现了二维未知物体的重构(Fhager A,Persson M.Comparison of two image reconstruction algorithms for microwavetomography[J].Radio Science,2005,40(3):1-15.)。然而，以上论述的方案虽然能够实现未知区域、物体的成像，但是收发分置同步扫描方案效率低下，该方案对于大的场景成像时则望而却步，并且成像质量较为粗糙。因此如何提高雷达获取探测区域回波的效率，如何得到高质量的未知场景图像仍需大量工作去验证。

发明内容

针对目前未知场景成像技术中收发分置同步扫描速率低下，成像效果不佳等缺陷，本发明提出一种基于窄带窄波束双基雷达的新型未知场景成像方案，在“感知”区域内部时能够有效提高效率，并得到高质量的场景图像。首先，设计了一种基于地基雷达和无人机的地空联合扫描方式来获取未知场景的数据，该方式具有执行效率高、环境适应性强的优点；接着，利用WKB近似方法建立起接收信号与未知场景之间的联系；最后，根据本发明提出的中值滤波下的总变分最小化代数迭代算法(TV-MF-ART)反演得到高精度场景成像结果。

本发明的技术方案为，一种基于双基雷达的三维未知场景成像方法，该方法包括：

步骤1：获取场景信息

对未知区域进行扫描，令一窄波束雷达在未知场景的外侧发射窄带连续波信号，在场景的另一侧，搭载了N_a个接收阵元的无人机按照规定好的路线移动，并接收穿透未知区域后的透射信号，在无人机移动期间，发射天线保证在测量时间内正对着相应的接收阵元，确保接收信号的衰减是直达波穿透场景所引起；

令n_p表示无人机的位置，n_v为视角；则n＝(n_a,n_p,n_v)表示为位于n_v视角下n_p位置处的第n_a个接收阵元的位置矢量；当完成N＝N_a×N_p×N_v次测量时，N_a表示接收阵元的总数量，N_p表示在固定视角下无人机移动的次数，N_v表示视角数，测量矩阵表示为：

P＝[p₁,p₂,...,p_N]^T (1)

其中，p_n表示在n位置处的接收功率。

测量值为信号功率的衰减量，它会随着不同电介质的尺寸、位置和介电常数的不同而不同，根据Wentzel-Kramers-Brillouin近似，测量的衰减量与传播路径的关系可以表示为：

σ_n∝exp(j2πf_c∫_T→Rα(r_n)dr) (2)

其中，f_c表示中心频率，σ_n为衰减值，α_n(r)为电磁波在位置r处的衰减速率，∫_T→R为收发阵元之间的线积分；将成像区域离散成M个单元，记作成像矩阵O；每个单元处的值取决于电磁波在该位置处的衰减速率，则成像矩阵表示为不同衰减速率的集合，即O＝[α(r₁),α(r₂),...,α(r_M)]，因此，测量向量与成像向量之间的联系为：

P＝A·O+b (3)

其中，A∈R^N×M为映射矩阵，表示成像向量与测量向量之间的映射；当j^th单元位于i^th次测量，则A(i,j)＝1，否则为0；b为测量误差，包括定位误差和环境噪声；至此，完成步骤1获取未知场景信息；

步骤2：对未知场景成像

执行步骤1得到包含未知场景信息的测量向量；因为测量的数据远少于成像单元数，并根据测量向量在稀疏采样数据下进行场景成像；

步骤2-1：代数迭代重建

将欠定方程组式(3)视为N个超平面，首先设定一个初始解，将初始解依次正交投影到每个超平面上，实现对初始解的迭代更新，逼近真实解；迭代方程为：

其中，λ为收敛因子，可以设置为1，O_q表示在第q次迭代时的成像向量，A_q,+表示为映射矩阵的第q行向量，N_ART为测量点数，初始的成像向量设置为全0，通过该步骤迭代得到一个初始的成像向量；

步骤2-2：对步骤2-1得到初始成像向量进行正约束；

对步骤2-1的初始成像向量中的负数置为0；表达式为：

O＝max(O,0) (5)

步骤2-3：总变分最小化步骤

将步骤2-2约束后的成像矩阵进行总变分最小化迭代，具体迭代方程为：

其中：α为人为设定的收敛因子，N_TV为迭代轮数，ΔO为步骤2-1与步骤2-2成像矩阵的差量，||·||₂为矩阵二范数，||·||_TV为图像的梯度，三维图像的梯度计算公式为：

其中：O_i,j,k为三维成像矩阵，i,j,k为其行列高索引，ρ为一小正数，防止分母为0；

步骤2-4：三维中值滤波步骤

对总变分最小化后的成像向量转换为三维矩阵执行三维中值滤波操作，具体表达式为：

其中为循环卷积，W为滑窗，窗大小视成像矩阵单元数而定；通过该步骤后得到最终的场景成像结果。

本发明提供了一种适用于收发分置雷达的场景成像方案，具有扫描速率快、成像精度高的优点。首先基于双基雷达对未知区域进行多位置、多角度的扫描，使得测量的数据能够包含场景的全部信息。该扫描方案仅记录接收信号的衰减值，因此也适用于窄带雷达，能够有效降低***成本。然后通过本专利提出的TV-MF-ART稀疏重建算法，对测量向量分别执行代数迭代重建、正约束、总变分最小化约束和中值滤波操作后，最终得到了高精度的三维场景图像。因此，本发明具有扫描速率快、计算量小和重建精度高的优点。

附图说明

图1为本发明双基雷达扫描示意图；

图2为重建算法框图；

图3为三维中值滤波示意图；

图4为待重建场景；

图5为双基雷达扫描场景的一个实例；

图6为传统ART算法重建的结果；

图7为本专利提出的重建算法的结果；

图8为RADON变换下测量曲线的比较。

具体实施方式

下面结合一个仿真给出本发明具体实施步骤。

建筑如图4所示，未知场景面积为3m×3m×1m。将该区域划分为0.01m×0.01m×0.01m的网格单元，也就是说该未知区域共有9×10⁶个单元。区域外设置一个地基雷达和无人机对该区域进行扫描，地基雷达发射的信号为2Ghz的正弦波信号，无人机搭载了100个等间隔阵元的线阵，用于接收透射场景后的信号并记录功率的衰减值。

步骤1：令无人机在场景另一侧移动，每个接收阵元移动间隔0.01m时记录一次信号衰减值，记录次数为300次。记录完毕后完成一个视角测量。该步骤记录的测量值包含直达路径的衰减量、穿透墙体时的多径分量以及环境噪声。为了抑制多径现象对本模型的影响，我们采用的窄波束天线，在增大发射功率的同时可以保证接收信号由直达波占主导地位。

步骤2：为了使得到的场景信息不会冗余，即测量时的张角要足够大。因此要完成多视角的测量。不同视角的选取会影响成像结果。本次仿真视角为4个，分别为0°，45°，90°和135°。完成本轮测量后共得到1200个测量点，记为P。记录每个测量点的位置，根据位置计算映射矩阵A，维度为1200×9000000。

步骤3：执行代数迭代重建步骤，初始成像向量设置为全0，迭代次数为1200次，根据式(4)迭代得到成像矩阵O^ART。

步骤4：执行正约束式(5)得到O^POS，并计算得到步骤3和步骤4成像矩阵的差量ΔO。

步骤5：根据式(6)和式(7)对总变分进行约束，迭代次数为20次，通过梯度下降方法降低图像的总能量。得到的结果记为O^MF。

步骤6：执行中值滤波操作式(8)，进一步平滑图像。

步骤7：若图像未收敛(收敛判定为ΔO的总值小于某个门限值)，继续执行步骤3-6。收敛后得到最终的未知场景图像。

本领域的普通技术人员将会意识到，这里所述的实施例是为了帮助读者理解本发明的原理，应被理解为本发明的保护范围并不局限于这样的特别陈述和实施例。对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的权利要求范围之内。

Claims (1)

1.一种基于双基雷达的三维未知场景成像方法，该方法包括：

步骤1：获取场景信息

对未知区域进行扫描，令一窄波束雷达在未知场景的外侧发射窄带连续波信号，在场景的另一侧，搭载了N_a个接收阵元的无人机按照规定好的路线移动，并接收穿透未知区域后的透射信号，在无人机移动期间，发射天线保证在测量时间内正对着相应的接收阵元，确保接收信号的衰减是直达波穿透场景所引起；

令n_p表示无人机的位置，n_v为视角；则n＝(n_a,n_p,n_v)表示为位于n_v视角下n_p位置处的第n_a个接收阵元的位置矢量；当完成N＝N_a×N_p×N_v次测量时，N_a表示接收阵元的总数量，N_p表示在固定视角下无人机移动的次数，N_v表示视角数，测量矩阵表示为：

P＝[p₁,p₂,...,p_N]^T (1)

其中，p_n表示在n位置处的接收功率。

测量值为信号功率的衰减量，它会随着不同电介质的尺寸、位置和介电常数的不同而不同，根据Wentzel-Kramers-Brillouin近似，测量的衰减量与传播路径的关系可以表示为：

σ_n∝exp(j2πf_c∫_T→Rα(r_n)dr) (2)

其中，f_c表示中心频率，σ_n为衰减值，α_n(r)为电磁波在位置r处的衰减速率，∫_T→R为收发阵元之间的线积分；将成像区域离散成M个单元，记作成像矩阵O；每个单元处的值取决于电磁波在该位置处的衰减速率，则成像矩阵表示为不同衰减速率的集合，即O＝[α(r₁),α(r₂),...,α(r_M)]，因此，测量向量与成像向量之间的联系为：

P＝A·O+b (3)

其中，A∈R^N×M为映射矩阵，表示成像向量与测量向量之间的映射；当j^th单元位于i^th次测量，则A(i,j)＝1，否则为0；b为测量误差，包括定位误差和环境噪声；至此，完成步骤1获取未知场景信息；

步骤2：对未知场景成像

执行步骤1得到包含未知场景信息的测量向量；因为测量的数据远少于成像单元数，并根据测量向量在稀疏采样数据下进行场景成像；

步骤2-1：代数迭代重建

将欠定方程组式(3)视为N个超平面，首先设定一个初始解，将初始解依次正交投影到每个超平面上，实现对初始解的迭代更新，逼近真实解；迭代方程为：

其中，λ为收敛因子，可以设置为1，O_q表示在第q次迭代时的成像向量，A_q,+表示为映射矩阵的第q行向量，N_ART为测量点数，初始的成像向量设置为全0，通过该步骤迭代得到一个初始的成像向量；

步骤2-2：对步骤2-1得到初始成像向量进行正约束；

对步骤2-1的初始成像向量中的负数置为0；表达式为：

O＝max(O,0) (5)

步骤2-3：总变分最小化步骤

将步骤2-2约束后的成像矩阵进行总变分最小化迭代，具体迭代方程为：

其中：α为人为设定的收敛因子，N_TV为迭代轮数，ΔO为步骤2-1与步骤2-2成像矩阵的差量，||·||₂为矩阵二范数，||·||_TV为图像的梯度，三维图像的梯度计算公式为：

其中：O_i,j,k为三维成像矩阵，i,j,k为其行列高索引，ρ为一小正数，防止分母为0；

步骤2-4：三维中值滤波步骤

对总变分最小化后的成像向量转换为三维矩阵执行三维中值滤波操作，具体表达式为：

其中为循环卷积，W为滑窗，窗大小视成像矩阵单元数而定；通过该步骤后得到最终的场景成像结果。