CN104346787A - 一种非局部均值图像去噪算法 - Google Patents

Abstract

本发明所述的属于图形图像处理领域，具体涉及一种非局部均值去噪算法。本发明的目的在于：针对现有技术存在的问题，提供一种非局部均值图像去噪算法。利用Pseudo-Zernike矩的旋转不变性特征不仅可以匹配到平移相似的patch，还可以匹配到旋转和对称而潜在相似的patch，从而为滤波提供更多的相似patch。进一步的使用二次高斯型核函数替代原算法中的指数型核函数来计算加权系数，从而有效解决图像去噪问题。本发明通过首先计算像素点i、p距离公式，结合像素点i的权重值，计算得到计算出i点的灰度估计值。

Description

一种非局部均值图像去噪算法

技术领域

本发明所述的属于图形图像处理领域，具体涉及一种非局部均值去噪算法。

背景技术

目前在图像形成、传输和存储过程中，不可避免地受到各种噪声的影响，图像去噪也成为图像处理领域中最基础和广泛研究的一个热门课题，也是图像处理中一个及其重要的预处理步骤，有效地抑制噪声可为图像的后续处理提供更为精准的信息。

在众多图像去噪算法中，Buades等人提出的原始非局部均值去噪算法被证明其性能优于其它经典局部去噪算法，如各向异性滤波、高斯平滑滤波、基于小波的滤波等。其基本思想是利用图像中存在的大量冗余信息，采用方形区域(patch)来描述像素间的相似性，并在整幅图像中搜索相似patch，当前像素的灰度估计值由图像中所有与它结构相似的像素加权平均得到。

由于原始非局部均值图像去噪算法匹配相似像素只利用patch的平移，并没有考虑到patch的旋转和对称变换，从而没有考虑经变换而潜在相似的patch；其次，原始算法法中使用指数型加权核函数计算权重系数，导致图像细节过度平滑而变得模糊。

发明内容

本发明的目的在于：针对现有技术存在的问题，提供一种非局部均值图像去噪算法。利用Pseudo-Zernike矩的旋转不变性特征不仅可以匹配到平移相似的patch，还可以匹配到旋转和对称而潜在相似的patch，从而为滤波提供更多的相似patch。进一步的使用二次高斯型核函数替代原算法中的指数型核函数来计算加权系数，从而有效解决图像去噪问题。

本发明目的通过下述技术方案来实现：

一种非局部均值图像去噪算法包括：

步骤1：在一幅数字图像中，计算以像素点i为中心的patch N_i，其n阶m重的pseudo-zernike矩Z_nm(i)；同理得到以像素点p为中心的patch N_p，其n阶m重的pseudo-zernike矩Z_nm(p)；

步骤2：根据步骤1中，像素点i的Z_nm(i)得到Z(i)；像素点p的Z_nm(p)得到Z(p)，进而得到像素点i、p距离公式变为：

d(i,p)＝||Z(i)-Z(p)||²   (1)；

步骤3：根据步骤2得到权重w(i,p)为：

$w(i,p)=\frac{1}{C\left(i\right)}\exp \{-{\left[\frac{{\left|\right|Z\left(i\right)-Z\left(p\right)\left|\right|}^2}{h^2}\right]}^2\}---\left(2\right);$

其中 $C\left(i\right)=\underset{p}{\mathrm{\Σ}}\exp \{-{\left[\frac{{\left|\right|Z\left(i\right)-Z\left(p\right)\left|\right|}^2}{h^2}\right]}^2\};$

步骤4：对整幅图像I内的所有点求和，计算出i点的灰度估计值：

$\mathrm{NL}\left[v\right]\left(i\right)=\underset{p\∈I}{\mathrm{\Σ}}w(i,p)v\left(i\right)---\left(3\right);$

其中v(i)为噪声图像I中被估计的像素点i的灰度值。

所述步骤1具体步骤包括：

