1　引言

地基光学望远镜是天文观测的重要工具，然而大气湍流使来自目标的光波产生波前畸变，导致光学望远镜的成像分辨力严重下降。自适应光学(Adaptive optics, AO)技术是减弱大气湍流影响最有效的手段之一，被广泛应用于地基光学望远镜中。AO技术能够实时探测并校正波前畸变，可有效提高望远镜的成像分辨力。然而受硬件性能等限制，AO系统波前探测与校正过程存在误差，只能实现部分校正，观测图像仍然存在残余模糊，需要结合图像事后重建才能获得更高分辨力的图像^[1]。

目前，国内外的大口径太阳望远镜几乎都配备了AO系统，采集到的太阳AO图像可通过盲解卷积^[2]、相位差法^[3-5]、斑点重建^[6]或深度学习^[7-9]等方法进行重建，进一步提高图像质量。其中，相比于其它几种图像事后重建方法，盲解卷积不需要系统设置额外的成像通道，所需的图像帧数较少，不依赖于大量数据训练，最具灵活性。

早期，由于空间观测的需求，盲解卷积主要应用于天文图像。Ayers和Dainty^[10]提出的迭代盲解卷积(iterative blind deconvolution, IBD)算法属于较早的盲解卷积研究工作，其利用简单的点扩散函数(point spread function, PSF)正性约束和能量守恒约束，在空域和频域进行交替迭代估计图像和PSF，然而IBD算法易放大噪声影响。随后，Davey^[11]等人和Fish^[12]等人分别将维纳滤波和Richardson-Lucy算法引入IBD中，提出Weiner-IBD和RL-IBD，提高了图像复原效果。Kundur和Hatzinakos^[13]等人提出非负约束和支持域约束的递归逆滤波(NAS-RIF)算法，模型为凸函数，具有较好的收敛性，但NAS-RIF算法对噪声敏感。Chan和Wong^[14]首次将总变分(total variation, TV)引入盲解卷积，并采用交替最小化方法优化目标函数，既提高了算法对噪声的鲁棒性，也更好地保持了图像的边缘。由于天文光学望远镜通常是连续曝光的，能够获得图像序列，且多帧有利于提高盲解卷积的稳定性，复原更多的小尺度结构，因此利用多帧图像进行盲解卷积也有不少研究工作。Schulz^[15]提出了基于最大似然估计的多帧盲解卷积算法，并指出可以利用惩罚最大似然估计或参数化PSF方法避免平凡解。Harmeling^[16]等人提出了一种在线多帧盲解卷积(online blind deconvolution, OBD)算法，图像序列依次进行处理，直到序列中所有图像都参与盲解卷积。Hirsch^[17]等人结合超分辨和饱和校正进一步扩展了OBD，并利用随机梯度下降优化模型。

进入二十一世纪后，由于摄像设备的普及，产生了更多运动去模糊的需求，图像盲去模糊(blind deblurring)领域得到了迅速发展。图像盲去模糊实际上就是盲解卷积的不同说法。Fergus^[18]等人分别以混合高斯分布和混合指数分布拟合图像先验和PSF先验，是基于变分贝叶斯(variational Bayes, VB)方法实现图像盲去模糊的重要工作。Fergus之后出现了更多基于VB的图像盲去模糊研究工作^[19-21]。然而此类方法理论相对复杂，计算复杂度也较高。另一类更为流行的是基于最大后验概率(maximum a posteriori, MAP)的方法，此类方法的盲解卷积模型以更简单形式的数据保真项、图像正则项(先验项)和PSF正则项(先验项)构成，更具灵活性，便于计算求解。基于MAP的盲解卷积方法关键在于选择图像和PSF先验，因此大量图像和PSF先验被提出和应用，如L1/L2规范化稀疏先验^[22]、近似L0稀疏^[23]、块先验^[24]、内部块递归性^[25]、L0范数图像强度和梯度先验^[26]、低秩先验^[27-28]、暗通道先验^[29]。此外，各种提高盲去模糊效果的思想和技巧也被大量提出。Cho和Lee^[30]通过对中间图像双边滤波和冲击滤波来抑制平坦区域的噪声和增强边缘，提高了盲去模糊中PSF估计的精度。Xu和Jia^[31]指出图像边缘并不总是有利于PSF估计，其中比PSF尺度小的边缘会增大PSF估计误差。这些抑制中间图像噪声、增强中间图像显著边缘的思想对盲去模糊研究产生了重要影响，不仅使得各种技巧广泛显式地应用在盲去模糊中^{[24, 27]}，也促进研究者寻找隐式的非自然稀疏先验^{[22-23, 26, 29]}。

