Как обнаружить круговую эрозию / расширение

Question

Я хотел бы обнаружить круговые эрозии и расширения на линии. Для расширений я пытался рекурсивно размывать изображение, и на каждой рекурсии я проверял соотношение сторон ширина / высота. Если отношение меньше 4, я предполагаю, что контур является круглым, и для каждого такого контура я вычисляю центр и радиус круга из моментов и площади. Это функция, которая обнаруживает круговые расширения:

def detect_circular_dilations(img, contours):
    contours_current, hierarchy = cv2.findContours(img, cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE)
    if len(contours_current) == 0:
        return get_circles_from_contours(contours)
    for c in contours_current:
        x, y, w, h = cv2.boundingRect(c)
        if w > h:
            aspect_ratio = float(w) / h
        else:
            aspect_ratio = float(h) / w
        if aspect_ratio < 4 and w < 20 and h < 20 and w > 5 and h > 5:
            contours.append(c)
    return detect_circular_dilations(cv2.erode(img, None, iterations=1), contours)

Примером круговых дилатаций, которые я хочу обнаружить, являются следующие:

Другая проблема, которую я не решил, - это обнаружение круговых эрозий. Пример круговой эрозии следующий:

Здесь я отметил круговую эрозию, которую я хотел бы обнаружить красным прямоугольником. Могут быть некоторые меньшие круглые образцы (слева), которые не должны рассматриваться как фактическая круговая эрозия.

Кто-нибудь знает, как лучше всего обнаружить такие круглые формы? Для круговых расширений я был бы признателен за любые комментарии / предложения, чтобы потенциально сделать обнаружение более устойчивым

Спасибо!

Answer 1

Я бы попытался найти два края линии с cv2.Canny() и найти контуры. Если вы отсортируете свой контур по ширине ограничивающего прямоугольника, первые два контура будут вашими линиями. После этого вы можете рассчитать минимальное расстояние каждой точки на одном ребре до другого ребра. Затем вы можете вычислить медиану расстояний и сказать, что если точка имеет большее или меньшее расстояние, чем медиана (+ - допуск), чем эта точка, то это расширение или эрозия линии и добавьте ее в список. При необходимости вы можете сортировать шумы, просматривая списки, и удалять точки, если они не последовательны (по оси x).

Вот простой пример:

import cv2
import numpy as np
from scipy import spatial

def detect_dilation(median, mindist, tolerance):
    count = 0
    for i in mindist:
        if i > median + tolerance:
            dilate.append((reshape_e1[count][0], reshape_e1[count][1]))
        elif i < median - tolerance:
            erode.append((reshape_e1[count][0], reshape_e1[count][1]))
        else:
            pass
        count+=1

def other_axis(dilate, cnt):
    temp = []
    for i in dilate:
        temp.append(i[0])
    for i in cnt:
        if i[0] in temp:
            dilate.append((i[0],i[1]))

img = cv2.imread('1.jpg')
gray = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
edges = cv2.Canny(gray,100,200)
_, contours, hierarchy = cv2.findContours(edges,cv2.RETR_TREE,cv2.CHAIN_APPROX_NONE)
contours.sort(key= lambda cnt :cv2.boundingRect(cnt)[3])
edge_1 = contours[0]
edge_2 = contours[1]
reshape_e1 = np.reshape(edge_1, (-1,2))
reshape_e2 =np.reshape(edge_2, (-1,2))
tree = spatial.cKDTree(reshape_e2)
mindist, minid = tree.query(reshape_e1)
median = np.median(mindist)
dilate = []
erode = []
detect_dilation(median,mindist,5)
other_axis(dilate, reshape_e2)
other_axis(erode, reshape_e2)

dilate = np.array(dilate).reshape((-1,1,2)).astype(np.int32)
erode = np.array(erode).reshape((-1,1,2)).astype(np.int32)
x,y,w,h = cv2.boundingRect(dilate)
cv2.rectangle(img,(x,y),(x+w,y+h),(255,0,0),2)
x,y,w,h = cv2.boundingRect(erode)
cv2.rectangle(img,(x,y),(x+w,y+h),(0,0,255),2)

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

Результат:

Edit:

Если на рисунке есть ломаная линия (что означает больше контуров), вам придется рассматривать каждый контур как отдельную линию. Этого можно достичь, создав интересующий регион с помощью cv2.boundingRect(). Но, как я попробовал с новой загруженной картинкой, процесс не очень надежный, так как вам нужно изменить допуск, чтобы получить желаемый результат. Поскольку я не знаю, как выглядят другие изображения, вам может понадобиться лучший способ получить среднее расстояние и коэффициент допуска. Любой способ здесь является примером того, что я описал (с 15 для допуска):

import cv2
import numpy as np
from scipy import spatial