步骤11：在一幅数字图像中，i的坐标点为f(x_i,y_i)，计算一个以像素点i为中心的patch N_i的pseudo-zernike矩时，把该patch N_i映射到单位圆内，圆心为像素点i，像素由直角坐标变为角坐标，pseudo-zernike是一种正交复数矩，它由在单位圆x²+y²≤1内的一组完备正交基{V_nm(x_i,y_i)}得到，其定义为：

V_nm(x_i,y_i)＝V_nm(ρ,θ)＝R_nm(ρ)e^jmθ   (4)

其中，n为非负整数，m为整数，且满足|m|≤n；ρ和θ为极坐标下像素的半径和角度，S表示0到n-m的任意整数值；patch N_i表示以i为中心点的方形区域，其中R_nm(ρ)为径向多项式：

$R_{\mathrm{nm}}\left(\mathrm{\ρ}\right)=\underset{s=0}{\overset{n-\left|m\right|}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}---\left(5\right);$

步骤12：对于数字图像中一个以像素点_i为中心的patch N_i，其n阶m重的pseudo-zernike矩定义为：

$Z_{\mathrm{nm}}\left(i\right)=\frac{n+1}{\mathrm{\π}}\underset{x_i}{\mathrm{\Σ}}\underset{y_i}{\mathrm{\Σ}}{V_{\mathrm{nm}}}^{*}(x_i,y_i)f(x_i,y_i)=\frac{n+1}{\mathrm{\π}}\underset{x_i}{\mathrm{\Σ}}\underset{y_i}{\mathrm{\Σ}}\underset{s=0}{\overset{n-\left|m\right|}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}{e}^{\mathrm{jm\θ}}*f(x_i,y_i)---\left(6\right);$

步骤13：在一幅数字图像中，p的坐标点为f(x_p,y_p)，计算一个以像素点p为中心的patch N_p的pseudo-zernike矩时，把该patch N_p映射到单位圆内，圆心为像素点p，像素由直角坐标变为角坐标，pseudo-zernike是一种正交复数矩，它由在单位圆x²+y²≤1内的一组完备正交基{V_nm(x_p,y_p)}得到，其定义为：

V_nm(x_p,y_p)＝V_nm(ρ,θ)＝R_nm(ρ)e^jmθ   (7)

其中，n为非负整数，m为整数，且满足|m|≤n；ρ和θ为极坐标下像素的半径和角度，S表示从0到n-|m|的全部整数值；patch N_p表示以p为中心点的方形区域，其中R_nm(ρ)为径向多项式：

$R_{\mathrm{nm}}\left(\mathrm{\ρ}\right)=\underset{s=0}{\overset{n-\left|m\right|}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}---\left(8\right);$

步骤14：根据步骤13，对于数字图像中一个以像素点p为中心的patch N_p，其n阶m重的pseudo-zernike矩定义为：

$Z_{\mathrm{nm}}\left(\text{p}\right)=\frac{n+1}{\mathrm{\π}}\underset{x_p}{\mathrm{\Σ}}\underset{y_p}{\mathrm{\Σ}}{V_{\mathrm{nm}}}^{*}(x_p,y_p)f(x_p,y_p)=\frac{n+1}{\mathrm{\π}}\underset{x_p}{\mathrm{\Σ}}\underset{y_p}{\mathrm{\Σ}}\underset{s=0}{\overset{n-\left|m\right|}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}{e}^{\mathrm{jm\θ}}*f(x_p,y_p)---\left(9\right).$

所述步骤2具体包括：

步骤21：当n＝3时，计算以像素i为中心的patch N_i的6个pseudo-zernike矩：Z₀₀(i)、Z₁₁(i)、Z₂₀(i)、Z₂₂(i)、Z₃₁(i)、Z₃₃(i)，令像素i的特征向量Z(i)＝{Z₀₀(i)，Z₁₁(i)，Z₂₀(i)，Z₂₂(i)，Z₃₁(i)，Z₃₃(i)}；

步骤22：当n＝3时，计算以像素p为中心的patch N_i的6个pseudo-zernike矩：Z₀₀(p)、Z₁₁(p)、Z₂₀(p)、Z₂₂(p)、Z₃₁(p)、Z₃₃(p)，令像素p的特征向量Z(p)＝{Z₀₀(p)，Z₁₁(p)，Z₂₀(p)，Z₂₂(p)，Z₃₁(p)，Z₃₃(p)}；

