如何用OpenCV Python提取单轮廓中的多段圆弧半径？

问题描述 #

我希望用不同半径的圆弧段拟合轮廓。当轮廓仅包含单一半径的圆弧段时，使用cf.least_squares_circle()函数可以正常拟合出对应半径；但当轮廓包含多段不同半径的圆弧段时，该函数会得到平均化的半径结果。请问除了cf.least_squares_circle()外，还有其他可行方法吗？

现有代码 #

import cv2
import numpy as np
import circle_fit as cf

# load image as color to draw
img = cv2.imread(r'image.png')

##img Height:
print((img.shape[0]))
##Bildbreite:
print((img.shape[2]))

PY = []
PX = []

for j in range (img.shape[2]):
    for i in range (img.shape[0]):
        if img[i,j,0] != 0:
            PY.append(i)
            PX.append(j)

PXY = np.zeros((len(PX), 2))
for n in range(len(PX)):
    PXY[(n, 1)] = PX[n]
    PXY[(n, 0)] = PY[n]

xc, yc, r, _ = cf.least_squares_circle(PXY)
cv2.circle(img, (int(yc), int(xc)), int(r), (0, 0, 255), thickness=1, lineType=8, shift=0)
cv2.imshow("output", img)
cv2.waitKey(0)

相关示意图 #

可行解决方案 #

1. RANSAC随机抽样一致性拟合 #

RANSAC能从混合数据中筛选出符合单一圆弧模型的点集，反复迭代即可拟合出多段不同半径的圆弧，避免平均化问题。

示例代码：

import cv2
import numpy as np

def ransac_fit_circle(points, iterations=100, threshold=2.0):
    best_inliers = []
    best_params = None
    n_points = len(points)
    if n_points < 3:
        return None, []
    
    for _ in range(iterations):
        # 随机选取3个点用于拟合圆
        sample_indices = np.random.choice(n_points, 3, replace=False)
        sample = points[sample_indices]
        (x1,y1), (x2,y2), (x3,y3) = sample
        
        # 解方程组求圆心
        A = np.array([
            [2*(x2-x1), 2*(y2-y1)],
            [2*(x3-x1), 2*(y3-y1)]
        ])
        B = np.array([
            x2**2 + y2**2 - x1**2 - y1**2,
            x3**2 + y3**2 - x1**2 - y1**2
        ])
        try:
            xc, yc = np.linalg.solve(A, B)
            r = np.sqrt((x1 - xc)**2 + (y1 - yc)**2)
        except np.linalg.LinAlgError:
            continue
        
        # 计算所有点到圆的距离，筛选内点
        distances = np.abs(np.sqrt((points[:,0]-xc)**2 + (points[:,1]-yc)**2) - r)
        inliers = points[distances < threshold]
        
        if len(inliers) > len(best_inliers):
            best_inliers = inliers
            best_params = (xc, yc, r)
    return best_params, best_inliers

# 优化轮廓点提取逻辑
img = cv2.imread('image.png')
gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
_, binary = cv2.threshold(gray, 1, 255, cv2.THRESH_BINARY)
contours, _ = cv2.findContours(binary, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
points = np.vstack(contours).squeeze()

remaining_points = points.copy()
all_circles = []

# 迭代拟合多段圆弧，剩余点过少时停止
while len(remaining_points) > 10:
    circle_params, inliers = ransac_fit_circle(remaining_points)
    if not circle_params:
        break
    all_circles.append(circle_params)
    
    # 移除已拟合的内点
    mask = np.ones(len(remaining_points), dtype=bool)
    for idx, p in enumerate(remaining_points):
        if any((p == inlier).all() for inlier in inliers):
            mask[idx] = False
    remaining_points = remaining_points[mask]

# 绘制所有拟合结果
for xc, yc, r in all_circles:
    cv2.circle(img, (int(yc), int(xc)), int(r), (0,0,255), 1)
cv2.imshow('output', img)
cv2.waitKey(0)

2. 曲率分段后单独拟合 #

通过计算轮廓点的曲率变化找到分段边界，将轮廓拆分为多个圆弧段后，再分别用最小二乘法拟合。

示例代码：

import cv2
import numpy as np
import circle_fit as cf

# 计算轮廓点的曲率
def calculate_curvature(contour):
    curvatures = []
    n = len(contour)
    for i in range(n):
        p_prev = contour[(i-1)%n][0]
        p_curr = contour[i][0]
        p_next = contour[(i+1)%n][0]
        
        vec1 = p_curr - p_prev
        vec2 = p_next - p_curr
        cross = np.cross(vec1, vec2)
        norm1 = np.linalg.norm(vec1)
        norm2 = np.linalg.norm(vec2)
        
        if norm1 == 0 or norm2 == 0:
            curvatures.append(0)
        else:
            curvature = 2 * cross / (norm1 * norm2 * (norm1 + norm2))
            curvatures.append(np.abs(curvature))
    return curvatures

# 加载图像并提取轮廓
img = cv2.imread('image.png')
gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
_, binary = cv2.threshold(gray, 1, 255, cv2.THRESH_BINARY)
contours, _ = cv2.findContours(binary, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
contour = contours[0]

# 计算曲率并找到分段点
curvatures = calculate_curvature(contour)
threshold = np.mean(curvatures) * 2
split_indices = [i for i, c in enumerate(curvatures) if c > threshold]
split_indices.append(0)  # 闭合轮廓，补全首尾
split_indices = sorted(list(set(split_indices)))

# 对每个分段单独拟合圆弧
for i in range(len(split_indices)-1):
    start = split_indices[i]
    end = split_indices[i+1]
    segment = contour[start:end]
    if len(segment) < 3:
        continue
    segment_points = np.array([p[0] for p in segment])
    xc, yc, r, _ = cf.least_squares_circle(segment_points)
    cv2.circle(img, (int(yc), int(xc)), int(r), (0,255,0), 1)

cv2.imshow('output', img)
cv2.waitKey(0)

3. 利用专业几何库分割拟合 #

使用shapely库将轮廓转换为几何线段，通过近似分割得到圆弧段后再拟合：

import cv2
import numpy as np
import circle_fit as cf
from shapely.geometry import LineString

# 加载图像提取轮廓点
img = cv2.imread('image.png')
gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
_, binary = cv2.threshold(gray, 1, 255, cv2.THRESH_BINARY)
contours, _ = cv2.findContours(binary, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
points = np.vstack(contours).squeeze()

# 转换为LineString并分割近似
line = LineString(points)
# tolerance参数控制近似精度，值越大分段越少
simplified = line.simplify(tolerance=1.0, preserve_topology=False)

# 对每个子段拟合圆弧
for geom in simplified.geoms:
    seg_points = np.array(list(geom.coords))
    if len(seg_points) >=3:
        xc, yc, r, _ = cf.least_squares_circle(seg_points)
        cv2.circle(img, (int(yc), int(xc)), int(r), (255,0,0),1)

cv2.imshow('output', img)
cv2.waitKey(0)

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