CN101882301B - 基于截断窗的剪切波图像去噪方法 - Google Patents

Abstract

本发明公开了一种基于截断窗的剪切波图像去噪方法，主要解决现有剪切波基于第一基函数算法时间复杂度大的问题与现有剪切波基于传统的拉普拉斯金字塔算法对频带划分不灵活且去噪效果差的问题。其中包括对含噪图像的离散傅里叶变换，矩形截断窗或圆形截断窗的构造；用矩形中心截断窗或圆形中心截断窗得到近似图像，并用剪切波对矩环截断窗或圆环截断窗对含噪图像频域作用的频带进行方向划分得到系数矩阵；继而对系数矩阵进行硬阈值去噪；而后用矩形截断窗或圆形截断窗结合剪切波再对去噪后的系数矩阵与近似图像进行重构，得到去噪后的图像。本发明构建的截断窗在图像去噪上能取得更加理想的效果，可用于图像分析。

Description

基于截断窗的剪切波图像去噪方法

技术领域

本发明属于图像处理领域，涉及图像的去噪方法，可用于基于全变差，硬阈值，以及非局部均值的图像去噪。

背景技术

剪切波分析是一种继承曲线波和方波的优点的新型多尺度几何分析工具，通过对基函数的缩放，剪切和平移等放射变换生成具有不同特征的剪切波函数，对于包含C²奇异曲线或曲线的高维信号具有最优特性。对于二维信号，它不仅可以检测到所有的奇异点，还可以自适应的跟踪奇异曲线的方向，随着尺度参数的变化，可以精确的描述函数的奇异性特征，实现以经典多尺度分析描述高维信号中的奇异性，同时也为方波建立了数学理论基础。

对于二维分段光滑信号，剪切波可以对其达到理论性的最稀疏表示，由于剪切波还是一个很新的多尺度几何分析工具，现在的应用领域还不是很广泛，但它有着其它多尺度集合工具所没有的一些优秀的特性，如：(1)定义在卡特森域，可以通过剪切变换获得各种各样的方向，因此他对剪切操作的方向不受限制，可以在更多的方向上分析图像；(2)逆变换仅需要剪切波滤波器的简单合成而不是方向滤波器的逆变换，具有简单快速的离散化实现形式等。这些特点使其在图像处中有所作为，如图像降噪，压缩，增强和水印等方面。

剪切波在去噪方面的应用较为广泛，比如：全变差，拉普拉斯金字塔分解，以及非局部均值都可以结合剪切波，进行去噪应用。

传统的剪切波基于拉普拉斯金字塔分解的阈值去噪算法，利用了拉普拉斯金字塔分解，将图像分为不同的频带，相应的得到细节图像与近似系数。而一般结合剪切波，将图像分三层，最外层所对应的频率最高，其次为次外层，最内层的是低频的近似图像。但是，传统的基于拉普拉斯金字塔分解结合剪切波，对图像作用得到的系数矩阵稀疏性不太高，从而在应用阈值去噪时，不能有效地提取图像的边缘信息，以及摒除图像的噪声，从而使得去噪的图像边缘信息保留的不够好，直接导致了图像的峰值信噪比PSNR比较低，去噪效果不理想。

发明内容

本发明的目的在于克服上述已有技术的不足，提出一种基于剪切波的截断窗的图像去噪方法，提高图像去噪结果的峰值信噪比PSNR，更加有效的保留了图像的边缘信息，而且更加有效的去除了图像的噪声，使得PSNR更加理想。

本发明是这样实现的：

本发明是在剪切波基础上提出，主要是利用剪切波的基函数

得到对称基函数

而后将剪切波基函数与对称基函数进行三维立体拓展，构建剪切波的矩形或圆形截断窗。剪切波的截断窗，继承传统拉普拉斯金字塔对图像分层：中心截断窗截，与原图像频域作用得到原图像的低频成份，环形截断窗与原图像的频域作用，得到原图像的次高频与高频成份。再利用剪切波对图像高频成分与次高频成份进行作用，得到图像的各层，各个方向的系数矩阵，利用阈值，对图像进行去噪。

根据构建截断窗的两种不同形状，本发明给出以下两种实施方案：

技术方案1，基于矩形截断窗的剪切波图像去噪方法，包括以下步骤：

(1)对含噪图像进行离散傅里叶变换；

(2)利用剪切波的基函数

构建截断窗：

2a)构建对称基函数

如下：

$\hat{f}c\left(\mathrm{\&xi;}\right)=\left\{\begin{array}{cc}{\widehat{\mathrm{\&psi;}}}_2\left(\mathrm{\&xi;}\right)& 0\&le;\mathrm{\&xi;}<1\\ {\widehat{\mathrm{\&psi;}}}_2(-\mathrm{\&xi;})& -1\&le;\mathrm{\&xi;}<0\end{array};\right.$

其中ξ是函数的自变量，

2c)在对称基函数的基础上，构建矩形中心截断窗cut(0，sξ₁，sξ₂)为：