步骤23：根据d(i,p)＝||Z(i)-Z(p)||²计算得到像素点i、p距离公式。

所述步骤2具体包括：

步骤21：当n＝2时，计算以像素i为中心的patch N_i的4个pseudo-zernike矩：Z₀₀(i)、Z₁₁(i)、Z₂₀(i)、Z₂₂(i)，令像素i的特征向量Z(i)＝{Z₀₀(i)，Z₁₁(i)，Z₂₀(i)，Z₂₂(i)}；

步骤22：当n＝2时，计算以像素p为中心的patch N_i的6个pseudo-zernike矩：Z₀₀(p)、Z₁₁(p)、Z₂₀(p)、Z₂₂(p)，令像素p的特征向量Z(p)＝{Z₀₀(p)，Z₁₁(p)，Z₂₀(p)，Z₂₂(p)}；

步骤23：根据d(i,p)＝||Z(i)-Z(p)||²计算得到像素点i、p距离公式。

本发明的有益效果：使用Pseudo-Zernike矩代替高斯加权欧氏距离来度量patch之间的相似度，利用Pseudo-Zernike矩的旋转不变性特征不仅可以匹配到平移相似的patch，还可以匹配到旋转和对称而潜在相似的patch，从而为滤波提供更多的相似patch。另外，本发明使用二次高斯型核函数替代原算法中的指数型核函数来计算加权系数，这样可以使得相似度高的patch获得更大的权重，相似度低的patch获得较小的权重，使得本发明比原始算法去噪后的视觉效果更好，尤其是对图像细节和纹理的保持，而且去噪后的峰值信噪比也高于原始算法。

附图说明

图1是本发明实施例以像素i为中心的3X3的patchN_i；

图2是指数型、高斯型核二次高斯型核函数权重随距离变换的加权特性曲线。

具体实施方式

下列非限制性实施例用于说明本发明。

本发明目的是：从被污染的图像中有效的去掉噪声，恢复原始图像；

在计算过程中，某一个像素的灰度估计值应该利用图像中的所有点来恢复，但是离该点越远贡献越小，因此，把搜索范围限制在某个RXR的方形区域，比如在计算第11行11列像素点的灰度值时，把搜索区域限制在R*R＝21*21的方形区域，由该区域内的所有像素点来恢复该点的估计值。

实施例1：

1、原始非局部均值去噪算法

1)计算高斯加权欧式距离：d(i,j)＝||v(N_i)-v(N_j)||² _2,a，V(Ni)表示以像素点i为中心的patch里所有点的灰度值组成的向量，向量v(Nj)表示以任意一个点为中心的patch所组成的灰度值向量。

2)计算w(i,j)：利用公式计数像素点j对像素点i贡献的权值。权值w(i,j)取决于像素i和j之间的相似度，像素i和j越相似，像素j对待恢复像素i的贡献越大，满足0≤w(i,j)≤1且像素i和j之间的相似度由灰度值向量v(N_i)与v(N_j)之间的相似度决定，其中N_i表示以像素i为中心的固定大小(2K+1)×(2K+1)的方形patch。

3)计算i点的灰度估计值：利用对整幅图像I内的所有点求和。具体过程是：给定一幅噪声图像v＝{v(i)|i∈I}，I表示图像域，由图像中的全部像素构成，对于图像中任何一个像素i，NLM算法利用该图像中所有像素灰度值的加权平均来得到i点的灰度估计值NL[v](i)，即：

$\mathrm{NL}\left[v\right]\left(i\right)=\underset{j\∈I}{\mathrm{\Σ}}w(i,j)v\left(i\right).$

二、本发明去噪详细过程：

实施例一：步骤1：步骤1：在一幅数字图像中，i的坐标点为f(x_i,y_i)，计算一个以像素点i为中心的patch N_i的pseudo-zernike矩时，把该patch N_i映射到单位圆内，圆心为像素点i，像素由直角坐标变为角坐标，pseudo-zernike是一种正交复数矩，它由在单位圆x²+y²≤1内的一组完备正交基{V_nm(x,y)}得到，其定义为：