def detect_dilation(median, mindist, tolerance):
    count = 0
    for i in mindist:
        if i > median + tolerance:
            dilate.append((reshape_e1[count][0], reshape_e1[count][1]))
        elif i < median - tolerance:
            erode.append((reshape_e1[count][0], reshape_e1[count][1]))
        else:
            pass
        count+=1

def other_axis(dilate, cnt):
    temp = []
    for i in dilate:
        temp.append(i[0])
    for i in cnt:
        if i[0] in temp:
            dilate.append((i[0],i[1]))

img = cv2.imread('2.jpg')
gray_original = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
_, thresh_original = cv2.threshold(gray_original, 0, 255, cv2.THRESH_BINARY+cv2.THRESH_OTSU)

# Filling holes
_, contours, hierarchy = cv2.findContours(thresh_original,cv2.RETR_CCOMP,cv2.CHAIN_APPROX_SIMPLE)
for cnt in contours:
    cv2.drawContours(thresh_original,[cnt],0,255,-1)
_, contours, hierarchy = cv2.findContours(thresh_original,cv2.RETR_EXTERNAL,cv2.CHAIN_APPROX_NONE)

for cnt in contours:
    x2,y,w2,h = cv2.boundingRect(cnt)
    thresh = thresh_original[0:img.shape[:2][1], x2+20:x2+w2-20] # Region of interest for every "line"
    edges = cv2.Canny(thresh,100,200)
    _, contours, hierarchy = cv2.findContours(edges,cv2.RETR_TREE,cv2.CHAIN_APPROX_NONE)
    contours.sort(key= lambda cnt: cv2.boundingRect(cnt)[3])
    edge_1 = contours[0]
    edge_2 = contours[1]
    reshape_e1 = np.reshape(edge_1, (-1,2))
    reshape_e2 =np.reshape(edge_2, (-1,2))
    tree = spatial.cKDTree(reshape_e2)
    mindist, minid = tree.query(reshape_e1)
    median = np.median(mindist)
    dilate = []
    erode = []
    detect_dilation(median,mindist,15)
    other_axis(dilate, reshape_e2)
    other_axis(erode, reshape_e2)
    dilate = np.array(dilate).reshape((-1,1,2)).astype(np.int32)
    erode = np.array(erode).reshape((-1,1,2)).astype(np.int32)
    x,y,w,h = cv2.boundingRect(dilate)
    if len(dilate) > 0:
        cv2.rectangle(img[0:img.shape[:2][1], x2+20:x2+w2-20],(x,y),(x+w,y+h),(255,0,0),2)
    x,y,w,h = cv2.boundingRect(erode)
    if len(erode) > 0:
        cv2.rectangle(img[0:img.shape[:2][1], x2+20:x2+w2-20],(x,y),(x+w,y+h),(0,0,255),2)

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

Результат:

Answer 2

Подобные проблемы часто решаются с помощью преобразования расстояния и преобразования средней оси.Они связаны в некоторой степени, так как медиальная ось проходит вдоль гребня трансформации расстояния.Общая идея такова:

1. Вычислить преобразование расстояния изображения (для каждого пикселя переднего плана возвращает расстояние до ближайшего пикселя фона; некоторые библиотеки реализуют это другим способом, и в этом случаевам нужно вычислить преобразование расстояния инвертированного изображения).

2. Вычислить медиальную ось (или скелет).

3. ЗначенияРасстояние преобразования по медиальной оси являются соответствующими значениями, мы игнорируем все остальные пиксели.Здесь мы видим локальный радиус линии.

4. Локальные максимумы являются центроидами дилатаций.Используйте пороговое значение, чтобы определить, какие из них являются важными дилатациями, а какие - нет (зашумленный контур вызовет множество локальных максимумов).

5. Локальные минимумы являются центроидами эрозии.

Например, я получил следующий вывод, используя приведенный ниже код MATLAB.

Вот код, который я использовал,Он использует MATLAB с DIPimage 3 , просто как быстрое доказательство принципа.Это должно быть просто для перевода на Python с любой библиотекой обработки изображений, которую вы хотите использовать.

% Read in image and remove the red markup:
img = readim('https://i.stack.imgur.com/bNOTn.jpg');
img = img{3}>100;
img = closing(img,5);

% This is the algorithm described above:
img = fillholes(img);               % Get rid of holes
radius = dt(img);                   % Distance transform
m = bskeleton(img);                 % Medial axis
radius(~m) = 0;                     % Ignore all pixels outside the medial axis
detection = dilation(radius,25)==radius & radius>25; % Local maxima with radius > 25
pos = findcoord(detection);         % Coordinates of detections
radius = double(radius(detection)); % Radii of detections

% This is just to make the markup:
detection = newim(img,'bin');
for ii=1:numel(radius)
   detection = drawshape(detection,2*radius(ii),pos(ii,:),'disk');
end
overlay(img,detection)