Sparse Targets Angular Super-resolution Reconstruction Method under Unknown Antenna Pattern Errors for Scanning Radar
-
摘要: 扫描雷达角超分辨技术是基于目标与天线方向图的关系,采用解卷积方法获取超越实波束的角分辨能力。目前的角超分辨方法大都是基于理想的无畸变天线方向图,未考虑实际过程中方向图的变化。然而,由于雷达天线罩、天线测量误差与平台非理想运动等因素的影响,天线方向图在实际中往往存在未知的误差,会导致目标分辨能力下降,甚至产生虚假目标。针对此问题,该文提出一种机载扫描雷达未知天线方向图误差下的角超分辨成像方法。首先,基于总体最小二乘(TLS)准则,该文考虑了方向图误差矩阵的影响,导出了相应的目标函数;其次,基于交替迭代的求解思路,利用迭代重加权优化方法实现了目标函数求解;最后,针对算法超参数选取,引入了一种自适应参数选取方法。仿真与实测结果表明,该文方法能实现未知天线误差条件下的超分辨重建,进一步提升了超分辨算法的稳健性。Abstract: Scanning radar angular super-resolution technology is based on the relationship between the target and antenna pattern, and a deconvolution method is used to obtain angular resolution capabilities beyond the real beam. Most current angular super-resolution methods are based on ideal distortion-free antenna patterns and do not consider pattern changes in the actual process due to the influence of factors such as radar radome, antenna measurement errors, and non-ideal platform motion. In practice, an antenna pattern often has unknown errors, which can result in reduced target resolution and even false target generation. To address this problem, this paper proposes an angular super-resolution imaging method for airborne radar with unknown antenna errors. First, based on the Total Least Square (TLS) criterion, this paper considers the effect of the pattern error matrix and derive the corresponding objective function. Second, this paper employs the iterative reweighted optimization method to solve the objective function by adopting an alternative iteration solution idea. Finally, an adaptive parameter update method is introduced for algorithm hyperparameter selection. The simulation and experimental results demonstrate that the proposed method can achieve super-resolution reconstruction even in the presence of unknown antenna errors, promoting the robustness of the super-resolution algorithm.
-
1. 引言
扫描雷达具有快速搜索和广域成像的优点,且能兼容任意的成像构型,在海面舰船搜索、导航避障、精确制导等领域均有着重要的研究价值[1−4]。在距离向上,扫描雷达通过发射大带宽信号并采用匹配滤波技术实现距离高分辨成像[5,6];然而,受雷达实际孔径尺寸的限制,方位角度分辨率远低于距离分辨率,极易造成远距离多目标信息的误判。因此,如何突破天线物理孔径的限制,提升扫描雷达角分辨能力,已成为实际应用中亟待解决的难点问题。
扫描雷达方位向信号可以建模为天线方向图与目标散射系数的卷积,利用卷积反演方法进行信号处理,可以得到超越实孔径限制的角分辨力。针对这一问题,近年来国内外学者对扫描雷达超分辨成像方法展开了广泛的研究。目前,国内外研究的各种超分辨成像方法可大致归纳为逆滤波方法[7−9]、贝叶斯方法[10,11]、谱估计方法[12,13]、正则化方法[14,15]等4大类。
逆滤波方法是通过在频域合理设计逆滤波器,实现超分辨成像。经典的逆滤波方法是维纳滤波算法[16],它可以避免对逆滤波器的零点的相除运算,从而解决直接逆滤波的病态问题。但是该方法需要目标分布的先验统计信息。Richards等人[17]在非相干雷达成像应用中提出了卷积核频域因式分解逆滤波方法,可在高信噪比条件下提高角分辨率约3倍。逆滤波方法通常对环境信噪比要求较高,在实际中难以实现。
基于贝叶斯的超分辨方法是通过利用目标的先验分布以及噪声的统计特性,在最大似然或最大后验准则下将超分辨问题转化为最优参数估计问题。比如Li等人[18]将光学超分辨成像领域中Lucy-Richardson方法应用到扫描雷达前视超分辨成像中,与传统的实波束相比,其角分辨提高了约4倍。随后,Tan等人[19]提出一种基于马尔可夫随机场(Markov Random Field, MRF)模型的贝叶斯超分辨方法,利用MRF模型的N阶邻域系统很好地描述和利用成像场景的二维空间结构特征,实现了更好的成像场景伪影抑制和轮廓恢复。Chen, Li等人[20,21]利用广义高斯分布提出了一种稀疏贝叶斯超分辨方法,在实现高分辨率成像的同时可有效抑制噪声。
阵列信号处理方法是通过利用天线方向图与阵列信号空间导引矢量之间的映射关系实现超分辨成像。美国麻省理工学院Capon[22]提出一种Capon频谱估计方法,该方法基于形状变化的波束窗口,实现了信号频谱高分辨估计。但Capon, MUSIC (MUltiple SIgnal Classification)和ESPRIT (Estimation of Signal Parameters using Rotational Invariance Techniques)等谱估计方法往往需要大量快拍数的积累,在单快拍数据下成像结果较差。为了解决大快拍问题,Zhang等人[23]提出了一种迭代自适应方法(Iterative Adaptive Approach, IAA)应用于扫描雷达超分辨成像中,实现了单快拍条件下的多目标分辨。此外,Li等人[24]将非相干信号模型应用于前视雷达角超分辨成像中,通过利用辐射方向图代替天线功率方向图,并采用一种改进的频域迭代自适应方法(IAA)进行反卷积求解,实现了快速前视角超分辨成像。
正则化方法是通过在目标信号模型中引入正则化项,限制解空间范围,获取稳定的解。Migliaccio等人[25]提出了Tikhonov正则化方法,由于L2范数的平滑性使得分辨率提升受限。Huo等人[26]提出一种动态平衡正则化方法,其实质是将L2范数与全变差范数相结合,通过设置动态正则化参数,调节L2范数和全变差范数的比重,实现了成像场景轮廓恢复。此外,Tang等人[27]提出了一种结合低秩和稀疏先验约束模型的多通道阵列雷达前视成像方法,针对低秩联合稀疏双先验约束优化问题,采用交替方向乘子法(Alternating Direction Method of Multipliers, ADMM)框架下的增广拉格朗日乘子(Augmented Lagrangian Method, ALM)重构方法,实现了目标场景的高效重建。
然而,现有的超分辨方法,通常假定成像过程中其卷积核不会发生畸变,即没有考虑天线方向图的误差。然而,在实际应用中由于雷达天线罩、天线方向图测量误差以及平台运动等因素,容易造成天线方向图误差[28,29]。当天线方向图存在未知误差时,会引起模型失配,卷积反演过程的病态性进一步扩大,导致目标分辨能力下降,甚至产生虚假目标,不利于方位向高分辨成像。
针对该问题,本文提出一种未知天线方向图误差条件下的稀疏目标角超分辨重建方法。首先,基于总体最小二乘法准则,考虑了方向图误差的影响与目标稀疏先验,构建目标函数;其次,基于变量交替迭代的求解思路,借助迭代重加权(Iterative Reweighted Norm, IRN)优化方法实现了目标函数的最优求解;更进一步,针对算法参数自适应选取难题,本文基于范数比准则导出了自适应参数选取方法;最后,通过点仿真、面仿真以及实测数据证明了本文方法的有效性。
2. 信号模型与问题形成
2.1 扫描雷达信号模型
扫描雷达以固定的频率发射线性调频脉冲信号,同时天线以机械转动扫描方式对目标区域进行扫描, 图1表示机载扫描雷达的几何模型。初始时刻载机位于点A,以速度V沿某一水平方向平飞,飞行高度为H,天线的扫描角速度为ω。假设扫描场景中有一感兴趣的点目标P,其相对于雷达的初始斜距为R0,方位角为α0,俯仰角为β。飞机在经过t时间飞行后,到达B点,此时飞机与目标P的斜距为
R(t)=√R20+V2t2−2R0Vtcosθ0 (1) 其中,θ0为点目标与雷达的初始空间方位角。通过泰勒公式近似展开,式(1)可以简化为
R(t)≈R0−Vt (2) 扫描雷达发射的线性调频脉冲信号可以表示为
f(τ)=rect(τT)exp[j2π(fcτ+12Kτ2)] (3) 其中,rect(⋅)为矩形窗,T为脉冲宽度,τ为距离向快时间,fc为信号载频,K为调频率。点目标P的回波经过下变频等预处理后表示为
yP(τ,t)=xPh(t)rect(τ−τdT)exp[jπK(τ−τd)2]⋅exp(jπfcτd) (4) 其中,xP表示目标后向散射系数,h(t)为天线调制函数,τd=2R(t)/2R(t)cc为双程回波时延,c为光速。
在对接收的回波信号进行脉冲压缩和距离走动校正处理后,回波信号可以建模为卷积形式:
y(R,θ)=x(R,θ)⊗h(R,θ)+n(R,θ) (5) 其中
h(R,θ)=h(θ)sinc(2BcR) (6) 在式(5)和式(6)中,y(R,θ)为接收信号,h(R,θ)为天线方向图函数,x(R,θ)为目标散射系数,n(R,θ)为加性噪声,⊗表示卷积运算。
将式(6)离散化后写成如下的形式:
y=Hx+n (7) 其中,y∈CN×1与x∈CN×1分别表示为第m个距离单元的接收信号向量与目标散射系数向量,N为原始场景中某一距离单元方位向的采样点数,n∈CN×1为加性噪声向量,H∈CN×N为天线方向图调制矩阵,其可以表示为
H=[h(θL/L22)h(θL/L2−1(2−1))⋯h(θ1)h(θL/L2(2+1))h(θL/L22)⋱h(θ1)⋮⋮⋱⋱⋱h(θL)⋱⋱h(θ1)h(θL)⋱⋱h(θ2)⋱⋱⋱⋮h(θL)h(θL−1)⋯h(θL/L22)]N×N (8) 其中,[h(θ1), h(θ2), ⋯, h(θL)]表示天线方向图采样向量。
2.2 问题形成
2.2.1 问题分析
在扫描雷达系统中,通常为了保护雷达天线会引入雷达天线罩,但天线罩会引起信号的散射、反射或吸收,使得天线方向图出现不同程度的误差[28]。在实测处理中,天线方向图根据原有的设计参数计算或测得,由于测量误差或平台运动等影响,也会不可避免地引入未知的天线方向图误差。在扫描雷达成像中,天线方向图误差会造成卷积过程失配,导致反演求解过程出现巨大偏差,降低分辨性能甚至出现虚假目标,影响后续的目标检测、识别。
2.2.2 总体最小二乘(Total Least Square, TLS)法
假设理想的天线方向图为h(θ),由其构造的天线卷积矩阵为H。而实际未知的方向图为h(θ)+Δh(θ),由其构造的天线方向图卷积矩阵为⌢H,其中Δh(θ)为未知的误差。对于天线方向图未知的误差值可以用总体最小二乘法[30]的思想来修正。
令A0和b0分别代表观测的无误差数据矩阵和无误差数据向量,实际观测的数据矩阵和数据向量为
A=A0+ΔA (9) b=b0+Δb (10) 其中,ΔA表示误差数据矩阵,Δb表示受ΔA影响造成的数据向量误差。
总体最小二乘法的基本思想是:不仅用校正向量去干扰数据向量,同时用校正矩阵去干扰数据矩阵,以便对二者内存在的误差或噪声进行联合补偿,抑制观测误差或者噪声对矩阵方程求解的影响,从而实现有误差的矩阵方程求解向精确矩阵方程的求解转换:
(A0+ΔA)x=b0+Δb⇒A0x=b0 (11) 在本文中,由式(9)和式(10)类比,实际天线卷积矩阵可以表示为⌢H=H+E,其中H表示无误差的天线方向图矩阵,E表示天线方向图误差矩阵,⌢y和y分别表示有误差和无误差的回波数据。根据2.1中的信号模型,结合式(11)中的校正方式,可以实现有误差的天线卷积矩阵的校正:
⌢Hx=(H+E)x=⌢y⇒Hx=y (12) 3. 自适应迭代重加权稀疏重建方法
本节基于总体最小二乘(TLS)准则,考虑了天线卷积矩阵的误差和目标稀疏先验,构建相应的目标函数;其次,基于变量交替迭代的求解思路,利用迭代重加权(IRN)优化方法推导了方向图误差条件下的稀疏重建方法;更进一步,针对算法的超参数选取难题,本文基于范数比准则导出了自适应参数选取方法。
3.1 自适应迭代重加权方法
根据2.2节的问题分析,本文基于总体最小二乘(TLS)理论,当观测数据y和天线卷积矩阵H均受到扰动时,此时的数学模型转化为
{ˆxTLS,ˆETLS,ˆnTLS}:=argminx,E,n∥[En]∥2Fs.t.y=(H+E)x+n (13) 其中,F代表矩阵的Frobenius范数。通常假设E和n服从复高斯分布且彼此独立,该方法的最优解应该满足由扰动E和噪声n这两项组成的矩阵[En]的F范数最小。
运用二次惩罚的思想,式(13)等价于式(14)的无约束优化问题:
J(ˆx,ˆE)=argminx,E‖ (14) 在雷达成像中,感兴趣的目标往往稀疏地分布在观测场景中,比如海面的舰船编队目标以及机场跑道的异物,此时可利用目标的稀疏先验增加对目标的惩罚项,可构造目标函数为
\begin{split} J(\hat {\boldsymbol{x}},\hat {\boldsymbol{E}}) = \,& \mathop {{\rm{arg}} {\rm{min}}}\limits_{{\boldsymbol{x}},{\boldsymbol{E}}} \left\| {{\boldsymbol{y}} - ({\boldsymbol{H}} + {\boldsymbol{E}}){\boldsymbol{x}}} \right\|_2^2 \\ & + \alpha {\left\| {\boldsymbol{x}} \right\|_1} + \beta \left\| {\boldsymbol{E}} \right\|_{\mathrm{F}}^2 \end{split} (15) 目标函数的求解,可采用变量分离的思想,每次求解时,固定一个变量,剩余的变量视为常数:
(1) 求解E问题:
J(\hat {\boldsymbol{x}},{\boldsymbol{E}}) = \mathop {{\mathrm{arg}} {\mathrm{min}}}\limits_{\boldsymbol{E}} \left\| {{\boldsymbol{y}} - ({\boldsymbol{H}} + {\boldsymbol{E}}){\boldsymbol{x}}} \right\|_2^2 + \beta \left\| {\boldsymbol{E}} \right\|_{\mathrm{F}}^2 (16) 对于E问题的求解,式(16)的最优解是对E的微分等于零,即:
{\nabla _{\boldsymbol{E}}}J(\hat {\boldsymbol{x}},{\boldsymbol{E}}) = {\boldsymbol{Ex}}{{\boldsymbol{x}}^{\mathrm{T}}} - ({\boldsymbol{y}} - {\boldsymbol{Hx}}){{\boldsymbol{x}}^{\mathrm{T}}} + \beta {\boldsymbol{E}} = 0 (17) 整理式(17)可得,对扰动矩阵的估计为
\hat {\boldsymbol{E}} = \frac{{({\boldsymbol{y}} - {\boldsymbol{Hx}}){{\boldsymbol{x}}^{\mathrm{T}}}}}{{{\boldsymbol{x}}{{\boldsymbol{x}}^{\mathrm{T}}} + \beta }} (18) (2) 求解x问题:
对于x问题的求解,可以将扰动矩阵E固定,建立以下目标函数,采用迭代重加权优化算法进行求解。
J({\boldsymbol{x}},\hat {\boldsymbol{E}}) = \mathop {{\mathrm{arg}} {\mathrm{min}}}\limits_{\boldsymbol{x}} \left\| {{\boldsymbol{y}} - ({\boldsymbol{H}} + {\boldsymbol{E}}){\boldsymbol{x}}} \right\|_2^2 + \alpha {\left\| {\boldsymbol{x}} \right\|_1} (19) 由于式(19)中 {l_1} 范数的不可导特性,迭代重加权范数方法可以通过用 {l_2} 范数逼近 {l_1} 范数以获得稀疏解。因此,本节定义一个对角矩阵W为
{\boldsymbol{W}} = {\mathrm{diag}}\left( {{{\left| {\boldsymbol{x}} \right|}^{ - 1}}} \right) (20) 其中,W中的对角元素值可以由目标散射系数迭代求得。IRN方法的求解步骤总结如下:
步骤1 初始化:选择初始系数估计 {{\boldsymbol{x}}^{\left( 0 \right)}} 和权重向量 {{\boldsymbol{W}}^{\left( 0 \right)}} 。
步骤2 迭代:
① 计算权重向量:
{{\boldsymbol{W}}^{\left( k \right)}} = {\mathrm{diag}}\left( {{{\left| {{{\boldsymbol{x}}^{\left( k \right)}}} \right|}^{ - 1}}} \right) (21) ② 求解带有权重的正则化问题:
{\hat {\boldsymbol{x}}^{\left( {k + 1} \right)}} = {\left( {{{({\boldsymbol{H}} + {\boldsymbol{E}})}^{\mathrm{T}}}({\boldsymbol{H}} + {\boldsymbol{E}}) + \alpha {{\boldsymbol{W}}^{\left( k \right)}}} \right)^{ - 1}}{ ({\boldsymbol{H}} + {\boldsymbol{E}})^{\mathrm{T}}}{\boldsymbol{y}} (22) ③ 更新迭代次数: k = k + 1 。
3.2 自适应参数选取方法
在正则化方法中,正则化参数影响着算法的最终性能,适当的正则化参数可以在误差达到最小值的同时,最大化地利用解的先验信息。目前被广泛使用的正则化参数选取方法主要有广义交叉验证(Generalized Cross-Validation, GCV)方法[31]、L-curve方法[32]、误差极小化准则[33]等。对于本文的稀疏正则化问题,在3.1节求解时,涉及两个参数的选取,即 \alpha 和 \beta ,其中通过调整 \alpha 的大小以控制 {l_1} 范数惩罚项的强度,以获取稀疏的重建结果。通过调整参数 \beta 能控制误差项的强度,以更好地校正天线方向图误差。
正则化参数的选取是目前稀疏正则化问题中的难点,参数的选取是否合适很大程度上影响超分辨成像的性能,近年来大多数正则化方法的参数依赖人工经验选取,算法泛化能力差;为解决这一问题,本文使用范数比正则化方法来自适应获取参数值[34]。范数比法的核心思想是基于正则化项和数据拟合项之间的范数比来确定最佳的正则化参数值,更好地权衡模型的拟合程度和正则化项的强度,在迭代过程中动态地调整正则化参数,以使正则化项和数据拟合项之间的范数比接近某个理想的目标值,帮助防止过拟合和提高模型的泛化能力。
本文中对于参数 \beta 选取方法为
\beta = \frac{{\left\| {{\boldsymbol{y}} - ({\boldsymbol{H}} + {\boldsymbol{E}}){\boldsymbol{x}}} \right\|_2^2}}{{\left\| {\boldsymbol{E}} \right\|_{\mathrm{F}}^2}} (23) 而对于正则化参数 \alpha ,其本质与噪声方差有关,从贝叶斯估计理论出发,其选取方法为
\alpha = \frac{{\left\| {{\boldsymbol{y}} - ({\boldsymbol{H}} + {\boldsymbol{E}}){\boldsymbol{x}}} \right\|_2^2}}{N} (24) 3.3 算法计算复杂度分析
由2.1节可知,第m个距离单元回波y的矩阵维度为{{\bf{C}}^{N \times 1}},假设采集的雷达回波数据共有M个距离单元,在雷达成像中为节省硬件资源,一般采用沿距离向逐行处理的方式。下面分别分析传统IRN-L1和本文所提方法的运算复杂度。
由3.1节的推导可知,假设IRN-L1方法的迭代次数为K。计算W的复杂度为 O(N) ;计算 {{\boldsymbol{H}}^{\mathrm{T}}} 的复杂度为 O({N^2}) ;计算 {{\boldsymbol{H}}^{\mathrm{T}}}{\boldsymbol{H}} 的复杂度为 O({N^3}) ;计算 {{\boldsymbol{H}}^{\mathrm{T}}}{\boldsymbol{H}} + \alpha {{\boldsymbol{W}}^{\left( k \right)}} 的复杂度为 O(2{N^2}) ;计算 {{\boldsymbol{H}}^{\mathrm{T}}}{\boldsymbol{y}} 的复杂度为 O(2{N^2} - N) ;计算 {({{\boldsymbol{H}}^{\mathrm{T}}}{\boldsymbol{H}} + \alpha {{\boldsymbol{W}}^{\left( k \right)}})^{ - 1}} 的复杂度为 O({N^3}) ;计算 {({{\boldsymbol{H}}^{\mathrm{T}}}{\boldsymbol{H}} + \alpha {{\boldsymbol{W}}^{\left( k \right)}})^{ - 1}}{{\boldsymbol{H}}^{\mathrm{T}}}{\boldsymbol{y}} 的复杂度为 O(2{N^2} - N) ;因此单次求解x需要的总计算复杂度为 O(2{N^3} + 7{N^2} - N) 。综上,IRN-L1方法的总计算复杂度约为 O(M({K_1}(2{N^3} + 7{N^2} - N))) 。
从3.1节关于算法推导可知,本文所提方法有内外两层循环,假设内外循环的迭代次数分别为 {K_1} 与 {K_2} 。对于内循环求解误差矩阵E:首先要计算x的转置 {{\boldsymbol{x}}^{\mathrm{T}}} ,其计算复杂度为 O(N) ;计算 {\boldsymbol{x}}{{\boldsymbol{x}}^{\mathrm{T}}} 的复杂度为 O({N^2}) ;计算{\boldsymbol{Hx}}的复杂度为 O(2{N^2} - N) ;计算 {\boldsymbol{y}} - {\boldsymbol{Hx}} 的复杂度为 O(N) ;计算 ({\boldsymbol{y}} - {\boldsymbol{Hx}}){{\boldsymbol{x}}^{\mathrm{T}}} 的复杂度为 O({N^2}) ;计算 {\boldsymbol{x}}{{\boldsymbol{x}}^{\mathrm{T}}} + \beta 的复杂度为 O({N^2}) ; ({\boldsymbol{y}} - {\boldsymbol{Hx}}){{\boldsymbol{x}}^{\mathrm{T}}}/({\boldsymbol{x}}{{\boldsymbol{x}}^{\mathrm{T}}} + \beta ) 的计算,可以看作 ({\boldsymbol{y}} - {\boldsymbol{Hx}}){{\boldsymbol{x}}^{\mathrm{T}}} 与 {\boldsymbol{x}}{{\boldsymbol{x}}^{\mathrm{T}}} + \beta 的逆矩阵相乘,计算 ({\boldsymbol{y }} - {\boldsymbol{Hx}}){{\boldsymbol{x}}^{\mathrm{T}}}/ ({\boldsymbol{x}}{{\boldsymbol{x}}^{\mathrm{T}}} + \beta ) 的复杂度为 O(2{N^3}) ;因此单次求解E需要的计算复杂度为 O(2{N^3} + 5{N^2} + N) 。
对于外循环求解目标散射系数矩阵x:计算W的复杂度为 O(N) ;计算 {({\boldsymbol{H}} + {\boldsymbol{E}})^{\mathrm{T}}} 的复杂度为 O(2{N^2}) ;计算 {({\boldsymbol{H}} + {\boldsymbol{E}})^{\mathrm{T}}}({\boldsymbol{H}} + {\boldsymbol{E}}) 的复杂度为 O({N^3}) ;计算 {({\boldsymbol{H}} + {\boldsymbol{E}})^{\mathrm{T}}}({\boldsymbol{H}} + {\boldsymbol{E}}) + \alpha {{\boldsymbol{W}}^{\left( k \right)}} 的复杂度为 O(2{N^2}) ;计算 {({\boldsymbol{H }}+ {\boldsymbol{E}})^{\mathrm{T}}}{\boldsymbol{y}} 的复杂度为 O(2{N^2} - N) ;计算 {({({\boldsymbol{H}} + {\boldsymbol{E}})^{\mathrm{T}}}({\boldsymbol{H}} + {\boldsymbol{E}}) + \alpha {{\boldsymbol{W}}^{\left( k \right)}})^{ - 1}} 的复杂度为 O({N^3}) ;计算 {({({\boldsymbol{H}} + {\boldsymbol{E}})^{\mathrm{T}}}({\boldsymbol{H}} + {\boldsymbol{E}}) + \alpha {{\boldsymbol{W}}^{\left( k \right)}})^{ - 1}} {({\boldsymbol{H}} + {\boldsymbol{E}})^{\mathrm{T}}}{\boldsymbol{y}} 的复杂度为 O(2{N^2} - N) ;计算正则化参数 \alpha , \beta 的复杂度为 O(3{N^2} + N) ;因此单次求解x需要的总计算复杂度为 O(2{N^3} + 11{N^2}) 。
综上,本文所提方法的计算复杂度为 O(M({K_1} \cdot (2{N^3} + 5{N^2} + N) + {K_1}{K_2}(2{N^3} + 11{N^2}))) 。
4. 仿真实验
在仿真实验中,本文采用一维点目标仿真以及二维面目标仿真来验证本文提出方法的成像性能。为了验证本文方法在未知方向图误差条件下的稳健性和适应性,本节结合实际情况,在仿真中主要考虑了3种天线方向图误差类型:展宽误差、随机误差、展宽误差+随机误差。仿真实验所使用的系统参数与仿真环境分别见表1和表2。
表 1 仿真参数Table 1. Simulation parameters参数 数值 载频(GHz) 10.75 信号带宽(MHz) 40 脉冲重复频率(Hz) 1000 天线扫描速度(°/s) 60 波束宽度(°) 3 扫描范围(°) –10~10 表 2 仿真环境Table 2. Simulation conditions硬件/软件 参数值 CPU Intel(R) Core(TM)i7-9700K RAM 64 GB 仿真软件 Matlab 2022a 首先定义天线展宽系数如下:
\gamma = \frac{{{\theta _{{\mathrm{Broaden}}}}}}{{{\theta _{{\mathrm{Ideal}}}}}} (25) 其中, {\theta _{{\mathrm{Broaden}}}} , {\theta _{{\mathrm{Ideal}}}} 分别表示展宽的天线3 dB主瓣宽度与理想的天线3 dB主瓣宽度。
随机误差定义为理想的天线方向图矩阵上附加一个微小的随机值构成随机误差下的天线方向图,其中,随机误差系数 \xi 分别取0.1, 0.2以模拟实际情况下存在的随机误差。
{{\boldsymbol{h}}_{{\mathrm{Random}}\; {\mathrm{errors}}}} = {{\boldsymbol{h}}_{{\mathrm{Ideal}}}} + \xi \times{\mathrm{rand}}(1,N) (26) 4.1 点目标仿真
4.1.1 不同误差类型下的点目标仿真
在本节中,通过设置一维点目标来验证提出方法的成像性能。一维点目标场景如图2所示,在某一距离单元–1.4°, 1.4°的角度位置分别设置了两个幅度为1的点目标。另外,回波信噪比设置为20 dB。
4.1.1.1 展宽误差下的点目标仿真
展宽系数设置为1.2时,天线方向图展宽误差下的一维点目标处理结果如图3所示。
图3(a)表示受展宽误差影响的天线方向图;图3(b)表示受到噪声污染的实波束回波,可以看出,受到方向图调制影响,目标的实波束回波混叠在一起。图3(c)和图3(d)分别表示展宽系数为1.2条件下维纳滤波方法与L2正则化方法的处理结果,由结果可知,维纳逆滤波方法与L2正则化方法的分辨率提升有限且旁瓣较高。图3(e)和图3(f)分别表示IRN-L1方法与本文所提方法的处理结果,展宽误差下IRN-L1方法的处理结果中目标幅度有所损失,且分辨率低于本文方法的处理结果。
为进一步验证本文所提算法在展宽误差较大时的性能,将展宽系数设置为1.4,一维点目标处理结果如图4所示。图4(a)表示受展宽误差影响的天线方向图;图4(b)表示受到噪声污染的实波束回波,可以看出,受到展宽方向图调制影响,目标的实波束回波混叠更为严重。图4(c)和图4(d)分别表示展宽系数为1.4条件下维纳滤波方法与L2正则化方法的处理结果,由结果可知,维纳逆滤波方法与L2正则化方法的分辨率下降明显,已几乎无法分辨两个点目标。图4(e)表示IRN-L1方法的处理结果,展宽系数为1.4下的IRN-L1方法的处理结果分辨率有所下降,与展宽系数为1.2条件下相比,目标有所展宽。图4(f)表示本文提出方法的处理结果,可以看出,本文提出的方法仍能有效地分辨出两个点目标,且目标几乎无展宽。
4.1.1.2 随机误差下的点目标仿真
设置随机误差系数为0.1,随机误差下的一维点目标处理结果如图5所示。图5(a)表示受随机误差影响的天线方向图;图5(b)表示受到噪声污染的实波束回波,可以看出在随机误差情况下,在实波束回波中仍无法对目标进行区分。图5(c)和图5(d)分别表示随机误差下的维纳逆滤波方法与L2正则化方法的处理结果。虽然上述两种方法可以分辨两个点目标,但存在较高的旁瓣。图5(e)表示IRN-L1方法的处理结果,由结果可知,在随机误差系数0.1条件下,IRN-L1方法分辨两个点目标;图5(f)表示本文提出方法的处理结果,本文提出的方法可有效地分辨出两个点目标,且分辨率要优于IRN-L1方法的处理结果。
为进一步验证本文所提算法的性能,将随机误差系数设置为0.2,随机误差下的一维点目标处理结果如图6所示。图6(a)表示受随机误差影响的天线方向图;图6(b)表示受到噪声污染的实波束回波,可以看出在随机误差情况下,在实波束回波中仍无法对目标进行区分。图6(c)和图6(d)分别表示随机误差下的维纳逆滤波方法与L2正则化方法的处理结果。虽然上述两种方法可以分辨两个点目标,分辨率有所下降,依旧存在较高的旁瓣。图6(e)表示IRN-L1方法的处理结果,由结果可知,在随机误差系数为0.2条件下,分辨率下降,目标幅度也出现失真;图6(f)表示本文提出方法的处理结果,本文提出的方法分辨率几乎无损失,优于 IRN-L1方法。
4.1.1.3 展宽误差+随机误差下的点目标仿真
本节将展宽系数设置为1.2,随机误差系数设置为0.1,天线方向图展宽误差与随机误差结合下的一维点目标处理结果如图7所示。图7(a)表示同时受随机误差与展宽误差影响的天线方向图,图7(b)表示受到噪声污染的实波束回波。可以看出,在展宽误差和随机误差结合情况下,图7(c)和图7(d)所示的维纳逆滤波方法与L2正则化方法的处理结果对目标的分辨能力明显下降,超分辨结果与目标原始分布偏差较大。图7(e)表示IRN-L1方法的处理结果,由结果可知,在展宽误差和随机误差结合下,成像结果目标分辨率进一步下降,目标展宽程度进一步加剧。图7(f)表示本文提出方法的处理结果,本文提出的方法分辨率优于IRN-L1方法,证明了本文方法的有效性。
为进一步验证本文所提算法的性能,将展宽系数设置为1.4,随机误差系数设置为0.2,探究天线方向图展宽误差与随机误差加剧时的一维点目标处理结果。图8(a)表示同时受随机误差与展宽误差影响的天线方向图,图8(b)表示受到噪声污染的实波束回波。可以看出,在展宽误差和随机误差进一步变大的情况下,图8(c)和图8(d)所示的维纳逆滤波方法与L2正则化方法的处理结果对两个点目标已几乎无法分辨,旁瓣更为严重。图8(e)表示IRN-L1方法的处理结果,成像结果目标展宽程度进一步加剧。图8(f)表示本文提出方法的处理结果,本文提出的方法虽然分辨率也出现一定的下降,但仍优于IRN-L1方法。
4.1.2 不同信噪比下的算法性能验证
本节通过设置不同的信噪比来验证本文方法的稳健性,仿真参数与仿真场景均与4.1.1节中的点目标仿真相同。为了定量评价点目标的超分辨重建结果,通过计算随机误差与展宽误差结合下的均方根误差(Root Mean Square Error, RMSE)指标作为衡量标准。
均方根误差定义为
{\mathrm{RMSE}} = \sqrt {\frac{1}{N}\sum\limits_{i = 1}^N {{{\left( {{{\hat {\boldsymbol{x}}}_i} - {{\boldsymbol{x}}_i}} \right)}^2}} } (27) 其中, {{\boldsymbol{x}}_i} 和 {\hat {\boldsymbol{x}}_i} 分别表示第i个真实的目标散射系数与求解的目标散射系数。均方根误差越小,表示求解的目标散射系数越接近真实的目标散射系数,目标重建性能越好。
由图9的均方根误差曲线图可知,随着信噪比的提高,超分辨方法的目标重建性会逐步提升。在信噪比为5 dB时,本文方法性能会明显下降,当信噪比≥10 dB时,本文方法能获得优于传统方法的均方根误差。
4.2 面目标仿真
4.2.1 不同误差类型下的面目标仿真
在本节中,通过设置二维面目标来验证提出方法在未知天线方向图误差下的成像性能。原始场景如图10(a)和图10(b)所示,这是一个由11艘舰艇组成的舰队。与整个成像场景相比,舰船编队可视为稀疏目标。仿真实验所使用的雷达系统参数与表1相同,信噪比设置为20 dB。
4.2.1.1 展宽误差下的面目标仿真
展宽系数为1.2条件下的仿真结果如图11所示,其中图11(a)表示加噪实波束回波,实波束回波分辨率较低,无法区分舰船编队。图11(b)和图11(c)分别表示维纳逆滤波方法和L2正则化方法在超分辨处理后的成像结果。上述两种方法的分辨率提升有限,难以获取有效的编队数量信息。图11(d)表示IRN-L1方法处理后的成像结果,由于受到方向图展宽误差的影响,11艘舰艇目标出现了一定的拓展。图11(e)表示本文提出方法的成像结果图,本文方法的分辨率优于IRN-L1方法,锐化效果更强。
展宽系数为1.4下的仿真结果如图12所示,其中图12(a)表示加噪实波束回波。图12(b)和图12(c)分别表示维纳逆滤波方法和L2正则化方法在超分辨处理后的成像结果。上述两种方法的分辨率进一步下降,目标展宽严重。图12(d)表示IRN-L1方法处理后的成像结果, 11艘舰艇目标展宽程度明显加剧。图12(e)表示本文提出方法的成像结果图,本文提出的方法几乎没有展宽,分辨率优于IRN-L1方法,但存在部分目标幅度失真。
4.2.1.2 随机误差下的面目标仿真
随机误差系数为0.1条件下的面目标仿真结果如图13所示,其中图13(a)表示加噪实波束回波,实波束回波中目标混叠,无法区分目标数量。图13(b)和图13(c)分别表示维纳逆滤波方法和L2正则化方法在超分辨处理后的成像结果。由结果可知,在随机误差情况下,上述两种方法处理后的分辨率略微提升,11艘舰艇仅能大致分辨,但存在较强的噪声和杂点。图13(d)表示IRN-L1方法处理后的成像结果,11艘舰艇可分辨,但目标出现了一定程度的展宽。图13(e)表示本文提出方法的成像结果图,11艘舰艇能清晰分辨。
随机误差系数为0.2条件下的面目标仿真结果如图14所示,其中图14(a)表示加噪实波束回波。图14(b)和图14(c)分别表示维纳逆滤波方法和L2正则化方法在超分辨处理后的成像结果。由结果可知,在随机误差情况下,上述两种方法处理后的11艘舰艇目标的分辨率有所下降,出现了较强的噪声和伪影。图14(d)表示IRN-L1方法处理后的成像结果,11艘舰艇可分辨,但出现了一定程度的伪影和虚假点。图14(e)表示本文提出方法的成像结果图,11艘舰艇仍可以分辨,但部分目标幅度存在失真且有微弱的伪影。
4.2.1.3 展宽误差+随机误差下的面目标仿真
本节将展宽系数设置为1.2,随机误差系数设置为0.1,天线方向图展宽误差与随机误差结合下的面目标处理结果如图15所示,其中图15(a)表示实波束回波,由于展宽误差和随机误差的双重影响,实波束回波中目标混叠严重,几乎无法区分被覆盖目标的数量。图15(b)和图15(c)分别表示维纳逆滤波方法和L2正则化方法在超分辨处理后的成像结果图。由结果图可知,在展宽误差和随机误差的双重影响下,上述两种方法处理后的结果恶化,成像场景充满了噪声。图15(d)表示IRN-L1方法处理后的成像结果,受到展宽误差和随机误差的影响,目标出现了一定程度的展宽,出现了部分的伪影。图15(e)表示本文提出方法的成像结果图,11艘舰艇可分辨且几乎不存在伪影。
为进一步验证本文所提算法的性能,将展宽系数设置为1.4,随机误差系数设置为0.2,探究天线方向图展宽误差与随机误差加剧时的面目标处理结果,展宽误差和随机误差下的面目标仿真结果如图16所示,其中图16(a)表示加噪实波束回波,由于展宽误差和随机误差的进一步加剧,实波束回波中目标混叠更为严重,完全无法区分编队目标的数量。图16(b)和图16(c)分别表示维纳逆滤波方法和L2正则化方法在超分辨处理后的成像结果图。由结果图可知,在展宽误差和随机误差的双重影响下,上述两种方法处理后的结果中11艘舰艇几乎无法分辨,且充满了噪声。图16(d)表示IRN-L1方法处理后的成像结果,受到展宽误差和随机误差的影响,编队目标展宽程度加剧,且出现了较严重的伪影。图16(e)表示本文提出方法的成像结果图,11艘舰艇仍可有效分辨,但也出现了部分的伪影与目标幅度失真。
4.2.2 不同信噪比下的面目标重建性能验证
在本节中,通过设置不同的信噪比来验证提出方法的性能,其仿真参数与仿真场景均与4.2.1节的面目标仿真相同。为了定量评价面目标的重建结果,我们通过计算随机误差与展宽误差结合下的图像熵(Entropy)指标作为衡量标准,图像熵越小,表明超分辨重建效果越好。
由图17图像熵曲线图可知,随着信噪比的提高,两种方法的目标重建性能都会呈现逐步提升的趋势。相比于传统的维纳逆滤波方法、L2正则化方法、IRN-L1方法,在信噪比≥10 dB时,本文所提出的方法在图像熵方面优于传统的超分辨方法,证明了其稳健性。
5. 实测验证
为验证所提算法在实际应用中的有效性,我们利用项目组研制的扫描雷达系统,录取实测数据并进行处理。图18为角反射器光学图,两个角反的间距为10 m。扫描雷达系统参数如表3所示。
表 3 扫描雷达系统参数Table 3. Scanning radar system parameters参数 数值 载频(GHz) 30.75 信号带宽(MHz) 200 脉冲重复频率(Hz) 4000 主瓣宽度(°) 4 天线扫描速度(°/s) 60 信号时宽(μs) 1 扫描范围(°) –35~35 图19(a)为经过距离维预处理后的实波束回波,我们可以看到两个角反的实波束回波混叠在一起,无法分辨。图19(b)和图19(c)分别是维纳逆滤波与L2正则化方法的处理结果,可以看到两种方法虽然平滑了噪声,但是由于分辨率提升有限,两个角反仍然无法区分。图19(d)表示IRN-L1方法的处理结果,虽然能大致辨别出两个角反,但是受到误差影响,角反目标在方位维存在扩展,仍然相连;图19(e)为本文提出方法的处理结果,本文方法可以完全分离2个角反目标,具有更优的分辨率。综上,实测数据处理结果证明了本文所提方法在实际应用中能获得优于传统方法的处理效果。
6. 结语
本文提出一种未知天线方向图误差条件下的扫描雷达角超分辨方法。基于总体最小二乘法准则,本文所提方法考虑了方向图误差的影响,导出了相应的目标函数,并采用迭代重加权方法实现了目标函数求解。点目标、面目标与实测数据处理结果表明,本文方法在一定的天线方向图误差条件下能实现稀疏目标的有效重建。相比于传统的超分辨方法,在信噪比≥10 dB时,能一定程度缓解天线方向图误差引起的展宽、虚假伪影和杂点等,获得更优的处理效果。
-
表 1 仿真参数
Table 1. Simulation parameters
参数 数值 载频(GHz) 10.75 信号带宽(MHz) 40 脉冲重复频率(Hz) 1000 天线扫描速度(°/s) 60 波束宽度(°) 3 扫描范围(°) –10~10 表 2 仿真环境
Table 2. Simulation conditions
硬件/软件 参数值 CPU Intel(R) Core(TM)i7-9700K RAM 64 GB 仿真软件 Matlab 2022a 表 3 扫描雷达系统参数
Table 3. Scanning radar system parameters
参数 数值 载频(GHz) 30.75 信号带宽(MHz) 200 脉冲重复频率(Hz) 4000 主瓣宽度(°) 4 天线扫描速度(°/s) 60 信号时宽(μs) 1 扫描范围(°) –35~35 -
[1] Rigelsford J. Introduction to airborne radar[J]. Sensor Review, 2002, 22(3): 265–266. [2] 张良, 祝欢, 吴涛. 机载预警雷达系统架构发展路径研究[J]. 现代雷达, 2015, 37(12): 11–18. doi: 10.16592/j.cnki.1004-7859.2015.12.003.ZHANG Liang, ZHU Huan, and WU Tao. A study on the evolution way of the system architecture of AEW radar[J]. Modern Radar, 2015, 37(12): 11–18. doi: 10.16592/j.cnki.1004-7859.2015.12.003. [3] CLARKE J, 徐映和, 译. 英国预警雷达的发展概况[J]. 现代雷达, 1987, 9(2): 1–6. doi: 10.16592/j.cnki.1004-7859.1987.02.001.CLARKE J, XU Yinghe, translation. Overview of the development of early warning radar in the UK[J]. Modern Radar, 1987, 9(2): 1–6. doi: 10.16592/j.cnki.1004-7859.1987.02.001. [4] ZHANG Qiping, ZHANG Yin, HUANG Yulin, et al. TV-sparse super-resolution method for radar forward-looking imaging[J]. IEEE Transactions on Geoscience and Remote Sensing, 2020, 58(9): 6534–6549. doi: 10.1109/TGRS.2020.2977719. [5] LING Hao. Novel radar techniques and applications[J]. IEEE Antennas and Propagation Magazine, 2018, 60(1): 132–134. doi: 10.1109/MAP.2017.2776153. [6] BEKKADAL F. Novel radar technology and applications[C]. 17th International Conference on Applied Electromagnetics and Communications, Dubrovnik, Croatia, 2003: 6–12. doi: 10.1109/ICECOM.2003.1290942. [7] ZHANG Yongchao, ZHANG Yin, LI Wenchao, et al. Super-resolution surface mapping for scanning radar: Inverse filtering based on the fast iterative adaptive approach[J]. IEEE Transactions on Geoscience and Remote Sensing, 2018, 56(1): 127–144. doi: 10.1109/TGRS.2017.2743263. [8] ZHANG Yongchao, JAKOBSSON A, ZHANG Yin, et al. Wideband sparse reconstruction for scanning radar[J]. IEEE Transactions on Geoscience and Remote Sensing, 2018, 56(10): 6055–6068. doi: 10.1109/TGRS.2018.2830100. [9] CHEN Rui, LI Wenchao, LI Kefeng, et al. A super-resolution scheme for multichannel radar forward-looking imaging considering failure channels and motion error[J]. IEEE Geoscience and Remote Sensing Letters, 2023, 20: 3501305. doi: 10.1109/LGRS.2023.3234264. [10] KANG Yao, ZHANG Yin, MAO Deqing, et al. Bayesian azimuth super-resolution of sea-surface target in forward-looking imaging[C]. 2020 IEEE Radar Conference (RadarConf20), Florence, Italy, 2020: 1–5. doi: 10.1109/RadarConf2043947.2020.9266692. [11] CHEN Hongmeng, WANG Zeyu, ZHANG Yingjie, et al. Data-driven airborne Bayesian forward-looking superresolution imaging based on generalized Gaussian distribution[J]. Frontiers in Signal Processing, 2023, 3: 1093203. doi: 10.3389/frsip.2023.1093203. [12] MAO Deqing, ZHANG Yin, ZHANG Yongchao, et al. Super-resolution Doppler beam sharpening method using fast iterative adaptive approach-based spectral estimation[J]. Journal of Applied Remote Sensing, 2018, 12(1): 015020. doi: 10.1117/1.JRS.12.015020. [13] LIU Sijia and PAN Minghai. Research on a forward-looking scanning imaging algorithm for a high-speed radar platform[J]. IET Signal Processing, 2023, 17(6): e12221. doi: 10.1049/sil2.12221. [14] MAO Deqing, YANG Jianyu, TUO Xingyu, et al. Angular superresolution of real aperture radar for target scale measurement using a generalized hybrid regularization approach[J]. IEEE Transactions on Geoscience and Remote Sensing, 2023, 61: 5109314. doi: 10.1109/TGRS.2023.3315310. [15] TUO Xingyu, MAO Deqing, ZHANG Yin, et al. Radar forward-looking super-resolution imaging using a two-step regularization strategy[J]. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2023, 16: 4218–4231. doi: 10.1109/JSTARS.2023.3270309. [16] YOUNG P. Alternative Recursive Approaches to Time-series Analysis[M]. YOUNG P. Recursive Estimation and Time-Series Analysis: An Introduction. Berlin, Heidelberg: Springer, 1984: 205–230. [17] RICHARDS M A. Iterative noncoherent angular superresolution (radar)[C]. 1988 IEEE National Radar Conference, Ann Arbor, USA, 1988: 100–105. doi: 10.1109/NRC.1988.10940. [18] LI Dongye, HUANG Yulin, and YANG Jianyu. Real beam radar imaging based on adaptive Lucy-Richardson algorithm[C]. 2011 IEEE CIE International Conference on Radar, Chengdu, China, 2011: 1437–1440. doi: 10.1109/CIE-Radar.2011.6159830. [19] TAN Ke, LU Xingyu, YANG Jianchao, et al. A novel Bayesian super-resolution method for radar forward-looking imaging based on Markov random field model[J]. Remote Sensing, 2021, 13(20): 4115. doi: 10.3390/rs13204115. [20] CHEN Hongmeng, LI Yachao, GAO Wenquan, et al. Bayesian forward-looking super-resolution imaging using Doppler deconvolution in expanded beam space for high-speed platform[J]. IEEE Transactions on Geoscience and Remote Sensing, 2022, 60: 5105113. doi: 10.1109/TGRS.2021.3107717. [21] LI Weixin, LI Ming, ZUO Lei et al. A computationally efficient airborne forward-looking super-resolution imaging method based on sparse Bayesian learning[J]. IEEE Transactions on Geoscience and Remote Sensing, 2023, 61: 5102613. doi: 10.1109/TGRS.2023.3260094. [22] CAPON J. High-resolution frequency-wavenumber spectrum analysis[J]. Proceedings of the IEEE, 1969, 57(8): 1408–1418. doi: 10.1109/PROC.1969.7278. [23] ZHANG Yongchao, LI Wenchao, ZHANG Yin, et al. A fast iterative adaptive approach for scanning radar angular superresolution[J]. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 2015, 8(11): 5336–5345. doi: 10.1109/JSTARS.2015.2449090. [24] LI Yueli, LIU Jianguo, JIANG Xiaoqing, et al. Angular superresolution for signal model in coherent scanning radars[J]. IEEE Transactions on Aerospace and Electronic Systems, 2019, 55(6): 3103–3116. doi: 10.1109/TAES.2019.2900133. [25] GAMBARDELLA A and MIGLIACCIO M. On the superresolution of microwave scanning radiometer measurements[J]. IEEE Geoscience and Remote Sensing Letters, 2008, 5(4): 796–800. doi: 10.1109/LGRS.2008.2006285. [26] HUO Weibo, TUO Xingyu, ZHANG Yin, et al. Balanced tikhonov and total variation deconvolution approach for radar forward-looking super-resolution imaging[J]. IEEE Geoscience and Remote Sensing Letters, 2022, 19: 3505805. doi: 10.1109/LGRS.2021.3072389. [27] TANG Junkui, LIU Zheng, RAN Lei, et al. Enhancing forward-looking image resolution: Combining low-rank and sparsity priors[J]. IEEE Transactions on Geoscience and Remote Sensing, 2023, 61: 5100812. doi: 10.1109/TGRS.2023.3237332. [28] KOZAKOFF D J. Analysis of radome-enclosed antennas[J]. Artech House, 1997 (685): 217–227. [29] PERSSON K and GUSTAFSSON M. Reconstruction of equivalent currents using a near-field data transformation-with radome applications[J]. Progress in Electromagnetics Research, 2005, 54: 179–198. doi: 10.2528/PIER04111602. [30] Fierro R D, Golub G H, Hansen P C, et al. Regularization by truncated total least squares[J]. SIAM Journal on Scientific Computing, 1997, 18(4): 1223–1241. doi: 10.1137/S106482759426383. [31] GOLUB G H, HEATH M, and WAHBA G. Generalized cross-validation as a method for choosing a good ridge parameter[J]. Technometrics, 1979, 21(2): 215–223. doi: 10.1080/00401706.1979.10489751. [32] JOHNSTON P R and GULRAJANI R M. Selecting the corner in the L-curve approach to Tikhonov regularization[J]. IEEE Transactions on Biomedical Engineering, 2000, 47(9): 1293–1296. doi: 10.1109/10.867966. [33] ENGL H W. Discrepancy principles for Tikhonov regularization of ill-posed problems leading to optimal convergence rates[J]. Journal of optimization theory and applications, 1987, 52: 209–215. doi: 10.1007/BF00941281. [34] RUDIN L I, OSHER S, and FATEMI E. Nonlinear total variation based noise removal algorithms[J]. Physica D: Nonlinear Phenomena, 1992, 60(1/4): 259–268. doi: 10.1016/0167-2789(92)90242-F. -