使用像素图查找对象的测地中心

在各种图像分析情况下，将图像的像素或图像区域视为网络或图是有用的，其中每个像素与其邻居相连（有或没有对角线）。一种这样的情况是找到对象的测地中心，即最接近所有其他点的点，如果您只允许在对象的像素上移动，而不是直线移动。这个点是网络中具有最大接近中心性的点 [1]。

在此示例中，我们创建一个骨架的这种像素图，并找到该骨架的中心像素。这证明了它与质心（也称为质心）相比的实用性，质心实际上可能落在对象之外。

参考文献

[1]

Linton C. Freeman: Centrality in networks: I. Conceptual clarification. Social Networks 1:215-239, 1979. DOI:10.1016/0378-8733(78)90021-7

import matplotlib.pyplot as plt
import numpy as np
from scipy import ndimage as ndi
from skimage import color, data, filters, graph, measure, morphology

我们首先加载数据：人视网膜的图像。

retina_source = data.retina()

_, ax = plt.subplots()
ax.imshow(retina_source)
ax.set_axis_off()
_ = ax.set_title('Human retina')

我们将图像转换为灰度，然后使用 Sato 血管性 滤波器 来更好地识别图像中的主要血管。

retina = color.rgb2gray(retina_source)
t0, t1 = filters.threshold_multiotsu(retina, classes=3)
mask = retina > t0
vessels = filters.sato(retina, sigmas=range(1, 10)) * mask

_, axes = plt.subplots(nrows=1, ncols=2)
axes[0].imshow(retina, cmap='gray')
axes[0].set_axis_off()
axes[0].set_title('grayscale')
axes[1].imshow(vessels, cmap='magma')
axes[1].set_axis_off()
_ = axes[1].set_title('Sato vesselness')

基于观察到的血管性值，我们使用 滞后 阈值处理 来定义主要血管。

thresholded = filters.apply_hysteresis_threshold(vessels, 0.01, 0.03)
labeled = ndi.label(thresholded)[0]

_, ax = plt.subplots()
ax.imshow(color.label2rgb(labeled, retina))
ax.set_axis_off()
_ = ax.set_title('Thresholded vesselness')

最后，我们可以 骨骼化 此标签图像，并将其作为找到该骨架中 中心 像素 的基础。将其与质心的位置进行比较！

largest_nonzero_label = np.argmax(np.bincount(labeled[labeled > 0]))
binary = labeled == largest_nonzero_label
skeleton = morphology.skeletonize(binary)
g, nodes = graph.pixel_graph(skeleton, connectivity=2)
px, distances = graph.central_pixel(
    g, nodes=nodes, shape=skeleton.shape, partition_size=100
)

centroid = measure.centroid(labeled > 0)

_, ax = plt.subplots()
ax.imshow(color.label2rgb(skeleton, retina))
ax.scatter(px[1], px[0], label='graph center')
ax.scatter(centroid[1], centroid[0], label='centroid')
ax.legend()
ax.set_axis_off()
ax.set_title('Vessel graph center vs centroid')

plt.show()