From b29ad8ef6d9086e8908c6dced1bf5784aa1b3697 Mon Sep 17 00:00:00 2001
From: Robin Wilson <robin@rtwilson.com>
Date: Fri, 28 Aug 2015 22:30:42 +0100
Subject: [PATCH] Added quantile_threshold option to canny edge detection

This allows you to specify the high and low thresholds as quantiles of the edge magnitude
image, rather than as absolute edge magnitude values
---
 skimage/feature/_canny.py           | 14 +++++++++++++-
 skimage/feature/tests/test_canny.py | 27 +++++++++++++++++++++++++++
 2 files changed, 40 insertions(+), 1 deletion(-)

diff --git a/skimage/feature/_canny.py b/skimage/feature/_canny.py
index 73566d00..91598ce0 100644
--- a/skimage/feature/_canny.py
+++ b/skimage/feature/_canny.py
@@ -50,7 +50,8 @@ def smooth_with_function_and_mask(image, function, mask):
     return output_image
 
 
-def canny(image, sigma=1., low_threshold=None, high_threshold=None, mask=None):
+def canny(image, sigma=1., low_threshold=None, high_threshold=None, mask=None,
+    quantile_threshold=False):
     """Edge filter an image using the Canny algorithm.
 
     Parameters
@@ -67,6 +68,9 @@ def canny(image, sigma=1., low_threshold=None, high_threshold=None, mask=None):
         If None, high_threshold is set to 20% of dtype's max.
     mask : array, dtype=bool, optional
         Mask to limit the application of Canny to a certain area.
+    quantile_threshold : bool, optional
+        If True then treat low_threshold and high_threshold as quantiles of the
+        edge magnitude image, rather than absolute edge magnitude values.
 
     Returns
     -------
@@ -246,6 +250,14 @@ def canny(image, sigma=1., low_threshold=None, high_threshold=None, mask=None):
     c2 = magnitude[1:, :-1][pts[:-1, 1:]]
     c_minus = c2 * w + c1 * (1 - w) <= m
     local_maxima[pts] = c_plus & c_minus
+
+    #
+    #---- If quantile_threshold is set then calculate the thresholds to use
+    #
+    if quantile_threshold:
+        high_threshold = np.percentile(magnitude, 100.0 * high_threshold)
+        low_threshold = np.percentile(magnitude, 100.0 * low_threshold)
+
     #
     #---- Create two masks at the two thresholds.
     #
diff --git a/skimage/feature/tests/test_canny.py b/skimage/feature/tests/test_canny.py
index 400a7f70..0bfb1fcf 100644
--- a/skimage/feature/tests/test_canny.py
+++ b/skimage/feature/tests/test_canny.py
@@ -1,7 +1,10 @@
 import unittest
 import numpy as np
+from numpy.testing import assert_equal
 from scipy.ndimage import binary_dilation, binary_erosion
 import skimage.feature as F
+from skimage import filters, data
+from skimage import img_as_float
 
 
 class TestCanny(unittest.TestCase):
@@ -66,3 +69,27 @@ class TestCanny(unittest.TestCase):
         result1 = F.canny(np.zeros((20, 20)), 4, 0, 0, np.ones((20, 20), bool))
         result2 = F.canny(np.zeros((20, 20)), 4, 0, 0)
         self.assertTrue(np.all(result1 == result2))
+
+
+def test_quantile_threshold():
+    image = img_as_float(data.camera()[::50,::50])
+
+    # Correct output produced manually with quantiles
+    # of 0.8 and 0.6 for high and low respectively
+    correct_output = np.array([[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
+       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],
+       [0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0],
+       [0, 0, 1, 1, 0, 0, 1, 0, 0, 0, 0],
+       [0, 0, 1, 0, 0, 0, 1, 0, 0, 0, 0],
+       [0, 1, 1, 0, 0, 0, 1, 0, 0, 0, 0],
+       [0, 1, 0, 0, 0, 1, 0, 0, 0, 0, 0],
+       [0, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0],
+       [0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0],
+       [0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0],
+       [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]], dtype=bool)
+
+    result = F.canny(image, low_threshold=0.6, high_threshold=0.8, quantile_threshold=True)
+
+    assert_equal(result, correct_output)
+
+