V_nm(x_i,y_i)＝V_nm(ρ,θ)＝R_nm(ρ)e^jmθ   (1)

其中，n为非负整数，m为整数，且满足|m|≤n；ρ和θ为极坐标下像素的半径和角度，S表示0到n-m的任意数值；其中R_nm(ρ)为径向多项式：

$R_{\mathrm{nm}}\left(\mathrm{\ρ}\right)=\underset{s=0}{\overset{n-\left|m\right|}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}---\left(2\right);$

步骤2：对于数字图像中一个以像素点i为中心的patch N_i，其n阶m重的pseudo-zernike矩定义为：

$Z_{\mathrm{nm}}\left(i\right)=\frac{n+1}{\mathrm{\π}}\underset{x_i}{\mathrm{\Σ}}\underset{y_i}{\mathrm{\Σ}}{V_{\mathrm{nm}}}^{*}(x_i,y_i)f(x_i,y_i)---\left(3\right);$

步骤3：以像素点i为中心的patchN_i的pseudo-zernike伪泽尼克矩计算，本算法计算了伪泽尼克的n＝3阶矩。如图1所述是一个以像素i为中心的3X3的patch，把该patchN_i映射到单位圆内，在圆上的像素点到圆心的角度为45度，135度，225度，315度，半径都是1，而在圆内的像素点半径都为二分之根号二，角度分别为0度，90度，180度，270度。计算出以像素_i为中心的patch N_i的6个pseudo-zernike矩：Z₀₀(i)、Z₁₁(i)、Z₂₀(i)、Z₂₂(i)、Z₃₁(i)、Z₃₃(i)，令像素i的特征向量Z(i)＝{Z₀₀(i)，Z₁₁(i)，Z₂₀(i)，Z₂₂(i)，Z₃₁(i)，Z₃₃(i)}；

步骤4：重复步骤1-3，以像素p为中心的patch N_p，当n＝3时，计算以像素p为中心的patch N_i的6个pseudo-zernike矩：Z₀₀(p)、Z₁₁(p)、Z₂₀(p)、Z₂₂(p)、Z₃₁(p)、Z₃₃(p)，令像素p的特征向量Z(p)＝{Z₀₀(p)，Z₁₁(p)，Z₂₀(p)，Z₂₂(p)，Z₃₁(p)，Z₃₃(p)}。

步骤5：把原始算法中灰度值向量v(N_i)与v(N_p)替换为特征向量Z(i)与Z(p)，则N_i与N_p之间的相似度由向量Z(i)与Z(p)间去掉高斯加权的欧式距离来描述，则像素点i、p距离公式变为：

d(i,p)＝||Z(i)-Z(p)||²   (4)；

步骤6：权重系数是影响去噪性能最直接的因素，而权重又主要由核函数的选取和相似度的计算方法决定。相似度高的patch应赋予较大的权值，而不相似的patch应得到较小的权值，甚至为0。因此，如何选取合适的核函数就成了影响去噪的关键，本发明对多种形式的加权核函数进行了分析和对比，在高斯型核函数的基础上，提出了二次高斯型核函数：

$w(i,p)=e^{-{\left[\frac{d^2(i,p)}{h^2}\right]}^2}---\left(5\right);$

其中h是常数，表示滤波参数，取值在30～90；

由图2可以看出，二次高斯型核函数在距离小的区间更为平坦，随着距离的增大而迅速衰减，即相似性大的patch可以得到更大的权值，不相似的patch被赋予更小的权值，从而有效地弥补了原始算法中指数型核函数加权不足和过度加权的缺点，减轻了不相似patch对去噪性能的影响。

综上，由式(4)和(5)权重w(i,j)可重新定义为：

$w(i,j)=\frac{1}{C\left(i\right)}\exp \{-{\left[\frac{{\left|\right|Z\left(i\right)-Z\left(j\right)\left|\right|}^2}{h^2}\right]}^2\}---\left(6\right);$

