Heuristic Canny edge detection for circular objects

    Two Canny-based edge detections with different filter sizes are
    performed to find the outmost contour of an object in a phase image
    while keeping artifacts at a minimum.

    Parameters
    ----------
    image: 2d ndarray
        Image containing an approximately spherically symmetric object
    radius: float
        The approximate object radius in pixels (required for filtering)
    mult_coarse: float
        The coarse edge detection has a filter size of
        ``sigma = mult_coarse * radius``
    mult_fine: float
        The fine edge detection has a filter size of
        ``sigma = mult_fine * radius``
    clip_rmin: float
        Removes edge points that are closer than `clip_rmin` times the
        average radial edge position from the center of the image.
    clip_rmax: float
        Removes edge points that are further than `clip_rmin` times the
        average radial edge position from the center of the image.
    maxiter: int
        Maximum number iterations for coarse edge detection, see Notes
    verbose: bool
        If set to `True`, issues EdgeDetectionWarning where applicable

    Returns
    -------
    edge : 2d boolean ndarray
        The detected edge positions of the object.

    Notes
    -----
    If no edge is found using the filter size defined by `mult_coarse`,
    then the coarse filter size is reduced by a factor of 2 until an
    edge is found or until `maxiter` is reached.

    The edge found using the filter size defined by `mult_fine` is
    heuristically filtered (parts at the center and at the edge of the
    image are removed). This heuristic filtering assumes that the
    circular object is centered in the image.

    See Also
    --------
    skimage.feature.canny: Canny edge detection algorithm used

def contour_canny(image, radius, mult_coarse=.40, mult_fine=.1,
                  clip_rmin=.9, clip_rmax=1.1, maxiter=20,
                  verbose=True):
    """Heuristic Canny edge detection for circular objects

    Two Canny-based edge detections with different filter sizes are
    performed to find the outmost contour of an object in a phase image
    while keeping artifacts at a minimum.

    Parameters
    ----------
    image: 2d ndarray
        Image containing an approximately spherically symmetric object
    radius: float
        The approximate object radius in pixels (required for filtering)
    mult_coarse: float
        The coarse edge detection has a filter size of
        ``sigma = mult_coarse * radius``
    mult_fine: float
        The fine edge detection has a filter size of
        ``sigma = mult_fine * radius``
    clip_rmin: float
        Removes edge points that are closer than `clip_rmin` times the
        average radial edge position from the center of the image.
    clip_rmax: float
        Removes edge points that are further than `clip_rmin` times the
        average radial edge position from the center of the image.
    maxiter: int
        Maximum number iterations for coarse edge detection, see Notes
    verbose: bool
        If set to `True`, issues EdgeDetectionWarning where applicable

    Returns
    -------
    edge : 2d boolean ndarray
        The detected edge positions of the object.

    Notes
    -----
    If no edge is found using the filter size defined by `mult_coarse`,
    then the coarse filter size is reduced by a factor of 2 until an
    edge is found or until `maxiter` is reached.

    The edge found using the filter size defined by `mult_fine` is
    heuristically filtered (parts at the center and at the edge of the
    image are removed). This heuristic filtering assumes that the
    circular object is centered in the image.

    See Also
    --------
    skimage.feature.canny: Canny edge detection algorithm used
    """
    image = (image - image.min()) / (image.max() - image.min())

    if radius > image.shape[0] / 2:
        msg = "`radius` in pixels exceeds image size: {}".format(radius)
        raise RadiusExceedsImageSizeError(msg)
    # 1. Perform a coarse Canny edge detection. If the edge found is empty,
    # the coarse filter size is reduced by a factor of 2.
    for ii in range(maxiter):
        fact_coarse = .5**ii
        sigma_coarse = radius * mult_coarse * fact_coarse
        edge_coarse = feature.canny(image=image,
                                    sigma=sigma_coarse)
        if np.sum(edge_coarse) != 0:
            break
    else:
        msg = "Could not find edge! Try to reducing `mult_coarse` " \
              + "or increasing `maxiter`."
        raise EdgeDetectionError(msg)

    fact_fine = .7**ii

    if fact_fine != 1 and verbose:
        msg = "The keyword argument `mult_coarse` is too large. " \
              + "If errors occur, adjust `mult_fine` as well.\n" \
              + "Given `mult_coarse`: {}\n".format(mult_coarse) \
              + "New `mult_coarse`: {}\n".format(mult_coarse * fact_coarse) \
              + "Given `mult_fine`: {}\n".format(mult_fine) \
              + "New `mult_fine`: {}".format(mult_fine * fact_fine)
        warnings.warn(msg, EdgeDetectionWarning)

    # 2. Perform a fine Canny edge detection.
    sigma_fine = radius * mult_fine * fact_fine
    edge_fine = feature.canny(image, sigma_fine)

    # 3. Remove parts from the fine edge
    # Assume that the object is centered.
    sx, sy = image.shape
    x = np.linspace(-sx / 2, sx / 2, sx, endpoint=True).reshape(-1, 1)
    y = np.linspace(-sy / 2, sy / 2, sy, endpoint=True).reshape(1, -1)
    # 3.a. Remove detected edge parts from the corners of the image
    # Radius of this disk is approximately
    ellipse = (x / sx)**2 + (y / sy)**2 < .25
    edge_fine *= ellipse
    edge_coarse *= ellipse

    if np.sum(edge_fine):
        # 3.b Also filter inside of edge
        rad = np.sqrt(x**2 + y**2)
        # Filter coarse edge with `clip_rmin` and `clip_rmax`
        rad_coarse = rad * edge_coarse
        avg_coarse = np.sum(rad_coarse) / np.sum(edge_coarse)
        rad_coarse[rad_coarse < avg_coarse * clip_rmin] = 0
        rad_coarse[rad_coarse > avg_coarse * clip_rmax] = 0
        # Filter inside of fine edge with smallest radius of
        # coarse edge, i.e. `rad_coarse.min()`.
        rad_fine = rad * edge_fine
        if np.sum(rad_coarse):
            edge_fine[rad_fine < rad_coarse[rad_coarse != 0].min()] = 0
        # Filter outside of fine edge with `clip_rmax` twice
        for __ in range(2):
            rad_fine = rad * edge_fine
            avg_fine = np.sum(rad_fine) / np.sum(edge_fine)
            edge_fine[rad_fine > avg_fine * clip_rmax] = 0
    elif np.sum(edge_coarse):
        # No fine edge detected.
        edge_fine = edge_coarse
    else:
        msg = "Could not find edge! Try reducing `mult_coarse` " \
              + "and `mult_fine`."
        raise EdgeDetectionError(msg)
    # make sure there are more than 4 points
    if np.sum(edge_fine) < 4:
        msg = "Detected edge too small! Try increasing `maxiter`, " \
              + "modifying `radius`, or reducing `mult_coarse`."
        raise EdgeDetectionError(msg)
    return edge_fine