¹新疆师范大学数学科学学院，新疆 乌鲁木齐

²新疆克拉玛依南湖小学，新疆 克拉玛依

收稿日期：2022年1月15日；录用日期：2022年2月17日；发布日期：2022年2月24日

摘要

传统的图像复原的方法会造成图像细节的丢失，且去噪效果一般，结合二进小波变换、高斯滤波和阈值函数去噪的优点，我们提出了一种基于B-样条二进小波变换的图像恢复方法。本文利用小波变换，用新构造的B-样条二进小波滤波器将方差为0.005的高斯噪声图像分解三次，分解得到的每一层高频系数分别使用传统的软阈值模型进行阈值去噪，仅对第一层分解得到的低频系数使用二维高斯滤波器进行去噪，接着，将处理后的高频低频系数利用小波逆变换从第三层重构到第二层，从第二层重构到第一层，从第一层重构到第零层，最后得到复原图像，结果显示本文得到的去噪图像充分保留了原图像的细节，人物的边缘，物体的边缘都能很好地被人眼观察到，具有很好的实用性。

关键词

B-样条二进小波变换，阈值去噪，高斯滤波，图像复原

Image Restoration Method Based on B-Spline Dyadic Wavelet Transform

Ruirui Ma¹, Gang Wang^1*, Jing Zhang¹, Yongxiu Dong²

¹School of Mathematical Sciences, Xinjiang Normal University, Urumqi Xinjiang

²Nanhu Primary School, Karamay Xinjiang

Received: Jan. 15^th, 2022; accepted: Feb. 17^th, 2022; published: Feb. 24^th, 2022

ABSTRACT

The traditional image restoration method will cause the loss of image details, and the denoising effect is general. Combining the advantages of binary wavelet transform, Gaussian filtering and threshold function denoising, we propose a B-spline binary wavelet transform based method—Image restoration method. In this paper, the wavelet transform is used to decompose the Gaussian noise image with a variance of 0.005 three times with the newly constructed B-spline binary wavelet filter. Only the low-frequency coefficients decomposed in the first layer are denoised using a two-dimensional Gaussian filter, and then the processed high-frequency and low-frequency coefficients are reconstructed from the third layer to the second layer using inverse wavelet transform, and the second layer is repeated from the second layer. The first layer is constructed, reconstructed from the first layer to the zeroth layer, and finally the restored image is obtained. The results show that the denoised image obtained in this paper fully retains the details of the original image, and the edges of characters and objects can be well observed by the human eye, it has good practicality.

Keywords:B-Spline Binary Wavelet Transform, Threshold Denoising, Gaussian Filtering, Image Restoration

Copyright © 2022 by author(s) and Hans Publishers Inc.

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

http://creativecommons.org/licenses/by/4.0/

1. 引言

在图像收集、传输、转换的过程中，增加了模糊、失真、噪声等噪声。在实际应用中，需要清晰而高品质的图像。图像复原的最基本的任务是在不丢失图像的详细信息的情况下去除图像中的噪声 [1] [2]。然而，噪声抑制和细节保持往往是矛盾的，这也是图像处理中一直没有解决的问题。

由于小波变换保留了图像的频率信息和空间信息特性使之成为图像去噪领域研究的一个热点。目前常用的滤波方法有高斯滤波、中值滤波、维纳滤波、小波软硬阈值滤波法等。为了获得模糊图像的清晰图像，传统的图像恢复算法通常对整个图像进行处理，信息全部被还原，使得需要的信息不够突出。文献 [3] 提出一种新的自适应高斯滤波图像去噪算法，算法具有一定的自适应效果，但图像细节还是会丢失。而随着小波理论的发展，不同的小波阈值去噪方法被许多学者提出，最早的小波阈值去噪方法是由Donoh提出的软阈值去噪和硬阈值去噪，但硬阈值函数处理后的图像会有振铃、伪吉布斯效应，而软阈值函数相对灵活，因此很多人对它进行了改进，尽量避免软阈值的弊端，文献 [4] 采用阈值去噪，去噪效果一般，并且文中使用的是使用的正交小波变换缺乏平移不变性，会容易使去噪后的图像出现Gibbs效应。

本文引入具有平移不变性的二进小波变换，我们提出了基于B-样条二进小波变换的图像恢复方法，并且结合了高斯滤波，阈值函数的特点，不仅改进了高斯滤波使图像细节丢失的缺点，而且改进了阈值函数去噪和高斯滤波器去噪效果一般的缺点。实验结果表明，基于B-样条二进小波变换的图像恢复方法不仅能够达到很好的去噪效果，并且可以很好地保留图像的细节。

