joining dotted line by interpolation in an image

Question

I have this image as shown below. It is a binary mask

I created this image using the below code. Basically I got the x_idx, y_idx for just those white pixels, and I know the actual image size, so I first create an empty array and filled those lines with the help of x_idx and y_idx

image = np.empty((x_shape, y_shape))

def line_to_img(linedf, image):
    
    x_idx = linedf.x
    y_idx = linedf.y

    for i,j in zip(x_idx, y_idx):
        image[i, j] = 1
            
    return image

as you can see all pixels are the same except 2 lines, one on the left and one in right.

As you can see that the right line is not continuous and I want to make this line continuous by some interpolation method

I've tried to work on two different methods to achieve this, but no luck so far

1st Method using skimage

new_image = skimage.morphology.remove_small_holes(old_image, 40, connectivity=2, in_place=False)

Output:

Interpretation of output: Same image without any interpolation

2nd Method using cv2

kernel = cv2.getStructuringElement(cv2.MORPH_ELLIPSE,(3,3))
new_image = cv2.morphologyEx(old_image,cv2.MORPH_OPEN,kernel)

Output:

Interpretation of output: Line got removed for some reason

Please help me on how to achieve this task and interpolate the line in the image to get a continuous line

EDIT (Use-Case): Basically I got the x_idx, y_idx for just those white pixels, and I know the actual image size, so I first create an empty array and filled those lines with the help of x_idx and y_idx I don't have control of the data, this is what it is, now I want to join the line on the right size. Basically, I've to create segmentation labels, in which above the line is one label and below the line is one label, left side is fine, I can divide the image into two labels based on that line, while the middle portion will remain for class 1 i.e., Upper part, while I know for sure that right side is a single line, it's just that the data I got is degraded, so I want this interpolation to come into the picture

Answer 1

Since you're just looking to connect the region gaps in data, the Bresenham algorithm (which many know as a common line-drawing algorithm) should perform well in this case.

https://en.wikipedia.org/wiki/Bresenham%27s_line_algorithm

Pseudo Algorithm:

1. Get all (x, y) coordinate pairs from the binary mask
2. Sort from left to right by X (this assumes you'll have a horizontal line, and if not, you may choose to sort by a different means for various segmentation masks)
3. Iterate over each coordinate pair and use Bresenham to connect them.

Implementation:

import numpy as np
import cv2
from matplotlib import pyplot as plt


def _bresenham(x0: int, y0: int, x1: int, y1: int):
    dx = x1 - x0
    dy = y1 - y0

    xsign = 1 if dx > 0 else -1
    ysign = 1 if dy > 0 else -1

    dx = abs(dx)
    dy = abs(dy)

    if dx > dy:
        xx, xy, yx, yy = xsign, 0, 0, ysign
    else:
        dx, dy = dy, dx
        xx, xy, yx, yy = 0, ysign, xsign, 0

    D = 2 * dy - dx
    y = 0

    for x in range(dx + 1):
        yield x0 + x * xx + y * yx, y0 + x * xy + y * yy
        if D >= 0:
            y += 1
            D -= 2 * dx
        D += 2 * dy

# Read in image and convert to binary mask
img = cv2.imread("C:\\Test\\so1.png", 0)
ret, thresh = cv2.threshold(img, 1, 255, cv2.THRESH_BINARY)

# Get xy coordinate list of points
pairs = []
points = np.nonzero(thresh)
points_row = points[0]
points_col = points[1]

for row, col in zip(points_row, points_col):
    pairs.append((col, row))

# Sort coordinates by X
coords_sorted = sorted(pairs, key=lambda x: x[0])

# Apply bresenham algorithm
result_coords = []
for n in range(len(coords_sorted) - 1):
    for p in _bresenham(coords_sorted[n][0], coords_sorted[n][1], coords_sorted[n + 1][0], coords_sorted[n + 1][1]):
        result_coords.append(p)

# Update the binary mask with the connected lines
for x, y in result_coords:
    thresh[y][x] = 255

plt.imshow(thresh, 'gray', vmin=0, vmax=255)
plt.show()

Output Mask:

Answer 2

Here is a simple approach using the extreme points of contours.

Advantage?

This approach has a slight advantage. For every contour obtained, there are 4 extreme points; these are the topmost, bottommost, rightmost and leftmost points of a contour. We only iterate these 4 points for every contour. Unlike the approach using Bresenham algorithm which iterates through every non-zero point in the image.

Flow:

• Obtain binary image
• Perform morphological operations to join nearby dots
• Find contours
• Iterate through the extreme points found for each contour and draw a line between the closest among them.

Code:

img = cv2.imread('image_path', 0)
img1 = cv2.imread(f, 1)

# binary image
th = cv2.threshold(img, 150, 255, cv2.THRESH_BINARY+cv2.THRESH_OTSU)[1]

# morphological operations
k1 = cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (5, 5))
dilate = cv2.morphologyEx(th, cv2.MORPH_DILATE, k1)
k2 = cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (3, 3))
erode = cv2.morphologyEx(dilate, cv2.MORPH_ERODE, k2)

# find contours
cnts1 = cv2.findContours(erode, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
#cnts = cnts[0] if len(cnts1) == 2 else cnts[1]
cnts = cnts1[0]

# For each contour, find the closest distance between their extreme points and join them
for i in range(len(cnts)):
    min_dist = max(img.shape[0], img.shape[1])
    cl = []
    
    ci = cnts[i]
    ci_left = tuple(ci[ci[:, :, 0].argmin()][0])
    ci_right = tuple(ci[ci[:, :, 0].argmax()][0])
    ci_top = tuple(ci[ci[:, :, 1].argmin()][0])
    ci_bottom = tuple(ci[ci[:, :, 1].argmax()][0])
    ci_list = [ci_bottom, ci_left, ci_right, ci_top]
    
    for j in range(i + 1, len(cnts)):
        cj = cnts[j]
        cj_left = tuple(cj[cj[:, :, 0].argmin()][0])
        cj_right = tuple(cj[cj[:, :, 0].argmax()][0])
        cj_top = tuple(cj[cj[:, :, 1].argmin()][0])
        cj_bottom = tuple(cj[cj[:, :, 1].argmax()][0])
        cj_list = [cj_bottom, cj_left, cj_right, cj_top]
        
        for pt1 in ci_list:
            for pt2 in cj_list:
                dist = int(np.linalg.norm(np.array(pt1) - np.array(pt2)))     #dist = sqrt( (x2 - x1)**2 + (y2 - y1)**2 )
                if dist < min_dist:
                    min_dist = dist             
                    cl = []
                    cl.append([pt1, pt2, min_dist])
    if len(cl) > 0:
        cv2.line(erode, cl[0][0], cl[0][1], 255, thickness = 2)

The final result:

I talked about extreme points, but where are they located? The following snippet shows that:

# visualize extreme points for each contour
for c in cnts:
    left = tuple(c[c[:, :, 0].argmin()][0])
    right = tuple(c[c[:, :, 0].argmax()][0])
    top = tuple(c[c[:, :, 1].argmin()][0])
    bottom = tuple(c[c[:, :, 1].argmax()][0])
    
    # Draw dots onto image
    #cv2.drawContours(img1, [c], -1, (36, 255, 12), 2)
    cv2.circle(img1, left, 2, (0, 50, 255), -1)
    cv2.circle(img1, right, 2, (0, 255, 255), -1)
    cv2.circle(img1, top, 2, (255, 50, 0), -1)
    cv2.circle(img1, bottom, 2, (255, 255, 0), -1)

The extreme points shown above are obtained from contours of image erode.