第40卷第2期
2021年4月红外与毫米波学报
J.Infrared Millim.Waves Vol.40,No.2 April,2021
基于L1范数优化模型的遥感图像条纹去除方法
李凯1,2,3,李文力1,2,3,韩昌佩1,2*
(1.中国科学院上海技术物理研究所,上海200083;
2.中国科学院红外探测与成像技术重点实验室,上海200083;
3.中国科学院大学,北京100049)
摘要:从条纹噪声的结构属性进行分析,通过分离出条纹成分来实现去条纹的目的。在优化模型中,基于L1范数的正则化表示条纹的全局稀疏特性;基于差分的约束条件用于描述条纹方向上的平滑度和条纹垂直方向上的不连续性。为了更好地保护图像的细节信息,在条纹垂直方向的约束上引入了边缘权重因子,最后通过交替方向乘子法(ADMM)对所提模型进行求解和优化。用多通道扫描辐射计(AGRI)获取的在
轨数据对算法进行了验证并与典型方法进行了比较,结果表明,消除条纹噪声的同时更好地保留了细节信息,并且呈现出较好的定性和定量结果。
关键词:L1稀疏优化模型;图像去条纹;边缘权重因子;交替方向乘子法;AGRI图像
中图分类号:TP751.1文献标识码:A
The method based on L1norm optimization model
for stripe noise removal of remote sensing image
LI Kai1,2,3,LI Wen-Li1,2,3,HAN Chang-Pei1,2*
(1.Shanghai Institute of Technical Physics,Chinese Academy of Sciences,Shanghai200083,China;
2.Key Laboratory of Infrared Detection and Imaging Technology,Chinese Academy of Sciences,Shanghai
200083,China;
3.University of Chinese Academy of Sciences,Beijing100049,China)
Abstract:Structural properties of stripe noise are analyzed and the purpose of destriping is achieved by separating the stripe components.In the proposed optimization model,the L1-norm-based is used to describe global sparse property of stripes.In addition,difference-based constraints are adopted to describe the smoothness and discontinuity in the along-stripe and across-stripe directions,respectively.In order to better protect the detailed information of an image,an edge weighting factor is introduced in the constraints of across-stripe direction.Finally,the proposed model is solved and op‐timized by the alternating direction method of multipliers(ADMM).The algorithm is verified by the in-orbit images ob‐tained by Advanced Geosynchronous Radiation Imager(AGRI)in comparison with typical destriping methods.Experi‐mental results show that the proposed algorithm completely eliminates the stripe noise and preserves more details,which shows better qualitative and quantitative result.
Key words:L1sparse optimization model,image destriping,edge weighting factor,alternating direction method of multipliers,AGRI image
PACS:42.79.Pw,05.40.Ca,42.68.Sq
引言
条纹噪声主要是由多元探测器响应不一致而引起的[1],通常出现在摆扫式和推扫式遥感成像系统中[2-3]。条纹噪声不仅影响成像质量,而且还会严重影响后续的遥感图像应用,例如分类、目标检测、定量应用等。因此,为了改善图像质量,在保留细节结构信息的同时,消除条纹噪声至关重要。
文章编号:1001-9014(2021)02-0272-12DOI:10.11972/j.issn.1001-9014.2021.02.018收稿日期:2020⁃04⁃26,修回日期:2020⁃06⁃08Received date:2020⁃04⁃26,Revised date:2020⁃06⁃08
基金项目:中国科学院上海技术物理研究所创新专项(CX-208)
Foundation items:Supported by Innovative Special Foundation of Shanghai Institute of Technical Physics(CX-208)
作者简介(Biography):李凯(1992-),男,河南周口人,博士研究生,主要从事遥感图像处理,模式识别方面的研究.E-mail:likai_sitp@163 *通讯作者(Corresponding author):E-mail:changpei_han@mail.sitp.ac
2期李凯等:基于L1范数优化模型的遥感图像条纹去除方法
目前,去除条纹噪声的方法大致可分为三类:基于滤波的方法、基于统计的方法和基于优化的方法,本文所提方法属于基于优化的方法。
基于滤波的方法通常可通过各种滤波器来获得去条纹结果,例如傅里叶滤波器[4]、小波分析[5]和小波-傅里叶组合滤波器[6]。这些方法通过分析和滤除变换域中的特定条纹成分来实现去条纹的目的。由于条纹的周期性特性可以在功率谱分析中轻松识别,并可以从变换域的信号中提取出来,因此基于滤波的方法通常对周期性条纹噪声能够输出较好的结果。受截止频率的影响,这些方法常会滤除比条纹噪声更多的信息,使滤除结果出现明显的模糊和振铃伪像。
基于统计的方法主要依赖于各传感器具有相同的灰度统计特性,例如矩匹配[7]和直方图匹配[8]。矩匹配假设每个传感器的均值和标准差是一致的,而直方图匹配是通过将未校准信号的直方图与参考信号进行匹配来消除条纹噪声。尽管这些方法效率相对较高,但是当条纹为非线性或不规则条纹时,并不能获得好的去除效果。
基于优化的方法通过极小化能量泛函,估计去条纹图像,因此这类方法的关键在于能量泛函的构造。文献[9]中假设正常像素值(DN值)符合Huber-
Markov分布,利用最大后验框架去除条纹。文献[10]通过考虑条纹噪声的方向性结构信息,设计了一种较为复杂的去条纹方法。为了能够同时去除随机噪声和条纹噪声,文献[11]提出一种结合了单向总
变分和稀疏表示的联合模型。由于优化模型设计的灵活性,这类方法是目前较为有效的去条纹方法。
但是,现有大多数优化方法都是通过去除噪声的角度来实现的,而没有考虑条纹的典型特性例如方向和结构特性,这将导致许多图像细节可能与条纹一起被去除[12]。文献[13]中使用基于L0范数的正则化来描述条纹的全局稀疏性来估计出条纹。但是L0范数是非凸的,其算法通常不具有收敛性[14]。因此本文设计了一个基于L1的正则化模型来描述条纹的稀疏特征,由于条纹噪声具有明确的方向特征,本文优化模型中还引入了沿条纹方向和垂直条纹方向的差分先验约束。同时为了进一步保护结构信息和防止过度平滑,本文在垂直条纹方向的约束上引入了边缘权重因子,使得在细节丰富区域上的惩罚因子较小,在平坦区域上的惩罚因子较大。为了有效地求解该凸优化模型,本文提出了一种基于交替方向乘子法(ADMM)的计算方法。本文使用多通道扫描辐射计(AGRI)在轨数据对所提方法进行实验验证,结果表明,与典型去条纹方法相比,所提方法在细节信息保护和条纹去除上呈现出较好的定性和定量结果。
1条纹去除算法
1.1问题描述
假设遥感图像f∈L2(Ω),条纹噪声可看作为加性噪声[10],则图像的退化模型为:
f(x,y)=u(x,y)+s(x,y),(1)这里,f(x,y)为遥感仪器所获取的实际观测数据,u(x,y)为探测器响应一致情况下的理想数据,s(x,y)是条纹噪声。x轴和y轴是分别沿遥感图像的水平方向和垂直方向,坐标原点位于图像的左上角。在本文中,算法主要部分是根据条纹的属性从退化图像f 中提取条纹分量s,其所提模型求解框图如图1所示。为了便于讨论算法的数值计算,可将式(1)写成矩阵向量形式,如式(2)所示:
f=u+s.(2)
1.2条纹成分的估计模型
这里需要充分利用条纹噪声的属性,并以适当的正则化约束条件来进行描述,用以构成s的估计模型。
1.2.1沿条纹方向的平滑度
遥感图像的条纹都是以逐行(x轴)或逐列(y轴)出现的,即沿着扫描方向的[2-3]。不失一般性,本文将条纹视作逐行情况来建立优化模型。通常,从沿着条纹的方向上看,单个条纹内的相邻像素具有相对较小的变化,这意味着条纹成分(图像)沿条纹方向上具有良好的平滑度,
因此条纹图像的水平梯度非图1本文算法框图
Fig.1The framework of the proposed model
273
40卷
红外与毫米波学报常稀疏。为了在水平方向上能够较好地保留条纹梯度,本文在梯度域上采用稀疏正则化对其进行约束,这里使用了具有凸函数性质的L1范数作为稀疏正则化,如下式所示:
R 1(s )=
∇x s
1
,(3)
式中∇x 表示水平方向上的一阶差分算子。1.2.2
全局稀疏性
近年来,稀疏表示模型广泛用于图像恢复中,例如去模糊[15]、去噪[16]和重建[17],并取得了很好的效果。在很多遥感图像中,条纹噪声的比例相对较低。基于这样一个事实,条纹图像应具有全局稀疏的属性。因此,本文采用条纹图像的L1范数来描述条纹的全局稀疏性,如式(4)所示:R 2(s )=  s 1
.(4)1.2.3
理想图像的局部连续性
理想图像u 通常在y 轴方向(垂直条纹方向)上可看作为连续的,在u 逐行添加s 时,会破坏u 的连续性。因此,为了保证其连续性应使∇y u 为一个较小的值。根据式(2)的假设关系,可以用式(5)的基于L1范数的正则化来描述理想图像的局部连续性:
R 3(s )=
∇y f -∇y s
1
,(5)
式中,
∇y u =∇y f -∇y s ,∇y 表示垂直方向上的一阶差分算子。值得注意的是R 3(s )正是单向总变分模型(UTV )用于平滑条纹的正则约束项[10]。
综上所述,结合三个正则化项可得出用于遥感图像去条纹的L1稀疏优化模型:s =arg min s
{
正则化描述正确的是
∇x
s 1
+λ1  s 1+λ2
∇y f -∇y s
1
}
,(6)
式中,
λ1与λ2是两个正的正则化参数用于平衡上述三项。1.2.4
边缘权重因子
最小化能量方程(6)还不能得到令人满意的去
条纹效果。考虑到图像中存在不同的特征,为整个图像设置相同的正则化参数λ2是不合理的,即分配给平坦区域与图像边缘区域的权重不应该相同,否则去条纹效果会变得过于平滑,并且细节信息将会丢失。为了克服这个问题,本文引入边缘权重因子分配给R 3(s )项,其中,高权重值用于条纹噪声很重的区域(条纹非稀疏区域),对于包含图像纹理和细节的区域则选择较小的权值。因此,引入边缘权重因子后的优化模型为式(7)所示:
s =arg min s
{
∇x
s 1
+λ1  s 1+λ2W f
∇y f -∇y s
1
}
,(7)
式中,
W f 为边缘权重因子。根据梯度域优化[18]和边缘感知加权理论[19],本
文定义了描述图像边缘的显式加权因子Φf (x ,y ),它
是由原始图像f 中所有像素的3*3窗口和r *r 窗口的局部方差来计算的:
Φf (x ,y )=σf ,3(x ,y )σf ,r (x ,y )
,(8)
式中,σf ,3(x ,y )和σf ,r (x ,y )分别是图像f 以像素(x ,y )为中心所计算的3*3窗口标准差和r *r 窗口标准差;M*N 是图像f 的尺寸大小,即像素总数。
根据Φf (x ,y )的定义可知,它的作用是衡量给定
像素(x ,
y )相对于整个图像f 的重要性。因为同时使用了较大的尺度窗口r 和较小的尺度窗口3,因此它
能够有效地将边缘分离出来,增强了加权因子的性能。图2(b )显示了权重因子Φf 的图像边缘和细节描述能力,但是受条纹噪声的影响,权重因子错误地将条纹噪声看作为图像边缘。为了克服这个问题,本文首先使用垂直滤波器[20]将图像f 进行平滑得到平滑部分f g (如图2(c )所示)和高频细节部分f d (如图2(d )所示),其中f d =f -f g 。利用这两个部分重新设置权重因子可得式(9),图2(e )显示了新权重因子,可见新权重因子能够准确描述原始遥感图像的边缘和细节并不受条纹噪声的干扰。
Φf (x ,y )=σf g
,3(x ,y )σf d
,r (x ,y )
.(9)
综上,可将优化模型中的边缘权重因子W f 定义如下:
W f (x ,y )={
1,if --Φf ()x ,y <S
δ,if --Φf ()x ,y ≥S
,(10)
式中W f (x ,y )表示在R 3(s )中像素(x ,y )的边缘权重,S 是将图像边缘与平滑区域分开的阈值,δ为一个较小的正值,--Φf 为Φf 的归一化值。边缘权重图像如图2
(f )所示。
2ADMM 最优化
由于L1范数具有不可导和不可分离性质,因此
很难直接从方程(7)中解出s 。为了克服这个问题,需要到一种有效的解决方法来处理这种非平滑的凸优化模型。ADMM 算法是一种高效地用于解决凸模型的优化方法,例如在图像处理中解决基于L1范数的最小化和总变分模型[12]。因此,本文选择该方
274
2期
李凯等:基于L1范数优化模型的遥感图像条纹去除方法
法来处理方程(7)。
首先,通过引入三个辅助变量Z =∇x s ,V =s 和H =∇y f -∇y s 将无约束最小化问题(7)转换为有约
束最小化问题,如式(11)所示:
arg min s ,Z ,V ,H
{
Z 1
+λ1  V 1+λ2W f  H 1
}
Z 1+λ1  V 1+λ2W f  H 1+p T 1(∇x s -Z )+p T
2
(s -V )+p T 3
(∇y f -∇y s -H )+
β12
∇x s -Z 22+β22
s -V 2
2+β32
∇y f -∇y s -H 22,(12)式中p 1,
p 2,p 3为拉格朗日系数;β1,β2,β3为正的惩罚参数;ADMM 的每次迭代都可以分解为四个较为简单的子问题。
2.1Z -子问题
关于Z 的子问题由式(13)所示:
Z =arg min Z {
Z 1+p T 1(∇x s -Z )+β1
2
∇x s -Z 22}
.(13)方程(13)可以通过以下软阈值收缩算子计算[21]Z k +1=shrink (∇x s k
+p k 1β1,1β1
)
,(14)
式中
shrink (r ,θ)=
r
|
|r *max (|r |-θ,0).(15)
2.2V -子问题
类似于Z 的子问题,V 子问题可由(16)式求解
V =arg min V
{
λ1  V 1+p T 2(s -V )+β22
s -V 2
2}
.(16)因此
V
k +1
=shrink (s k
+p k 2
β2,
λ1β2
)
.(17)
2.3H -子问题
同理,H 子问题可由(18)式求解
H =arg min H
{
λ2W f  H 1+p T 3(∇y f -∇y s -H )+β3
2
∇y f -∇y s -H 22
}
.(18)
因此H
k +1
=shrink (
∇y f -∇y s k
+p k 3β3,λ2β3
∙W f )
.(19)
2.4s -子问题
关于s 的子问题,可由(20)式所示:
s =arg min s
{
p T 1(∇x s -Z )+p T 2(s -V )+p T
3(∇y f -∇y s -H )+
β1
2
∇x s -Z 22
+
β22
s -V 2
2+β3
2
∇y f -∇y s -H 22
}
.(20)
这是一个二次最小化问题,等效于求解以下线性方程(21)(β1∇T x ∇x +β2+β3∇T y ∇y )s k +1=β1∇T x (
Z k +1
-p k 1β1)
+β2(
V k +1-p k 2β2)
+β3∇T y (
∇y f -H k +1
+p k 3β3
)
.(21)对方程(21)可以通过快速傅里叶(FFT )有效
求解。
最后,在每次迭代中对拉格朗日系数p 1,
p 2,p 3按照式(22)进行更新:ìíîïïïïp k +11=p k 1+β1()∇x s k +1-Z k +1p k +1
2=p k 2+β2()s k +1-V k +1p k +13=p k 3+β3(∇y f -∇y s k +1
-H k +1).(22)综上所述,本文利用ADMM
算法将复杂的优化
图2(a )原始遥感图像,(b )式(8)权重因子,(c )平滑部分,(d )高频部分,(e )式(9)权重因子,(f )式(10)边缘权重Fig.2
(a)The original remote sensing image,(b)weighting
factor image in Eq.8,(c)the smooth part,(d)the high frequen ‐cy part,(e)weighting factor image in Eq.9,(f)edge weighting image in Eq.10
275
40卷
红外与毫米波学报模型(7)分解为四个简单的子问题。具体来说,Z ,V 和H 子问题可以通过软阈值算子解决。对于s 子问题,可以选择FFT 进行高效求解。此外,拉格朗日系数p 1,p 2,p 3可以并行更新。因此,本文提出的去条纹
算法可以总结如表1所示。
3实验与分析
为了验证所提方法的有效性,本文进行了实际
遥感图像数据实验,并与3种典型方法:基于滤波的
方法(WFAF )[6]、基于统计的方法(SLD )[22]
和基于优化的方法(UTV )[10]进行比较。本文方法将输入参数
设置如下:正则化参数λ1=0.001,
λ2∈[0.005,0.01],惩罚参数β1=β2=β3=0.1,尺度窗口r =33,阈值S =0.1,
因子常数δ=0.2。为了便于定量评价,本文将所有实验数据进行了归一化。本文实验的计算机硬件运行平台为:Intel (R )Core (TM )i5-7300HQ CPU @2.50GHz 和8GB RAM ;软件运行平台为:MATLAB
(R2014a )。3.1
数据源
风云四号A 星(FY-4A )于2016年12月11日发射,是新一代中国地球静止气象卫星系列风云四号的第一颗卫星。其中,多通道扫描辐射计(AGRI )是FY-4A 的主要载荷之一[23]。AGRI 成像方式是多元
线列探测器并扫结合二维扫描,因此AGRI 遥感图像
条纹噪声呈现在东西扫描方向上,如图2(a )所示。AGRI 共有14个波段,光谱范围0.45~13.8µm ,虽然经过发射前辐射定标和在轨辐射定标,但某些波段仍受到条纹噪声影响严重,如水汽波段9和波段10以及长波红外波段14,其参数如表2所示。本文所用数据来源为FY-4A 上AGRI 于2018年12月9日06:00:00(Universal Time Coordinated ,UTC )生成。3.2
周期条纹噪声去除
在扫描成像中,同一个波段内,光谱一致的像元一般不受随机条纹噪声的影响,只有较浅的周期条纹。以AGRI 波段11为例,行内没有条纹(四个像元为一个扫描行),行间有周期条纹出现,如图3(a )所示。因此,本小节对波段11图像数据进行去条纹处理。
图3(b )~(e )为分别使用WFAF 、SLD 、UTV 和本文方法对图3(a )周期条纹噪声进行去除的结果。从图中可以看出,几种去条纹方法处理后,周期条纹噪声得到了不同程度的抑制。如图3(b )所示,WFAF 方法虽然能够一定程度上的去除条纹,但容易丢失细节信息,造成图像模糊。SLD 方法也仅
在条纹强度非常一致的地方去除的比较干净,其他区域仍有条纹残留。相对来说,UTV 方法呈现出了较好的去条纹效果,但是由于方向约束性太强导致图像沿垂直于
表1
去条纹算法
Table 1
The proposed destriping algorithm
1:Input :Stripe image f ,parameters λ1,
λ2,β1,β2,β3,δand S .2:Initialize :Set s 0=0,
Z 0=V 0=0,H 0=∇y f ,p 1=0,p 2=0,p 3=0,and ε=10-4.
3:Solve W f by (10)
4:While
()f -s k
-()f -s k -1
/  f -s k
>εand k <N
max
do
5:Solve Z k +1,V k +1,H k +1using a thresholding method by (14),
(17),(19)
6:Solve s k +1
using FFT by (21)
7:Update p k +11,p k +12,and p k +1
3
by (22)8:End while
9:Output :
u k +1=f -s k +1.表2
光谱参数
Table 2
Spectral parameters
No.9
101114Central Band /µm
6.25
7.10
8.50
13.50Spectral Band /µm 5.80∼6.706.90∼7.308.00∼9.00
13.20∼13.80Spatial Resolution
4km 4km 4km 4km
Number of pixels
4*14*14*14*1
Main Application
upper-level water vapor mid-level water vapor
integrated water vapor ,cloud
cloud ,water
vapor 图3AGRI 波段11子图像去条纹结果(a )原始图像,(b )WFAF 方法,(c )SLD 方法,d )UTV 方法,(e )本文方法Fig.3
Destriped results of AGRI band 11subimage (a)origi ‐
nal image,(b)WFAF,(c)SLD,(d)UTV,(e)proposed method
276

版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198,我们将在24小时内删除。