Extracting the (size of the grid)features from an image

Question

I am working with a image

Actualimage:

I have used techniques like Canny, Hough transform to detect the line, and got this

Output:

Now I would like to extract the features like thickness size on each side in the grid cell indicated in the image. I need calculate area of specific

Part:

Can you enlighten me on this?

# Dilate the edge image to make the edges thicker
   dialted_edge_image = cv2.dilate(edge_canny_image,kernel,iterations = 2)  

  # Perform a Hough Transform to detect lines
  lines = probabilistic_hough_line(edge_canny_image, threshold=40, line_length=1, line_gap=2)

  # Create a separate image for line detection
  line_detected_image = np.dstack([edge_canny_image] * 3)  # Convert to RGB for colored lines
   for line in lines:
    p0, p1 = line
    cv2.line(line_detected_image, p0, p1, (255, 255, 255), 2)  `

Answer 1

Okay, i am not sure with all the calculation you need but i can give you an approach to find the thicknesses of the sides.

To be able to provide accuracy and stability for future calculations:

1. You need to make sure the light levels of the image stays relatively constant with each different sample so that the inRange method gives proper filtering.
2. I used spesific roi (region of interst) for thickness calculation for each edge, you need to make sure each brightedge stays inside corresponding roi for different images.

I crop the image near the bright edges and since the thicknes changes accross the edge i measure it from different locations and calculate min max and average thickness. You may use the one suits your needs best.

and this is my code:

import cv2

def estimateThickness(img):
    #Determine if it is a vertical or a horizantal edege
    height,width = img.shape
    if height<=width:
        img = cv2.rotate(img,cv2.ROTATE_90_CLOCKWISE)
    height,width = img.shape

    #Estimate the thickness of sides from various locations
    #and extract min max average thickness
    thicknesess = []
    for nh in range(height//10):
        first,last = None,None
        for nw in range(width):
            # print(nl,ns)
            #Find the first white pixel on the direction 
            if img[10*nh][nw] == 255 and first is None:
                first = nw
            #Find the last white pixel on the direction 
            if img[10*nh][width-nw-1] == 255 and last is None:
                last = width-nw-1
            
            if first is not None and last is not None:
                thicknesess.append(last-first)

    return max(thicknesess),min(thicknesess),sum(thicknesess)/len(thicknesess)


#Read the image
src_image = cv2.imread('img\grid.png')
gray = cv2.cvtColor(src_image,cv2.COLOR_BGR2GRAY)


#Extract the bright part in the image and filter the rest to measure thickness
bright_part = cv2.inRange(gray,110,255)
bright_part = cv2.morphologyEx(bright_part,cv2.MORPH_OPEN,cv2.getStructuringElement(cv2.MORPH_RECT,(3,3)))
bright_part = cv2.morphologyEx(bright_part,cv2.MORPH_CLOSE,cv2.getStructuringElement(cv2.MORPH_RECT,(15,15)))

#Crop top left bot and right edges from the filtered image
left_edge = bright_part[200:500,180:280]
right_edge = bright_part[200:500,750:850]
top_edge = bright_part[20:120,400:700]
bot_edge = bright_part[580:680,400:700]

#Use the defined function with cropped image
minL,maxL,avgL = estimateThickness(left_edge)
minR,maxR,avgR = estimateThickness(right_edge)
minT,maxT,avgT = estimateThickness(top_edge)
minB,maxB,avgB = estimateThickness(bot_edge)

print('L',minL,maxL,avgL)
print('R',minR,maxR,avgR)
print('T',minT,maxT,avgT)
print('B',minB,maxB,avgB)
    
cv2.imshow('L',left_edge)
cv2.imshow('R',right_edge)
cv2.imshow('T',top_edge)
cv2.imshow('B',bot_edge)

cv2.imshow('Bright Part',bright_part)
cv2.imshow('Source',src_image)
key = cv2.waitKey(0)

For your other calculations if you can explain them further i can try to help you.

Answer 2

Image segmentation

The first sub-problem is to isolate the three different structures: the grid, the seam, and the negative space within each grid cell. As the pixel intensities fall into three barely overlapping distributions, this can be done with a multi-otsu threshold:

Grid cell width

The second sub-problem is determining the grid cell width. The approach taken here is to find the area of the large negative space in the middle, and then computing the square root to get the width in pixel: 519.1 px

Grid line width

The third sub-problem is finding the width of the grid lines. The basic idea is to create a clean mask of all grid lines, skeletonize that mask, and then divide the number of non-zero mask pixels by the number of skeleton pixels. However, the skeleton is far from perfect, as specifically the ends of each grid line are not well captured. Instead we process each straight identifiable grid line segment separately. The specific steps taken here are:

1. Isolate and clean up the grid mask.
2. Skeletonize the grid mask.
3. Convert non-zero skeleton mask pixels to a connectivity graph, in which nodes are connected, if the corresponding non-zero pixels are adjacent to each other.
4. Find all chains in the connectivity graph.
5. Isolate the grid segment corresponding to each chain.
6. Count the number of non-zero pixels in the segment and divide by the extent of the chain to determine the average width of the grid segment: 50.5 px

Seam width

The last sub-problem is finding the width of the seams. As the seam width is variable, it seemed most useful to determine the distribution of seam widths. The approach taken here is to use profile lines that span from the centroid of the largest grid cell to its perimeter, and determine the number of non-zero pixels along each profile line.

Code

import numpy as np
from matplotlib import pyplot as plt
import networkx as nx

from itertools import product
from skimage.io import imread
from skimage.filters import threshold_multiotsu
from skimage.morphology import (
    binary_opening,
    binary_closing,
    skeletonize,
)
from skimage.measure import (
    label,
    find_contours,
    centroid,
    profile_line,
)
from skimage.draw import polygon2mask


def to_graph(skeleton_image, connectivity=8):
    """Convert a skeleton image to a networkx graph object."""
    perpendicular_steps = set([(-1, 0), (1, 0), (0, -1), (0, 1)])
    diagonal_steps = set([(-1, 1), (-1, -1), (1, -1), (1, 1)])
    if connectivity == 8:
        allowed_steps = perpendicular_steps.union(diagonal_steps)
    elif connectivity == 4:
        allowed_steps = perpendicular_steps
    else:
        raise ValueError(f"The parameter connectivity is either 4 or 8 not {connectivity}.")

    nodes = list(zip(*np.where(skeleton_image)))
    edges = [((x, y), (x+dx, y+dy)) for (x, y), (dx, dy) in product(nodes, allowed_steps) if (x+dx, y+dy) in nodes]

    return nx.Graph(edges)


def get_orthogonal_unit_vector(v):
    """Determine the orthogonal unit vector to a given vector.

    Parameters
    ----------
    v : numpy.array
        The input vector.

    Returns
    -------
    w : numpy.array
        The output vector.

    Notes
    -----
    Adapted from https://stackoverflow.com/a/16890776/2912349

    """
    v = np.atleast_2d(v)
    if not np.all(np.isclose(v, 0)):
        v = v / np.linalg.norm(v, axis=-1)[:, None] # unit vector
        w = np.c_[-v[:,1], v[:,0]]                  # orthogonal vector
        w = w / np.linalg.norm(w, axis=-1)[:, None] # orthogonal unit vector
        return np.squeeze(w)
    else:
        raise ValueError("Cannot determine the orthogonal vector. Input vector has zero length.")


if __name__ == "__main__":

    img = imread("~/wdir/tmp/grid.jpg")
    bw = img.mean(axis=-1)

    # --------------------------------------------------------------------------------
    # image decomposition

    fig, axes = plt.subplots(1, 3, figsize=(10, 5))
    axes = axes.ravel()
    axes[0].imshow(bw, cmap="gray")
    axes[1].hist(bw.ravel(), bins=100)*2
    thresholds = threshold_multiotsu(bw)
    for threshold in thresholds:
        axes[1].axvline(threshold, color="k", linestyle="--")
    axes[1].set_xlabel("Pixel intensity")
    axes[1].set_ylabel("Count")
    regions = np.digitize(bw, bins=thresholds)
    axes[2].imshow(regions, cmap="bwr")
    axes[1].set_aspect("auto")
    fig.tight_layout()

    # --------------------------------------------------------------------------------
    # determine grid cell width by finding the area of the largest negative space (LNS)

    fig, axes = plt.subplots(1, 3, figsize=(10, 5))
    axes = axes.ravel()

    negative_space = regions < regions.max()
    cleaned = binary_opening(negative_space, np.ones((25, 25), dtype=bool))
    axes[0].imshow(cleaned, cmap="gray")

    labelled = label(cleaned)
    largest_negative_space_label = np.argmax(np.bincount(labelled.ravel())[1:]) + 1
    largest_negative_space_mask = labelled == largest_negative_space_label
    axes[1].imshow(largest_negative_space_mask, cmap="gray")

    lns_contour = sorted(find_contours(largest_negative_space_mask), key=lambda x: len(x))[-1]
    axes[2].imshow(bw, cmap="gray")
    axes[2].plot(lns_contour[:, 1], lns_contour[:, 0], color="blue")

    area = largest_negative_space_mask.sum()
    interior_width = np.sqrt(area)
    print(f"Grid cell interior width: {interior_width:.1f} pixels")

    lns_centroid = centroid(largest_negative_space_mask)
    row, col = lns_centroid
    x = [col - interior_width/2, col + interior_width/2]
    y = [row, row]
    axes[2].plot(x, y, color="yellow")

    fig.tight_layout()

    # --------------------------------------------------------------------------------
    # determine grid line width

    fig, axes = plt.subplots(2, 2)
    axes = axes.ravel()
    axes[0].imshow(bw, cmap="gray")

    mask_gridline = regions == regions.max()
    axes[1].imshow(mask_gridline, cmap="gray")

    # morphological cleaning
    d = 20
    selem = np.ones((d, d), dtype=bool)
    clean_gridline = binary_opening(binary_closing(mask_gridline, selem), selem)

    # skeletonize
    skeleton = skeletonize(clean_gridline)

    axes[2].imshow(clean_gridline.astype(int) + skeleton.astype(int), cmap="gray")
    axes[3].imshow(clean_gridline.astype(int) + skeleton.astype(int), cmap="gray")

    # get a first estimate of the width
    area = clean_gridline.sum()
    length = skeleton.sum()
    estimate = area / length

    # decompose skeleton into chains and estimate area of each chain
    g = to_graph(skeleton)
    chains = list(nx.connected_components(nx.subgraph(g, [node for node, degree in g.degree() if degree == 2])))
    chain_widths = []
    chain_lengths = []
    for chain in chains:
        if len(chain) > 100: # exclude short chains
            # ensure nodes in chain are correctly ordered
            start, stop = [node for node, degree in nx.subgraph(g, chain).degree() if degree == 1]
            chain = nx.shortest_path(g, start, stop)
            chain = np.array(list(chain))

            # find a polygon enclosing the grid segment corresponding to the chain
            start = chain[0]
            stop = chain[-1]
            delta = stop - start
            v = get_orthogonal_unit_vector(np.atleast_2d(delta)) * (1.5 * estimate / 2)
            polygon = np.array([start - v, start + v, stop + v, stop - v], dtype=int)
            polygon_mask = polygon2mask(bw.shape, polygon)

            # determine segment area and average width
            area = clean_gridline[polygon_mask].sum()
            length = np.linalg.norm(delta)
            width = area / length
            chain_widths.append(width)
            chain_lengths.append(length)

            # plot individual estimates
            color = np.random.rand(3)
            axes[3].plot(chain[:, 1], chain[:, 0], color=color)
            axes[3].plot(np.r_[polygon[:, 1], polygon[0, 1]], np.r_[polygon[:, 0], polygon[0, 0]], color=color)
            r, c = chain[int(len(chain) / 2)]
            dr, dc = get_orthogonal_unit_vector(np.atleast_2d(delta)) * width / 2
            axes[3].plot([c - dc, c + dc], [r - dr, r + dr], color=color)

    # estimate mean width weighted by chain length
    gridline_width = np.sum(np.array(chain_widths) * np.array(chain_lengths) / np.sum(chain_lengths).astype(float))
    print(f"Grid line width: {gridline_width:.1f} pixels")

        # --------------------------------------------------------------------------------
    # determine seam widths

    seam = regions == 1

    fig, axes = plt.subplots(2, 2)
    axes = axes.ravel()
    axes[0].imshow(seam, cmap="gray")

    # isolate seam within largest grid cell
    clean_seam = seam.copy()
    clean_seam[clean_gridline] = 0
    clean_seam[~largest_negative_space_mask] = 0
    axes[1].imshow(clean_seam, cmap="gray")

    # walk around contour; determine seam widths
    seam_thickness = []
    src = lns_centroid
    for ii, idx in enumerate(np.linspace(0, len(lns_contour), 72)[:-1]): # i.e. every 5 degrees
        dst = lns_contour[int(idx)]
        color = np.random.rand(3)
        if (ii % 2) == 0: # i.e. every 10 degrees
            axes[1].plot([src[1], dst[1]], [src[0], dst[0]], color=color)
        y = profile_line(clean_seam, src, dst)
        seam_thickness.append(np.sum(y))

    axes[2].plot(y, color=color)
    axes[2].set_ylabel("Intensity")
    axes[2].set_xlabel("Line length")
    axes[2].set_title("Example profile line")
    axes[3].hist(seam_thickness)
    axes[3].set_ylabel("Count")
    axes[3].set_xlabel("Seam thickness [pixel]")
    fig.tight_layout()

    plt.show()