diff --git a/skimage/filters/LICENSE.txt b/skimage/filters/LICENSE.txt new file mode 100644 index 00000000..3507299a --- /dev/null +++ b/skimage/filters/LICENSE.txt @@ -0,0 +1,28 @@ +threshold_li based on imagej code: + +Copyright (c) 2009 - 2015, Board of Regents of the University of +Wisconsin-Madison, Broad Institute of MIT and Harvard, and Max Planck +Institute of Molecular Cell Biology and Genetics. +All rights reserved. + +Redistribution and use in source and binary forms, with or without modification, +are permitted provided that the following conditions are met: + +1. Redistributions of source code must retain the above copyright notice, this + list of conditions and the following disclaimer. + +2. Redistributions in binary form must reproduce the above copyright notice, + this list of conditions and the following disclaimer in the documentation + and/or other materials provided with the distribution. + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE +ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE +LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +POSSIBILITY OF SUCH DAMAGE. diff --git a/skimage/filters/__init__.py b/skimage/filters/__init__.py index 795819c6..7cffb140 100644 --- a/skimage/filters/__init__.py +++ b/skimage/filters/__init__.py @@ -9,7 +9,7 @@ from .edges import (sobel, hsobel, vsobel, sobel_h, sobel_v, from ._rank_order import rank_order from ._gabor import gabor_kernel, gabor_filter from .thresholding import (threshold_adaptive, threshold_otsu, threshold_yen, - threshold_isodata) + threshold_isodata, threshold_li) from . import rank from .rank import median @@ -66,4 +66,5 @@ __all__ = ['inverse', 'threshold_otsu', 'threshold_yen', 'threshold_isodata', + 'threshold_li', 'rank'] diff --git a/skimage/filters/tests/test_thresholding.py b/skimage/filters/tests/test_thresholding.py index 5de758e8..a47b6f20 100644 --- a/skimage/filters/tests/test_thresholding.py +++ b/skimage/filters/tests/test_thresholding.py @@ -5,6 +5,7 @@ import skimage from skimage import data from skimage.filters.thresholding import (threshold_adaptive, threshold_otsu, + threshold_li, threshold_yen, threshold_isodata) @@ -28,6 +29,17 @@ class TestSimpleImage(): image = np.float64(self.image) assert 2 <= threshold_otsu(image) < 3 + def test_li(self): + assert int(threshold_li(self.image)) == 2 + + def test_li_negative_int(self): + image = self.image - 2 + assert int(threshold_li(image)) == 0 + + def test_li_float_image(self): + image = np.float64(self.image) + assert 2 <= threshold_li(image) < 3 + def test_yen(self): assert threshold_yen(self.image) == 2 @@ -152,6 +164,25 @@ def test_otsu_astro_image(): img = skimage.img_as_ubyte(data.astronaut()) assert 109 < threshold_otsu(img) < 111 +def test_li_camera_image(): + camera = skimage.img_as_ubyte(data.camera()) + assert 63 < threshold_li(camera) < 65 + + +def test_li_coins_image(): + coins = skimage.img_as_ubyte(data.coins()) + assert 95 < threshold_li(coins) < 97 + + +def test_li_coins_image_as_float(): + coins = skimage.img_as_float(data.coins()) + assert 0.37 < threshold_li(coins) < 0.38 + + +def test_li_astro_image(): + img = skimage.img_as_ubyte(data.astronaut()) + assert 66 < threshold_li(img) < 68 + def test_yen_camera_image(): camera = skimage.img_as_ubyte(data.camera()) assert 197 < threshold_yen(camera) < 199 diff --git a/skimage/filters/thresholding.py b/skimage/filters/thresholding.py index e10bca5c..2921017b 100644 --- a/skimage/filters/thresholding.py +++ b/skimage/filters/thresholding.py @@ -1,7 +1,8 @@ __all__ = ['threshold_adaptive', 'threshold_otsu', 'threshold_yen', - 'threshold_isodata'] + 'threshold_isodata', + 'threshold_li', ] import numpy as np import scipy.ndimage @@ -302,3 +303,71 @@ def threshold_isodata(image, nbins=256, return_all=False): return thresholds else: return thresholds[0] + + +def threshold_li(image): + """Return threshold value based on adaptation of Li's Minimum Cross Entropy method. + + Parameters + ---------- + image : array + Input image. + + Returns + ------- + threshold : float + Upper threshold value. All pixels intensities that less or equal of + this value assumed as foreground. + + References + ---------- + .. [1] Li C.H. and Lee C.K. (1993) "Minimum Cross Entropy Thresholding" + Pattern Recognition, 26(4): 617-625 + .. [2] Li C.H. and Tam P.K.S. (1998) "An Iterative Algorithm for Minimum + Cross Entropy Thresholding" Pattern Recognition Letters, 18(8): 771-776 + .. [3] Sezgin M. and Sankur B. (2004) "Survey over Image Thresholding + Techniques and Quantitative Performance Evaluation" Journal of + Electronic Imaging, 13(1): 146-165 + http://citeseer.ist.psu.edu/sezgin04survey.html + .. [4] ImageJ AutoThresholder code, http://fiji.sc/wiki/index.php/Auto_Threshold + + Adapted for skimage by J. Metz from ImageJ plugin by G.Landini + + Examples + -------- + >>> from skimage.data import camera + >>> image = camera() + >>> thresh = threshold_li(image) + >>> binary = image <= thresh + """ + # Requires positive image (because of log(mean)) + offset = image.min() + # Can not use fixed tolerance for float image + imrange = image.max() - offset + image -= offset + + tolerance = 0.5 * imrange / 256.0 + # Calculate the mean gray-level + mean = image.mean() + + # Initial estimate + new_thresh = mean + old_thresh = new_thresh + 2 * tolerance + + # Stop the iterations when the difference between the + # new and old threshold values is less than the tolerance + while abs(new_thresh - old_thresh) > tolerance: + old_thresh = new_thresh + threshold = old_thresh + tolerance # range + # Calculate the means of background and object pixels + mean_back = image[image <= threshold].mean() + mean_obj = image[image > threshold].mean() + + temp = (mean_back - mean_obj) / (np.log(mean_back) - np.log(mean_obj)) + + if temp < 0: + new_thresh = temp - tolerance + else: + new_thresh = temp + tolerance + + return threshold + offset