$\mathrm{cut}(0,s{\mathrm{\&xi;}}_1,s{\mathrm{\&xi;}}_2)$

$=\hat{f}c\left(s{\mathrm{\&xi;}}_1\right){\mathrm{\&chi;}}_{D_0}({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)+\hat{f}c\left(s{\mathrm{\&xi;}}_2\right){\mathrm{\&chi;}}_{D_1}({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)$

其中S是大于等于1的正整数，d＝0，1，

定义为： ${\mathrm{\&chi;}}_{D_d}({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)=\left\{\begin{array}{cc}1& ({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)\&Element;D_d\\ 0& ({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)\&NotElement;D_d\end{array}\right.,$

D_d定义为：

其中

表示二维实数集；

2d)在基函数

的基础上构建矩形环截断窗cut(k，sξ₁，sξ₂)为：

$\mathrm{cut}(k,s{\mathrm{\&xi;}}_1,s{\mathrm{\&xi;}}_2)$

$=[{\widehat{\mathrm{\&psi;}}}_2(s{\mathrm{\&xi;}}_1-k)+{\widehat{\mathrm{\&psi;}}}_2(-(s{\mathrm{\&xi;}}_1+k))]{\mathrm{\&chi;}}_{D_0}({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)$

$+[{\widehat{\mathrm{\&psi;}}}_2(s{\mathrm{\&xi;}}_2-k)+{\widehat{\mathrm{\&psi;}}}_2(-(s{\mathrm{\&xi;}}_2+k))]{\mathrm{\&chi;}}_{D_1}({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)$

其中k＝1，2，......s；

2e)令s＝2，将矩形中心截断窗与矩形环截断窗统一表示为截断窗cut：

cut(1)＝cut(O，2ξ₁，2ξ₂)，cut(2)＝cut(1，2ξ₁，2ξ₂)，cut(3)＝cut(2，2ξ₁，2ξ₂)；

(3)将步骤(1)得到的含噪图像离散傅里叶变换与cut²(1)点乘，再进行逆离散傅里叶变换得到近似图像a₀，再用步骤(1)含噪图像进行离散傅里叶变换再分别与cut(2)和cut(3)点乘，得到带通域矩阵F₁和高通域矩阵F₂；

(4)将剪切波

按

重新规划为W_j，l，

其中l＝1，2，...4j；n＝2，3，.......2j，

L＝[2j+1]，[2j+1]表示对2j+1取四舍五入整，

j为大于等于1的整数，m取-j到j的整数，

(5)用重新规划的剪切波W_j，l对步骤(3)得到的带通域F₁与高通域F₂进行去噪；

(6)对个系数矩阵f′_j，l与近似图像a₀进行如下重构，得到去噪图像：

6a)对去噪后的系数矩阵f′_j，l进行离散傅里叶变换，再与cut(j+1)进行点乘运算，将其结果再与W_j，l做点乘运算；

6b)对6a)进行逆离散傅里叶变换，得到rf′_i，j，将所有的rf′_i，j相加，再加上近似图像a₀，得到重构图像。

技术方案2，基于圆形截断窗的剪切波图像去噪方法，包括以下步骤：

1)对含噪图像进行离散傅里叶变换；

2)利用剪切波的基函数

构建截断窗：

2.1)构建对称基函数

如下：

$\hat{f}c\left(\mathrm{\&xi;}\right)=\left\{\begin{array}{cc}{\widehat{\mathrm{\&psi;}}}_2\left(\mathrm{\&xi;}\right)& 0\&le;\mathrm{\&xi;}<1\\ {\widehat{\mathrm{\&psi;}}}_2(-\mathrm{\&xi;})& -1\&le;\mathrm{\&xi;}<0\end{array};\right.$

其中ξ是函数的自变量，

2.2)在对称基函数的基础上，构建圆形中心截断窗cut(0，sξ₁，sξ₂)为：

$\mathrm{cut}(0,s{\mathrm{\&xi;}}_1,s{\mathrm{\&xi;}}_2)$

$=\hat{f}c\left(s\sqrt{{{\mathrm{\&xi;}}_1}^2+{{\mathrm{\&xi;}}_2}^2}\right)$

其中s是大于等于1的正整数，d＝0，1，

