2.6.6 测量对象属性: ndimage.measurements

合成数据:

In [59]:

  1. n = 10
  2. l = 256
  3. im = np.zeros((l, l))
  4. points = l*np.random.random((2, n**2))
  5. im[(points[0]).astype(np.int), (points[1]).astype(np.int)] = 1
  6. im = ndimage.gaussian_filter(im, sigma=l/(4.*n))
  7. mask = im > im.mean()
  • 联通分支分析

标记联通分支:ndimage.label:

In [60]:

  1. label_im, nb_labels = ndimage.label(mask)
  2. nb_labels # 多少区域?
  3. plt.imshow(label_im)

Out[60]:

  1. <matplotlib.image.AxesImage at 0x11543ab10>

2.6.6 测量对象属性: ndimage.measurements - 图1

2.6.6 测量对象属性: ndimage.measurements - 图2

[Python source code]

计算每个区域的大小、mean_value等等:

In [61]:

  1. sizes = ndimage.sum(mask, label_im, range(nb_labels + 1))
  2. mean_vals = ndimage.sum(im, label_im, range(1, nb_labels + 1))

清理小的联通分支:

In [62]:

  1. mask_size = sizes < 1000
  2. remove_pixel = mask_size[label_im]
  3. remove_pixel.shape
  4. label_im[remove_pixel] = 0
  5. plt.imshow(label_im)

Out[62]:

  1. <matplotlib.image.AxesImage at 0x114644a90>

2.6.6 测量对象属性: ndimage.measurements - 图3

现在用np.searchsorted重新分配标签:

In [63]:

  1. labels = np.unique(label_im)
  2. label_im = np.searchsorted(labels, label_im)

2.6.6 测量对象属性: ndimage.measurements - 图4

[Python source code]

找到包含感兴趣对象的区域:

In [65]:

  1. slice_x, slice_y = ndimage.find_objects(label_im==4)[0]
  2. roi = im[slice_x, slice_y]
  3. plt.imshow(roi)

Out[65]:

  1. <matplotlib.image.AxesImage at 0x1147fa8d0>

2.6.6 测量对象属性: ndimage.measurements - 图5

[Python source code]

其他空间测量: ndimage.center_of_mass, ndimage.maximum_position, 等等,可以在分割应用这个有限的部分之外使用。

实例: 块平均数:

In [66]:

  1. from scipy import misc
  2. f = misc.face(gray=True)
  3. sx, sy = f.shape
  4. X, Y = np.ogrid[0:sx, 0:sy]
  5. regions = (sy//6) * (X//4) + (Y//6) # 注意我们使用广播
  6. block_mean = ndimage.mean(f, labels=regions, index=np.arange(1, regions.max() +1))
  7. block_mean.shape = (sx // 4, sy // 6)

2.6.6 测量对象属性: ndimage.measurements - 图6

[Python source code]

当区域是标准的块时,使用步长刻度更加高效 (实例: 使用刻度的虚假维度)。

非标准空间块: radial平均:

In [67]:

  1. sx, sy = f.shape
  2. X, Y = np.ogrid[0:sx, 0:sy]
  3. r = np.hypot(X - sx/2, Y - sy/2)
  4. rbin = (20* r/r.max()).astype(np.int)
  5. radial_mean = ndimage.mean(f, labels=rbin, index=np.arange(1, rbin.max() +1))

2.6.6 测量对象属性: ndimage.measurements - 图7

[Python source code]

  • 其他测量

相关函数、Fourier/wavelet频谱、等㩐。

使用数学形态学的一个实例: granulometry

In [69]:

  1. def disk_structure(n):
  2. struct = np.zeros((2 * n + 1, 2 * n + 1))
  3. x, y = np.indices((2 * n + 1, 2 * n + 1))
  4. mask = (x - n)**2 + (y - n)**2 <= n**2
  5. struct[mask] = 1
  6. return struct.astype(np.bool)
  7. def granulometry(data, sizes=None):
  8. s = max(data.shape)
  9. if sizes == None:
  10. sizes = range(1, s/2, 2)
  11. granulo = [ndimage.binary_opening(data, \
  12. structure=disk_structure(n)).sum() for n in sizes]
  13. return granulo
  14. np.random.seed(1)
  15. n = 10
  16. l = 256
  17. im = np.zeros((l, l))
  18. points = l*np.random.random((2, n**2))
  19. im[(points[0]).astype(np.int), (points[1]).astype(np.int)] = 1
  20. im = ndimage.gaussian_filter(im, sigma=l/(4.*n))
  21. mask = im > im.mean()
  22. granulo = granulometry(mask, sizes=np.arange(2, 19, 4))
  1. /Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/site-packages/ipykernel/__main__.py:11: FutureWarning: comparison to `None` will result in an elementwise object comparison in the future.

2.6.6 测量对象属性: ndimage.measurements - 图8

[Python source code]

也看一下: 关于图像处理的更多内容: