一种高速超分辨率MIMO阵列成像方法与流程

本发明涉及合成孔径成像算法，尤其涉及一种高速超分辨率mimo阵列成像方法。

背景技术：

多输入多输出(multiinputmultioutput，mimo)阵列是一种具有广阔发展前景的阵列成像方案。利用分集技术，大大减少天线和通道数量，并保持系统的分辨率性能。合成孔径(syntheticaperture，sa)成像算法可以应用于mimo阵列成像系统，该技术的特点是以脉冲压缩技术为基础，相对于传统的电磁成像方法，计算量更小，内存需求更低，是实时成像系统的首选方法。sa采用一种born近似的电场信号模型，为保证图像的精度，不考虑任何远场近似。mimo阵列成像最经典的两种sa成像方法是反投影成像算法(backprojectionalgorithm,bpa)和距离偏移成像算法(rangemigrationalgorithm，rma)。两种方法各有特点。bpa的重要优势在于对天线的空间采样方式没有要求，允许非均匀采样。这一点对于非均匀排布的阵列成像和一些复杂的双基成像雷达具有重要意义。另一方面，bpa的抗噪性能也更强。它的主要缺点是计算量非常大，难以应用在实时系统中。为克服这一问题，有人提出了快速bpa(fastbpa，fbpa)和快速分解bpa(fastfactorizedbpa，ffbpa)，这一类算法是以牺牲图像精度为代价换取提高速度，因此需要均衡考虑。相对来说，rma得益于快速傅立叶变换(fft)，其速度优势是与生俱来的，虽然精度比bpa稍低，但是整体性能仍然可接受。mimo-rma的最大限制是发射和接收天线的空间排布必须是均匀的。

基于rma的超分辨率成像算的主要原理是，将rma的最后一步快速傅里叶逆变换(ifft)用功率谱分析方法代替，因为某些功率谱分析方法得到的谱密度函数比信号的频谱具有更高的频谱分辨率，如capon谱分析和子空间类谱分析方法。但是谱分析方法的计算量和内存资源占用都很高，因此这一类方法多数局限于二维rma成像算法中，且难以应用于实时成像系统中。

随着压缩感知(compressivesensing,cs)技术的快速发展，雷达成像的一个重要的技术应用就是压缩感知成像。在成像应用中，cs技术一般与其他正则化方法相结合，有利于提高图像质量，达到提高分辨率的效果。除了在3d大场景应用中内存需求过高的问题，计算效率也是这一类方法难以用于实时系统的巨大障碍。

由以上分析可见，目前的提出的超分辨率成像算法普遍特点是，计算量大，内存占用率高。对于实时成像系统，尚不存在可用的超分辨率mimo阵列成像算法。

技术实现要素：

针对以上所述超分辨成像方法的技术缺陷，本发明提出一种高速超分辨率mimo阵列成像方法，相干因子技术的提出是用于压制稀疏阵列雷达系统的旁瓣和栅瓣，与bpa相结合，应用于空间域成像算法中。本发明将空间域相干因子扩展到波数域，并与mimo-rma相结合，具有实时成像能力；本方法在压低图像旁瓣和基底噪声的同时，进一步提高了分辨率。

本发明的目的就是为了解决上述问题，提供：

一种高速超分辨率mimo阵列成像方法，包括如下步骤：

步骤(1)：对mimo阵列进行数据采样，得到回波数据；

步骤(2)：利用mimo-rma对回波数据进行处理，得到反射率函数，即待重构的初始三维图像；所述步骤(2)在处理的过程中得到频域响应信号；

步骤(3)：利用mimo-rma对步骤(2)求得的频域响应信号加速求取相干因子；

步骤(4)：利用步骤(3)得到的相干因子校正步骤(2)得到的反射率函数的初始三维图像。

所述步骤(1)的步骤为：在mimo阵列中，发射天线阵元与接收天线阵元是分开的，发射天线阵元将电磁信号辐射到空间中，遇到被测目标体被散射点反射，接收天线阵元记录回波数据。

所述步骤(2)的步骤为：

步骤(2.1)：根据驻定相位原理，对回波信号在垂直距离向做四维傅立叶变换，得到垂直距离向波数域信号；

步骤(2.2)：对垂直距离向波数域信号进行高通滤波，实现信号脉冲补偿，得到频域响应信号；