定义为： ${\mathrm{\&chi;}}_{D_d}({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)=\left\{\begin{array}{cc}1& ({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)\&Element;D_d\\ 0& ({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)\&NotElement;D_d\end{array}\right.,$

D_d定义为：

其中

表示二维实数集；

2.3)在基函数

的基础上构建圆形环截断窗cut(k，sξ₁，sξ₂)为：

$\mathrm{cut}(k,s{\mathrm{\&xi;}}_1,s{\mathrm{\&xi;}}_2)$

$={\widehat{[\mathrm{\&psi;}}}_2(s\sqrt{({{\mathrm{\&xi;}}_1}^2+{{\mathrm{\&xi;}}_2}^2)}-k)$

其中

2.4)令s＝2，将圆形中心截断窗与圆形环截断窗统一表示为截断窗cut：

cut(1)＝cut(0，2ξ₁，2ξ₂)，cut(2)＝cut(1，2ξ₁，2ξ₂)

$\mathrm{cut}\left(3\right)=\sqrt{\mathrm{cut}(\mathrm{2,2}{\mathrm{\&xi;}}_1,2{\mathrm{\&xi;}}_2)^2+\mathrm{cut}(\mathrm{3,2}{\mathrm{\&xi;}}_1,2{\mathrm{\&xi;}}_2)^2}$

3)将步骤1)得到的含噪图像离散傅里叶变换与cut²(1)点乘，再进行逆离散傅里叶变换得到近似图像a₀，再用步骤1)含噪图像进行离散傅里叶变换，再分别与cut(2)和cut(3)点乘，得到带通域矩阵F₁和高通域矩阵F₂；

4)将剪切波

按

重新规划为W_j，l，

其中l＝1，2，...4j；n＝2，3，.......2j；

L＝[2j+1]，[2j+1]表示对2j+1取四舍五入整；

j为大于等于1的整数，m取-j到j的整数，

5)用重新规划的剪切波W_j，l对步骤3)得到的带通域F₁与高通域F₂进行去噪；

6)对个系数矩阵f′_j，l与近似图像a₀进行重构，得到去噪图像：

6.1)对去噪后的系数矩阵f′_j，l进行离散傅里叶变换，再与cut(j+1)进行点乘运算，将其结果再与W_j，l做点乘运算；

6.2)对6.1)进行逆离散傅里叶变换，得到rf′_i，j，将所有的rf′_i，j相加，再加上近似图像a₀，得到重构图像。

本发明具有如下优点：

(1)本发明所构建的截断窗可以任意的对图像频带进行划分层，与传统的拉普拉斯金字塔对图像频带只能固定的划分相比，具有很强的灵活性，针对性更强；

(2)剪切波基于截断窗的硬阈值去噪与基于拉普拉斯金字塔的剪切波硬阈值去噪方法相比，去噪效果更好；

(3)本发明构建的矩形截断窗或圆形截断窗与基于剪切波基函数

构建的频带划分窗相比，在去噪等其他应用上，速度有明显的提高，复杂度有明显的减小；

(4)仿真结果表明，在单纯的硬阈值去噪上，本发明较其他几种方法，圆形截断窗与圆形截断窗在硬阈值去噪上，具有很大的优越性；

(5)本发明构建的截断窗还可以作为图像的边缘检测等的图像分析应用的基础框架。

附图说明

图1是本发明的流程图；

图2是剪切波的各个方向与区域划分图；

图3是本发明的矩形截断窗示意图；

图4是本发明的圆形截断窗示意图；

图5是对lena512.bmp图像加标准差为30的高斯白噪声之后，利用现有拉普拉斯的剪切波分解与本发明基于截断窗分解的剪切波去噪后局部区域对比图。

具体实施方式

实施例一：

参照图1，本发明基于矩形截断窗的剪切波图像去噪方法，包括以下步骤：

步骤1、对含噪图像im进行离散傅里叶变换DFT(im)。

步骤2、在频域构造现有的剪切波基函数

(2.1)构造肿块函数h₁，该h₁要满足在(-2，2)区间内0到无穷阶可微，且在(-2，2)内0≤h₁≤1，在[-1，1]内h₁＝1，构造的h₁表达式如下：

$h_1\left(\mathrm{\&xi;}\right)=e^{-\frac{28{(\mathrm{\&xi;}-1)}^{14}}{2^{14}}}$

其中ξ是函数的自变量

(2.2)构造辅助函数h₂该h₂表达式如下：

$h_2\left(\mathrm{\&xi;}\right)=\sqrt{1-e^{\frac{1}{(\mathrm{\&xi;}-1)}}}$

(2.3)在(2.1)，(2.2)的基础上，得到左边支函数h，该h表达式如下：

h(ξ)＝h₁(ξ)h₂(ξ)

(2.4)构造函数右边支函数g，该函数g表达式如下：

$g\left(\mathrm{\&xi;}\right)=\sqrt{1-\frac{(1-e^{-\frac{1}{\mathrm{\&xi;}-2}})}{e^{{\left(\frac{28{(\mathrm{\&xi;}-2)}^{14}}{2^{14}}\right)}^2}}}$

(2.5)利用前面得到的左边支函数h和右边支函数g，构造剪切波基函数其表达式如下：

