文章目录
前言
之前提到图像映射变换必须要用的插值法，其作用是使得图像像素更加平滑些。
今天将提到的算法也具有平滑像素的作用，主要用于图像增强、图像降噪。那么，为什么要降噪呢？
首先，我们需要知道图像的噪声是什么。在图像传感器成像过程中，光电转换及数模放大时，不可避免地会产生噪声，在图像传输过程中，也将引入二次噪声。比如我们平时看到的视频中的“雪花”，就是一种噪声。
本次实验将采用滤波算法对图像进行降噪操作，本质上是使用一个类似于滑动窗口的滤波器核（滤波器核大小势必大于单个噪点所占像素地大小）对图像进行滤波运算。因此，我们需要先看看什么是空间滤波运算。
一、滤波运算——相关运算与卷积运算
相关运算，是利用卷积核对图像进行邻近操作的，将卷积核的重心移动到待处理地像素点，对卷积核对应地像素点加权求和。
卷积运算，也是利用卷积核对图像进行邻域操作的，只是把相关运算的卷积核翻转了180度。
#空间相关运算与卷积运算(根据上述几张图,此实验例子的滤波器核的高和宽正好相等)
import numpy as np
#创建一个图像矩阵，整型
image = np.zeros((5,5))
image[2,2] = 1
# print(image)
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.matrix([[1,2,3],[4,5,6],[7,8,9]],dtype=float)
# print(filter_correlation)
h_f,w_f = filter.shape
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image(image,h_f,w_f):
h_img,w_img = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new))
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2)] = image[i,j]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间相关运算
def correlation_operation(image,filter):
h_img,w_img = image.shape
h_f,w_f = filter.shape
image_out = np.zeros((h_img,w_img))
image_padding = padding_image(image,h_f,w_f)
for i in range(h_img):
for j in range(w_img):
image_out[i,j] = np.sum(np.multiply(filter[:,:],image_padding[i:i+h_f,j:j+w_f]))
return image_out
image_correlation = correlation_operation(image,filter)
# print(image_correlation)
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img))
image_padding = padding_image(image,h_f,w_f)
for i in range(h_img):
for j in range(w_img):
image_out[i,j] = np.sum(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f]))
return image_out
image_convolution = convolution_operation(image,filter)
# print(image_convolution)
[[0. 0. 0. 0. 0.]
[0. 1. 2. 3. 0.]
[0. 4. 5. 6. 0.]
[0. 7. 8. 9. 0.]
[0. 0. 0. 0. 0.]]
二、盒式滤波器（均值滤波算法）
盒式滤波器，又称作算术均值滤波。盒式核的所有像素权重都相同，它是以一个像素值全是1的滤波器核除以滤波器总像素个数为原型的。
f
(
x
,
y
)
=
1
m
∗
n
∗
∑
(
r
,
c
)
∈
S
x
y
g
(
r
,
c
)
f(x,y) =\frac{1}{m*n} * \sum_ {(r,c) \in S_{xy}} g(r,c)
f(x,y)=m∗n1∗(r,c)∈Sxy∑g(r,c)
其中，
S
x
y
S_{xy}
Sxy表示中心为
(
x
,
y
)
(x,y)
(x,y)、大小为
m
×
n
m \times n
m×n的矩形子图像窗口（邻域）的一组坐标。
f
(
x
,
y
)
f(x,y)
f(x,y)表示经过滤波算法后的输出图像像素值，
g
(
x
,
y
)
g(x,y)
g(x,y)表示原图像素值,r和c是邻域
S
x
y
S_{xy}
Sxy中包含的像素的行坐标和列坐标。
#盒式滤波器——实验一：创建图像
import numpy as np
#创建一个图像矩阵，整型
image = np.zeros((5,5))
image[2,2] = 1
# print(image)
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.matrix([[1,1,1],[1,1,1],[1,1,1]],dtype=float)
# print(filter_correlation)
h_f,w_f = filter.shape
filter[:,:] = filter[:,:] / (h_f * w_f)
print(filter)
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image(image,h_f,w_f):
h_img,w_img = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new))
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2)] = image[i,j]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img))
image_padding = padding_image(image,h_f,w_f)
for i in range(h_img):
for j in range(w_img):
image_out[i,j] = np.sum(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f]))
return image_out
image_convolution = convolution_operation(image,filter)
print(image_convolution)
[[0.11111111 0.11111111 0.11111111]
[0.11111111 0.11111111 0.11111111]
[0.11111111 0.11111111 0.11111111]]
[[0. 0. 0. 0. 0. ]
[0. 0.11111111 0.11111111 0.11111111 0. ]
[0. 0.11111111 0.11111111 0.11111111 0. ]
[0. 0.11111111 0.11111111 0.11111111 0. ]
[0. 0. 0. 0. 0. ]]
#盒式滤波器——实验二：选择图像(对彩色的滤波运算/卷积运算，必须考虑其颜色通道，需要对三个通道分别进行卷积),3*3的盒式核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.matrix([[1,1,1],[1,1,1],[1,1,1]],dtype=float)
# print(filter_correlation)
h_f,w_f = filter.shape
filter[:,:] = filter[:,:] / (h_f * w_f)
# print(filter)
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = np.sum(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]))
return image_out
image_convolution = convolution_operation(image,filter)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
#盒式滤波器——实验三：选择图像(对彩色的滤波运算/卷积运算，必须考虑其颜色通道，需要对三个通道分别进行卷积),7*7的盒式核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.ones((7,7))
# print(filter_correlation)
h_f,w_f = filter.shape
filter[:,:] = filter[:,:] / (h_f * w_f)
# print(filter)
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = np.sum(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]))
return image_out
image_convolution = convolution_operation(image,filter)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
我们会发现，随着盒式核的大小在增大，图像开始模糊了。但另一方面，这个过程中图像像素之间也更平滑了。我们接下来再试一个
11
×
11
11 \times 11
11×11的盒式核。
#盒式滤波器——实验四：选择图像(对彩色的滤波运算/卷积运算，必须考虑其颜色通道，需要对三个通道分别进行卷积),11*11的盒式核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.ones((11,11))
# print(filter_correlation)
h_f,w_f = filter.shape
filter[:,:] = filter[:,:] / (h_f * w_f)
# print(filter)
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = np.sum(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]))
return image_out
image_convolution = convolution_operation(image,filter)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
三、其他均值滤波
1.几何均值滤波
几何均值滤波实际上是将图像在某一邻域内（假设该邻域大小为
m
×
n
m \times n
m×n）的诸多像素值求乘积，再求这个乘积的
1
m
∗
n
\frac{1}{m*n}
m∗n1次方。大致公式如下：
f
(
x
,
y
)
=
[
∏
(
r
,
c
)
∈
S
x
y
g
(
r
,
c
)
]
1
m
∗
n
f(x,y) = [\prod_ {(r,c) \in S_{xy}} g(r,c)]^\frac{1}{m*n}
f(x,y)=[(r,c)∈Sxy∏g(r,c)]m∗n1
其中，
∏
\prod_{}
∏表示相乘，
S
x
y
S_{xy}
Sxy表示中心为
(
x
,
y
)
(x,y)
(x,y)、大小为
m
×
n
m \times n
m×n的矩形子图像窗口（邻域）的一组坐标。
f
(
x
,
y
)
f(x,y)
f(x,y)表示经过滤波算法后的输出图像像素值，
g
(
x
,
y
)
g(x,y)
g(x,y)表示原图像素值,r和c是邻域
S
x
y
S_{xy}
Sxy中包含的像素的行坐标和列坐标。
#几何均值滤波算法：3*3的滤波器核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.ones((3,3))
# print(filter_correlation)
h_f,w_f = filter.shape
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = np.power(np.prod(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k])),1/(h_f * w_f))
return image_out
image_convolution = convolution_operation(image,filter)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
#几何均值滤波算法：7*7的滤波器核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.ones((7,7))
# print(filter_correlation)
h_f,w_f = filter.shape
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = np.power(np.prod(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k])),1/(h_f * w_f))
return image_out
image_convolution = convolution_operation(image,filter)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
2、谐波均值滤波
f ( x , y ) = m ∗ n ∑ ( r , c ) ∈ S x y ∗ 1 g ( r , c ) f(x,y) = \frac {m*n}{\sum_ {(r,c) \in S_{xy}} * \frac {1}{g(r,c)}} f(x,y)=∑(r,c)∈Sxy∗g(r,c)1m∗n
#谐波均值滤波算法：3*3的滤波器核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.ones((3,3))
# print(filter_correlation)
h_f,w_f = filter.shape
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = (h_f * w_f) / np.sum(1 / np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]))
return image_out
image_convolution = convolution_operation(image,filter)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
#谐波均值滤波算法：11*11的滤波器核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.ones((11,11))
# print(filter_correlation)
h_f,w_f = filter.shape
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = (h_f * w_f) / np.sum(1 / np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]))
return image_out
image_convolution = convolution_operation(image,filter)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
3、反谐波均值滤波
f
(
x
,
y
)
=
∑
(
r
,
c
)
∈
S
x
y
g
(
r
,
c
)
Q
+
1
∑
(
r
,
c
)
∈
S
x
y
g
(
r
,
c
)
Q
f(x,y) = \frac {\sum_{(r,c) \in S_{xy}} g(r,c)^{Q+1}}{\sum_{(r,c) \in S_{xy}} g(r,c)^{Q}}
f(x,y)=∑(r,c)∈Sxyg(r,c)Q∑(r,c)∈Sxyg(r,c)Q+1
其中Q称为滤波器的阶数。这种滤波器适用于降低或消除椒盐噪声。Q值为正时，该滤波器消除胡椒噪声；Q值为负时，该滤波器消除盐粒噪声。然而，该滤波器不能同时消除这两种噪声。注意，当
Q
=
0
Q=0
Q=0时，反谐波均值滤波器简化为算术均值滤波器；
Q
=
−
1
Q=-1
Q=−1时，简化为谐波均值滤波。
#反谐波均值滤波算法：3*3的滤波器核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.ones((3,3))
# print(filter_correlation)
h_f,w_f = filter.shape
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter,Q):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = np.sum(np.power(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]),Q+1)) / np.sum(np.power(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]),Q))
return image_out
image_convolution = convolution_operation(image,filter,1)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
注意：numpy.power([2,3,4],2)是分别求2、3、4的2次方，numpy.power(2,[2,3,4])是分别求2的2、3、4次方。
四、统计排序滤波器
统计排序滤波器是空间滤波器，其响应基于滤波器所围邻域中的像素值的顺序，排序结果决定了滤波器的响应。一般包括：中值滤波器、最大值滤波器、最小值滤波器、重点滤波器和修正阿尔法值均值滤波器。
1、中值滤波器
中值滤波器属于非线性滤波器。它能有效降低某些随机噪声，且模糊较小，但运算时间很长。
f
(
x
,
y
)
=
m
e
d
i
a
n
(
r
,
c
)
∈
S
x
y
{
g
(
r
,
c
)
}
f(x,y) = {median}_{(r,c) \in S_{xy}}\{g(r,c)\}
f(x,y)=median(r,c)∈Sxy{g(r,c)}
#中值滤波,3*3滤波器核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.ones((3,3))
# print(filter_correlation)
h_f,w_f = filter.shape
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image(image,h_f,w_f):
h_img,w_img = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new))
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2)] = image[i,j]
return image_padding
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = np.median(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]))
return image_out
image_convolution = convolution_operation(image,filter)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
2、最大值滤波器
f ( x , y ) = max ( r , c ) ∈ S x y { g ( r , c ) } f(x,y) = \max_{(r,c) \in S_{xy}}\{g(r,c)\} f(x,y)=(r,c)∈Sxymax{g(r,c)}
#最大值滤波,3*3滤波器核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.ones((3,3))
# print(filter_correlation)
h_f,w_f = filter.shape
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image(image,h_f,w_f):
h_img,w_img = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new))
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2)] = image[i,j]
return image_padding
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = np.max(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]))
return image_out
image_convolution = convolution_operation(image,filter)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
3、最小值滤波器
f ( x , y ) = min ( r , c ) ∈ S x y { g ( r , c ) } f(x,y) = \min_{(r,c) \in S_{xy}}\{g(r,c)\} f(x,y)=(r,c)∈Sxymin{g(r,c)}
#最小值滤波,3*3滤波器核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.ones((3,3))
# print(filter_correlation)
h_f,w_f = filter.shape
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image(image,h_f,w_f):
h_img,w_img = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new))
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2)] = image[i,j]
return image_padding
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = np.min(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]))
return image_out
image_convolution = convolution_operation(image,filter)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
4、中点滤波器
f
(
x
,
y
)
=
(
min
(
r
,
c
)
∈
S
x
y
{
g
(
r
,
c
)
}
+
min
(
r
,
c
)
∈
S
x
y
{
g
(
r
,
c
)
}
)
2
f(x,y) = \frac {(\min_{(r,c) \in S_{xy}}\{g(r,c)\} + \min_{(r,c) \in S_{xy}}\{g(r,c)\})}{2}
f(x,y)=2(min(r,c)∈Sxy{g(r,c)}+min(r,c)∈Sxy{g(r,c)})
中点滤波器是统计排序滤波器和均值滤波器的结合，适合处理随机分布的噪声，如高斯噪声和均匀噪声。注意中点的值median一般是不同的。这里的中点指的是图像像素在某一邻域内最大值和最小值的均值。
#中点滤波,3*3滤波器核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.ones((3,3))
# print(filter_correlation)
h_f,w_f = filter.shape
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image(image,h_f,w_f):
h_img,w_img = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new))
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2)] = image[i,j]
return image_padding
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算
def convolution_operation(image,filter):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = (np.max(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k])) + np.min(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]))) / 2
return image_out
image_convolution = convolution_operation(image,filter)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
5、修正阿尔法均值滤波器
假设我们要在邻域
S
x
y
S_{xy}
Sxy内删除
g
(
r
,
c
)
g(r,c)
g(r,c)的
d
2
\frac{d}{2}
2d个最低灰度值和
d
2
\frac{d}{2}
2d个最高灰度值。令
g
R
(
r
,
c
)
g_{R}(r,c)
gR(r,c)表示
S
x
y
S_{xy}
Sxy中剩下的
m
×
n
−
d
m \times n-d
m×n−d个像素。通过平均这些剩余像素所形成的滤波器，成为修正阿尔法均值滤波器。其形式为：
f
(
x
,
y
)
=
1
m
∗
n
−
d
∗
∑
(
r
,
c
)
∈
S
x
y
g
R
(
r
,
c
)
f(x,y) = \frac {1}{m*n-d} * \sum_{(r,c) \in S_{xy}} g_{R}(r,c)
f(x,y)=m∗n−d1∗(r,c)∈Sxy∑gR(r,c)
其中，d的取值范围是从0到
m
×
n
−
1
m \times n-1
m×n−1。当
d
=
0
d=0
d=0时修正阿尔法均值滤波器简化为算术均值滤波器。如果
d
=
m
×
n
−
1
d=m \times n-1
d=m×n−1，那么该滤波器变为中值滤波器。d去其他值时，修正阿尔法均值滤波器适合于处理多种混合噪声，如高斯噪声和椒盐噪声。（d一定是偶数）
#修正阿尔法均值滤波器，3*3滤波器核
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
filter = np.ones((3,3))
# print(filter_correlation)
h_f,w_f = filter.shape
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image(image,h_f,w_f):
h_img,w_img = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new))
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2)] = image[i,j]
return image_padding
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算,注意：这里的d一定是偶数
def convolution_operation(image,filter,d):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = np.sum(np.sort(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]),axis=None)[int(d/2):h_f*w_f-int(d/2)]) / (h_f * w_f - d)
return image_out
image_convolution = convolution_operation(image,filter,2)
# print(image_convolution)
image_out = np.uint8(image_convolution)
plt.title("image_convolution")
plt.imshow(image_out,cmap="brg")
plt.show()
以上提到的均值滤波器以及统计排序滤波器都是低通滤波器，又称为平滑空间滤波器，它主要用于降低灰度的急剧过渡。由于随机噪声通常由灰度的急剧过渡组成，因此平滑的一个明显应用就是降噪。它的另一个应用时平滑灰度级数量不足导致的图像中的伪轮廓。
换言之，低通滤波器衰减或去除高频而通过低频，这样会模糊图像的边缘信息。
还有一种低通滤波器叫做高斯滤波器。
五、高斯滤波器
均值滤波和统计排序滤波，对于滤波窗口内每个像素的权重是一样的。但是噪声在图像当中往往比较突兀，即和原图有较大的差异，那么它必然是不是平均分布的。
我们先看看正态分布是什么。
正态分布是概率分布中的一种，由于历史的原因，正态分布也成为高斯分布。在正态分布里，认为中间状态是常态，过高和过低都属于少数，因此正态分布具有很高的普遍性，比如我们的身高、寿命、血压、成绩、测量误差等都遵循正要分布。
当然，自然界也有部呈正态分布的例子，如财富分布，最初可能是正态分布的，但是优质的资源总是掌握在少数人手里，是大部分人共享了少量的资源，结果导致富者越富，贫者越贫，这就是所谓的“马太效应（两极分化现象）”。
我们要处理的噪声包括高斯噪声和椒盐噪声。高斯噪声顾名思义，就是服从正态分布（高斯分布）的噪声，高斯噪声是画面上每个点都存在着不同程度的，与当前像素成高斯分布的噪声。
由中心极限定理可知，如果样本量足够大，样本的均值的分布慢慢变成正态分布（中心极限定理是指给定一个任意分布的总体，首先，每次从这些总体中随机抽取n个样本，一共抽取m次，然后对着m组样本分别求出平均值，这些平均值的分布接近正态分布）。那么噪声的分布也应该符合正态分布（高斯白噪声），所以就有了高斯滤波算法。
接下来我们看看两个高斯分布函数。
一维高斯分布函数：
G
(
r
)
=
1
2
π
σ
e
−
r
2
2
σ
2
G(r) = \frac {1}{\sqrt{2\pi}\sigma}e^{-\frac {r^2}{2\sigma^2}}
G(r)=2πσ1e−2σ2r2
二维高斯分布函数：
G
(
r
)
=
1
2
π
σ
2
e
−
r
2
2
σ
2
G(r) = \frac {1}{2\pi\sigma^2}e^{-\frac {r^2}{2\sigma^2}}
G(r)=2πσ21e−2σ2r2
其中r表示高斯滤波器核的中心点到核上任意一点的距离，
σ
\sigma
σ时高斯核的标准差（尺度因子），
σ
\sigma
σ值越大，数据越分散，反之数据月向中心集中。
这里我们需要知道：1，高斯卷积核是圆对称的（各向同性），中心点权重最大，离中心点越远权重就越小。2，高斯卷积核是可分离卷积核，可以通过水平卷积核和垂直卷积核实现对图像的卷积。3，高斯卷积核的有效尺寸是
(
6
σ
+
1
)
(
6
σ
+
1
)
(6\sigma+1)(6\sigma+1)
(6σ+1)(6σ+1)，尺寸越大，平滑程度就越高。
我们可以通过代码生成三种不同的
σ
\sigma
σ值对比三个以为高斯分布函数的曲线图。
#一维高斯分布函数图
import numpy as np
import matplotlib.pyplot as plt
r = np.linspace(-10,10,200)
fig, ax = plt.subplots(figsize=(5, 5), layout='constrained')
sigma1,sigma2,sigma3 = 1,2,3
y1 = np.exp(r**2 * (-1)/(sigma1 ** 2 * 2)) / (np.sqrt(2 * np.pi) * sigma1)
y2 = np.exp(r**2 * (-1)/(sigma2 ** 2 * 2)) / (np.sqrt(2 * np.pi) * sigma2)
y3 = np.exp(r**2 * (-1)/(sigma3 ** 2 * 2)) / (np.sqrt(2 * np.pi) * sigma3)
plt.grid(color='gray',linestyle='--', linewidth=1)
ax.plot(r,y1,label="sigma1")
ax.plot(r,y2,label="sigma2")
ax.plot(r,y3,label="sigma3")
ax.set_xlabel('r') # Add an x-label to the axes.
ax.set_ylabel('G(r)') # Add a y-label to the axes.
plt.title("1-dimension-Gaussian")
plt.legend()
plt.show()
#一维高斯分布模拟
import cv2
import numpy as np
a = cv2.getGaussianKernel(5,sigma=0,ktype=cv2.CV_64F)
a = np.array(a,dtype=float)
a1 = a.T
a2 = a * a1
sigma = 1.1
b = np.exp(-4/((sigma ** 2) * 2)) / (np.pi * 2 * sigma ** 2)
c = np.exp(-1/((sigma ** 2) * 2)) / (np.pi * 2 * sigma ** 2)
d = np.exp(-0/((sigma ** 2) * 2)) / (np.pi * 2 * sigma ** 2)
sum = b*2 + c*2 + d
b = np.float64(b / sum)
c = np.float64(c / sum)
d = np.float64(d / sum)
a,b,c,d,cv2.getGaussianKernel(5,sigma=1.1,ktype=cv2.CV_64F)
(array([[0.0625],
[0.25 ],
[0.375 ],
[0.25 ],
[0.0625]]),
0.07076637133154648,
0.24446039891162388,
0.3695464595136593,
array([[0.07076637],
[0.2444604 ],
[0.36954646],
[0.2444604 ],
[0.07076637]]))
在利用numpy和opencv分别对一维高斯分布进行模拟时，我发现两个非常奇怪的地方：
第一，当
σ
\sigma
σ值为0时，也就是高斯核元素标准差为0时，高斯分布函数就会出现分母有0的情况，这很显然有问题。不过经过极限运算后，我发现其结果是无穷大，这也恰恰和前面提到的三种
σ
\sigma
σ值的对比图对应起来了，当
σ
\sigma
σ越接近0，分布函数的值就越接近无穷大。
很显然，我们单纯地利用numpy实现公式，是不全面的。于是我利用OpenCV地生成一维高斯分布的函数，发现其中
σ
\sigma
σ值可以为0，而且该方法还做了特别说明：当
σ
\sigma
σ为非正数时，
σ
\sigma
σ由ksize得到：
σ
=
0.3
(
(
k
s
i
z
e
−
1
)
/
2
−
1
)
+
0.8
\sigma = 0.3((ksize -1) / 2 - 1) + 0.8
σ=0.3((ksize−1)/2−1)+0.8。于是我尝试了将numpy方法中的
σ
\sigma
σ值改为1.1，结果发现这与OpenCV方法中的
σ
\sigma
σ为0时的值并不同。
于是我又尝试在Opencv的方法里将
σ
\sigma
σ设为1.1，我发现此时一维高斯分布函数的值和利用numpy实现的一维高斯分布函数值接近了。
这里面就涉及到Opencv在进行高斯滤波运算时的特殊情况，就是当
σ
\sigma
σ为0时，OpenCV会对高斯核进行特殊处理，尽可能使得高斯核各元素值服从正态分布。
第二，在进行卷积运算时，尤其是涉及到高斯滤波器核的卷积运算，我们尽可能选择服从离散分布的高斯核，一般来说，在numpy实现高斯滤波算法时，我们会选择numpy.linespace获取离散的数值。这些数值并不是随机分布的，而是等间隔分布的，数值范围可大可小。但由于高斯核尺寸是有要求的
(
6
σ
+
1
,
6
σ
+
1
)
(6\sigma+1,6\sigma+1)
(6σ+1,6σ+1)，因此等间隔分布的离散数值在设置时要注意范围。
接下来，我们利用服从离散分布的高斯核，进行高斯分布函数的运算以及滤波运算吧。
#离散分布高斯核，一维高斯分布函数
import numpy as np
import matplotlib.pyplot as plt
def continuous_gaussian(x, sigma):
"""
连续高斯分布函数
"""
return 1/(np.sqrt(2*np.pi*sigma**2)) * np.exp(-x**2/(2*sigma**2))
def discretize_gaussian(num_points, sigma):
"""
离散化高斯分布函数
"""
# 选择离散采样点
x_values = np.linspace(-3*sigma, 3*sigma, num_points)
# 计算连续高斯分布在采样点上的值
continuous_values = continuous_gaussian(x_values, sigma)
# 将连续函数值映射到离散序列中
# discrete_values = np.round(continuous_values * num_points).astype(int)
discrete_values = continuous_values
return discrete_values
# 参数
num_points = 21
sigma = 1.0
# 离散化高斯分布
discrete_values = discretize_gaussian(num_points, sigma)
# 可视化
plt.plot(discrete_values, marker='o', linestyle='-', label='Discretized Gaussian')
plt.title('Discretized Gaussian Distribution')
plt.xlabel('Index')
plt.ylabel('Value')
plt.legend()
plt.show()
#一维高斯分布函数图（Opencv方法）
import numpy as np
import matplotlib.pyplot as plt
import cv2
fig, ax = plt.subplots(figsize=(5, 5), layout='constrained')
sigma1,sigma2,sigma3 = 1,2,3
y1 = cv2.getGaussianKernel(61,sigma=sigma1,ktype=cv2.CV_64F).tolist()
y2 = cv2.getGaussianKernel(121,sigma=sigma2,ktype=cv2.CV_64F).tolist()
y3 = cv2.getGaussianKernel(181,sigma=sigma3,ktype=cv2.CV_64F).tolist()
plt.grid(color='gray',linestyle='--', linewidth=1)
ax.plot([(i-30)/10 for i in range(61)],y1,label="sigma1")
ax.plot([(i-60)/10 for i in range(121)],y2,label="sigma2")
ax.plot([(i-90)/10 for i in range(181)],y3,label="sigma3")
ax.set_xlabel('r') # Add an x-label to the axes.
ax.set_ylabel('G(r)') # Add a y-label to the axes.
plt.title("1-dimension-Gaussian")
plt.legend()
plt.show()
此处，我们会发现，条用Opencv库中的生成高斯核方法后，调整了高斯核的大小，严格按照
r
∈
(
−
3
σ
,
3
σ
)
r \in (-3\sigma,3\sigma)
r∈(−3σ,3σ),即高斯核大小最大范围是
(
6
σ
+
1
,
6
σ
+
1
)
(6\sigma + 1,6\sigma + 1)
(6σ+1,6σ+1)。但是对比之前的实验结果，我们又会发现，当
σ
\sigma
σ不变而高斯核越大时，高斯分布函数值的函数图中间部分越集中，也就是图中所示的越来越细的曲线图。
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import cm
from mpl_toolkits.mplot3d import Axes3D #这里只是用Axes3D函数，所以只导入了Axes3D
fig = plt.figure()
ax = Axes3D(fig)
x = np.linspace(-3,3,61) # 生成数据
y = np.linspace(-3,3,61)
x,y=np.meshgrid(x,y) #生成x,y轴数据
sigma1 = 1
z = np.exp((x ** 2 + y ** 2) * (-1)/(sigma1 ** 2 * 2)) / (2 * np.pi * sigma1 ** 2)
surf = ax.plot_wireframe(x,y,z,rstride=1,cstride=1,cmap=cm.viridis)
# ax.plot_surface(x,y,z,rstride=1,cstride=1)
# ax.contourf(x,y,z,zdir='z',offset=-2)
plt.show()
#高斯滤波器
import numpy as np
import matplotlib.pyplot as plt
from PIL import Image
#所选图片为彩色图，维度为3，通道数为3
image = Image.open("erha.jpg")
image = np.array(image)
h_img,w_img,ch = image.shape
# print(image)
plt.title("image_origin")
plt.imshow(image,cmap="brg")
plt.show()
#创建相关运算二维核，而且一般来说，滤波器核的大小都是奇数。
def creat_filter(ksize,random_type,num1,num2):
if random_type == "normal":
mean,std = num1,num2
filter = np.random.normal(mean,std,size=(ksize,ksize))
elif random_type == "int":
low,high = num1,num2
filter = np.random.randint(num1,num2,size=(ksize,ksize))
elif random_type == "std_normal":
filter = np.random.randn(ksize,ksize) * num1 + num2
elif random_type == "float_0_1":
filter = np.random.random(size=(ksize,ksize))
elif random_type == "full_one":
filter = np.ones((ksize,ksize))
# sigma = 0,此时标准差sigma需按照公式取值
# sigma = 0.3 * ((ksize - 1) / 2 - 1) + 0.8
if np.std(filter) == 0:
sigma = 0.3 * ((ksize - 1) / 2 - 1) + 0.8
elif random_type == "normal":
sigma = num2
else:
sigma = np.std(filter)
for x in range(ksize):
for y in range(ksize):
filter[x,y] = filter[x,y] * np.exp(((x - (ksize - 1)/2) ** 2 + (y - (ksize - 1)/2) ** 2) * (-1)/(sigma ** 2 * 2)) / (2 * np.pi * sigma ** 2)
return filter
#边界扩充：填充，相关运算和卷积运算对图像边界点进行特殊处理时，需要适当扩充图像边界。默认填充0，代表黑色。其中padding为填充0的圈数。
#为了使经过滤波器核输出的图像大小和原图一致，本次实验选择高和宽都填充两层。即：padding = h_f - 1。（padding = w_f - 1）
#注意此处滤波器核大小是奇数。
def padding_image(image,h_f,w_f):
h_img,w_img = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new))
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2)] = image[i,j]
return image_padding
def padding_image_3d(image,h_f,w_f):
h_img,w_img,ch = image.shape
h_new,w_new = h_img + h_f - 1,w_img + w_f - 1
image_padding = np.zeros((h_new,w_new,ch))
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_padding[i + int((h_f-1)/2),j + int((w_f-1)/2),k] = image[i,j,k]
return image_padding
# image_padding = padding_image(image,h_f,w_f)
# print(image_padding)
#无论是相关运算，还是卷积运算，经过滤波器核之后的图像输出结果
#空间卷积运算,注意：这里的d一定是偶数
def convolution_operation(image,filter):
h_img,w_img,ch = image.shape
h_f,w_f = filter.shape
filter_convolution = np.flip(filter)#卷积运算所需二维核是对相关运算二维核的翻转。
image_out = np.zeros((h_img,w_img,ch))
image_padding = padding_image_3d(image,h_f,w_f)
for k in range(ch):
for i in range(h_img):
for j in range(w_img):
image_out[i,j,k] = np.sum(np.multiply(filter_convolution[:,:],image_padding[i:i+h_f,j:j+w_f,k]))
return image_out
def draw_convolution(random_type,ksize,num1,num2):
image_convolution = convolution_operation(image,creat_filter(ksize,random_type,num1,num2))
image_out = np.uint8(image_convolution)
plt.title("image_convolution_" + random_type)
plt.imshow(image_out,cmap="brg")
plt.show()
if __name__ == "__main__":
# draw_convolution("std_normal",11,1,1)
# draw_convolution("normal",11,2,2)
# draw_convolution("normal",43,2,2)
# draw_convolution("float_0_1",11,1,1)
draw_convolution("full_one",11,1,1)
draw_convolution("full_one",43,1,1)
上述二维高斯核滤波算法所使用的滤波器核均是随机生成的数据。
#使用OpenCV库进行二维高斯滤波运算
import cv2
import matplotlib.pyplot as plt
img = cv2.imread("erha.jpg")
img = cv2.cvtColor(img,cv2.COLOR_BGR2RGB)
ksize_11 = (11,11)
Gauss_filter11 = cv2.GaussianBlur(img,ksize_11,0) #simga为0时，sigma值由ksize公式计算得出
ksize_43 = (43,43)
Gauss_filter43 = cv2.GaussianBlur(img,ksize_43,0)
plt.title("img_origin")
plt.imshow(img,cmap="brg",vmin=0,vmax=255)
plt.show()
plt.title("img_Gauss_filter11")
plt.imshow(Gauss_filter11,cmap="brg",vmin=0,vmax=255)
plt.show()
plt.title("img_Gauss_filter43")
plt.imshow(Gauss_filter43,cmap="brg",vmin=0,vmax=255)
plt.show()
结语
后续补上不太常见的低通滤波。
926
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
