2.6.6 测量对象属性: ndimage.measurements
合成数据:
In [59]:
n = 10
l = 256
im = np.zeros((l, l))
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()
- 联通分支分析
标记联通分支:ndimage.label
:
In [60]:
label_im, nb_labels = ndimage.label(mask)
nb_labels # 多少区域?
plt.imshow(label_im)
Out[60]:
<matplotlib.image.AxesImage at 0x11543ab10>
计算每个区域的大小、mean_value等等:
In [61]:
sizes = ndimage.sum(mask, label_im, range(nb_labels + 1))
mean_vals = ndimage.sum(im, label_im, range(1, nb_labels + 1))
清理小的联通分支:
In [62]:
mask_size = sizes < 1000
remove_pixel = mask_size[label_im]
remove_pixel.shape
label_im[remove_pixel] = 0
plt.imshow(label_im)
Out[62]:
<matplotlib.image.AxesImage at 0x114644a90>
现在用np.searchsorted
重新分配标签:
In [63]:
labels = np.unique(label_im)
label_im = np.searchsorted(labels, label_im)
找到包含感兴趣对象的区域:
In [65]:
slice_x, slice_y = ndimage.find_objects(label_im==4)[0]
roi = im[slice_x, slice_y]
plt.imshow(roi)
Out[65]:
<matplotlib.image.AxesImage at 0x1147fa8d0>
其他空间测量: ndimage.center_of_mass
, ndimage.maximum_position
, 等等,可以在分割应用这个有限的部分之外使用。
实例: 块平均数:
In [66]:
from scipy import misc
f = misc.face(gray=True)
sx, sy = f.shape
X, Y = np.ogrid[0:sx, 0:sy]
regions = (sy//6) * (X//4) + (Y//6) # 注意我们使用广播
block_mean = ndimage.mean(f, labels=regions, index=np.arange(1, regions.max() +1))
block_mean.shape = (sx // 4, sy // 6)
当区域是标准的块时,使用步长刻度更加高效 (实例: 使用刻度的虚假维度)。
非标准空间块: radial平均:
In [67]:
sx, sy = f.shape
X, Y = np.ogrid[0:sx, 0:sy]
r = np.hypot(X - sx/2, Y - sy/2)
rbin = (20* r/r.max()).astype(np.int)
radial_mean = ndimage.mean(f, labels=rbin, index=np.arange(1, rbin.max() +1))
- 其他测量
相关函数、Fourier/wavelet频谱、等㩐。
使用数学形态学的一个实例: granulometry
In [69]:
def disk_structure(n):
struct = np.zeros((2 * n + 1, 2 * n + 1))
x, y = np.indices((2 * n + 1, 2 * n + 1))
mask = (x - n)**2 + (y - n)**2 <= n**2
struct[mask] = 1
return struct.astype(np.bool)
def granulometry(data, sizes=None):
s = max(data.shape)
if sizes == None:
sizes = range(1, s/2, 2)
granulo = [ndimage.binary_opening(data, \
structure=disk_structure(n)).sum() for n in sizes]
return granulo
np.random.seed(1)
n = 10
l = 256
im = np.zeros((l, l))
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()
granulo = granulometry(mask, sizes=np.arange(2, 19, 4))
/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.
也看一下: 关于图像处理的更多内容:
- 关于Scikit-image的章节
- 其他更有力更完整的模块: OpenCV (Python绑定), CellProfiler, ITK 带有Python绑定