${\widehat{\mathrm{\&psi;}}}_2\left(\mathrm{\&xi;}\right)=\left\{\begin{array}{cc}h(2\mathrm{\&xi;}+1)& -1\&le;\mathrm{\&xi;}<0\\ g\left(2\mathrm{\&xi;}\right)& 0\&le;\mathrm{\&xi;}<1\\ 0& \end{array}\right.$

步骤3、在步骤2构建的剪切波基函数

的基础之上，构建矩形截断窗

(3.1)先构建对称基函数

构建公式如下：

$\hat{f}c\left(\mathrm{\&xi;}\right)=\left\{\begin{array}{cc}{\widehat{\mathrm{\&psi;}}}_2\left(\mathrm{\&xi;}\right)& 0\&le;\mathrm{\&xi;}<1\\ {\widehat{\mathrm{\&psi;}}}_2(-\mathrm{\&xi;})& -1\&le;\mathrm{\&xi;}<0\end{array}\right.$

(3.2)在(3.1)的基础上，构建矩形窗的中心截断窗cut(0，sξ₁，sξ₂)，构建公式为

$\mathrm{cut}(0,s{\mathrm{\&xi;}}_1,s{\mathrm{\&xi;}}_2)$

$=\hat{f}c\left(s{\mathrm{\&xi;}}_1\right){\mathrm{\&chi;}}_{D_0}({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)+\hat{f}c\left(s{\mathrm{\&xi;}}_2\right){\mathrm{\&chi;}}_{D_1}({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)$

其中s为大于等于1的整数

定义为： ${\mathrm{\&chi;}}_{D_d}({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)=\left\{\begin{array}{cc}1& ({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)\&Element;D_d\\ 0& ({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)\&NotElement;D_d\end{array}\right.,$

D_d定义为：

D₀，D₁的区域表示如图4所示，

其中

表示二维实数集；

3.3)构建矩形截断窗cut(k，sξ₁，sξ₂)

$\mathrm{cut}(k,s{\mathrm{\&xi;}}_1,s{\mathrm{\&xi;}}_2)$

$=[{\widehat{\mathrm{\&psi;}}}_2(s{\mathrm{\&xi;}}_1-k)+{\widehat{\mathrm{\&psi;}}}_2(-(s{\mathrm{\&xi;}}_1+k))]{\mathrm{\&chi;}}_{D_0}({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)$

$+[{\widehat{\mathrm{\&psi;}}}_2(s{\mathrm{\&xi;}}_2-k)+{\widehat{\mathrm{\&psi;}}}_2(-(s{\mathrm{\&xi;}}_2+k))]{\mathrm{\&chi;}}_{D_1}({\mathrm{\&xi;}}_1,{\mathrm{\&xi;}}_2)$

其中k＝1，2，......s

3.4)取s＝2，重新规划矩形截断窗：

cut(1)＝cut(0，2ξ₁，2ξ₂)

cut(2)＝cut(1，2ξ₁，2ξ₂)

cut(3)＝cut(2，2ξ₁，2ξ₂)

当s＝2时，重新规划的矩形窗截断窗如图3所示，图a为矩形中心截断窗cut(1)的平面图，图b为矩形中心截断窗cut(1)的立体图，图c为矩形环截断窗cut(2)的平面图，图e为矩形环截断窗cut(2)的立体图，图e为矩形环截断窗cut(3)的平面图，图f为矩形环截断窗cut(3)的立体图；

步骤4、将步骤1得到的含噪图像离散傅里叶变换与cut²(1)点乘，再进行逆离散傅里叶变换得到近似图像a₀，再用步骤(1)得到的含噪图像离散傅里叶变换分别与cut(2)和cut(3)点乘，对其结果再进行逆离散傅里叶变换，得到带通域矩阵F₁和高通域矩阵F₂。

步骤5、构建剪切波，并对其进行重行规划

(5.1)利用步骤2构建的剪切波基函数

构建剪切波

$\left.\begin{array}{c}{W}_{j,m}^{\left(0\right)}\left(\mathrm{\&xi;}\right)=s{\widehat{\mathrm{\&psi;}}}_2(j\frac{{\mathrm{\&xi;}}_2}{{\mathrm{\&xi;}}_1}-m){\mathrm{\&chi;}}_{D_0}\\ W_{j,m}^{\left(1\right)}\left(\mathrm{\&xi;}\right)=s{\widehat{\mathrm{\&psi;}}}_2(j\frac{{\mathrm{\&xi;}}_1}{{\mathrm{\&xi;}}_2}-m){\mathrm{\&chi;}}_{D_1}\end{array}\right\}$

j为大于等于1的整数，m取-j到j的整数，

(5.2)对上一步构建的剪切波

按

重新规划，得到剪切波W_j，l；

其中n＝2，2，....2j，l＝1，2，...4j；L＝[2j+1]，其中[2j+1]表示对2j+1取四舍五入整，

这样，根据j的变化，剪切波W_j，l的个数将以4，8，12.......4n(n≥1)规律出现。

步骤6、用重新规划的剪切波W_j，l对步骤4得到的带通域F₁与高通域F₂进行去噪：

(6.1)用步骤4得到F₁与W_1，l点乘，得到F_1，l，并对其进行逆离散傅里叶变换IDFT(F_1，l)，得到四个截断剪切波方向窗的系数f_1，l，同理用步骤4得到F₂与W_2，l点乘，得到F_2，l，并对其进行逆离散傅里叶变换IDFT(F_2，l)，得到八个截断剪切波方向窗系数f_2,l；

(6.2)对上面得到的系数矩阵f_j，l进行硬阈值去噪，令阈值

其中，

为对应系数矩阵f_j，l的方差，

为对应系数矩阵f_j，l的噪声估计方差，噪声方差σ_j，l，n用蒙特卡罗表达式

进行估计得到；

其中E()表示求期望，δ表示噪声标准差；阈值选取规则是：对于系数矩阵f_j，l元素的绝对值大于等于τ_j，l的保留原值，对于系数矩阵f_j，l元素的绝对值小于τ_j，l的置零，这样就去除了系数矩阵f_j，l所含的大部分非边缘的噪声；

步骤7、对个系数矩阵f′_j，l与近似图像a₀进行如下重构，得到去噪图像：

(7.1)对去噪后的系数矩阵f′_j，l进行离散傅里叶变换，再与cut(j+1)进行点乘运算，将其结果再与W_j，l做点乘运算；

(7.2)对(7.1)得到的结果进行逆离散傅里叶变换，得到rf′_i，j，将所有的rf′_i，j相加，再加上近似图像a₀，得到重构图像，即去噪图像xr，用公式表示如下：

$\mathrm{xr}=a_0+\underset{i,l}{\mathrm{\&Sigma;}}{\mathrm{rf}}_{j,l}^{\&prime;\left(d\right)}.$

实施例二：

参照图1，本发明基于圆形截断窗的剪切波图像去噪方法，包括以下步骤：

步骤A，对含噪图像进行离散傅里叶变换。

步骤B，在频域构造现有的剪切波基函数

(B1)构造肿块函数h₁，该h₁要满足在(-2，2)区间内0到无穷阶可微，且在(-2，2)内0≤h₁≤1，在[-1，1]内h₁＝1，构造的h₁表达式如下：

$h_1\left(\mathrm{\&xi;}\right)=e^{-\frac{28{(\mathrm{\&xi;}-1)}^{14}}{2^{14}}}$

其中ξ是函数的自变量；

(B2)构造辅助函数h₂该h₂表达式如下：

$h_2\left(\mathrm{\&xi;}\right)=\sqrt{1-e^{\frac{1}{(\mathrm{\&xi;}-1)}}};$

(B3)在步骤(B1)和(B2)基础上，得到左边支函数h，该h表达式如下：

h(ξ)＝h₁(ξ)h₂(ξ)；

(B4)构造右边支函数g，该函数g表达式如下：

$g\left(\mathrm{\&xi;}\right)=\sqrt{1-\frac{(1-e^{-\frac{1}{\mathrm{\&xi;}-2}})}{e^{{\left(\frac{28{(\mathrm{\&xi;}-2)}^{14}}{2^{14}}\right)}^2}}};$

(B5)利用步骤(B3)得到的左边支函数h和步骤(B4)得到的右边支函数g，构造剪切波基函数

其表达式如下：

${\widehat{\mathrm{\&psi;}}}_2\left(\mathrm{\&xi;}\right)=\left\{\begin{array}{cc}h(2\mathrm{\&xi;}+1)& -1\&le;\mathrm{\&xi;}<0\\ g\left(2\mathrm{\&xi;}\right)& 0\&le;\mathrm{\&xi;}<1\\ 0& \end{array}\right..$

步骤C，利用剪切波的基函数

构建圆形截断窗

C1)构建对称基函数如下：

$\hat{f}c\left(\mathrm{\&xi;}\right)=\left\{\begin{array}{cc}{\widehat{\mathrm{\&psi;}}}_2\left(\mathrm{\&xi;}\right)& 0\&le;\mathrm{\&xi;}<1\\ {\widehat{\mathrm{\&psi;}}}_2(-\mathrm{\&xi;})& -1\&le;\mathrm{\&xi;}<0\end{array}\right.;$

C2)在对称基函数的基础上，构建圆形中心截断窗cut(0，sξ₁，sξ₂)为：

$\mathrm{cut}(0,s{\mathrm{\&xi;}}_1,s{\mathrm{\&xi;}}_2)$

$=\hat{f}c\left(s\sqrt{{{\mathrm{\&xi;}}_1}^2+{{\mathrm{\&xi;}}_2}^2}\right)$

s为大于等于1的整数；

C3)在基函数的基础上构建圆形环截断窗cut(k，sξ₁，sξ₂)为：

$\mathrm{cut}(k,s{\mathrm{\&xi;}}_1,s{\mathrm{\&xi;}}_2)$

$=[{\widehat{\mathrm{\&psi;}}}_2(s\sqrt{({{\mathrm{\&xi;}}_1}^2+{{\mathrm{\&xi;}}_2}^2)}-k)$

其中

C4)令s＝2，将圆形中心截断窗与圆形环截断窗统一表示为截断窗cut：

cut(1)＝cut(0，2ξ₁，2ξ₂)，cut(2)＝cut(1，2ξ₁，2ξ₂)

$\mathrm{cut}\left(3\right)=\sqrt{\mathrm{cut}(\mathrm{2,2}{\mathrm{\&xi;}}_1,2{\mathrm{\&xi;}}_2)^2+\mathrm{cut}(\mathrm{3,2}{\mathrm{\&xi;}}_1,2{\mathrm{\&xi;}}_2)^2},$

当s＝2时，重新规划的矩形窗截断窗如图4所示，其中图4a为圆形中心截断窗cut(1)的平面图，图4b为圆形中心截断窗cut(1)的立体图，图4c为圆形环截断窗cut(2)的平面图，图4d为圆形环截断窗cut(2)的立体图，图4e为圆形环截断窗cut(3)的平面图，图4f为圆形环截断窗cut(3)的立体图。

步骤D，将步骤A得到的含噪图像离散傅里叶变换与cut²(1)点乘，再进行逆离散傅里叶变换得到近似图像a₀，再用步骤A含噪图像进行离散傅里叶变换，再分别与cut(2)和cut(3)点乘，得到带通域矩阵F₁和高通域矩阵F₂。

步骤E，构建剪切波，并对其进行重行规划

(E1)利用步骤B构建的剪切波基函数

构建剪切波

$\left.\begin{array}{c}{W}_{j,m}^{\left(0\right)}\left(\mathrm{\&xi;}\right)=s{\widehat{\mathrm{\&psi;}}}_2(j\frac{{\mathrm{\&xi;}}_2}{{\mathrm{\&xi;}}_1}-m){\mathrm{\&chi;}}_{D_0}\\ W_{j,m}^{\left(1\right)}\left(\mathrm{\&xi;}\right)=s{\widehat{\mathrm{\&psi;}}}_2(j\frac{{\mathrm{\&xi;}}_1}{{\mathrm{\&xi;}}_2}-m){\mathrm{\&chi;}}_{D_1}\end{array}\right\}$

j为大于等于1的整数，m取-j到j的整数，

(E2)对(E1)构建的剪切波

按

重新规划，得到剪切波W_j，l；

其中n＝2，2，....2j，l＝1，2，...4j；

L＝[2j+1]，其中[2j+1]表示对2j+1取四舍五入整，

这样，根据j的变化，剪切波W_j，l的个数将以4，8，12.......4n(n≥1)规律出现。

步骤F，用重新规划的剪切波W_j，l对步骤D得到的带通域F₁与高通域F₂进行去噪：

(F1)用步骤D得到F₁与W_1，l点乘，得到F_1，l，并对其进行逆离散傅里叶变换IDFT(F_1，l)，得到四个截断剪切波方向窗的系数f_1，l，同理用步骤D得到F₂与W_2，l点乘，得到F_2，l，并对其进行逆离散傅里叶变换IDFT(F_2，l)，得到八个截断剪切波方向窗系数f_2，l；

(F2)对(F1)得到的系数矩阵f_j，l进行硬阈值去噪，令阈值其中，

为对应系数矩阵f_j，l的方差，

为对应系数矩阵f_j，l的噪声估计方差，噪声方差σ_j，l，n用蒙特卡罗表达式

进行估计得到；

其中E()表示求期望，δ表示噪声标准差；阈值选取规则是：对于系数矩阵f_j，l元素的绝对值大于等于τ_j，l的保留原值，对于系数矩阵f_j，l元素的绝对值小于τ_j，l的置零，这样就去除了系数矩阵f_j，l所含的大部分非边缘的噪声。

步骤G，对各个系数矩阵f′_j，l与近似图像a₀进行如下重构，得到去噪图像。

(G1)对去噪后的系数矩阵f′_j，l进行离散傅里叶变换，再与cut(j+1)进行点乘运算，将其结果再与W_j，l做点乘运算；

(G2)对(G1)得到的结果进行逆离散傅里叶变换，得到rf′_i，j，将所有的rf′_i，j相加，再加上近似图像a₀，得到重构图像，即去噪图像xr，用公式表示如下：

$\mathrm{xr}=a_0+\underset{i,l}{\mathrm{\&Sigma;}}{\mathrm{rf}}_{j,l}^{\&prime;\left(d\right)}$

本发明的效果可以通过下面的仿真进一步说明：

1、仿真条件

将一幅256×256的自然图像lena_256.bmp，三幅512×512的自然图像lena512.bmp，barbara512.bmp，flower.bmp在matalab7.01下运行。

2、仿真内容

A：按照说明书具体实施例一中的步骤1，步骤2，......步骤5，步骤6中的(6.1)与具体实施案例二中的步骤A，步骤B，......步骤E，步骤F中的(F1)，对仿真条件中的一幅256×256lena 256的自然图像进行分解，得到系数矩阵，并得并对系数的元素实部就近取整，到它们的灰度直方图，并计算出各直方图的最大值与陡度，并与用现有剪切波第一剪切波基函数得到的灰度直方图的指标比较，结果如表1所示：

表1截断窗，psi1截断窗的灰度直方图表

B：用说明书具体实施例二方式中的步骤，对仿真条件所述的三幅512×512的自然图像lena512.bmp，barbara512.bmp，flower.bmp加高斯白噪声之后，进行硬阈值去噪，在对去噪图像中将小于零的元素全部置为0，大于255的元素全部置为255，再对其它元素进行四舍五入的整数化，得到最后去噪后的峰值信噪比结果如表2，其中用上述具体实施案例二中的步骤对lena512.bmp的去噪结果如图5，其中图5的(a)为lena512的原局部图像，(b)为加上标准差为30图像，(c)为现有剪切波基于拉普拉斯金字塔分解的阈值去噪结果图，(d)为剪切波基于圆形截断窗的阈值去噪结果图，；

表2图像的阈值去噪结果对比表

3、仿真结果分析

从表1我们可以观察出，这些值所表现出的一个共同规律就是：系数直方图所呈现的特性是高尖峰，长拖尾。从系数矩阵要求高的稀疏性来看，在内层四方向，本发明构造的截断窗比psi1截断窗不论在陡度上，还是最大值上，都有着明显的优越性，而在外层八个方向上，截断窗则在最大值上优点突出，psi1截断窗的优点则体现在陡度上。截断窗可以任意的去调节，这样就能任意的去改变系数矩阵，从而随意的改变系数直方图，这一点上，psi1截断窗不具备这样的能力。

再看参照表2与图5，表2中对三幅自然图像加不同噪声进行去噪，用基于截断窗的剪切波的阈值去噪结果明显要优于其余各方法，从峰值信噪比能明显的看到这一点；从图5c与图5d的对比可以看出，基于截断窗的剪切波去噪比传统的拉普拉斯分解的阈值去噪的效果无论从对噪声的去除还是对边缘的保持，都有明显的改进。

Claims (3)

1.一种基于矩形截断窗的剪切波图像去噪方法，包括以下步骤：

(1)对含噪图像进行离散傅里叶变换；

(2)采用如下方法构造剪切波的基函数 

2a)构造肿块函数h₁，该h₁要满足在(-2，2)区间内0到无穷阶可微，且在(-2，2)内0≤h₁≤1，在[-1，1]内h₁＝1，构造的h₁表达式如下：

其中ξ是函数的自变量

2b)构造辅助函数h_2，该h₂表达式如下：

2c)在2a)，2b)的基础上，得到左边支函数h，该h表达式如下：

