2.6.4 图像过滤
局部过滤器: 通过临近像素值的函数来替换像素的值。
邻居: 方块 (选择大小)、 硬盘、或者更复杂的结构化元素。
2.6.4.1 模糊 / 光滑
高斯过滤器 来自 scipy.ndimage
:
In [29]:
from scipy import misc
face = misc.face(gray=True)
blurred_face = ndimage.gaussian_filter(face, sigma=3)
very_blurred = ndimage.gaussian_filter(face, sigma=5)
均匀过滤器
In [30]:
local_mean = ndimage.uniform_filter(face, size=11)
2.6.4.2 锐化
锐化模糊的图像:
In [31]:
from scipy import misc
face = misc.face(gray=True).astype(float)
blurred_f = ndimage.gaussian_filter(face, 3)
通过添加Laplacian的近似值来增加边缘的权重:
In [32]:
filter_blurred_f = ndimage.gaussian_filter(blurred_f, 1)
alpha = 30
sharpened = blurred_f + alpha * (blurred_f - filter_blurred_f)
2.6.4.3 降噪
有噪音的脸:
In [33]:
from scipy import misc
f = misc.face(gray=True)
f = f[230:290, 220:320]
noisy = f + 0.4 * f.std() * np.random.random(f.shape)
高斯过滤器光滑了噪音… 以及边缘:
In [34]:
gauss_denoised = ndimage.gaussian_filter(noisy, 2)
绝大多数局部线性各向同性过滤器模糊图像 (ndimage.uniform_filter
)
中位数过滤器更好的保留的边缘:
In [35]:
med_denoised = ndimage.median_filter(noisy, 3)
中位数过滤器: 对直边缘的结果更好 (低曲度):
In [37]:
im = np.zeros((20, 20))
im[5:-5, 5:-5] = 1
im = ndimage.distance_transform_bf(im)
im_noise = im + 0.2 * np.random.randn(*im.shape)
im_med = ndimage.median_filter(im_noise, 3)
其他排名过滤器:ndimage.maximum_filter
、ndimage.percentile_filter
其他的局部非线性过滤器: Wiener (scipy.signal.wiener
) 等等.
非局部过滤器
练习: 降噪
- 创建一个带有一些对象 (圆形、椭圆、正方形或随机形状) 的二元图像 (0 和 1) .
- 添加一些噪音 (例如, 20% 的噪音)
- 用两种不同的降噪方法来降噪这个图像: 高斯过滤器和中位数过滤器。
- 比较两种不同降噪方法的直方图。哪一个与原始图像(无噪音)的直方图最接近?
也看一下: 在skimage.denoising
中有更多的的降噪过滤器可用,见教程Scikit-image: 图像处理.
2.6.4.4 数学形态学
看一下wikipedia上的数学形态学定义。
探索一个简单形状 (结构化的元素) 的图像,然后根据这个形状是如何局部适应或不适应这个图像来修改这个图像。
结构化元素:
In [38]:
el = ndimage.generate_binary_structure(2, 1)
el
Out[38]:
array([[False, True, False],
[ True, True, True],
[False, True, False]], dtype=bool)
In [39]:
el.astype(np.int)
Out[39]:
array([[0, 1, 0],
[1, 1, 1],
[0, 1, 0]])
风化(Erosion) = 最小值过滤器。用结构化元素所覆盖的最小值来替换像素的值。:
In [41]:
a = np.zeros((7,7), dtype=np.int)
a[1:6, 2:5] = 1
a
Out[41]:
array([[0, 0, 0, 0, 0, 0, 0],
[0, 0, 1, 1, 1, 0, 0],
[0, 0, 1, 1, 1, 0, 0],
[0, 0, 1, 1, 1, 0, 0],
[0, 0, 1, 1, 1, 0, 0],
[0, 0, 1, 1, 1, 0, 0],
[0, 0, 0, 0, 0, 0, 0]])
In [42]:
ndimage.binary_erosion(a).astype(a.dtype)
Out[42]:
array([[0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 1, 0, 0, 0],
[0, 0, 0, 1, 0, 0, 0],
[0, 0, 0, 1, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0]])
In [43]:
#风化删除了比结构小的对象
ndimage.binary_erosion(a, structure=np.ones((5,5))).astype(a.dtype)
Out[43]:
array([[0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0]])
扩大(Dilation):最大值过滤器:
In [44]:
a = np.zeros((5, 5))
a[2, 2] = 1
a
Out[44]:
array([[ 0., 0., 0., 0., 0.],
[ 0., 0., 0., 0., 0.],
[ 0., 0., 1., 0., 0.],
[ 0., 0., 0., 0., 0.],
[ 0., 0., 0., 0., 0.]])
In [45]:
ndimage.binary_dilation(a).astype(a.dtype)
Out[45]:
array([[ 0., 0., 0., 0., 0.],
[ 0., 0., 1., 0., 0.],
[ 0., 1., 1., 1., 0.],
[ 0., 0., 1., 0., 0.],
[ 0., 0., 0., 0., 0.]])
对于灰度值图像同样适用:
In [46]:
np.random.seed(2)
im = np.zeros((64, 64))
x, y = (63*np.random.random((2, 8))).astype(np.int)
im[x, y] = np.arange(8)
bigger_points = ndimage.grey_dilation(im, size=(5, 5), structure=np.ones((5, 5)))
square = np.zeros((16, 16))
square[4:-4, 4:-4] = 1
dist = ndimage.distance_transform_bf(square)
dilate_dist = ndimage.grey_dilation(dist, size=(3, 3), structure=np.ones((3, 3)))
Opening: erosion + dilation:
In [48]:
a = np.zeros((5,5), dtype=np.int)
a[1:4, 1:4] = 1; a[4, 4] = 1
a
Out[48]:
array([[0, 0, 0, 0, 0],
[0, 1, 1, 1, 0],
[0, 1, 1, 1, 0],
[0, 1, 1, 1, 0],
[0, 0, 0, 0, 1]])
In [49]:
# Opening 删除小对象
ndimage.binary_opening(a, structure=np.ones((3,3))).astype(np.int)
Out[49]:
array([[0, 0, 0, 0, 0],
[0, 1, 1, 1, 0],
[0, 1, 1, 1, 0],
[0, 1, 1, 1, 0],
[0, 0, 0, 0, 0]])
In [50]:
# Opening 也可以光滑转角
ndimage.binary_opening(a).astype(np.int)
Out[50]:
array([[0, 0, 0, 0, 0],
[0, 0, 1, 0, 0],
[0, 1, 1, 1, 0],
[0, 0, 1, 0, 0],
[0, 0, 0, 0, 0]])
应用: 去除噪音:
In [51]:
square = np.zeros((32, 32))
square[10:-10, 10:-10] = 1
np.random.seed(2)
x, y = (32*np.random.random((2, 20))).astype(np.int)
square[x, y] = 1
open_square = ndimage.binary_opening(square)
eroded_square = ndimage.binary_erosion(square)
reconstruction = ndimage.binary_propagation(eroded_square, mask=square)
Closing: dilation + erosion 许多其他的数学形态学操作: hit 和 miss 转换、tophat等等。
2.6.5.2 分割
- 直方图分割 (没有空间信息)
In [52]:
n = 10
l = 256
im = np.zeros((l, l))
np.random.seed(1)
points = l*np.random.random((2, n**2))
im[(points[0]).astype(np.int), (points[1]).astype(np.int)] = 1
im = ndimage.gaussian_filter(im, sigma=l/(4.*n))
mask = (im > im.mean()).astype(np.float)
mask += 0.1 * im
img = mask + 0.2*np.random.randn(*mask.shape)
hist, bin_edges = np.histogram(img, bins=60)
bin_centers = 0.5*(bin_edges[:-1] + bin_edges[1:])
binary_img = img > 0.5
用数学形态学来清理结果:
In [53]:
# 删除小白色区域
open_img = ndimage.binary_opening(binary_img)
# 删除小黑洞
close_img = ndimage.binary_closing(open_img)
练习
检查重建操作 (erosion + propagation) 生成一个比opening/closing更好的结果:
In [55]:
eroded_img = ndimage.binary_erosion(binary_img)
reconstruct_img = ndimage.binary_propagation(eroded_img, mask=binary_img)
tmp = np.logical_not(reconstruct_img)
eroded_tmp = ndimage.binary_erosion(tmp)
reconstruct_final = np.logical_not(ndimage.binary_propagation(eroded_tmp, mask=tmp))
np.abs(mask - close_img).mean()
Out[55]:
0.0072783654884356133
In [56]:
np.abs(mask - reconstruct_final).mean()
Out[56]:
0.00059502552803868841
练习
检查第一步降噪 (例如,中位数过滤器) 如何影响直方图,检查产生的基于直方图的分割是否更加准确。
也可以看一下:更高级分割算法可以在scikit-image
中找到: 见Scikit-image: 图像处理。
也可以看一下:其他科学计算包为图像处理提供的算法。在这个例子中,我们使用scikit-learn
中的谱聚类函数以便分割黏在一起的对象。
In [57]:
from sklearn.feature_extraction import image
from sklearn.cluster import spectral_clustering
l = 100
x, y = np.indices((l, l))
center1 = (28, 24)
center2 = (40, 50)
center3 = (67, 58)
center4 = (24, 70)
radius1, radius2, radius3, radius4 = 16, 14, 15, 14
circle1 = (x - center1[0])**2 + (y - center1[1])**2 < radius1**2
circle2 = (x - center2[0])**2 + (y - center2[1])**2 < radius2**2
circle3 = (x - center3[0])**2 + (y - center3[1])**2 < radius3**2
circle4 = (x - center4[0])**2 + (y - center4[1])**2 < radius4**2
# 4个圆形
img = circle1 + circle2 + circle3 + circle4
mask = img.astype(bool)
img = img.astype(float)
img += 1 + 0.2*np.random.randn(*img.shape)
# 将图像转化为边缘带有坡度值的图。
graph = image.img_to_graph(img, mask=mask)
# 用一个梯度递减函数: 我们用它轻度依赖于接近voronoi的梯度细分
graph.data = np.exp(-graph.data/graph.data.std())
labels = spectral_clustering(graph, n_clusters=4, eigen_solver='arpack')
label_im = -np.ones(mask.shape)
label_im[mask] = labels