虽然盲解卷积最具简单灵活性，但由于太阳图像具有精细复杂的纹理结构、较强的噪声以及非等晕等特点，太阳图像盲解卷积具有一定的挑战性。盲解卷积的性能依赖于图像中的大尺度边缘，虽然从整体上看太阳图像具有较多的大尺度边缘，但非等晕性使得太阳图像中的许多大尺度边缘结构超出了等晕区而难以有效发挥作用。而在太阳图像局部的等晕区内，小尺度的精细纹理经过模糊后淹没混合成为大量的平坦区域，虽然较大尺度边缘经过模糊后能够得到保留，但在等晕区内可能极少而不足以估计精确的PSF。另外，太阳图像中噪声的存在也会严重影响算法的收敛，加剧盲解卷积的困难。

本文将盲去模糊领域的图像和PSF先验引入太阳图像重建中，提出了基于二阶广义总变分的空变多帧盲解卷积算法。首先利用交替最小化和半二次分裂方法求解模型，形成多尺度框架下的空不变多帧盲解卷积算法，以提高重建过程的稳定性和对噪声的鲁棒性，复原更多的纹理细节。针对非等晕的大视场太阳图像的PSF空变特性，本文利用图像重叠分块、子块序列盲解卷积与子块加权拼接的方案，减小图像非等晕导致的重建误差。

2　基于二阶广义总变分的多帧盲解卷积

2.1　多帧盲解卷积

图像的线性空不变退化模型可表示为

1g=f*h+n,

其中：g为观测到的退化图像，f为原始目标图像，h为PSF，n为加性噪声，﹡为卷积符号。当观测目标保持不变，在时间上连续成像或在空间上多通道成像时，则可由式(1)得到多帧图像退化模型为

2gi=f*hi+ni,

其中：gi 、hi 和ni 分别为第i帧成像时的退化图像、PSF和噪声，1⩽i⩽K 。

由于未知量个数大于方程个数，解空间无穷大，因此盲解卷积具有病态性。对于盲解卷积这类病态性反问题，最流行且有效的手段是基于最大后验概率，建立最优化模型：

3minf,{hi}∑i=1K‖f*hi−gi‖22+αKJ1(f)+∑i=1KβJ2(hi),

其中：J1(f) 为图像正则项，J2(h) 为PSF正则项，α和β为权重参数。图像正则项对于盲解卷积起到关键的作用，因此大部分的盲解卷积算法都是围绕着图像正则项设计。PSF正则项往往较为简单，但对于算法收敛和避免平凡解(无模糊解)也起到重要作用。依据图像先验特征，合理引入并结合图像和PSF正则项，可以缩小解空间，提高盲解卷积过程的稳定性、收敛性。

2.2　模型与算法

图像正则项中应用较早的是L2范数梯度正则化，它的形式最为简单，然而其会对梯度较大的边缘过度惩罚，导致复原图像过度平滑。相比于L2正则化，总变分正则化即TV正则化能够较好地保持边缘，广泛应用于图像去噪、去模糊等领域^{[14, 32-33]}。TV正则化提高了图像复原的效果，然而TV允许分片常数解，容易导致复原图像出现阶梯效应。为了在保持图像的边缘和细节的同时抑制阶梯效应，Bredies^[34]等人提出了广义总变分(total generalized variation, TGV)。作为TV的推广，TGV能够有效地逼近任意阶的多项式函数，例如分片常数、分片仿射函数等^[35]。随着TGV阶数的提高，相比于计算复杂度的提高程度，图像复原精度提升较小，因此常用的是二阶TGV，即：

4TGVα2(f)=α0‖∇f−v‖1+α1‖ε(v)‖1,

其中：

∇为梯度算子,∇f=[∇xf;∇yf],ε(v)=[∇xvx;12(∇xvy+∇yvx);12(∇xvy+∇yvx);∇yvy].

TGVα2 基于自然图像特征，能复原更多图像细节，适用于非盲解卷积而并不适用于盲解卷积。为了将TGVα2 引入盲解卷积，避免平凡解，Shao^[36]等人结合L0范数的非自然稀疏特性，采取了一种简单有效的改进策略，将原始的L1-L1TGVα2 变为L0-L1TGVα2 ，即将式(4)中‖∇f−v‖1 变为‖∇f−v‖0 。