2. 高斯滤波

2.1. 二维高斯分布的定义

二维高斯分布函数是X和Y方向上一维高斯函数的乘积：

$f\left(x,y\right)=f\left(x\right)f\left(y\right)=\frac{1}{\sqrt{2\pi }{\sigma }_x}{\text{e}}^{-\frac{{\left(x-u_x\right)}^2}{2{\sigma }_x{}^2}}\frac{1}{\sqrt{2\pi }{\sigma }_y}{\text{e}}^{-\frac{{\left(x-u_y\right)}^2}{2{\sigma }_y{}^2}}$ (1)

根据二维高斯函数的旋转对称性： ${\sigma }_x={\sigma }_y$

此时，二维高斯函数如下：

$f\left(x,y\right)=f\left(x\right)f\left(y\right)=\frac{1}{2\pi {\sigma }^2}{\text{e}}^{-\frac{x^2+y^2}{2{\sigma }^2}}$ (2)

2.2. 二维高斯分布的曲线

二维高斯分布曲线如下图1所示：

图1. 二维高斯分布曲线

2.3. 高斯滤波器简介

高斯滤波是对整个图像进行加权平均的过程。每个像素的值通过自身和邻域中其他像素值的加权平均来获得。高斯滤波的具体操作是用模型(或卷积掩模)扫描图像中的每个像素，并用模型确定的附近的加权平均像素值替换模型的中心像素值 [5]。本文选取的是15 × 15的高斯滤波器模板，选取 $\sigma =1.6$，计算整数模板的公式，假设我们的模板放在矩阵M中，M的大小为 $\left(2k+1\right)\ast \left(2k+1\right)$，那么 $M\left(i,j\right)$ 的值为 [4]：

$M\left(i,j\right)=\frac{1}{2\pi {\sigma }^2}\exp \left(-\frac{{\left(i-k-1\right)}^2+{\left(j-k-1\right)}^2}{2{\sigma }^2}\right)$ (3)

3. B-样条二进小波滤波器

3.1. 二进小波快速算法—— ${}^{‵}aTrous$

在对图像进行处理时，大部分是在对象为离散状况下进行的，所以需要用到二进小波变换的离散形式计算方法。在1992年， $Mallat$ 创建了能够完成二进小波变换的计算方法，即 ${}^{‵}aTrous$ 算法，与 $Mallat$ 计算方法具有一样的优势，代表了不一样尺度的二进小波系数的联系。需要说明文中所提的 ${}^{‵}aTrous$ 算法是 $T.Abdukirim$ 在 $Mallat$ 创建的算法的基础上修改得到的。

设 $\phi ,\tilde{\phi }$ 是滤波器 $h,\tilde{h}$ 的尺度函数， $\psi ,\tilde{\psi }$ 是滤波器 $g,\tilde{g}$ 的小波函数，假设 $a_0\left[n\right]\ne f\left(t\right)$ 表示输入的样本值，当 $t=n$ 时， $f$ 的局部平均值就可以写为

$a_0\left[n\right]=\langle f\left(t\right),\phi \left(t-n\right)\rangle ={\displaystyle {\int }_{-\infty }^{+\infty }f\left(t\right){\phi }^{*}\left(t-n\right)\text{d}t}$ (4)

当 $j\ge 0$ 时，记

$a_j\left[n\right]=\langle f\left(t\right),{\phi }_{2^{-j}}\left(t-n\right)\rangle$ (5)

其中

${\phi }_{2^{-j}}\left(t-n\right)=2^{\frac{j}{2}}\phi \left(2^j\left(t-n\right)\right)$

若 $j>0$ 在所有点上，那么二进小波的系数计算方法为

$d_j\left[n\right]={\left(W_{\psi },f\right)}_{\left(2^{-j},n\right)}\left(t\right)=\langle f\left(t\right),{\psi }_{2^{-j}}\left(t-n\right)\rangle ,n\in Z$ (6)

也就是说任意尺度 $2^j>1,\left\{d_{1,}{d}_2,\cdots ,d_j,a_j\right\}$ 这一离散信号的序列便可以被叫做 $a_0$ 离散信号的二进小波变换。

定理2.1二维 ${}^{‵}aTrous$ 算法 [6] [7] [8]

二维DyWT分解公式