h(ξ)＝h₁(ξ)h₂(ξ)

2d)构造函数右边支函数g，该函数g表达式如下：

2e)利用前面得到的左边支函数h和右边支函数g，构造剪切波基函数 

其表达式如下：

(3)利用步骤(2)构造的剪切波的基函数 构建矩形截断窗： 

3a)构建对称基函数 

如下：

其中ξ是函数的自变量，

3b)在对称基函数的基础上，构建矩形中心截断窗cut(0，sξ₁，sξ₂)为：

其中s是大于等于1的正整数，d＝0，1，

定义为：

D_d定义为： 

其中 

表示二维实数集；

3c)在基函数 

的基础上构建矩环截断窗cut(k，sξ₁，sξ₂)为：

其中k＝1，2，......s；

3d)令s＝2，将矩形中心截断窗与矩环截断窗统一表示为截断窗cut：

cut(1)＝cut(0，2ξ₁，2ξ₂)，cut(2)＝cut(1，2ξ₁，2ξ₂)，cut(3)＝cut(2，2ξ₁，2ξ₂)；

(4)将步骤(1)得到的含噪图像离散傅里叶变换与cut²(1)点乘，再进行逆离散傅里叶变换得到近似图像a₀，再用步骤(1)含噪图像进行离散傅里叶变换再分别与cut(2)和cut(3)点乘，得到带通域矩阵F₁和高通域矩阵F₂； 

