技术领域
[0001] 本
发明属于雷达
信号处理技术领域，具体涉及太赫兹孔径编码成像体制下的面目标成像成像方法。
背景技术
[0002] 相对于传统雷达采用的
微波频段，太赫兹波具有更高
频率和更短
波长，使得太赫兹雷达能够提供更大的绝对带宽，更高的成像
分辨率以及比光波更好的穿透能
力，同时，系统收发链路相对简单，易实现小型化。而孔径编码技术可通过产生多样性的照射模式使得成像系统生成更加丰富的目标回波信息，从而有望在相同成像孔径条件下获得超出传统成像系统衍射极限的分辨率。
[0003] 太赫兹频段孔径编码的成像过程是求解欠定方程的过程，该成像体制目前主要借鉴了微波关联成像的思想[1]。由于该体制是针对DARPA于2014年提出的先进扫描
成像雷达技术需求而设计的一种成像技术[2]，因此该体制下的面目标成像方法鲜有人研究。在微波关联成像体制下的某些方法可以对面目标进行成像，如最小二乘方法，相关法和正则化方法等。但是最小二乘方法和相关法受
时空二维随机
辐射场非相关性影响较大，且成像效果并不理想；同时正则化方法是一种普适性的
算法，其针对面目标的成像效果一般，且成像速度较慢[3-5]。在太赫兹孔径编码实时成像过程中，矩阵方程大多数时候是欠定方程，稀疏重构算法是一种有效求解欠定方程的方法。但是稀疏重构算法以目标稀疏为先验条件，而面目标的稀疏性极差，稀疏重构算法的成像效果受到很大限制。针对太赫兹孔径编码成像体制下的面目标成像方法，仍是一个有待研究的问题，针对这一问题，当前并没有有效的手段。
现有技术的主要缺点在于尚未有在太赫兹孔径编码成像体制下针对面目标成像的有效方法，且现有成像方法对面目标成像，成像效果不理想，成像速度慢。现有技术中的相关参考文献如下：
[0004] [1]D.Z.Li,Xiang L.,and Y.L.Qin,et al.,"Radar Coincidence Imaging:An Instantaneous Imaging Technique with Stochastic Signals,"IEEE Transactions on Geoscience and Remote Sensing,vol.52,no.4,pp.2261-2277,2014.
[0005] [2]H.B.Wallace,"Advanced Scanning Technology for Imaging Radars(ASTIR),"DARPA,Approved for Public Release,Distribution Unlimited,DISTAR Case#23357,Aug.2014.
[0006] [3]李东泽.雷达关联成像技术研究[D].国防科技大学,2014.
[0007] [4]查国锋.运动目标微波关联成像技术研究[D].国防科技大学,2016.[0008] [5]徐浩.基于空间谱理论和时空两维随机辐射场的雷达成像研究[D].中国科学技术大学,2011.
发明内容
[0009] 针对上述技术问题，本发明旨在提供太赫兹孔径编码成像体制下针对面目标的有效、快速成像方法。基于
小波变换的方法能够将非稀疏的面目标转换到稀疏的小波域。基于稀疏重构算法的成像方法对小波域的稀疏目标具有良好的适应性，计算效率高。因此，通过小波变换和稀疏重构算法的联合处理能实现对面目标的有效、快速成像。具体技术方案如下：
[0010] 一种基于太赫兹孔径编码成像体制的面目标成像方法，包括以下步骤：
[0011] (S1)构造小波稀疏基矩阵，包括：
[0012] (S11)对待成像面目标区域进行网格划分；设待成像面目标区域大小为h×v，h为目标区域
水平宽度，v为目标区域竖直高度；设单个网格单元的尺寸为hcell×vcell，hcell为网格单元水平宽度，vcell为网格单元竖直高度；根据单个网格单元尺寸和待成像面目标区域大小，待成像面目标区域由P×Q个网格单元组成，其中P＝h/hcell，Q＝v/vcell，设h和v分别为hcell和vcell的整数倍，待成像面目标区域的网格总数目为K＝P×Q；
[0013] (S12)选择小波母函数和采用matlab
软件wmaxlev函数确定最大分解层数；
[0014] (S13)采用matlab软件的wavedec函数求解小波稀疏基矩阵Ψ：
[0015] (S2)接收天线接收面目标回波信号，对待成像面目标稀疏重构成像，包括：
[0016] (S21)计算参考信号矩阵 及接收信号Sr(t)；
[0017]
[0018] Wr(·)表示调频信号的包络，fc为信号中心频率，为快时间，tm为慢时间，k为线性调频系数；τij(x,y)表示信号经过电控次反射面的第i行j列
位置处阵元的辐射场时延，(x,y)表示待成像面目标区域上任意一点的位置矢量， 表示根据编码时间施加至阵元(i,j)的随机
相位调制；
[0019] 则接收信号为：Sr(t)＝∫IβrSI(t)dr，
[0020] 其中 I表示待成像面目标区域，点(xp,yq)表示网格上的任一点；r表示待成像面目标区域网格点的矢量位置；t为全时间，全时间与快时间和慢时间的关系为 βr表示连续的目标散射系数；
[0021] (S22)根据参考信号矩阵和小波稀疏基矩阵构造字典矩阵；
[0022] 建立太赫兹频段孔径编码成像的数学模型：Sr＝S·β+w，
[0023] 其中S表示参考信号矩阵，Sr是接收信号矢量，β表示目标散射系数矢量，β是βr的离散化向量表示，β为未知量，w为噪声矢量；根据t1,t2,…,tN个时间点上的离散
采样，太赫兹频段孔径编码成像的数学模型表示为：
[0024]
[0025] r1,r2,…,rK表示待成像面目标区域网格划分后的K个网格点的矢量位置；
[0026] 令θ＝Ψ·β，将非稀疏解β换成稀疏解θ，构造稀疏字典矩阵A：A＝S·Ψ-1，则将太赫兹频段孔径编码成像的数学模型表示为：
[0027] Sr＝S·Ψ-1·θ+w＝A·θ+w
[0028] (S23)忽略噪声矢量w不计，采用稀疏重构算法，进行最优化求解得到[0029]
[0030] 其中λ为稀疏度控制参数，||·||2表示l2范数，||·||0表示l0范数；argmin函数表示求取使目标函数取最小值时的变量值。
[0031] (S24)通过小波逆变换求解面目标散射系数
[0032] 进一步地，所述步骤(S23)中l0范数替换成lp范数，p取值不为0。
[0033] 进一步地，所述步骤中(S23)中所述lp范数为l1范数。
[0034] 进一步地，所述小波母函数为Daubechies小波系或Coiflet小波系或Biorthogonal小波系。
[0035] 为更好的理解本发明方法内容，现将有关理论及推导内容介绍如下：
[0036] 小波变换是信号和
图像处理领域提取有效提取信号和图像特征的有效手段。针对面目标，小波变换能将面目标散射分布转换成表征面目标特征信息的稀疏解。
[0037] 上述建立的太赫兹孔径编码数学模型表达式为Sr＝S·β+w，忽略噪声影响，Sr＝S·β+w可表示成
压缩感知的数学模型y＝A·θ，其中θ为待求稀疏系数解。对于稀疏目标，Sr＝S·β+w中的β即为待求稀疏解，同y＝A·θ中的θ一致；而对于面目标，β为非稀疏解。
[0038] 本文将采用小波
正交基构造稀疏矩阵Ψ，通过小波变换Ψ·β，将非稀疏解β换成稀疏解θ，即有β＝Ψ-1·θ。针对待成像面目标，Sr＝S·β+w可重新建模表述成：
[0039] Sr＝S·β+w＝S·Ψ-1·θ+w＝A·θ+w (1)
[0040] 其中A＝S·Ψ-1，为字典矩阵。针对新建立的数学模型，通过稀疏重构算法进行求解，可以得到小波稀疏解θ。最后通过β＝Ψ-1·θ小波逆变换，求得面目标散射系数，即太赫兹频段孔径编码成像体制下的面目标成像结果。如图1为太赫兹频段孔径编码成像技术原理
框图。太赫兹频段孔径编码成像系统主要包括发射天线、电控次反射面、目标平面、接收天线和系统控制主机。发射天线发射信号，经过电控次反射面，到目标平面，被目标平面反射后，接收天线接收传至系统控制主机，经过系统控制主机解码处理，输出成像结果。
[0041] 本发明中构造小波稀疏基矩阵的过程就是
离散小波变换的过程。本发明采用离散小波变换(DWT)基于Mallat(法国学者，音译为
马拉特)算法实现。首先选取小波母函数，选取小波母函数的主要原则有：(1)正交小波；(2)具有最大可分解层数。太赫兹频段孔径编码成像的数学模型Sr＝S·β+w中，成像分辨率与参考信号矩阵S正交性成正比，新引入小波稀疏基矩阵不能破坏参考信号矩阵S正交性，因此选取小波母函数的原则之一为正交小波。小波母函数的多分辨分析能力与最大可分解层数成正比，且小波稀疏基矩阵稀疏表征功能与其多分辨率分析能力成正比，因此选取小波母函数的另一个原则为具有最大可分解层数。常用小波母函数包括：Daubechies小波系、Coiflet小波系、Biorthogonal小波系和Haar小波等。
[0042] 采用本发明获得的有益效果：本发明方法解决了太赫兹频段孔径编码体制下面目标成像难得问题，相比稀疏重构算法，该方法成像模型一致性更高，成像效果和成像速度均更优；相比于最小二乘方法、相关法和正则化方法，本方法对面目标针对性更强，成像
质量更佳。本发明采用基于小波变换的小波稀疏基构造，将非稀疏的面目标有效转换到稀疏的小波域表示。本发明针对面目标提出新的数学模型，基于本方法数学模型的稀疏重构算法能够极大程度改善太赫兹频段孔径编码体制下的面目标成像质量。
附图说明
[0043] 图1为太赫兹频段孔径编码成像技术原理框图；
[0044] 图2为本发明方法流程示意图；
[0045] 图3为太赫兹孔径编码成像系统发射与扫描光路示意图；
[0046] 图4为具体
实施例中手枪面目标小波变换前后图像对比图，图(a)为待成像目标；图(b)为小波变换后图像表示；
[0047] 图5为不同
迭代次数下OMP法和WT-OMP法(本发明)成像结果对比图；
[0048] 图6为不同迭代次数下OMP法和WT-OMP法NMSE对比图。
具体实施方式
[0049] 下面，结合附图和实施例对本发明作进一步说明。
[0050] 如图2所示，本文提出了一种基于小波变换和稀疏重构算法的太赫兹频段孔径编码成像体制下的面目标成像方法，利用小波变换将面目标转换到适于稀疏重构算法模型的稀疏小波域，然后采用稀疏重构算法进行面目标成像。具体技术方案包括两大步，第一步为小波稀疏基构造，第二步为基于稀疏重构算法的面目标成像。
[0051] (S1)小波稀疏基构造；
[0052] (S11)对待成像面目标区域进行网格划分。设置成像基本参数，主要包括设置太赫兹频段孔径编码成像系统，包括发射天线、电控次反射面、目标平面、接收天线和系统控制主机的位置。设待成像面目标区域大小为h×v，h为目标区域水平宽度，v为目标区域竖直高度；设单个网格单元的尺寸为hcell×vcell，hcell为网格单元水平宽度，vcell为网格单元竖直高度；根据单个网格单元尺寸和待成像面目标区域大小，待成像面目标区域由P×Q个网格单元组成，其中P＝h/hcell，Q＝v/vcell，设h和v分别为hcell和vcell的整数倍，待成像面目标区域的网格总数目为K＝P×Q；
[0053] (S12)选择小波母函数和采用matlab软件中的wmaxlev函数确定最大分解层数；本实施例中选择小波母函数Daubechies10。
[0054] (S13)根据wmaxlev函数确定的最大分解层数和确定的目标平面网格数目K，采用matlab软件中的wavedec函数求解小波稀疏基矩阵Ψ，Ψ∈RK×K
[0055] (S2)接收天线接收面目标回波信号，对待成像面目标稀疏重构成像，包括：
[0056] 如图3所示，为太赫兹孔径编码成像系统发射与扫描光路示意图，下面简述“快拍”成像系统中的信号采集过程。
[0057] 成像系统采用宽带线性调频连续信号，信号
波形为：
[0058]
[0059] 其中，Wr(·)表示调频信号的包络，fc为信号中心频率，为快时间，tm为慢时间，t为全时间，全时间与快时间和慢时间的关系为 k为线性调频系数。
[0060] 成像系统工作过程中，对电控次反射面的第(i,j)个阵元用均匀分布随机数移相形成一次反射场快照，亦即形成一次孔径编码，生成所需的空间随机辐射场，假定电控次反射面相当于一个M行G列的矩阵，则第(i,j)阵元即是矩阵中第i行第j列位置处的元素，i＝1,2,…,M，j＝1,2,…,G。根据采样时间(快时间) 进行采样，然后根据编码时间(慢时间)tm切换至下一次随机移相方式，此时时间维信号从确定的波形变为随机信号波形，辐射场从空间随机场变为时空二维随机场，即如图1中时间随机指辐射场分布在不同的时间相互独立，空间随机指辐射场分布在不同的空间位置相互独立。最终形成的瞬时辐射场波形 表示为：
[0061]
[0062] 其中，τij(x,y)表示经过第(i,j)个阵元的辐射场时延，由两者之间的空间距离决定。 表示根据编码时间施加至各阵元的随机相位调制。
[0063] 在待成像目标区域I内，(x,y)为其中任意一点的位置矢量。待成像目标区域划分为P×Q的网格，在任一网格点(xp,yq)处的参考信号根据形成的瞬时辐射场波形即可表示为探测信号用SI(t)表示：
[0064]
[0065] p＝1,2,…,P；q＝1,2,…,Q。
[0066] 接收信号为成像区域内探测参考信号的
叠加：
[0067] Sr(t)＝∫IβrSI(t)dr (5)
[0068] βr是连续的目标散射系数表示，β是βr的离散化向量表示，
[0069] 基于式(5)所示的模型，获得孔径编码成像的离散“快拍”采集过程和参数
化成像模型：
[0070] Sr＝S·β+w
[0071] 其中，S是参考信号矩阵，Sr是接收信号矢量，β是未知的目标散射系数矢量，w是噪声系数向量；在t1,t2,…,tN时间点进行离散采样，上式表示为：
[0072]
[0073] r1,r2,…,rK表示待成像面目标区域网格划分后的所有网格点的矢量位置；
[0074] 对于非稀疏面目标，压缩感知/稀疏重构类算法在有限次迭代条件下成像效果受到限制。在此
基础上，根据第一步确定的小波稀疏基矩阵：Ψ∈RK×K。通过小波变换Ψ·β将非稀疏解β换成稀疏解θ，明显有β＝Ψ-1·θ。综上所述，针对面目标，Sr＝S·β+w可重新建模表述成：
[0075] Sr＝S·β+w＝S·Ψ-1·θ+w＝A·θ+w (7)
[0076] 现在已经有严格的数学证明，若θ稀疏，可将θ的求解过程建立成以下优化模型：
[0077]
[0078] 式(8)为稀疏重构算法的基本原理模型，其中λ为稀疏度控制参数。
[0079] 进一步，实施例中小波域稀疏重构成像过程如下：
[0080] 因为式(8)最小l0范数问题是一个NP-hard问题，因而无法求解，通过近似处理能够很好的解决该问题，最常见的一种方法就是用lp范数代替l0范数
[0081]
[0082] 当p为1时，为基追踪方法(Basis Pursuit，BP)法。此外，稀疏重构类方法还包括贪婪类方法、稀疏贝叶斯学习方法、非凸优化方法等。其中贪婪类算法中的正交
匹配追踪算法在太赫兹频段孔径编码成像中较为常用。选取稀疏重构算法求出θ的估计值
[0083] 对上一步中求出的小波稀疏基求逆，通过小波逆变换，求出非稀疏的面目标散射向量β：β＝Ψ-1·θ，求出β即完成了面目标成像。
[0084] 实施例中的目标平面基于太赫兹孔径编码成像体制，且其网格大小受成像系统分辨率限制。小波母函数具有正交性和多分辨性能，其正交性保证了字典矩阵的行和列的非相关性，多分辨性能保证了非稀疏面目标在小波域的稀疏表示。将参考信号矩阵和小波稀疏基矩阵相乘构造新的参考信号矩阵，而且模型中的待求解参数具有稀疏性，可采用稀疏重构算法进行求解。
[0085] 下面对本发明方法进行仿真实验，以手枪面目标作为对比研究的对象，本文方法中的稀疏重构算法采用正交匹配追踪算法，与单纯的正交匹配追踪算法之间的对比，说明了本方法的有效性。
[0086] 仿真实验场景设置如图3所示，具体参数如下：采用载频340GHz，带宽24GHz的馈源发射天线，
采样频率5MHz。电控次反射面阵元规模为30x30，对斜距为1.5m处的目标平面(待成像面目标区域)进行成像。目标平面(待成像面目标区域)处网格划分为32x32，网格单元尺寸设置为5mm，小波母函数采用Daubechies10小波。坐标原点到目标平面中心的连线为一方向矢量R′，θ′和 分别是该方向矢量的
俯仰角和方位角，俯仰角θ′为方向矢量和y轴正向的夹角，方位角 为方向矢量在xoz平面投影与z轴正方向的夹角。Sh和Sv分别为成像扫描过程中雷达波束在目标明面上的扫描光斑大小。B′(eh,ev,a+b)点是目标平面定点坐标，A′(x0,y0,a)为发射天线位置点坐标。具体参数设置如表1。
[0087] 表1仿真参数设置
[0088]
[0089]
[0090] 待成像手枪面目标，稀疏性极差，经小波变换后的小波域表示，其稀疏性有极大改善。图中WT-OMP法表示小波变换-正交匹配追踪方法，即本发明所提方法，OMP法表示正交匹配追踪算法。图5为对比不同迭代次数下WT-OMP法和OMP法之间的成像结果，其中图(a)为OMP法迭代250次，图(b)为WT-OMP法迭代250次；图(c)为OMP法迭代500次，图(d)为WT-OMP法迭代500次；图(e)为OMP法迭代750次，图(f)为WT-OMP法迭代750次；图(g)为OMP法迭代1000次，图(h)为WT-OMP法迭代1000次；从各图的显示效果来看，证明了本文所提方法的有效性。
[0091] 如图6所示，为量化成像结果，采用归一化均方误差(NMSE)衡量成像误差，对比了不同迭代次数下OMP和WT-OMP方法的NMSE，图(a)为OMP法和WT-OMP法的NMSE同尺度对比显示，图(b)为WT-OMP法的NMSE，进一步证明了该方法的有效性。
[0092] 本发明中的太赫兹孔径编码成像借鉴微波关联成像思想，通过电控次反射面阵列天线的实时编码调制形成时-空二维随机分布的辐射场，最后利用探测回波和辐射场参考矩阵通过矩阵方程求解的方式实现高分辨、前视和凝视成像，弥补了合成孔径高分辨成像依赖目标运动的不足。
[0093] 以上仅是实施例仅用于说明本发明的效果，本发明的保护范围并不仅局限于上述实施例，凡属于本发明思路下的技术方案均属于本发明的保护范围。应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理前提下的若干改进和润饰，应视为本发明的保护范围。