综上所述，本文建立基于L0-L1 TGVα2 的多帧盲解卷积模型：

5minf,v,{hi}∑i=1K‖f*hi−gi‖22+α0K‖∇f−v‖0+α1K‖ε(v)‖1+∑i=1Kβ0‖hi‖22+∑i=1Kβ1‖∇hi‖22,

其中：简化的ε(v)=[∇xvx;∇xvy+∇yvx;∇yvy] 。除了图像正则化项外，模型(5)中还包含PSF正则化项，利用‖hi‖22 和‖∇hi‖22 来分别保持PSF的空域稀疏性和平滑性，提高PSF估计的精确性^[37]。模型(5)的优化求解可采用交替最小化方法，将其分解为图像子模型和PSF子模型，进行交替迭代优化。对应的图像子模型和PSF子模型分别为

6minf,v∑i=1K‖f*hi−gi‖22+α0K‖∇f−v‖0+α1K‖ε(v)‖1,

7min{hi}∑i=1K‖f*hi−gi‖22+∑i=1Kβ0‖hi‖22+∑i=1Kβ1‖∇hi‖22.

2.2.1　图像子问题

图像子模 型(6)可利用半二次分裂法(half-quadratic splitting, HQS)求解^{[26, 29, 38]}。引入辅助变量A、B、C、D、E，令：

A=∇xf−vx,B=∇yf−vy,C=∇xvx,D=∇yvy,E=∇xvy+∇yvx,

则模型(6)转化为

8minf,v,A,B,C,D,E∑i=1K‖f*hi−gi‖22+α0K‖A‖0+α0K‖B‖0+α1K‖C‖1+α1K‖D‖1+α1K‖E‖1+λ1K‖∇xf−vx−A‖22+λ1K‖∇yf−vy−B‖22+λ2K‖∇xvx−C‖22+λ2K‖∇yvy−D‖22+λ2K‖∇xvy+∇yvx−E‖22,

其中：λ1 和λ2 是可更新的惩罚参数。当λ1 和λ2 足够大时，模型(8)的解将趋近模型(6)的解。为简化模型(8)，令λ=λ1=λ2 ，得到：

9minf,v,A,B,C,D,E∑i=1K‖f*hi−gi‖22+α0K‖A‖0+α0K‖B‖0+α1K‖C‖1+α1K‖D‖1+α1K‖E‖1+λK‖∇xf−vx−A‖22+λK‖∇yf−vy−B‖22+λK‖∇xvx−C‖22+λK‖∇yvy−D‖22+λK‖∇xvy+∇yvx−E‖22.

由模型(9)，固定vx 、vy 、A、B、C、D、E时，可得到关于f的子模型，即：

10minf∑i=1K‖f*hi−gi‖22+λK‖∇xf−vx−A‖22+λK‖∇yf−vy−B‖22.

模型(10)的求解可令导数为0，并结合快速傅里叶变换(FFT)得到：

11f=F−1⋅(∑i=1KF(hi)¯F(gi)+λKF(∇x)¯F(vx+A)+λKF(∇y)¯F(vy+B)∑i=1KF(hi)¯F(hi)+λKF(∇x)¯F(∇x)+λKF(∇y)¯F(∇y)),

其中：F 和F−1 分别表示傅里叶正变换和逆变换。

由模型(9)，固定f、A、B、C、D、E时，可得到关于vx 和vy 的子模型，即：

12minvx,vy‖∇xf−vx−A‖22+‖∇yf−vy−B‖22+‖∇xvx−C‖22+‖∇yvy−D‖22+‖∇xvy+∇yvx−E‖22.

将式(12)分别对vx 和vy 求导，由线性方程组可得模型的解为

13vx=F−1[F(℧)F(Δ1)−F(Λ)F(Δ2)F(℧)F(℧)−F(ΛT)F(Λ)],

14vy=F−1[F(℧)F(Δ2)−F(ΛT)F(Δ1)F(℧)F(℧)−F(ΛT)F(Λ)],

其中：

Λ=∇yT∇x�℧=I+∇xT∇x+∇yT∇y,Δ1=∇xTC+∇yTE+∇xf−A,Δ2=∇yTD+∇xTE+∇yf−B.

