CSDN搬家失敗，手動匯出markdown後再匯入部落格園
參考 Matlab 計算輪廓內切圓
初衷是為了求裂縫的最大寬度
![[output/attachments/5ecf17abcb54aaa4fb35b00c3f243f32_MD5.png]]
直接上程式碼
import random
import cv2
import math
import numpy as np
from numpy.ma import cos, sin
import matplotlib.pyplot as plt
def max_circle(f):
img = cv2.imread(f, cv2.IMREAD_COLOR)
img_gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
# _, img_gray = cv2.threshold(img_gray, 0, 255, cv2.THRESH_BINARY + cv2.THRESH_OTSU)
contous, hierarchy = cv2.findContours(img_gray, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
"""
第二個參數列示輪廓的檢索模式，有四種（本文介紹的都是新的cv2介面）：
cv2.RETR_EXTERNAL表示只檢測外輪廓
cv2.RETR_LIST檢測的輪廓不建立等級關係
cv2.RETR_CCOMP建立兩個等級的輪廓，上面的一層為外邊界，裡面的一層為內孔的邊界資訊。如果內孔內還有一個連通物體，這個物體的邊界也在頂層。
cv2.RETR_TREE建立一個等級樹結構的輪廓。
第三個引數method為輪廓的近似辦法
cv2.CHAIN_APPROX_NONE儲存所有的輪廓點，相鄰的兩個點的畫素位置差不超過1，即max（abs（x1-x2），abs（y2-y1））==1
cv2.CHAIN_APPROX_SIMPLE壓縮水平方向，垂直方向，對角線方向的元素，只保留該方向的終點座標，例如一個矩形輪廓只需4個點來儲存輪廓資訊
cv2.CHAIN_APPROX_TC89_L1，CV_CHAIN_APPROX_TC89_KCOS使用teh-Chinl chain 近似演算法
"""
for c in contous:
left_x = min(c[:, 0, 0])
right_x = max(c[:, 0, 0])
down_y = max(c[:, 0, 1])
up_y = min(c[:, 0, 1])
upper_r = min(right_x - left_x, down_y - up_y) / 2
# 定義相切二分精度
precision = math.sqrt((right_x - left_x) ** 2 + (down_y - up_y) ** 2) / (2 ** 13)
# 構造包含輪廓的矩形的所有畫素點
Nx = 2 ** 8
Ny = 2 ** 8
pixel_X = np.linspace(left_x, right_x, Nx)
pixel_Y = np.linspace(up_y, down_y, Ny)
# [pixel_X, pixel_Y] = ndgrid(pixel_X, pixel_Y);
# pixel_X = reshape(pixel_X, numel(pixel_X), 1);
# pixel_Y = reshape(pixel_Y, numel(pixel_Y), 1);
xx, yy = np.meshgrid(pixel_X, pixel_Y)
# % 篩選出輪廓內所有畫素點
in_list = []
for c in contous:
for i in range(pixel_X.shape[0]):
for j in range(pixel_X.shape[0]):
if cv2.pointPolygonTest(c, (xx[i][j], yy[i][j]), False) > 0:
in_list.append((xx[i][j], yy[i][j]))
in_point = np.array(in_list)
pixel_X = in_point[:, 0]
pixel_Y = in_point[:, 1]
# 隨機搜尋百分之一畫素提高內切圓半徑下限
N = len(in_point)
rand_index = random.sample(range(N), N // 100)
rand_index.sort()
radius = 0
big_r = upper_r
center = None
for id in rand_index:
tr = iterated_optimal_incircle_radius_get(c, in_point[id][0], in_point[id][1], radius, big_r, precision)
if tr > radius:
radius = tr
center = (in_point[id][0], in_point[id][1]) # 只有半徑變大才允許位置變更，否則保持之前位置不變
# 迴圈搜尋剩餘畫素對應內切圓半徑
loops_index = [i for i in range(N) if i not in rand_index]
for id in loops_index:
tr = iterated_optimal_incircle_radius_get(c, in_point[id][0], in_point[id][1], radius, big_r, precision)
if tr > radius:
radius = tr
center = (in_point[id][0], in_point[id][1]) # 只有半徑變大才允許位置變更，否則保持之前位置不變
# 效果測試
plot_x = np.linspace(0, 2 * math.pi, 100)
circle_X = center[0] + radius * cos(plot_x)
circle_Y = center[1] + radius * sin(plot_x)
print(radius * 2)
plt.figure()
plt.imshow(img_gray)
plt.plot(circle_X, circle_Y)
plt.show()
def iterated_optimal_incircle_radius_get(contous, pixelx, pixely, small_r, big_r, precision):
radius = small_r
L = np.linspace(0, 2 * math.pi, 360) # 確定圓散點剖分數360, 720
circle_X = pixelx + radius * cos(L)
circle_Y = pixely + radius * sin(L)
for i in range(len(circle_Y)):
if cv2.pointPolygonTest(contous, (circle_X[i], circle_Y[i]), False) < 0: # 如果圓散集有在輪廓之外的點
return 0
while big_r - small_r >= precision: # 二分法尋找最大半徑
half_r = (small_r + big_r) / 2
circle_X = pixelx + half_r * cos(L)
circle_Y = pixely + half_r * sin(L)
if_out = False
for i in range(len(circle_Y)):
if cv2.pointPolygonTest(contous, (circle_X[i], circle_Y[i]), False) < 0: # 如果圓散集有在輪廓之外的點
big_r = half_r
if_out = True
if not if_out:
small_r = half_r
radius = small_r
return radius
if __name__ == '__main__':
max_circle('thresh_crack.png')