DEV Community

David Haley
David Haley

Posted on

Algorithm edge cases vs notation

I just spent 2.5 hours debugging an infinite loop in this algorithm, caused by an edge case (center point of structuring element is excluded). I'd missed a single but crucial requirement…

Here is the looping code:

def opencv_reconstruct(
    marker: np.ndarray, mask: np.ndarray, kernel: np.ndarray, anchor: tuple = (-1, -1)
):
    # .item() converts the numpy scalar to a python scalar
    pad_value = np.min(marker).item()

    while True:
        expanded = cv2.dilate(
            src=marker,
            kernel=kernel,
            borderType=cv2.BORDER_CONSTANT,
            borderValue=pad_value,
            anchor=anchor,
        )
        expanded = np.fmin(expanded, mask)

        # Termination criterion: Expansion didn't change the image at all
        if (marker == expanded).all():
            return expanded

        marker = expanded
Enter fullscreen mode Exit fullscreen mode

See below for the input data I was using.

I was using the footprint:

100
000
101
Enter fullscreen mode Exit fullscreen mode

Crucially, the center point is false– which I gather is a bit unusual in these operations, every example I've ever seen has a true in the center.

Anyhow, I was confused, because Vincent'93 says the algorithm terminates but mine never converged. Let's say the above is an implementation of the following:

Screenshot of Vincent'93 parallel reconstruction algorithm

opencv for its part takes the footprint very literally. If you didn't ask for the center point, it won't give you the center point.

But Vincent'93 actually defines reconstruction as including the center point:

Vincent '93 dilation step of parallel reconstruction

In English this is saying,

  • for each pixel in the image,
  • the dilated pixel value is equal to the maximum of:
    • the pixel value of each point q, where q is in the neighborhood of p, or p itself

It's that last little that does the trick…

If we change the algorithm to select the center point (hackity hack)

kernel[1, 1] = True
Enter fullscreen mode Exit fullscreen mode

then … voilà, no more infinite loop (and we get the correct answer to boot).

In terms of my algorithm, my bug is here:

    for footprint_row_offset in range(-footprint_center_row, footprint_center_row + 1):
        for footprint_col_offset in range(
            -footprint_center_col, footprint_center_col + 1
        ):
            # Skip this point if not in the footprint.
            footprint_x = footprint_center_row + footprint_row_offset
            footprint_y = footprint_center_col + footprint_col_offset
            if not footprint[footprint_x, footprint_y]:
                continue
Enter fullscreen mode Exit fullscreen mode

I skip if its not in the footprint– I should also ensure the center is always kept.

That union symbol really mattered 😎 all because of an odd footprint that doesn't include the center.

Input data:

    seed = np.array(
        [
            [
                49,
                -1,
                -54,
                -16,
                55,
                -20,
                56,
                41,
                -5,
                -28,
            ],
            [
                -18,
                56,
                -16,
                56,
                -26,
                56,
                41,
                -6,
                7,
                7,
            ],
            [
                -17,
                37,
                56,
                -16,
                56,
                41,
                -35,
                -69,
                -15,
                95,
            ],
            [
                -31,
                -17,
                -18,
                -2,
                -19,
                18,
                41,
                45,
                95,
                -0,
            ],
            [
                31,
                7,
                42,
                -22,
                10,
                -52,
                45,
                95,
                4,
                28,
            ],
            [
                31,
                42,
                31,
                10,
                31,
                -50,
                -21,
                45,
                95,
                2,
            ],
            [
                42,
                31,
                31,
                44,
                86,
                -54,
                107,
                14,
                45,
                17,
            ],
            [
                -103,
                45,
                31,
                66,
                -63,
                66,
                40,
                -64,
                40,
                49,
            ],
            [
                -46,
                41,
                42,
                -28,
                107,
                45,
                -36,
                60,
                -26,
                18,
            ],
            [
                41,
                -29,
                -37,
                45,
                23,
                48,
                40,
                -26,
                4,
                -40,
            ],
        ],
        dtype=np.int16,
    )
    mask = np.array(
        [
            [
                49,
                -1,
                -54,
                35,
                55,
                -20,
                122,
                46,
                17,
                -28,
            ],
            [
                37,
                57,
                -16,
                77,
                -26,
                109,
                52,
                -6,
                7,
                7,
            ],
            [
                36,
                91,
                56,
                27,
                56,
                94,
                -35,
                -69,
                -15,
                126,
            ],
            [
                -31,
                -17,
                -18,
                -2,
                -19,
                18,
                41,
                99,
                112,
                -0,
            ],
            [
                31,
                7,
                44,
                -22,
                92,
                -52,
                102,
                123,
                4,
                28,
            ],
            [
                50,
                82,
                81,
                10,
                61,
                -50,
                -21,
                45,
                95,
                2,
            ],
            [
                42,
                105,
                31,
                52,
                86,
                -54,
                124,
                14,
                72,
                17,
            ],
            [
                -103,
                97,
                31,
                107,
                -63,
                126,
                60,
                -64,
                124,
                93,
            ],
            [
                -46,
                48,
                45,
                -28,
                107,
                86,
                -36,
                113,
                -0,
                18,
            ],
            [
                51,
                125,
                -37,
                66,
                23,
                48,
                45,
                -26,
                4,
                -40,
            ],
        ],
        dtype=np.int16,
    )
Enter fullscreen mode Exit fullscreen mode

Top comments (0)