与上面类似，固定不相关变量，可分别得到关于A、B、C、D、E的子模型，并利用硬阈值和软阈值方法，得到对应的解分别为

15A={∇xf−vx,|∇xf−vx|⩾α0λ0,otherwise,

16B={∇yf−vy,|∇yf−vy|⩾α0λ0,otherwise,

17C=sign(∇xvx)max(|∇xvx)−α12λ,0|,

18D=sign(∇yvy)max(|∇yvy)−α12λ,0|,

19E=sign(∇xvy+∇yvx)max(|∇xvy+∇yvx)−α12λ,0|.

2.2.2　PSF子问题

对于PSF子模型(7)，为提高PSF估计精度，常将其转化为梯度域PSF子模型^{[26, 29-31]}，即

20min{hi}∑i=1K‖∇f*hi−∇gi‖22+∑i=1Kβ0‖hi‖22+∑i=1Kβ1‖∇hi‖22,

令导数为0，并结合FFT，可得：

21hi=F−1[F(∇xf)¯F(∇xgi)+F(∇yf)¯F(∇ygi)F(∇xf)¯F(∇xf)+F(∇yf)¯F(∇yf)+β0+β1F(∇x)¯F(∇x)+β1F(∇y)¯F(∇y)].

由于hi 满足元素非负且‖hi‖1=1 ，所以求得hi 后需要将其中的负值置0以及规范化，即将hi 投影到约束集{hi⩾0,‖hi‖1=1} 。Perrone和Favaro^[33]分析了此步骤对盲解卷积算法收敛的重要性。

2.2.3　图像非盲解卷积

估计出PSF后，就可以对图像进行非盲解卷积求解。由于自然图像的梯度统计服从重尾分布，Krishnan和Fergus^[39]提出了基于超拉普拉斯先验的快速非盲解卷积模型：

22minfγ2‖f*h−g‖22+‖∇f‖pp,

其中：γ 为权重参数，0.5⩽p⩽0.8 。由于本文盲解卷积模型输入是多帧的退化图像序列，所以将模型(22)拓展成多帧版本，得到：

23minf∑i=1Kγ2‖f*hi−gi‖22+K‖∇f‖pp,

其中：i为帧号，K为总帧数。模型(23)的求解过程与文献[39]类似，此处不再做过多描述。

2.3　太阳图像重建流程

需要补充说明的是，观测得到的退化图像视野有限，即边界外的内容未知，然而退化图像视野内边界附近的像素生成与原始目标图像边界外的内容有关，所以观测所得的退化图像的周期延拓并不满足傅里叶变换求解图像的周期性条件假设，这将导致盲/非盲解卷积后的图像边界振铃效应^[40]。为了减弱边界振铃，需要预先对输入的退化图像进行边界处理。本文采用的是文献[40]中的图像边界处理算法。

综上所述，在多尺度框架下，基于二阶广义总变分的多帧盲解卷积算法如下：

输入：退化图像集合{gi} ，权重参数α0 、α1 、β0 、β1

输出：重建图像f

1) 设置参数：最大尺度层L、每层最大迭代次数T；

2) for l=L→1

{gi} 下采样到尺度层l：{gi(l)} ；

{gi(l)} 边界处理；

for t=1→T

根据图像子模型(6)计算得到f(l) ；

根据PSF子模型(7)计算得到{hi(l)} ；

将{hi(l)} 投影到约束集{hi(l)⩾0,‖hi(l)‖1=1} ；

end

　end

3) 求解非盲解卷积模型(23)得到重建图像f。

本文以上提出的多帧盲解卷积算法是基于PSF空不变的，然而由于大气湍流的非等晕效应，大视场图像的PSF通常是空变的，因此对大视场图像直接进行空不变盲解卷积必然导致不同程度的误差。一种常用的解决方案是对大视场图像进行重叠分块，从而使每一个子块看作等晕的或近似等晕的，即每一个子块内PSF是空不变的，这样对子块序列进行空不变盲解卷积就可以减小误差。依次对每个子块序列进行盲解卷积，全部子块序列完成盲解卷积后即可将复原子块加权拼接成最终的大视场重建图像^[41]。