(5)将剪切波 

按

重新规划为W_j，l，

其中l＝1，2，...4j；n＝2，3，......2j，

L＝[2j+1]，[2j+1]表示对2j+1取四舍五入整，

j为大于等于1的整数，m取-j到j的整数， 

(6)用重新规划的剪切波W_j，l对步骤(4)得到的带通域F₁与高通域F₂进行去噪，得到系数矩阵f′_j，l：

6a)将F₁与W_1，l点乘，得到四方向带通频谱F_1，l，并对其进行逆离散傅里叶变换IDFT(F_1，l)，得到四个截断剪切波方向窗的系数f_1，l，同理将F₂与W_2，l点乘，得到八个方向高通频谱F_2，l，并对其进行逆离散傅里叶变换IDFT(F_2，l)，得到八个截断剪切波方向窗系数f_2，l；

6b)确定阈值 

对系数矩阵进行硬阈值去噪，得到系数矩阵f′_j，l，其中 

为对应系数矩阵f_j，l的方差， 为对应系数矩阵f_j，l的噪声估计方差，这里σ_j，l，n用蒙特卡罗表示式 

进行估计，其中E()表示求期望，δ表示噪声标准差；

(7)对系数矩阵f′_j，l与近似图像a₀进行如下重构，得到去噪图像：

