基于优化后中轴的环状图像厚度计算及边界距离求解

一、优化后中轴上等距n个点的厚度计算步骤 #

假设你已通过NetworkX得到无分支的单环状中轴（二值图像optimized_medial_axis），按以下步骤计算等距点的厚度：

1. 提取中轴像素坐标 #

先把中轴的所有非零像素坐标提取出来：

import numpy as np
from skimage import measure

# 提取中轴的像素坐标（y, x）
y_coords, x_coords = np.where(optimized_medial_axis)
points = np.column_stack((x_coords, y_coords))  # 转为(x, y)格式，方便后续拟合

2. 拟合闭合环状曲线 #

由于中轴是闭合环状，用样条拟合生成连续曲线，确保首尾衔接：

from scipy.interpolate import splprep, splev

# 拟合闭合B样条曲线，per=True表示周期闭合
tck, u = splprep(points.T, s=0, per=True)

3. 生成n个等距点 #

在0-1的参数区间取n个等间隔点，得到对应坐标：

n = 10  # 你需要的点数
u_new = np.linspace(0, 1, n, endpoint=False)  # 不包含终点避免重复
x_new, y_new = splev(u_new, tck)
# 转为整数像素坐标（如果需要）
sample_points = np.round(np.column_stack((x_new, y_new))).astype(int)

4. 计算每个点到内外边界的距离 #

方法：预计算内外边界的距离变换图

• 先提取环状的内外轮廓：

# 假设你的原始二值图像为ring_img，前景（环状区域）值为1，背景为0
contours = measure.find_contours(ring_img, level=0.5)
# 通常轮廓中长度大的是外边界，小的是内边界
contours_sorted = sorted(contours, key=lambda c: len(c), reverse=True)
outer_contour = contours_sorted[0]
inner_contour = contours_sorted[1]

# 生成内外边界的二值图
outer_boundary = np.zeros_like(ring_img)
inner_boundary = np.zeros_like(ring_img)
# 给边界像素赋值为1（用round取整坐标）
for (x, y) in np.round(outer_contour).astype(int):
    if 0 <= y < ring_img.shape[0] and 0 <= x < ring_img.shape[1]:
        outer_boundary[y, x] = 1
for (x, y) in np.round(inner_contour).astype(int):
    if 0 <= y < ring_img.shape[0] and 0 <= x < ring_img.shape[1]:
        inner_boundary[y, x] = 1

• 计算距离变换：

from skimage import distance

# 计算环状区域内的点到外边界的距离
dist_outer = distance.distance_transform_edt(1 - outer_boundary)
# 计算环状区域内的点到内边界的距离
dist_inner = distance.distance_transform_edt(1 - inner_boundary)

• 对每个等距点取距离值，相加得到厚度：

thicknesses = []
for (x, y) in sample_points:
    # 确保坐标在图像范围内
    if 0 <= y < ring_img.shape[0] and 0 <= x < ring_img.shape[1]:
        thickness = dist_inner[y, x] + dist_outer[y, x]
        thicknesses.append(thickness)

二、更优的环状图像厚度计算方法 #

1. 基于轮廓的垂直对应点法（无需中轴） #

直接利用内外轮廓，对每个外轮廓等距点，找内轮廓上的最近点，计算距离：

# 对外轮廓生成n个等距点
tck_outer, u_outer = splprep(outer_contour.T, s=0, per=True)
u_outer_new = np.linspace(0, 1, n, endpoint=False)
x_outer, y_outer = splev(u_outer_new, tck_outer)
outer_sample = np.column_stack((x_outer, y_outer))

# 计算每个外轮廓点到内轮廓的最近距离作为厚度
from sklearn.neighbors import NearestNeighbors

nn = NearestNeighbors(n_neighbors=1)
nn.fit(inner_contour)
distances, _ = nn.kneighbors(outer_sample)
thicknesses_contour = distances.flatten()

这种方法跳过中轴提取，避免分支问题，适合规则或不规则环状。

2. 改进的中轴提取方法（避免分支） #

如果坚持用中轴法，可先提取骨架后去除分支：

from skimage.morphology import skeletonize, remove_small_objects

# 提取骨架
skeleton = skeletonize(ring_img)
# 去除小的连通分支（保留最大的环状骨架）
skeleton_cleaned = remove_small_objects(skeleton.astype(bool), min_size=100)  # min_size根据你的图像调整

再结合之前的等距点计算流程，无需额外NetworkX处理，适合简单分支情况。

3. 极坐标转换法（近似圆形环状） #

如果环状接近圆形，将图像转为极坐标，每个角度下的径向长度差即为厚度：

import cv2

# 找到环状的中心（比如外轮廓的质心）
M = cv2.moments(outer_contour.astype(np.float32))
cx = int(M['m10']/M['m00'])
cy = int(M['m01']/M['m00'])

# 转换为极坐标
height, width = ring_img.shape
max_radius = int(np.sqrt((cx)**2 + (cy)**2))
polar_img = cv2.warpPolar(ring_img, (max_radius, 360), (cx, cy), max_radius, cv2.WARP_POLAR_LINEAR)

# 对每个角度（0-359）计算厚度
thicknesses_polar = []
for angle in range(360):
    row = polar_img[:, angle]
    # 找到前景的起始和结束半径
    foreground_indices = np.where(row == 1)[0]
    if len(foreground_indices) > 0:
        thickness = foreground_indices[-1] - foreground_indices[0]
        thicknesses_polar.append(thickness)
# 取n个等角度的厚度值
step = 360 // n
selected_thicknesses = thicknesses_polar[::step]

内容的提问来源于stack exchange，提问作者Savindi