步骤(2.3)：场景中心补偿：将需要探测的成像区域事先划定，对频域响应信号进行成像区域的场景中心补偿；

步骤(2.4)：对场景中心补偿后的频域响应信号进行降维度变换和距离向插值，得到在垂直距离向和距离向均匀排布的三维数据，再进行三维傅立叶逆变换得到反射率函数。

所述步骤(3)的步骤为：

步骤(3.1)：通过数据重排化简反射率函数的非相干功率，得到重排后的回波数据；

步骤(3.2)：在垂直距离向对重排后的回波数据做四维傅立叶变换，得到波数域数据；

步骤(3.3)：与步骤(2.3)类似，对步骤(3.2)输出的波数域数据做场景中心补偿；

步骤(3.4)：在重排后的新维度空间中，对场景中心补偿后的信号进行降维度变换和距离向插值，得到的三维数据在各维度上均匀分布，再进行三维傅立叶逆变换后得到反射率函数的非相关功率；

步骤(3.5)：由步骤(2)得到的反射率函数和步骤(3.4)得到的非相干功率计算相干因子。

所述步骤(1)的步骤为：在mimo阵列中，发射天线阵元与接收天线阵元是分开的，假设所有天线阵元位于同一个平面，在笛卡尔坐标系中，天线平面与x-y平面相互平行；坐标为(xt,yt,zt)的发射阵元将电磁信号辐射到空间中，下标t代表发射阵元，遇到被测目标体被某一散射点(x,y,z)反射，位于(xr,yr,zr)的接收阵元记录回波数据为，下标r代表接收阵元；

其中，rt和rr分别是从散射点到发射阵元和接收阵元的单程距离；d(x,y,z)是目标域，σ(x,y,z)是成像目标的反射率函数；p(k)是探测的脉冲波形；f是探测信号的频率，k是f对应的波数，c是电磁波在自由空间中的传播速度；s(xt,yt,zt,xr,yr,zr；k)表示对应发射和接收阵元的回波数据。

所述步骤(2.1)的步骤为：

根据驻定相位原理，对回波信号s(xt,yt,zt,xr,yr,zr；k)的xt、yt、xr和yr做fft，取傅立叶变换对

xt→kxt

yt→kyt

xr→kxr

yr→kyr

其中→代表傅里叶变换中变量关系，得到变换结果为

其中，kzt是发射阵列在z方向上的波数，kzr是接收阵列在z方向上的波数，s(kxt,kyt,zt,kxr,kyr,zr,k)是垂直距离向波数域信号。

所述步骤(2.2)的步骤为：

假设mimo平面阵列与x-y平面重合，即zt＝zr＝0；做信号的高通滤波，补偿结果为

其中sc(kxt,kyt,kxr,kyr,k)是频域响应信号。

所述步骤(2.2)虽然合成孔径以相位聚焦为主，但是自由空间的衰减补偿有利于提高图像精度；信号调制脉冲波形体现了探测信号在各频率分量上的信噪比，大带宽信号有利于深度向分辨率的提高。信号脉冲补偿，从性能而言是做信号的高通滤波。

所述步骤(2.3)的步骤为：将需要探测的成像区域事先划定，设在z方向上场景中心与天线平面的距离为hc，则场景中心的补偿结果scs(kxt,kyt,kxr,kyr,k)是

scs(kxt,kyt,kxr,kyr,k)＝sc(kxt,kyt,kxr,kyr,k)exp[i(kzt+kzr)hc](6)

其中，hc表示在z方向上场景中心与天线平面的距离。

由于傅立叶变换的周期性特点，z方向的成像范围并不一定包含目标的真实距离位置，有可能周期性地反褶进入成像区，从而难以确定真实距离。因此更合理的办法是，将需要探测的成像区域事先划定，设在z方向上场景中心与天线平面的距离为hc，对场景中心进行补偿。

所述步骤(2.4)的步骤为：

根据脉冲压缩原理，由公式(2)变形得到

kx＝kxt+kxr(8)

ky＝kyt+kyr(9)

kz＝kzt+kzr(10)

公式(7)中的积分核形式满足傅立叶变换形式，但是积分限需要先做变换才能满足fft的使用条件。