7a)对去噪后的系数矩阵f′_j，l进行离散傅里叶变换，再与cut(j+1)进行点乘运算，将其结果再与W_j，l做点乘运算；

7b)对7a)的处理结果进行逆离散傅里叶变换，得到rf′_j，l，将所有的rf′_j，l相加，再加上近似图像a₀，得到重构图像。

2.根据权利要求1步骤(6)所述的方法，其中步骤6b)所述的用确定的阈值τ_j，l对系数矩阵f_j，l进行去噪，是用f_j，l元素的绝对值与选取的阈值τ_j，l进行比较，如果系数矩阵f_j，l元素的绝对值大于等于τ_j，l，则该元素被保留，如果系数矩阵τ_j，l元素的绝对值小于 τ_j，l，则将该元素置零，以去除系数矩阵f_j，l所含的噪声。

3.一种基于圆形截断窗的剪切波图像去噪方法，包括以下步骤：

1)对含噪图像进行离散傅里叶变换；

2)采用如下方法构造剪切波的基函数 

2.1)构造肿块函数h₁，该h₁要满足在(-2，2)区间内0到无穷阶可微，且在(-2，2)内0≤h₁≤1，在[-1，1]内h₁＝1，构造的h₁表达式如下：

其中ξ是函数的自变量

2.2)构造辅助函数h_2，该h₂表达式如下：