其中 $C\left(i\right)=\underset{j}{\mathrm{\Σ}}\exp \{-{\left[\frac{{\left|\right|Z\left(i\right)-Z\left(j\right)\left|\right|}^2}{h^2}\right]}^2\};$

步骤8：计算i点的灰度估计值：

$\mathrm{NL}\left[v\right]\left(i\right)=\underset{p\∈I}{\mathrm{\Σ}}w(i,p)v\left(i\right)---\left(7\right);$

对整幅图像I内的所有点求和，其中v(i)为噪声图像中被估计的像素点。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (4)

1.一种非局部均值图像去噪算法，其特征在于包括：

步骤1：在一幅数字图像中，计算以像素点i为中心的patch N_i，其n阶m重的pseudo-zernike矩Z_nm(i)；同理得到以像素点p为中心的patch N_p，其n阶m重的pseudo-zernike矩Z_nm(p)；patch N_i指的是以像素点i为中心的方形区域；patch N_p指的是以像素点p为中心的方形区域；

步骤2：根据步骤1中，像素点i的Z_nm(i)得到Z(i)；像素点p的Z_nm(p)得到Z(p)，进而得到像素点i、p距离公式为：

d(i,p)＝||Z(i)-Z(p)||²                    (1)；

步骤3：根据步骤2得到权重w(i,p)为：

$w(i,p)=\frac{1}{C\left(i\right)}\exp \{-{\left[\frac{{\left|\right|Z\left(i\right)-Z\left(p\right)\left|\right|}^2}{h^2}\right]}^2\}---\left(2\right);$

其中 $C\left(i\right)=\underset{p}{\mathrm{\Σ}}\exp \{-{\left[\frac{{\left|\right|Z\left(i\right)-Z\left(p\right)\left|\right|}^2}{h^2}\right]}^2\};$

步骤4：对整幅图像I内的所有点求和，计算出i点的灰度估计值：

$\mathrm{NL}\left[v\right]\left(i\right)=\underset{p\∈I}{\mathrm{\Σ}}w(i,p)v\left(i\right)---\left(3\right);$

其中v(i)为噪声图像I中被估计的像素点值i的灰度值。

2.根据权利要求1所述的一种非局部均值图像去噪算法，其特征在于所述步骤1具体步骤包括：

步骤11：在一幅数字图像中，i的坐标点为f(x_i,y_i)，计算一个以像素点i为中心的patch N_i的pseudo-zernike矩时，把该patch N_i映射到单位圆内，圆心为像素点i，像素由直角坐标变为角坐标，pseudo-zernike是一种正交复数矩，它由在单位圆x²+y²≤1内的一组完备正交基{V_nm(x_i,y_i)}得到，其定义为：

V_nm(x_i,y_i)＝V_nm(ρ,θ)＝R_nm(ρ)e^jmθ           (4)

其中，n为非负整数，m为整数，且满足|m|≤n；ρ和θ为极坐标下像素的半径和角度，S表示从0到n-|m|的全部整数值；patch N_i表示以i为中心点的方形区域，其中R_nm(ρ)为径向多项式：

$R_{\mathrm{nm}}\left(\mathrm{\ρ}\right)=\underset{s=0}{\overset{n-\left|m\right|}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}---\left(5\right);$

步骤12：对于数字图像中一个以像素点i为中心的patch N_i，其n阶m重的pseudo-zernike矩定义为：

$Z_{\mathrm{nm}}\left(i\right)=\frac{n+1}{\mathrm{\π}}\underset{x_i}{\mathrm{\Σ}}\underset{y_i}{\mathrm{\Σ}}{V_{\mathrm{nm}}}^{*}(x_i,y_i)f(x_i,y_i)=\frac{n+1}{\mathrm{\π}}\underset{x_i}{\mathrm{\Σ}}\underset{y_i}{\mathrm{\Σ}}\underset{s=0}{\overset{n-\left|m\right|}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}{e}^{\mathrm{jm\θ}}*f(x_i,y_i)---\left(6\right);$