公式(8)和(9)对应的变换称为降维度操作，降维度操作需要通过补零或插值保证kx和ky在变换后保持等间隔采样。

公式(10)的变换称为stolt插值，由于kz与kxt,kyt,kxr,kyr和k是非线性关系，因此需要做插值。

若降维度是通过补零完成，则stolt插值能够通过一维插值实现；

若降维度需要通过插值实现，则最佳方法是与stolt插值合并成三维插值实现；

公式(7)重新写成公式(11)形式

因为公式(11)的形式满足傅立叶变换格式，所以求取σ(x,y,z)通过ifft实现。

步骤(3.1)的步骤为：

反射率函数的非相干功率定义为

其中，a(xr,yr)表示接收阵元孔径，a(xt,yt)表示发射阵元孔径；

公式(12)重新整理成公式(13)：

公式(13)中，s(k)是回波信号s(xt,yt,zt,xr,yr,zr；k)的简写。*表示复数的共轭操作；波数k1和k2的取值范围与k相同，都是[kmin,kmax]。

通过数据重排将公式(13)化成能够使用快速运算的形式，数据重排的变量代换关系是

ks＝k1-k2(14)

其中，ks是重排后的波数。在物理意义上，回波信号s(k)应满足时间反对称要求，因此s(k1)s^*(k2)只与k1-k2相关，于是公式(13)中的k1和k2双重积分必定能够变换成ks的单重积分。

具体而言，波数域数据重排的过程是这样实现的：

假设采样的回波数据包含n个波数域等间隔采样数据，则ki(i＝1,2)对应的n个离散的波数为ki0,ki1,…,ki(n-1)。

以s(k1)的数据为行，s^*(k2)的数据为列，构成s(k1)s^*(k2)的二维矩阵。

所述二维矩阵中的每个对角线元素即对应一个重排后的波数(k1-k2)，即ks；

将各对角线的数据相加得到重排后的新数据e(ks)。

所述步骤(3.2)的步骤为：

经过回波数据重排后，公式(9)写成如下形式：

根据反投影成像算法bpa原理，σ(x,y,z)通过如下公式计算

对比公式(15)和(16)，发现其形式完全相同，因此采用rma快速计算ip(x,y,z)。

根据匹配滤波原理，将公式(15)变换成回波形式，并在两端同时做四维fft得到

其中，kzt′和kzr′分别是发射和接收阵列的z方向波数。

所述步骤(3.3)过程与步骤(2.3)相似，场景中心补偿结果为

ec(kxt,kyt,kxr,kyr；ks)＝e(kxt,kyt,kxr,kyr；ks)exp[i(k′zt+kzr′)hc](20)

其中，hc是mimo阵列平面与场景中心的z方向距离。

所述步骤(3.4)的步骤为：

由公式(17)变型得到

kx＝kxt+kxr(22)

ky＝kyt+kyr(23)

kz′＝kzt′+kzr′(24)

降维度公式(22)和(23)能够通过补零或插值实现。

经过stolt插值(24)并重新整理公式(21)，得到

公式(25)可以采用三维ifft快速实现。

所述步骤(3.5)的步骤为：

相干因子定义为

公式(26)中的参数α和β用于调节相干因子的性能，α∈[0.5,4.5],β∈[1/2,1/500]。

α和β的取值以经验为主，但是大量的仿真和实验发现，更为合理的取值方式是α∈[0.5,4.5],β∈[1/2,1/500]。不同系统中参数的取值不同，但是一旦确定，相干因子的性能会非常稳定，不需要进一步修改。

所述步骤(4)由步骤(3)得到的相干因子，计算最终图像为

σm(x,y,z)＝σ(x,y,z)cf(x,y,z)(27)。

本发明的有益效果：

1通过重排回波数据化简反射率非相干功率，利用mimo-rma加速计算相干因子，大幅度提高相干因子的计算速度，保证成像速度的实时性；

2利用相干因子校正mimo-rma重构的三维像，能够有效压制旁瓣和基底噪声，同时获得更高的分辨率；

3针对相干因子对弱散射点的压制缺陷，通过控制参数α和β平衡超分辨率性能。

附图说明

图1是mimo阵列成像原理图；

图2是基于波数域相干因子的高速超分辨率mimo阵列成像算法原理框图；

图3是十字型mimo阵列几何排布示意图；