$\begin{array}{l}{a}_{j+1}\left[n_1,n_2\right]={\displaystyle \underset{k_1,k_2}{\sum }h\left[k_1\right]h\left[k_2\right]a_j\left[n_1+2^{-j}{k}_1,n_2+2^{-j}{k}_2\right]}\\ d_{j+1}^1\left[n_1,n_2\right]={\displaystyle \underset{k_1,k_2}{\sum }h\left[k_1\right]g\left[k_2\right]a_j\left[n_1+2^{-j}{k}_1,n_2+2^{-j}{k}_2\right]}\\ d_{j+1}^2\left[n_1,n_2\right]={\displaystyle \underset{k_1,k_1}{\sum }g\left[k_1\right]h\left[k_2\right]a_j\left[n_1+2^{-j}{k}_1,n_2+2^{-j}{k}_2\right]}\\ d_{j+1}^3\left[n_1,n_2\right]={\displaystyle \underset{k_1,k_2}{\sum }g\left[k_1\right]g\left[k_2\right]a_j\left[n_1+2^{-j}{k}_1,n_2+2^{-j}{k}_2\right]}\end{array}$ (7)

二维DyWT重构公式

$\begin{array}{c}{a}_j\left[n_1,n_2\right]=\frac{1}{4}({\displaystyle \underset{k_1,k_2}{\sum }\tilde{h}\left[k_1\right]\tilde{h}\left[k_2\right]a_{j+1}\left[n_1-2^{-j}{k}_1,n_2-2^{-j}{k}_2\right]}\\ +{\displaystyle \underset{k_1,k_2}{\sum }\tilde{h}\left[k_1\right]\tilde{g}\left[k_2\right]d_{j+1}^1\left[n_1-2^{-j}{k}_1,n_2-2^{-j}{k}_2\right]}\\ +{\displaystyle \underset{k_1,k_2}{\sum }\tilde{g}\left[k_1\right]\tilde{h}\left[k_2\right]d_{j+1}^2\left[n_1-2^{-j}{k}_1,n_2-2^{-j}{k}_2\right]}\\ +{\displaystyle \underset{k_1,k_2}{\sum }\tilde{g}\left[k_1\right]\tilde{g}\left[k_2\right]d_{j+1}^3\left[n_1-2^{-j}{k}_1,n_2-2^{-j}{k}_2\right]})\end{array}$ (8)

其中， $j=0,1,2,\cdots$。

3.2. B-样条二进小波滤波器的构造

该方法是对 $Mallat$ 的B-样条二进小波的推广 [6] [7] [9]，这个过程是 $T.Abdukirim$ 完成的。具体过程如下， $m$ 次 $B$ -样条函数 $\varphi \left(t\right)$ 是特征函数 ${\chi }_{\left[0,1\right]}$ 与其自身的卷积的平移 [10]： $\varphi \left(t\right)=\chi \ast \chi \ast \chi \cdots {\chi }_{\left[0,1\right]}(\; t\; )$

对 $\varphi \left(t\right)$ 进行 $Fourier$ 变换，可以得到：