由于在大气湍流影响下，光学望远镜观测的图像从微观上来看像素存在不规则的偏移、抖动现象，从宏观上来看单帧图像表现出扭曲畸变，而序列图像则表现出抖动畸变^[42]。从波前的角度来说，在等晕区内大气湍流导致了波前倾斜像差，使目标成像产生刚性的整体偏移，而整个非等晕的大视场图像则呈现出复杂的非刚性扭曲畸变。本文上述算法的处理对象是多帧的子块，这时找到最佳匹配的子块对重建效果有较大的影响。本文通过限定子块大小为接近等晕区大小，并利用简单有效的最大归一化互相关图像配准方法实现找出最佳匹配子块序列，有效降低子块间的偏移效应(低阶像差)，从而提高重建效果。其它高阶像差等效应则被包含在需要估计的PSF中。

综上所述，太阳图像重建流程如图1所示。

图 1. 太阳图像重建流程图

Fig. 1. Flow chart of the solar image reconstruction

下载图片 查看所有图片

3　实验结果与分析

3.1　图像重建结果

本文利用中国云南抚仙湖太阳观测站的1 m新真空望远镜(new vacuum solar telescope, NVST)于2016年10月7日[UT] 5:35:55前后观测的TiO波段太阳图像进行重建实验，活动区域为NOAA 12599，曝光时间为1.5 ms，图像大小为1792 × 1792，像元比例尺为0.0345″。由于实验所用的图像是大视场图像，利用上述介绍的重叠分块和加权拼接方法对其实现重建，子块大小设为256 × 256，重叠度为50%。

经过大量现有的盲解卷积或盲去模糊算法的太阳图像重建实验，本文选择表现较好的文献[33]、[17]中的空不变盲解卷积算法作为对比，分别简记为TVBD、OBD。其中TVBD是单帧盲解卷积算法，以经典的总变分作为图像正则项；OBD是多帧盲解卷积算法，常应用于天文图像处理。将本文基于二阶广义总变分的多帧盲解卷积算法简记为M-TGV，对应的单帧版本简记为S-TGV。

在S-TGV、M-TGV中，多尺度层数均设为4层，每层图像和PSF交替迭代次数均设为10，PSF大小均设为85 × 85，α0 和α1 分别设为0.5和0.05，β0 和β1 均设为10。另外，根据对比算法的特点和实验的效果，TVBD、OBD的PSF大小分别设为65 × 65、80 × 80，OBD中PSF迭代次数设为1000，其余参数均为算法默认参数。OBD和M-TGV两种多帧盲解卷积算法的输入帧数均设为5帧。最后，除OBD外的其余算法，在非盲解卷中均采用本文2.2.3小节中介绍的算法，其中γ 均设为3 × 10³，p均设为2/3。

大视场图像的重建结果如图2所示，其中图2(a)为OBD和M-TGV两种多帧盲解卷积算法的输入图像序列中总体质量最好的一帧，也是TVBD、S-TGV两种单帧盲解卷积算法的输入图像。图2(a)中矩形框标出的米粒1、米粒2、半影、本影区域的重建细节放大对比分别如图3、图4、图5、图6所示。图2(b)~2(e)分别为TVBD、S-TGV、OBD和M-TGV重建图像。图2(f)为斑点重建图像，本文将其看作参考图像，对盲解卷积算法重建图像进行质量评价。米粒1、米粒2、半影、本影区域和整帧重建图像的峰值信噪比(PSNR)、结构相似性(SSIM)评价结果如表1所示。

图 2. 大视场图像重建结果。(a) 输入序列中总体质量最好帧；(b) TVBD；(c) S-TGV；(d) OBD；(e) M-TGV；(f) 斑点重建

Fig. 2. Reconstruction results of the wide field-of-view image. (a) The image with the best overall quality in the input sequence; (b) TVBD; (c) S-TGV; (d) OBD; (e) M-TGV; (f) Speckle reconstruction

下载图片 查看所有图片

图 3. 大视场图像子区域重建结果。(a) 输入序列中总体质量最好帧；(b) TVBD；(c) S-TGV；(d) OBD；(e) M-TGV；(f) 斑点重建

Fig. 3. Subregion reconstruction results of the wide field-of-view image. (a) The image with the best overall quality in the input sequence; (b) TVBD; (c) S-TGV; (d) OBD; (e) M-TGV; (f) Speckle reconstruction

下载图片 查看所有图片

图 4. 不同输入帧数的重建结果。(a) 输入序列的其中1帧；(b) ~ (g) 分别对应输入1、2、3、5、10、20帧的重建结果；(h) 斑点重建

