fix percentile autolevel

doc/examples/plot_local_autolevels.py

@@ -0,0 +1,45 @@
+"""
+=====================
+Local Autolevel
+=====================
+
+Local autolevel stretch local histogram between 0 and max_graylevel (e.g. 255 for 8 bit image).
+The following code shows the difference between autolevel and percentile auto_level where [min,max] interval
+is replaced by [p0,p1] percentiles interval
+
+"""
+import matplotlib.pyplot as plt
+
+from skimage import data
+
+from skimage.rank import percentile_autolevel,autolevel
+from skimage.morphology import disk
+
+
+image = data.camera()
+
+selem = disk(20)
+loc_autolevel = autolevel(image,selem=selem)
+loc_perc_autolevel = percentile_autolevel(image,selem=selem,p0=.0,p1=1.0)
+
+assert (loc_autolevel==loc_perc_autolevel).all()
+
+loc_perc_autolevel = percentile_autolevel(image,selem=selem,p0=.01,p1=.99)
+
+fig, axes = plt.subplots(nrows=3, figsize=(7, 8))
+ax0, ax1, ax2 = axes
+plt.gray()
+
+ax0.imshow(image)
+ax0.set_title('Image')
+
+ax1.imshow(loc_autolevel)
+ax1.set_title('Autolevel')
+
+ax2.imshow(loc_perc_autolevel,vmin=0,vmax=255)
+ax2.set_title('percentile autolevel')
+
+for ax in axes:
+    ax.axis('off')
+
+plt.show()

skimage/rank/_core16.pxd

@@ -14,10 +14,6 @@ import numpy as np
 cimport numpy as np
 from libc.stdlib cimport malloc, free
 
-# generic cdef functions
-cdef inline int int_max(int a, int b): return a if a >= b else b
-cdef inline int int_min(int a, int b): return a if a <= b else b
-
 #---------------------------------------------------------------------------
 # 16 bit core kernel receives extra information about data bitdepth
 #---------------------------------------------------------------------------

skimage/rank/_core16b.pxd

@@ -14,10 +14,6 @@ import numpy as np
 cimport numpy as np
 from libc.stdlib cimport malloc, free
 
-# generic cdef functions
-cdef inline int int_max(int a, int b): return a if a >= b else b
-cdef inline int int_min(int a, int b): return a if a <= b else b
-
 #---------------------------------------------------------------------------
 # 16 bit core kernel receives extra information about data bitdepth and bilateral interval
 #---------------------------------------------------------------------------

skimage/rank/_core8.pxd

@@ -14,10 +14,6 @@ import numpy as np
 cimport numpy as np
 from libc.stdlib cimport malloc, free
 
-# generic cdef functions
-cdef inline int int_max(int a, int b): return a if a >= b else b
-cdef inline int int_min(int a, int b): return a if a <= b else b
-
 #---------------------------------------------------------------------------
 # 8 bit core kernel
 #---------------------------------------------------------------------------

skimage/rank/_core8p.pxd

@@ -15,8 +15,8 @@ cimport numpy as np
 from libc.stdlib cimport malloc, free
 
 # generic cdef functions
-cdef inline int int_max(int a, int b): return a if a >= b else b
-cdef inline int int_min(int a, int b): return a if a <= b else b
+cdef inline np.uint8_t uint8_max(np.uint8_t a, np.uint8_t b): return a if a >= b else b
+cdef inline np.uint8_t uint8_min(np.uint8_t a, np.uint8_t b): return a if a <= b else b
 
 #---------------------------------------------------------------------------
 # 8 bit core kernel receives extra information about data inferior and superior percentiles

skimage/rank/_crank8.pyx

@@ -33,11 +33,13 @@ cdef inline np.uint8_t kernel_autolevel(int* histo, float pop, np.uint8_t g):
             if histo[i]:
                 imin = i
                 break
-    delta = imax-imin
-    if delta>0:
-        return <np.uint8_t>(255.*(g-imin)/delta)
+        delta = imax-imin
+        if delta>0:
+            return <np.uint8_t>(255.*(g-imin)/delta)
+        else:
+            return <np.uint8_t>(imax-imin)
     else:
-        return <np.uint8_t>(imax-imin)
+        return <np.uint8_t>(0)
 
 cdef inline np.uint8_t kernel_bottomhat(int* histo, float pop, np.uint8_t g):
     cdef int i

skimage/rank/_crank8_percentiles.pyx

@@ -15,7 +15,7 @@ import numpy as np
 cimport numpy as np
 
 # import main loop
-from _core8p cimport _core8p
+from _core8p cimport _core8p,uint8_max,uint8_min
 
 # -----------------------------------------------------------------
 # kernels uint8 (SOFT version using percentiles)
@@ -27,25 +27,29 @@ cdef inline np.uint8_t kernel_autolevel(int* histo, float pop, np.uint8_t g, flo
     if pop:
         sum = 0
         p1 = 1.0-p1
+        imin = 0
+        imax = 255
+
         for i in range(256):
             sum += histo[i]
-            if sum>=p0*pop:
+            if sum>(p0*pop):
                 imin = i
                 break
         sum = 0
         for i in range(255,-1,-1):
             sum += histo[i]
-            if sum>=p1*pop:
+            if sum>(p1*pop):
                 imax = i
                 break
-
         delta = imax-imin
         if delta>0:
-            return <np.uint8_t>(255.*(g-imin)/delta)
+#            return <np.uint8_t>(255.)
+#            return <np.uint8_t>(delta)
+            return <np.uint8_t>(255*(uint8_min(uint8_max(imin,g),imax)-imin)/delta)
         else:
-            return <np.uint8_t>(0)
+            return <np.uint8_t>(imax-imin)
     else:
-        return <np.uint8_t>(0)
+        return <np.uint8_t>(128)
 
 
 cdef inline np.uint8_t kernel_gradient(int* histo, float pop, np.uint8_t g, float p0, float p1):

skimage/rank/tests/test_suite.py

@@ -4,7 +4,10 @@ import numpy as np
 
 from skimage.rank import _crank8,_crank8_percentiles
 from skimage.rank import _crank16,_crank16_bilateral,_crank16_percentiles
-from skimage.morphology import cmorph
+from skimage.morphology import cmorph,disk
+from skimage import data
+from skimage import rank
+
 
 class TestSequenceFunctions(unittest.TestCase):
 
@@ -92,6 +95,15 @@ class TestSequenceFunctions(unittest.TestCase):
         elem = np.ones((3,3),dtype='uint8')
         f = _crank16_percentiles.mean(image=a16,selem = elem,shift_x=0,shift_y=0,p0=.1,p1=.9,bitdepth=4)
 
+    def test_compare_autolevels(self):
+        image = data.camera()
+
+        selem = disk(20)
+        loc_autolevel = rank.autolevel(image,selem=selem)
+        loc_perc_autolevel = rank.percentile_autolevel(image,selem=selem,p0=.0,p1=1.)
+
+        assert (loc_autolevel==loc_perc_autolevel).all()
+
 if __name__ == '__main__':
 
     suite = unittest.TestLoader().loadTestsFromTestCase(TestSequenceFunctions)