步骤13：在一幅数字图像中，p的坐标点为f(x_p,y_p)，计算一个以像素点p为中心的patch N_p的pseudo-zernike矩时，把该patch N_p映射到单位圆内，圆心为像素点p，像素由直角坐标变为角坐标，pseudo-zernike是一种正交复数矩，它由在单位圆x²+y²≤1内的一组完备正交基{V_nm(x_p,y_p)}得到，其定义为：

V_nm(x_p,y_p)＝V_nm(ρ,θ)＝R_nm(ρ)e^jmθ                     (7)

其中，n为非负整数，m为整数，且满足|m|≤n；ρ和θ为极坐标下像素的半径和角度，S表示从0到n-|m|的全部整数值；patch N_p表示以p为中心点的方形区域，其中R_nm(ρ)为径向多项式：

$R_{\mathrm{nm}}\left(\mathrm{\ρ}\right)=\underset{s=0}{\overset{n-\left|m\right|}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}---\left(8\right);$

步骤14：根据步骤13，对于数字图像中一个以像素点p为中心的patch N_p，其n阶m重的pseudo-zernike矩定义为：

$Z_{\mathrm{nm}}\left(p\right)=\frac{n+1}{\mathrm{\π}}\underset{x_p}{\mathrm{\Σ}}\underset{y_p}{\mathrm{\Σ}}{V_{\mathrm{nm}}}^{*}(x_p,y_p)f(x_p,y_p)=\frac{n+1}{\mathrm{\π}}\underset{x_p}{\mathrm{\Σ}}\underset{y_p}{\mathrm{\Σ}}\underset{s=0}{\overset{n-\left|m\right|}{\mathrm{\Σ}}}\frac{{(-1)}^s(2n+1)!{\mathrm{\ρ}}^{n-s}}{s!(n+|m|+1-s)!(n-|m|-s)!}{e}^{\mathrm{jm\θ}}*f(x_p,y_p)---\left(9\right).$

3.根据权利要求1所述的一种非局部均值图像去噪算法，其特征在于所述步骤2具体包括：

步骤21：当n＝3时，计算以像素i为中心的patch N_i的6个pseudo-zernike矩：Z₀₀(i)、Z₁₁(i)、Z₂₀(i)、Z₂₂(i)、Z₃₁(i)、Z₃₃(i)，令像素i的特征向量Z(i)＝{Z₀₀(i)，Z₁₁(i)，Z₂₀(i)，Z₂₂(i)，Z₃₁(i)，Z₃₃(i)}；

步骤22：当n＝3时，计算以像素p为中心的patch N_i的6个pseudo-zernike矩：Z₀₀(p)、Z₁₁(p)、Z₂₀(p)、Z₂₂(p)、Z₃₁(p)、Z₃₃(p)，令像素p的特征向量Z(p)＝{Z₀₀(p)，Z₁₁(p)，Z₂₀(p)，Z₂₂(p)，Z₃₁(p)，Z₃₃(p)}；

步骤23：根据d(i,p)＝||Z(i)-Z(p)||²计算得到像素点i、p距离公式。

4.根据权利要求1所述的一种非局部均值图像去噪算法，其特征在于所述步骤2具体包括：

步骤21：当n＝2时，计算以像素i为中心的patch N_i的4个pseudo-zernike矩：Z₀₀(i)、Z₁₁(i)、Z₂₀(i)、Z₂₂(i)，令像素i的特征向量Z(i)＝{Z₀₀(i)，Z₁₁(i)，Z₂₀(i)，Z₂₂(i)}；

步骤22：当n＝2时，计算以像素p为中心的patch N_i的6个pseudo-zernike矩：Z₀₀(p)、Z₁₁(p)、Z₂₀(p)、Z₂₂(p)，令像素p的特征向量Z(p)＝{Z₀₀(p)，Z₁₁(p)，Z₂₀(p)，Z₂₂(p)}；

步骤23：根据d(i,p)＝||Z(i)-Z(p)||²计算得到像素点i、p距离公式。