Fig. 4. Reconstruction results of different input frames. (a) One frame in the input sequence; (b) ~ (g) The reconstruction results corresponding to input 1, 2, 3, 5, 10 and 20 frames respectively (h) Speckle reconstruction

下载图片 查看所有图片

图 5. 图像正则项的有效性。(a) 输入序列的其中一帧；(b) 无正则项；(c) M-TGV；(d) 斑点重建

Fig. 5. The effectiveness of image regularization. (a) One frame in the input sequence; (b) Without regularization; (c) M-TGV; (d) Speckle reconstruction

下载图片 查看所有图片

图 6. 不同K和β₁的重建图像和估计的PSF。(a) K=1, β₁=0；(b) K=1, β₁=10；(c) K=5, β₁=0；(d) K=5, β₁=10

Fig. 6. Reconstructed images and estimated PSF with different K and β₁ values. (a) K=1, β₁=0; (b) K=1, β₁=10; (c) K=5, β₁=0; (d) K=5, β₁=10

下载图片 查看所有图片

在单帧盲解卷积算法中，由图2、图3可以看到TVBD和S-TGV对图像部分半影区域质量改善明显，然而同时也在米粒等区域导致了较强的伪影。从表1可以看到，相比于输入图像，TVBD重建图像总体PSNR和SSIM均有所下降，而S-TGV重建图像虽然总体PSNR下降，但SSIM却略有提升。实验结果表明本文提出的多帧盲解卷积算法即使是在单帧情况下也具有一定图像复原效果。

在多帧盲解卷积算法中，由图2、图3可以看到OBD重建图像仍然非常模糊，包含许多噪声，而M-TGV复原效果明显，米粒、半影等区域的边缘锐度提升，抑制了噪声，显现出了更多的细节。从表1可以看到，M-TGV重建图像整体和放大显示的四个子区域的PSNR和SSIM均高于OBD。此外，不论是重建图像的视觉效果，还是定量评价，M-TGV均好于S-TGV，表明了多帧的有效性。

3.2　分析与讨论

3.2.1　多帧的有效性

在本文3.1节已表明了多帧的有效性，这里再给出输入图像帧数对M-TGV算法重建结果的影响，如图4所示，其中图4(a)为输入序列中的其中1帧，图4(b)~4(g)分别为输入1、2、3、5、10、20帧时的重建图像，图4(h)为斑点重建图像。可以发现当输入帧数在10帧之内时，随着帧数增多，重建图像质量明显改善。而当输入帧数大于10帧时，随着帧数增多，重建图像质量提升较小。

本文认为产生上述现象可能有两方面的原因。第一方面是算法本身性能的限制。当输入图像序列大小在一定范围内时，随着帧数的增多，帧间高频信息得到互补，在算法迭代过程中减弱了输入图像中的噪声和PSF不精确的影响，能够有效抑制重建图像中的噪声和伪影，图像主要的结构、部分细节等得到更好地复原。然而当输入帧数继续增多后，此时影响重建图像质量的就不再是噪声和伪影了，而是图像的精细纹理，然而目前的盲解卷积算法难以对其实现复原。因此，当输入图像序列大小超过一定范围时，随着帧数增多，重建图像质量提升较小。第二方面，当输入图像序列增大时，其中部分强噪声、质量较差的观测图像也被包含在内，会对整体的图像重建造成一定程度的影响，也就是说重建图像的质量不仅与输入图像的帧数有一定关系，也与输入图像的质量有关系。

3.2.2　图像正则项的有效性

正则化是克服图像盲解卷积病态性的常用手段，本文的多帧盲解卷积模型以L0-L1 TGVα2 作为图像正则项，如图5显示了图像正则项的有效性。当无图像正则项时，重建图像完全失真，辨别不出任何信息，而当加入了L0-L1TGVα2 正则项后，重建图像复原出了许多纹理细节。

3.2.3　PSF正则项的有效性

PSF的空域正则化项‖hi‖22 在理论上能够避免平凡解，因此被普遍应用于盲解卷积模型中。而为了增强对噪声的鲁棒性，提高PSF估计精度，本文在提出的多帧盲解卷积模型中增加了PSF的梯度域正则化项‖∇hi‖22 。如图6所示为其它参数不变，在不同输入图像帧数K和PSF梯度域正则化项权重参数β1 下的重建图像和估计的PSF。