2.3)在2.1)，2.2)的基础上，得到左边支函数h，该h表达式如下：

h(ξ)＝h₁(ξ)h₂(ξ)

2.4)构造函数右边支函数g，该函数g表达式如下：

2.5)利用前面得到的左边支函数h和右边支函数g，构造剪切波基函数 

其表达式如下：

3)利用步骤2)中构造的剪切波的基函数 

构建圆形截断窗：

3.1)构建对称基函数 

如下： 

其中ξ是函数的自变量，

3.2)在对称基函数的基础上，构建圆形中心截断窗cut(0，sξ₁，sξ₂)为：

其中s是大于等于1的正整数，d＝0，1，

定义为：

D_d定义为： 

其中 

表示二维实数集；

3.3)在基函数 

的基础上构建圆形环截断窗cut(k，sξ₁，sξ₂)为：

其中 

3.4)令s＝2，将圆形中心截断窗与圆形截断窗统一表示为截断窗cut：

cut(1)＝cut(0，2ξ₁，2ξ₂)，cut(2)＝cut(1，2ξ₁，2ξ₂)

4)将步骤1)得到的含噪图像离散傅里叶变换与cut²(1)点乘，再进行逆离散傅里叶变换得到近似图像a₀，再用步骤1)含噪图像进行离散傅里叶变换，再分别与cut(2)和cut(3)点乘，得到带通域矩阵F₁和高通域矩阵F₂；

5)将剪切波 

按

重新规划为W_j，l，

其中l＝1，2，...4j；n＝2，3，.......2j；

L＝[2j+1]，[2j+1]表示对2j+1取四舍五入整；

j为大于等于1的整数，m取-j到j的整数， 

6)用重新规划的剪切波W_j，l对步骤4)得到的带通域F₁与高通域F₂进行去噪，得到系数矩阵f′_j，l：

6a)将F₁与W_1，l点乘，得到四方向带通频谱F_1，l，并对其进行逆离散傅里叶变换IDFT(F_1，l)，得到四个截断剪切波方向窗的系数地f_1，l，同理将F₂与W_2，l点乘，得到八个方向高通频谱F_2，l，并对其进行逆离散傅里叶变换IDFT(F_2，l)，得到八个截断剪切波方向窗系数f_2，l；

6b)确定阈值 

对系数矩阵进行硬阈值去噪，得到系数矩阵f′_j，l，其中 

为对应系数矩阵f_j，l的方差， 

为对应系数矩阵f_j，l的噪声估计方差，这里σ_j，l，n用蒙特卡罗表示式 

进行估计，其中E()表示求期望，δ表示噪声标准差；

7)对系数矩阵f′_j，l与近似图像a₀进行重构，得到去噪图像：

7.1)对去噪后的系数矩阵f′_j，l进行离散傅里叶变换，再与cut(j+1)进行点乘运算，将其结果再与W_j，l做点乘运算；

7.2)对7.1)的处理结果进行逆离散傅里叶变换，得到rf′_j，l，将所有的rf′_j，l相加，再加上近似图像a₀，得到重构图像。 

Images