图4是x轴方向的点扩展函数；

图5是y轴方向的点扩展函数；

图6是用于散射场计算的八叶金属扇形目标的示意图；

图7是mimo-rma电磁仿真成像结果；

图8是相干因子mimo-rma电磁仿真成像结果；

图中：1、mimo阵列的发射天线阵元，2、mimo阵列的接收天线阵元，3、目标域。

具体实施方式

下面结合附图与实施例对本发明作进一步说明。

mimo阵列成像系统的工作方式如图1所示。mimo阵列的发射天线阵元(1)向自由空间辐射电磁信号，mimo阵列的接收天线阵元(2)接收并记录回波信号，所有天线阵元位于同一个平面内。发射和接收天线的二维孔径分别表示为a(xt,yt)和a(xr,yr)。mimo阵列辐射区域应覆盖目标域(3)，任何待成像的目标都应位于目标域d(x,y,z)的范围内。

实施例一：点扩展函数

本实施例分别计算了mimo-rma和相干因子mimo-rma两种成像算法的点扩展函数。采用的mimo阵列排布方式如图3所示。阵列中包含82个发射阵元和82个接收阵元，相邻阵元的间距为6.2mm，均匀地排布于x轴和y轴上，形成发射和接收孔径的尺寸均为500mm。用于计算的其他参数在表1中列出。

表1点扩展函数计算采用的主要参数

点扩展函数即是采用理想的点目标，由sa的信号模型产生回波信号，并利用图2所示的成像算法得到的点目标的图像。

分别采用mimo-rma和相干因子mimo-rma两种方法计算得到的点扩展函数结果如图4和图5所示。其中图4是x方向上的点扩展函数，图5是y方向上的点扩展函数，各图中十字曲线是mimo-rma方法的处理结果，点状曲线是相干因子mimo-rma方法的处理结果。由于图4和图5得到的曲线结果完全相同，只采用图4结果进行进一步分析，结果如表2所示。

表2中的几种参量的计算方法是：

半峰全宽(fwhm)：靠近0db处的两个负峰值被认为是主瓣的两端，主瓣幅值下降到3db处对应的宽度是fwhm；

背景噪声电平：主瓣采样点之外，所有采样点的幅值平均值；

最大旁瓣幅值：主瓣外区域中采样点的最大幅值。

表2点扩展函数的分析结果

对表2的结果进行总结。采用相干因子后，相对于mimo-rma：

1)点扩展函数主瓣的fmhw减少了1.496mm，相当于分辨率提高了大于20％；

2)背景噪声的平均电平下降了21.265db。该参量对于稀疏阵列成像格外有意义；

3)最大旁瓣幅值下降了11.000db。

实施例二：电磁仿真成像

该实施例主要用于验证相干因子mimo-rma对连续目标体的成像性能。仿真中采用的主要参数如表1所示，唯一的不同之处是α和β参数取值不同，即α＝1.5，β＝0.00375。

仿真包括两步：

步骤一：利用矩量法产生目标的散射场作为回波数据；

成像的目标体如图6所示，是八叶扇形理想电导体的二维模型。电磁仿真采用赫兹电偶极子源，源的位置相当于发射阵元的位置；点探测器相当于接收阵元。发射源按空间顺序激发，接收阵元平行接收，且只考虑同方向极化电场。

步骤二：分别采用mimo-rma和相干因子mimo-rma对回波数据进行处理，完成图像重构。成像结果如图7和图8所示。图6是mimo-rma的成像结果，图7是相干因子mimo-rma的成像结果。图像采用灰度图，动态范围是30db。图像的背景噪声在-30db以下，因此在图中被滤去。图7中可以明显看到在扇叶附近存在较多的旁瓣，使得边缘处变得模糊。对比于图8，采用相干因子后，扇叶附近的旁瓣被极大地压制，边缘处非常清晰，另外可以发现图8中的扇叶形状更为尖锐，目标中心处轮廓更清晰，这些都是相干因子mimo-rma的分辨率被提高的原因。

上述虽然结合附图对本发明的具体实施方式进行了描述，但并非对本发明保护范围的限制，所属领域技术人员应该明白，在本发明的技术方案的基础上，本领域技术人员不需要付出创造性劳动即可做出的各种修改或变形仍在本发明的保护范围以内。