从图6可以看到，在输入帧数少时，仅使用空域正则化项‖hi‖22 的模型对噪声比较敏感，重建图像完全失真，估计的PSF中也全是分布散乱的噪点。而增加了梯度域正则化项‖∇hi‖22 后，对噪声的鲁棒性有明显提高，虽然此时重建图像仍然包含大量伪影，但相比于不使用‖∇hi‖22 时图像质量有明显改善。当输入帧数较多时，增加‖∇hi‖22 前后的模型重建图像和估计的PSF非常接近，原因是多帧模型本身具有较高的噪声鲁棒性，即使没有将‖∇hi‖22 作为PSF正则项也能得到较好的重建结果。综上分析，本文在模型中增加的PSF梯度域正则项能够在输入帧数较少时发挥有效的作用。

3.2.4　算法收敛性

对图3所示的米粒2区域的子块序列进行重建，其中各尺度层的M-TGV算法目标函数的迭代曲线如图7所示。从图中可以看到，M-TGV具有较好的收敛性，开始保持着较快的收敛速度，迭代几次后收敛速度趋向于平缓。

图 7. 目标函数迭代曲线。(a) 尺度层4；(b) 尺度层3；(c) 尺度层2；(d) 尺度层1

Fig. 7. Iteration curves of the objective function. (a) Fourth scale; (b) Third scale; (c) Second scale; (d) First scale

下载图片 查看所有图片

3.2.5　算法复杂度

一般来说，大视场太阳图像的尺寸较大，且大视场太阳图像需要进行图像分块重建，如1792 × 1792的图像以256 × 256大小进行重叠分块后共有169个子块。若子块重建的空不变盲解卷积算法复杂度太高，则整个图像重建过程将会非常耗时。因此，太阳图像重建对空不变盲解卷积算法的复杂度有一定要求。

本文提出的算法复杂度为O(NlogN)。在Intel(R) Core(TM) i5-1135G7 @ 2.40 GHz 2.40 GHz处理器，16 GB 内存，Windows10 64位操作系统，MATLAB 2016运行环境下，对TVBD、S-TGV两种单帧盲解卷积算法(输入1帧256 × 256大小图像)和OBD、M-TGV两种多帧盲解卷积算法(输入5帧256 × 256大小图像)的运行时间进行测试，各参数均与前面保持一致，结果如表2所示。

从表2可以看到，TVBD十分耗时，OBD也较为耗时，而本文算法不仅在重建精度上超过了对比算法，在重建速率上也取得了较好的效果。

3.2.6　盲解卷积算法的局限性

虽然本文基于二阶广义总变分的多帧盲解卷积算法在太阳图像重建上取得了较好的效果，然而算法仍有其局限性。首先，大视场太阳AO图像具有等晕区大小不一致的非等晕性，本文算法在重建时使用了一个经验的统一大小对其进行分块；若子块相比等晕区过大，将导致图像重建误差；若子块相比等晕区过小，则在子块内图像的大尺度边缘等结构特征非常少，这并不利于PSF估计，并且子块太小也会在后续图像拼接时导致明显的边缘效应；后续研究可使用可变分块进一步提高重建效果。此外，太阳图像的结构纹理十分复杂，包含了大量的尺度小于PSF的精细纹理，经过退化后这些精细纹理被淹没，若估计的PSF存在误差，则这些精细纹理将难以被复原。并且在通常情况下，退化的太阳图像中也包含较强的噪声，而现有的绝大多数盲解卷积算法对噪声较为敏感，不能估计得到精确的PSF。即使是在PSF精确已知的情况下，虽然可消除PSF误差导致的图像伪影，但观测图像中的噪声仍然存在，受噪声影响，通过非盲解卷积重建的图像质量也并不能达到完全理想的状态。

4　结论

本文针对常用的太阳AO图像事后重建方法中最具灵活性的盲解卷积展开研究。针对太阳图像复杂的纹理结构、较强的噪声以及非等晕等特点，本文首先建立了基于二阶广义总变分的空不变多帧盲解卷积算法，然后利用重叠分块与加权拼接将其拓展为适用于非等晕大视场太阳图像的空变盲解卷积算法。通过NVST观测获得的真实太阳图像重建实验，表明本文算法在主观视觉效果和客观指标上均具有较好的图像重建效果，二阶广义总变分正则化和多帧能够提高图像重建的稳定性和可靠性。然而，盲解卷积仍然难以复原太阳图像中的部分精细纹理，后续需要采取更多的手段提高盲解卷积的性能。