摘 要：为有效抑制图像噪声，提高电气设备红外诊断的准确性，采用基于小波系数尺度间相关性和双变量收缩函数的方法对电缆瓷套终端红外图像进行去噪.将图像进行小波分解，计算小波系数尺度间的相关系数，使用模糊c-均值聚类法对相关系数聚类，即将小波系数分为有效系数和无效系数两类.对无效小波系数直接进行置零处理，对有效小波系数使用双变量收缩函数进行处理，得到真实图像小波系数的估计值.最后，对处理得到的真实图像小波系数的估计值进行重构，便得到去噪后图像.含噪图像的去噪结果表明，运用文中方法能有效地去除红外图像中的噪声，且与使用传统软阈值方法去噪得到的图像对比，文中方法去噪后的图像信噪比更高，最小均方误差更小.
关键词：图像去噪；尺度间相关性；小波变换；相关系数；模糊c-均值；双变量收缩
红外热成像能够实现非接触性测量，具有不容易被电磁干扰、准确、安全等特点[1]，近年来已广泛应用于电缆瓷套终端状态检测与故障诊断.红外热成像仪接收电缆瓷套终端的热辐射信号，并将其转换成相应的电信号，再经过信号处理后，被测电缆瓷套终端的表面热分布以人眼可见的图像呈现出来.当电缆终端存在某些缺陷时(如瓷套式终端进水)，电缆终端瓷套表面温度分布有所变化.因此，可以通过分析红外图像来发现电缆瓷套终端的某些缺陷，从而减少事故的发生.
红外图像在生成过程中由于受周围环境噪声和仪器本身噪声影响而具有高噪声、低对比度的特点[2].因此，为了准确地对电缆瓷套终端进行缺陷及故障诊断，有必要对其红外图像进行去噪处理.
在图像去噪的研究中，比较常用的方法有自适应加权均值滤波、形态学去噪、模糊去噪、小波去噪[3-6]等.其中小波变换由于其良好的时频局部化能力，具有较好的去噪效果[7-8].在小波去噪方法中，基于小波系数阈值的滤波方法由于实现简单、重构信号误差小而得到广泛应用.Donoho在文献[9]中提出了软阈值收缩函数，该函数通过将小波系数和选定阈值比较进行收缩，其去噪效果和阈值选取密切相关.常用的阈值选取规则有VisuShrink阈值法[10]、SureShrink阈值法[11]、BayesShrink阈值法[12].软阈值收缩函数简单有效，但该函数实质是在对小波系数进行单变量Laplase分布建模时得到的，忽略了小波系数尺度间相关性.文献[13]中提出了一种基于小波系数统计数据和隐式马尔可夫模型的信号处理方法，该方法考虑了小波系数尺度间相关性，采用多个状态的混合高斯模型来描绘小波系数的分布，并取得了比传统软阈值方法更好的去噪效果.文献[14-16]中提出了一种基于双变量收缩函数的图像去噪方法，该方法采用非高斯的双变量分布模型来描述小波系数的分布，基于贝叶斯估计理论得到相应的非线性阈值函数，并取得了比传统软阈值方法和文献[13]中提出的隐式马尔可夫树法(HMT)更好的去噪效果.
文中在文献[14-16]已有研究的基础上，结合小波系数尺度间相关性，将电缆瓷套终端红外图像进行小波分解，定义和计算小波系数尺度间相关系数，并使用模糊c-均值法对相关系数进行聚类，即将小波系数划分成有效小波系数和无效小波系数两类，对无效小波系数进行置零处理，对有效小波系数采用双变量收缩函数进行处理，对处理后的小波系数进行二维小波反变换即可达到去除红外图像噪声的目的.
1 含噪图像信号小波系数的模糊分类
含噪图像信号经过二维小波变换后得到各分解尺度不同方向上的小波系数和最高分解尺度上的尺度系数.其中，有用图像信号所对应的小波系数幅值大、数量少，噪声所对应的小波系数幅值小、数量多.而尺度系数则主要集中了有用图像信号的能量.文献[17]中指出，经过小波变换之后，信号所对应的小波系数在各尺度上具有比较强的相关性，而噪声信号所对应的小波系数在尺度间缺乏这种明显的相关性.因而文中利用各尺度上小波系数的尺度间相关性来对小波系数进行分类.根据信号处理需要可以将小波系数划分为两类：有效系数和无效系数，前者对应于信号，后者对应于噪声.具体实现如下.
1.1 信号小波系数相邻尺度的相关系数
设小波分解的最大尺度为
f(m,n)表示尺度j、kd方向、位置(m,n)处含噪信号的离散小波变换.位置(m,n)处、kd方向含噪信号小波系数在尺度j与j+l间的相关系数
f(m,n)定义如下：
(1)
式中:kd=H,V,D，分别表示水平、垂直和对角线方向;l为参与相关计算的尺度数.随着分解尺度的增大，信号突变部分的宽度也会增大，而且在粗尺度上相邻的突变点会相互影响[18]，因此文中取l=2，即式(1)可表示为
(2)
1.2 含噪图像信号小波系数的模糊c-均值分类
如前所述，图像信号的小波系数在各尺度上具有较强的相关性，即 f(m,n)绝对值较大,而噪声信号的 f(m,n)绝对值较小.根据图像信号和噪声信号的小波系数在相关系数方面的差异性，可以通过在每个尺度上设定一个阈值，将相应尺度计算出来的相关系数与之比较.经过对相关系数的阈值处理后，就可以把每个尺度上的小波系数分为有效系数和无效系数两类.为了使相关系数的阈值提取更符合人们的思维方式，文中使用模糊c-均值阈值法提取小波系数相关系数的阈值.该法通过比较每个数据点对所有聚类中心隶属度进行分类，具有自适应和无监督学习特点.
模糊c-均值聚类是一种基于划分的聚类算法，是普通c均值算法的改进，其实质上就是求解使得目标函数J最小时的隶属度矩阵[uij]D×N，以及类别中心V.
(3)
式中，m∈(1,∞),为加权指数，dij表示第j(j=1,2,…,N)个数据xj与第i(i=1,2,…,D)类的聚类中心vi之间的距离，即
.模糊c-均值聚类算法具体如下：
步骤1 根据式(3)，首先确定聚类数目D，然后初始化聚类中心V及加权指数m；
步骤2 对第t次迭代，需要更新隶属度和D个聚类中心，其计算公式分别为
(4)
(5)
若目标函数J连续两次迭代的误差满足|J(t)-J(t-1)|≤ε，则迭代过程停止；若目标函数的误差不满足精度，则返回步骤2继续迭代.
迭代过程结束后，得到各数据对不同聚类中心的隶属度矩阵U=[uij]D×N.在决定数据点分类时，文中采用最大隶属度原则[19]，则如果数据隶属于某聚类中心的程度最大，其就隶属于该聚类中心代表的类别.最大隶属度原则表示为
Ai0(u0)=max{A1(u0),A2(u0),…,AD(u0)}
(6)
则认为数据点u0相应隶属于类别Ai0.
文中将相关系数划分为两类，经聚类后阈值由下式决定：
(7)
式中,
为第2类即无效小波系数不同尺度间的相关系数最大值.实际应用中，也可根据小波系数相关系数的分布规律，使用第1类即有效系数不同尺度间的相关系数最小值作为阈值.通过模糊c-均值阈值处理后，可以将各尺度不同方向上的小波系数划分为有效小波系数和无效小波系数.对无效小波系数进行置零处理，得到初步去噪后的信号，其中仍含有少量噪声.在初步去噪的基础上，再对各尺度上的有效小波系数，采用双变量收缩函数和层内局域方差估计进行处理，即可恢复出真实信号的小波系数.
2 双变量收缩和局部方差估计
2.1 单变量收缩函数
假设图像被加进高斯白噪声污染，噪声的方差为
，对含噪图像进行小波变换，则有
y=w+n
(8)
式中，y为含噪图像的小波变换系数，w为原始图像的小波变换系数，n为噪声的小波变换系数.
则原始图像的小波变换系数的最大后验估计(MAP)为
(9)
式中，pw|y为条件概率密度函数，应用贝叶斯定律，由式(9)可得
(10)
式中，pn表示噪声的概率密度函数，pw表示原始图像小波系数的概率密度函数.其中，pn为零均值的高斯白噪声分布，即
(11)
式(10)可以变换为
(12)
记
，则由式(11)和(12)得
(13)
式(13)相当于解方程
(14)
在图像去噪方面，有许多统计模型被用来描述小波系数的分布，其中拉普拉斯模型应用比较广泛[20].拉普拉斯模型的数学表达式只有标准差σ一个参数，结合MAP估计和贝叶斯定律通常可以得到解析解.因此，文中对尺度内小波系数的分布特性进行拉普拉斯建模，即
(15)
那么式(14)的解为
(16)
式中，令
的定义如下：
(17)
把
看做阈值，式(16)实际上就是软阈值收缩函数.
2.2 双变量收缩函数
和一维小波变换的小波系数一样，二维小波变换的小波系数在纵向上具有传播性，即不同尺度间同一位置上的小波系数幅值之间存在传播性，为考虑图像小波系数的这一分布特性，需要对小波系数的分布模型作如下修改:
(18)
式中，
k=(w1k,w2k)，小波系数w2k为小波系数w1k的父系数，即w2k为下一个尺度上与第k个小波系数w1k在同一位置的小波系数.为方便叙述，均去掉下标k.σ为w1和w2的方差.
因为白高斯噪声
的概率密度函数为
(19)
即得出w1的MAP估计为
y1
(20)
式(20)表示的双参数收缩函数如图1所示，y1、y2均为含噪父系数.从图1可以看出，
1和y2有关.软阈值收缩函数如图2所示，可以看出1和y2无关.
图1 双参数收缩函数图
Fig.1 Figure of bivariate shrinkage function
图2 软阈值收缩函数图
Fig.2 Figure of soft threshold shrinkage function
2.3 局部方差估计
从式(20)可以看出，为了得到w1的估计，需要求出每个小波系数的局部方差σ2和噪声方差.对于局部方差σ2，文中采用局域法求取，即利用以第k个小波系数为中心的矩形窗区域内的小波系数来估计第k个小波系数的局部方差σ2.矩形窗区域N(k)的定义如图3所示.噪声方差可以通过图像小波分解各尺度上的对角线细节系数的中值估计得到，即
HH
(21)
式中，HH为对角线方向.
从观测模型可以得到
为观测值y1和y2的局部方差.由于y1和y2的均值为0，则的经验估计为
(22)
式中，M为矩形窗区域N(k)的大小.由此得到σ的估计为
(23)
图3 领域N(k)示意图
Fig.3 Illustration of neighborhood N(k)
3 电缆瓷套终端红外检测图像去噪
3.1 电缆瓷套终端红外检测
为提高电缆线路的供电可靠性和稳定性，广州供电局每个季度对所有电缆线路的瓷套终端进行红外检测.当瓷套终端存在发热缺陷时，会在其对应的红外图像里面反映出来，通过分析红外图像的特征，就可以判断出终端发热的位置和类型.图4为对某瓷套终端进行红外检测得到的红外图像.
图4 红外图像
Fig.4 Infrared image
从图4可以看出，瓷套终端的应力锥部位存在明显的局部发热现象，其中发热区域的最高温度达到27.1 ℃，而瓷套终端底部最高温度为25.0 ℃，两者相差为2.1 ℃，远远超过了正常范围，说明瓷套终端存在发热缺陷，需要查找和排除.
通过对该终端的缺陷查找发现，该终端由于密封不好存在进水现象，其内部进水水位已经达到应力锥的顶部.当终端内部进水时，由于水的密度比里面填充硅油的密度大，所以水分会集中在终端下部，使接近水的那部分硅油的介质损耗角正切变大.由介质损耗功率公式P=U2ωCtan δ(U为施加电压，ω为施加电压角频率，C为介质等效电容，δ为介质损耗角)可知，硅油的介质损耗功率与其介质损耗角正切和施加电压平方成正比.由于瓷套终端应力锥部分的电压分布比较集中，故该区域的发热比终端其他部位明显.此外，由于终端进水，可能导致泄露电流增大，进而使发热量增加.3.2 去噪算法
依据上述原理，基于双变量收缩函数和局部方差估计去噪算法的具体步骤如下：
步骤1 将红外图像信号进行二维小波分解，得到最高分解尺度上的尺度系数和各分解尺度不同方向上的小波系数.
步骤2 用式(2)计算每个小波系数的相关系数.
步骤3 使用模糊c-均值聚类法将小波系数分为有效小波系数和无效小波系数两类.
步骤4 将每一个无效小波系数置零.对每一个有效小波系数进行如下处理:
①用式(21)计算得到噪声方差的估计值
，用式(22)和式(23)计算得到局部方差的估计值2.
②将参数估计值和2代入式(20),得到小波系数的估计值
1.
步骤5 对最高分解尺度上尺度系数和由(4)得到的小波系数估计值进行二维小波重构，就可以得到去噪后的红外图像信号.
3.3 电缆瓷套终端红外图像去噪的仿真
为方便对图像进行处理，笔者在对红外图像进行去噪处理之前先对图像进行预处理，即将红外图像的真彩图转换为灰度图.将如图4所示的电缆瓷套式终端红外图像(大小为240像素×240像素)转换为对应的灰度图，加入不同噪声水平(噪声标准方差分别为5、10、15、20)的高斯白噪声，即可得到不同信噪比的含噪图像.如图5所示是标准方差为15的含噪图像.
图5 含噪图像
Fig.5 Noisy image
文中采用前述基于小波系数尺度间相关性和双变量收缩函数的自适应去噪方法，对不同高斯噪声水平的含噪电缆瓷套终端红外图像进行去噪.试验中采用sym4小波对图像进行小波分解，分解层数为3，矩形窗区域N(k)的大小取为3×3.
此外，作为与试验的对比，分别采用Rigrsure规则阈值法、Sqtwolog规则法、Minimaxi规则阈值法，并结合软阈值函数对红外图像进行去噪.
文中采用文献[21]中定义的最小均方误差MSE和信噪比SNR作为评估去噪效果好坏的指标，MSE越小，SNR越大，去噪效果越好.
(24)
式中，
(i,j)为去噪后图像的像素灰度值， f(i,j)为原始图像的像素灰度值.
(25)
式中，为去噪后图像像素灰度值的方差.
3.4 结果与分析
表1给出了采用不同去噪算法去噪后图像质量评价参数的计算结果.从表1可以看出，在噪声方差相同的情况下，与Rigrsure阈值法、Minimaxi阈值法以及Sqtwolog阈值法相比，利用文中方法对图像进行去噪后信噪比更大，最小均方误差更小，说明文中方法的去噪性能要比上述3种阈值方法更好.当噪声水平高时，采用文中方法去噪后图像的信噪比和最小均方误差较3种阈值法有明显的提高，说明本方法适用于信噪比小、噪声水平高的红外图像去噪.
表1 文中方法和传统阈值法对不同噪声水平的含噪红外图像去噪后的评价参数
Table 1 evaluating parameters of denoising noisy image contaminated by various noises using the proposed method and traditional threshold methods
图6的4个小图分别是采用Rigrsure阈值法、Minimaxi阈值法、Sqtwolog阈值法、文中方法去噪后的图像.对比发现，在保持完整图像细节的情况下，文中提出的方法去噪后的图像更加清晰，具有更好的去噪效果.
图6 不同去噪方法去噪后图像
Fig.6 Denoised images obtained by different denoising methods
4 结论
文中针对电缆瓷套终端红外图像去噪和小波系数尺度间存在相关性问题，提出了一种基于小波系数尺度间相关性和双变量收缩函数的局域自适应去噪方法.该方法的特点是：
(1)定义小波系数的尺度间相关系数，将相关系数的大小模糊化，使用模糊c-均值阈值法将相关系数划分为两类.对相关系数绝对值小于阈值的小波系数进行置零处理，对相关系数绝对值大于阈值的小波系数进行收缩处理.
(2)考虑小波系数尺度间相关性，引入双变量概率分布模型描述小波系数，进而由MAP估计和贝叶斯定律得到双变量收缩函数；再基于层内局域方差估计得到每个有效小波系数的估计值，对估计小波系数进行重构即可得到去噪后的图像.
文中分别采用3种阈值法和文中提出的基于小波系数尺度间相关性和双变量收缩函数去噪方法，对含有不同噪声水平的电缆瓷套终端红外图像进行仿真去噪处理，并使用信噪比和最小均方误差对去噪效果进行评价.对比不同去噪算法去噪后的图像和评价参数可以得出以下结论：文中提出的去噪方法能够在有效去除噪声的同时，完整保留图像信息；文中方法去噪后图像的信噪比更高，最小均方误差更小，具有更好的去噪效果，尤其是在强噪声下，去噪效果更佳.