$\hat{\phi }\left(\omega \right)={\left(\frac{\sin \left(\frac{\omega }{2}\right)}{\frac{\omega }{2}}\right)}^{m+1}{\text{e}}^{-i\xi \frac{\omega }{2}},\xi =\{\begin{array}{cc}1,& m=2n\\ 0,& m=2n+1\end{array}\left(m\in Z,m\ge 0\right)$ (9)

根据二尺度关系有

$\hat{h}\left(\omega \right)=\sqrt{2}\frac{\hat{\phi }\left(2\omega \right)}{\hat{\phi }\left(\omega \right)}\sqrt{2}{\left(\cos \left(\frac{\omega }{2}\right)\right)}^{m+1}{\text{e}}^{-i\xi \frac{\omega }{2}}$ (10)

选取滤波器 $\tilde{h}:\widehat{\stackrel{˜}{h}}\left(\omega \right)=\tilde{h}\left(\omega \right)$，选取滤波器 $g$ ：

$\hat{g}\left(\omega \right)={\left(-i\right)}^{\tau }\sqrt{2}{\text{e}}^{-i\frac{2-\tau }{2}\omega }{\left(\sin \left(\frac{\omega }{2}\right)\right)}^r,r=1,2;\tau =\{\begin{array}{cc}1,& r=1\\ 0,& r=2\end{array}$ (11)

由二进完全重构有：

$\widehat{\stackrel{˜}{g}}\left(\omega \right)={\left(-i\right)}^{\tau }\sqrt{2}{\text{e}}^{-i\frac{2-\tau }{2}\omega }{\left(\sin \left(\frac{\omega }{2}\right)\right)}^{2-\tau }{\displaystyle \underset{l=0}{\overset{m}{\sum }}{\left(\cos \left(\frac{\omega }{2}\right)\right)}^{2l}}$ (12)

因此由二尺度关系有

$\hat{\psi }\left(\omega \right)={\left(-\frac{\omega }{4}\right)}^r{\text{e}}^{-i\left(\frac{2+\xi -\tau }{4}\omega -\frac{\pi r}{2}\right)}{\left(\frac{\sin \left(\frac{\omega }{4}\right)}{\frac{\omega }{4}}\right)}^{m+r+1}$ (13)

上述 $\left\{h,g,\tilde{h},\tilde{g}\right\}$ 是二进小波滤波器组，如表1是本文构造的滤波器系数。

表1. r = 1，m = 0时的B-样条二进小波滤波器的系数

4. 阈值函数去噪

在阈值模型去噪中，阈值模型反映了阈值上下波动系数模块的不同处理策略和不同估计方法。如果阈值太小，表示的信号仍然有噪声，如果阈值太大，信号的重要特征将被过滤 [10]，估计小波系数的方法有两种：

1) 软阈值模型：

${\tilde{w}}_{j,k}=\mathrm{sgn}{w}_{j,k}\max \left\{\left|w_{j,k}\right|-\lambda ,0\right\}=\{\begin{array}{l}{w}_{j,k}-\lambda ,w_{j,k}\ge \lambda \\ 0,\left|w_{j,k}\right|\prec \lambda \\ w_{j,k}+\lambda ,w_{j,k}\le \lambda \end{array}$ (14)

2) 硬阈值模型

${\tilde{w}}_{j,k}=\{\begin{array}{l}{w}_{j,k},\left|w_{j,k}\right|\ge \lambda \\ 0,\left|w_{j,k}\right|\prec \lambda \end{array}$ (15)

硬阈值函数能在一定程度上保留图像边缘细节等局部信息，但其不连续性会时常导致图像失真现象，软阈值函数相对连续，重构后图像失真较少，本文选取软阈值模型对图像进行去噪。

5. 实验过程与分析

本文实验是以MATLAB R2016a为实验平台，以标准的“cameraman.tif”图像加入方差为0.005的高斯噪声作为分析对象。本文的实验过程是首先利用表1构造的二进小波滤波器对含噪图像分解三次。

1) 小波变换，利用新构造的B-样条二进小波滤波器对含噪图像进行3层小波分解。

2) 对图像三层分解得到的高频系数使用传统的软阈值模型进行阈值去噪。

3) 对图像第一层分解得到的低频系数使用二维高斯滤波器进行去噪。

4) 接着小波逆变换，将图像第三层到第二层，第二层到第一层，第一层到第零层进行重构。

5) 图像完成重构后，得到最终去噪后的图像。

分析：从图2和表2的实验结果图中可以看到，单独使用高斯滤波器去噪得到的图像去噪效果一般，图像中的人物边缘线条模糊，远处的建筑物几乎与背景融为一体，照相机支架边缘也模糊，不利于人眼观察，单独使用软硬阈值的去噪方法，也达不到很好的去噪效果，而本文方法去噪的图像人物线条清晰，人眼，照相机，包括远处的建筑物都能够很好的观察到，并且，达到了很好的去噪效果，分解流程图如图3所示。

图2. “Cameraman.tif”实验结果图

表2. 峰值信噪比

图3. 本文方法流程图

6. 结语

本文结合B-样条二进小波变换的特征和高斯滤波的优点，提出了基于B-样条二进小波变换的图像复原方法，人物图像的细节很多，边缘容易模糊，保持所有有用信息，特别重要的是达到去除噪音的效果。实验结果表明，结合B-样条二进小波滤波器去噪效果大于单纯的使用高斯滤波及软硬阈值传统方法的去噪效果，本文提到的方法复原后的图像没有模糊的现象，轮廓更完整，边缘更清晰，去噪效果更好，是一种有效可行的图像复原方法，在后续的工作中，本文中去噪方法还需要进一步提升，让去噪效果更佳，本文也证明了B-样条二进小波小波在图像去噪中的良好性能，为下一步研究图像去噪奠定了基础。

基金项目

新疆师范大学优秀青年教师科研启动基金资助项目(批准号：XJNU202014)。

文章引用

马瑞瑞,王 刚,张 静,董永秀. 基于B-样条二进小波变换的一种改进图像复原方法
Image Restoration Method Based on B-Spline Dyadic Wavelet Transform[J]. 理论数学, 2022, 12(02): 309-315. https://doi.org/10.12677/PM.2022.122035

参考文献