From 196205956e131f073e158facd6101096c0f06ab8 Mon Sep 17 00:00:00 2001 From: blink1073 Date: Sat, 8 Feb 2014 20:05:47 -0600 Subject: [PATCH 01/13] Set the unused portion of the image to 0 --- skimage/exposure/_adapthist.py | 13 +++++++++---- 1 file changed, 9 insertions(+), 4 deletions(-) diff --git a/skimage/exposure/_adapthist.py b/skimage/exposure/_adapthist.py index 699e9acf..aa5552f5 100644 --- a/skimage/exposure/_adapthist.py +++ b/skimage/exposure/_adapthist.py @@ -51,8 +51,7 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, Notes ----- * The algorithm relies on an image whose rows and columns are even - multiples of the number of tiles, so the extra rows and columns are left - at their original values, thus preserving the input image shape. + multiples of the number of tiles, so the extra rows and columns are to zero, thus preserving the input image shape. * For color images, the following steps are performed: - The image is converted to LAB color space - The CLAHE algorithm is run on the L channel @@ -73,14 +72,20 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, args[0] = rescale_intensity(l_chan, out_range=(0, NR_OF_GREY - 1)) new_l = _clahe(*args).astype(float) new_l = rescale_intensity(new_l, out_range=(0, 100)) - lab_img[:new_l.shape[0], :new_l.shape[1], 0] = new_l + col, row = new_l.shape + lab_img[:col, :row, 0] = new_l + lab_img[col:, :, 0] = 0 + lab_img[:, row:, 0] = 0 image = color.lab2rgb(lab_img) image = rescale_intensity(image, out_range=(0, 1)) else: image = skimage.img_as_uint(image) args[0] = rescale_intensity(image, out_range=(0, NR_OF_GREY - 1)) out = _clahe(*args) - image[:out.shape[0], :out.shape[1]] = out + col, row = out.shape + image[:col, :row] = out + image[col:, :] = 0 + image[:, row:] = 0 image = rescale_intensity(image) return image From 61fb831d31403d6725b0d786e5fb96a3c8302333 Mon Sep 17 00:00:00 2001 From: blink1073 Date: Sat, 8 Feb 2014 20:15:19 -0600 Subject: [PATCH 02/13] Implement a mode enum --- skimage/exposure/_adapthist.py | 28 +++++++++++++++++++--------- 1 file changed, 19 insertions(+), 9 deletions(-) diff --git a/skimage/exposure/_adapthist.py b/skimage/exposure/_adapthist.py index aa5552f5..241b8c8e 100644 --- a/skimage/exposure/_adapthist.py +++ b/skimage/exposure/_adapthist.py @@ -26,7 +26,7 @@ NR_OF_GREY = 16384 # number of grayscale levels to use in CLAHE algorithm def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, - nbins=256): + nbins=256, mode='unchanged'): """Contrast Limited Adaptive Histogram Equalization. Parameters @@ -42,6 +42,8 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, contrast). nbins : int, optional Number of gray bins for histogram ("dynamic range"). + mode : string, one of {'unchanged', 'zero', 'trim'}, optional + How to treat any pixels falling outside of the tiles. Returns ------- @@ -72,20 +74,28 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, args[0] = rescale_intensity(l_chan, out_range=(0, NR_OF_GREY - 1)) new_l = _clahe(*args).astype(float) new_l = rescale_intensity(new_l, out_range=(0, 100)) - col, row = new_l.shape - lab_img[:col, :row, 0] = new_l - lab_img[col:, :, 0] = 0 - lab_img[:, row:, 0] = 0 + if mode == 'trim': + lab_img = new_l + else: + col, row = new_l.shape + lab_img[:col, :row, 0] = new_l + if mode == 'zero': + lab_img[col:, :, 0] = 0 + lab_img[:, row:, 0] = 0 image = color.lab2rgb(lab_img) image = rescale_intensity(image, out_range=(0, 1)) else: image = skimage.img_as_uint(image) args[0] = rescale_intensity(image, out_range=(0, NR_OF_GREY - 1)) out = _clahe(*args) - col, row = out.shape - image[:col, :row] = out - image[col:, :] = 0 - image[:, row:] = 0 + if mode == 'trim': + image = out + else: + col, row = out.shape + image[:col, :row] = out + if mode == 'zero': + image[col:, :] = 0 + image[:, row:] = 0 image = rescale_intensity(image) return image From 9fa39e564a25aa706f4968f5b02871bb57db3240 Mon Sep 17 00:00:00 2001 From: blink1073 Date: Sat, 8 Feb 2014 20:15:36 -0600 Subject: [PATCH 03/13] Update docstring --- skimage/exposure/_adapthist.py | 9 +++++++-- 1 file changed, 7 insertions(+), 2 deletions(-) diff --git a/skimage/exposure/_adapthist.py b/skimage/exposure/_adapthist.py index 241b8c8e..552591f4 100644 --- a/skimage/exposure/_adapthist.py +++ b/skimage/exposure/_adapthist.py @@ -43,7 +43,7 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, nbins : int, optional Number of gray bins for histogram ("dynamic range"). mode : string, one of {'unchanged', 'zero', 'trim'}, optional - How to treat any pixels falling outside of the tiles. + How to treat any pixels falling outside of the tiles. See the notes. Returns ------- @@ -53,7 +53,12 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, Notes ----- * The algorithm relies on an image whose rows and columns are even - multiples of the number of tiles, so the extra rows and columns are to zero, thus preserving the input image shape. + multiples of the number of tiles, so the extra rows and columns are not + affected by the algorithm. The handling of those outlier pixels is + determined by the "mode" paramter. If mode is 'unchanged', the values + the same as the input values. If mode is 'zero', they are set to zero. + If mode is 'trim', only the portion of the image that was equalized will + be returned. * For color images, the following steps are performed: - The image is converted to LAB color space - The CLAHE algorithm is run on the L channel From 9f7119b1c0da17af7f7e9e4f671d9fcdf66e2ea0 Mon Sep 17 00:00:00 2001 From: blink1073 Date: Sun, 16 Feb 2014 19:08:14 -0600 Subject: [PATCH 04/13] Change 'trim' to 'crop' --- skimage/exposure/_adapthist.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/skimage/exposure/_adapthist.py b/skimage/exposure/_adapthist.py index 552591f4..f5a6d0a1 100644 --- a/skimage/exposure/_adapthist.py +++ b/skimage/exposure/_adapthist.py @@ -26,7 +26,7 @@ NR_OF_GREY = 16384 # number of grayscale levels to use in CLAHE algorithm def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, - nbins=256, mode='unchanged'): + nbins=256, mode='ignore'): """Contrast Limited Adaptive Histogram Equalization. Parameters @@ -42,7 +42,7 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, contrast). nbins : int, optional Number of gray bins for histogram ("dynamic range"). - mode : string, one of {'unchanged', 'zero', 'trim'}, optional + mode : string, one of {'unchanged', 'zero', 'crop'}, optional How to treat any pixels falling outside of the tiles. See the notes. Returns @@ -79,7 +79,7 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, args[0] = rescale_intensity(l_chan, out_range=(0, NR_OF_GREY - 1)) new_l = _clahe(*args).astype(float) new_l = rescale_intensity(new_l, out_range=(0, 100)) - if mode == 'trim': + if mode == 'crop': lab_img = new_l else: col, row = new_l.shape @@ -93,7 +93,7 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, image = skimage.img_as_uint(image) args[0] = rescale_intensity(image, out_range=(0, NR_OF_GREY - 1)) out = _clahe(*args) - if mode == 'trim': + if mode == 'crop': image = out else: col, row = out.shape From 95a208eec5fe5d7b7978c5a5cbfb62b0a5155408 Mon Sep 17 00:00:00 2001 From: blink1073 Date: Sun, 16 Feb 2014 21:10:51 -0600 Subject: [PATCH 05/13] Add test for adapthist modes --- skimage/exposure/tests/test_exposure.py | 35 +++++++++++++++++++++---- 1 file changed, 30 insertions(+), 5 deletions(-) diff --git a/skimage/exposure/tests/test_exposure.py b/skimage/exposure/tests/test_exposure.py index dcaf6837..fb3ba78a 100644 --- a/skimage/exposure/tests/test_exposure.py +++ b/skimage/exposure/tests/test_exposure.py @@ -194,6 +194,31 @@ def test_adapthist_color(): return data, adapted +def test_adapthist_modes(): + '''Test adaptist `mode` parameter values. + ''' + img = skimage.img_as_float(data.moon()) + full_scale = skimage.exposure.rescale_intensity(skimage.img_as_uint(img)) + ignore = exposure.equalize_adapthist(img.copy(), ntiles_x=9, ntiles_y=9) + zero = exposure.equalize_adapthist(img.copy(), ntiles_x=9, ntiles_y=9, + mode='zero') + crop = exposure.equalize_adapthist(img.copy(), ntiles_x=9, ntiles_y=9, + mode='crop') + assert ignore.shape == zero.shape + assert ignore.shape != crop.shape + assert_array_equal(zero[crop.shape[0]:, :], 0) + assert_array_equal(zero[:, crop.shape[1]:], 0) + + assert_almost_equal = np.testing.assert_almost_equal + full_cropped = full_scale[:crop.shape[0], :crop.shape[1]] + assert_almost_equal(peak_snr(full_scale, ignore), 56.175, 3) + assert_almost_equal(peak_snr(full_scale, zero), 94.7988, 3) + assert_almost_equal(peak_snr(full_cropped, crop), 120.4598, 3) + assert_almost_equal(norm_brightness_err(full_scale, ignore), 0.223, 3) + assert_almost_equal(norm_brightness_err(full_scale, zero), 0.01589, 3) + assert_almost_equal(norm_brightness_err(full_cropped, crop), 0.02878, 3) + + def peak_snr(img1, img2): '''Peak signal to noise ratio of two images @@ -236,11 +261,6 @@ def norm_brightness_err(img1, img2): return nbe -if __name__ == '__main__': - from numpy import testing - testing.run_module_suite() - - # Test Gamma Correction # ===================== @@ -409,3 +429,8 @@ def test_adjust_inv_sigmoid_cutoff_half(): def test_negative(): image = np.arange(-10, 245, 4).reshape(8, 8).astype(np.double) assert_raises(ValueError, exposure.adjust_gamma, image) + + +if __name__ == '__main__': + from numpy import testing + testing.run_module_suite() From 909f23baab96551af9614f6f21c30370f9bd60ba Mon Sep 17 00:00:00 2001 From: blink1073 Date: Mon, 17 Feb 2014 16:11:43 -0600 Subject: [PATCH 06/13] Switch to HSV, refactor equalize_adapthist, fix #757, cleanup docstring --- skimage/exposure/_adapthist.py | 77 ++++++++++++++++++---------------- 1 file changed, 40 insertions(+), 37 deletions(-) diff --git a/skimage/exposure/_adapthist.py b/skimage/exposure/_adapthist.py index f5a6d0a1..0887224f 100644 --- a/skimage/exposure/_adapthist.py +++ b/skimage/exposure/_adapthist.py @@ -42,7 +42,7 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, contrast). nbins : int, optional Number of gray bins for histogram ("dynamic range"). - mode : string, one of {'unchanged', 'zero', 'crop'}, optional + mode : {'unchanged', 'zero', 'crop'}, optional How to treat any pixels falling outside of the tiles. See the notes. Returns @@ -55,13 +55,13 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, * The algorithm relies on an image whose rows and columns are even multiples of the number of tiles, so the extra rows and columns are not affected by the algorithm. The handling of those outlier pixels is - determined by the "mode" paramter. If mode is 'unchanged', the values + determined by the `mode` parameter. If mode is 'unchanged', the values the same as the input values. If mode is 'zero', they are set to zero. If mode is 'trim', only the portion of the image that was equalized will be returned. * For color images, the following steps are performed: - - The image is converted to LAB color space - - The CLAHE algorithm is run on the L channel + - The image is converted to HSV color space + - The CLAHE algorithm is run on the V (Value) channel - The image is converted back to RGB space and returned * For RGBA images, the original alpha channel is removed. @@ -70,38 +70,36 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, .. [1] http://tog.acm.org/resources/GraphicsGems/gems.html#gemsvi .. [2] https://en.wikipedia.org/wiki/CLAHE#CLAHE """ - args = [None, ntiles_x, ntiles_y, clip_limit * nbins, nbins] - if image.ndim > 2: - lab_img = color.rgb2lab(skimage.img_as_float(image)) - l_chan = lab_img[:, :, 0] - l_chan /= np.max(np.abs(l_chan)) - l_chan = skimage.img_as_uint(l_chan) - args[0] = rescale_intensity(l_chan, out_range=(0, NR_OF_GREY - 1)) - new_l = _clahe(*args).astype(float) - new_l = rescale_intensity(new_l, out_range=(0, 100)) - if mode == 'crop': - lab_img = new_l - else: - col, row = new_l.shape - lab_img[:col, :row, 0] = new_l - if mode == 'zero': - lab_img[col:, :, 0] = 0 - lab_img[:, row:, 0] = 0 - image = color.lab2rgb(lab_img) - image = rescale_intensity(image, out_range=(0, 1)) - else: - image = skimage.img_as_uint(image) - args[0] = rescale_intensity(image, out_range=(0, NR_OF_GREY - 1)) - out = _clahe(*args) - if mode == 'crop': - image = out - else: - col, row = out.shape - image[:col, :row] = out - if mode == 'zero': - image[col:, :] = 0 - image[:, row:] = 0 + ndim = image.ndim + if ndim == 3: + if image.shape[2] == 4: + image = image[:, :, :3] + image = skimage.img_as_float(image) image = rescale_intensity(image) + hsv_img = color.rgb2hsv(image) + image = hsv_img[:, :, 2].copy() + image = skimage.img_as_uint(image) + image = rescale_intensity(image, out_range=(0, NR_OF_GREY - 1)) + out = _clahe(image, ntiles_x, ntiles_y, clip_limit * nbins, nbins) + if mode == 'crop': + image = image[:out.shape[0], :out.shape[1]] + if ndim == 3: + hsv_img = hsv_img[:out.shape[0], :out.shape[1], :] + image[:out.shape[0], :out.shape[1]] = out + image = skimage.img_as_float(image) + if ndim == 3: + hsv_img[:, :, 2] = rescale_intensity(image) + image = color.hsv2rgb(hsv_img) + else: + image = rescale_intensity(image) + if mode == 'zero': + mask = np.zeros(image.shape) + if ndim == 3: + for i in range(3): + mask[:out.shape[0], :out.shape[1], i] = 1 + else: + mask[:out.shape[0], :out.shape[1]] = 1 + image[mask == 0] = 0 return image @@ -134,8 +132,8 @@ def _clahe(image, ntiles_x, ntiles_y, clip_limit, nbins=128): """ ntiles_x = min(ntiles_x, MAX_REG_X) ntiles_y = min(ntiles_y, MAX_REG_Y) - ntiles_y = max(ntiles_x, 2) - ntiles_x = max(ntiles_y, 2) + ntiles_y = max(ntiles_y, 2) + ntiles_x = max(ntiles_x, 2) if clip_limit == 1.0: return image # is OK, immediately returns original image. @@ -144,6 +142,11 @@ def _clahe(image, ntiles_x, ntiles_y, clip_limit, nbins=128): y_res = image.shape[0] - image.shape[0] % ntiles_y x_res = image.shape[1] - image.shape[1] % ntiles_x + # make the tile size divisible by 2 + while y_res % (2 * ntiles_y): + y_res -= 1 + while x_res % (2 * ntiles_x): + x_res -= 1 image = image[: y_res, : x_res] x_size = image.shape[1] // ntiles_x # Actual size of contextual regions From ba010a3191ad03ec65d97e32f190d699c165255f Mon Sep 17 00:00:00 2001 From: blink1073 Date: Mon, 17 Feb 2014 16:13:52 -0600 Subject: [PATCH 07/13] Update adapthist tests, add test for alpha channel and rgb mode handling. --- skimage/exposure/tests/test_exposure.py | 86 +++++++++++++++++++------ 1 file changed, 67 insertions(+), 19 deletions(-) diff --git a/skimage/exposure/tests/test_exposure.py b/skimage/exposure/tests/test_exposure.py index fb3ba78a..01d10219 100644 --- a/skimage/exposure/tests/test_exposure.py +++ b/skimage/exposure/tests/test_exposure.py @@ -11,6 +11,7 @@ from skimage.exposure.exposure import intensity_range from skimage.color import rgb2gray from skimage.util.dtype import dtype_range +import matplotlib.pyplot as plt # Test histogram equalization # =========================== @@ -147,13 +148,13 @@ def test_adapthist_scalar(): ''' img = skimage.img_as_ubyte(data.moon()) adapted = exposure.equalize_adapthist(img, clip_limit=0.02) - assert adapted.min() == 0 - assert adapted.max() == (1 << 16) - 1 + assert adapted.min() == 0.0 + assert adapted.max() == 1.0 assert img.shape == adapted.shape - full_scale = skimage.exposure.rescale_intensity(skimage.img_as_uint(img)) + full_scale = skimage.exposure.rescale_intensity(skimage.img_as_float(img)) assert_almost_equal = np.testing.assert_almost_equal - assert_almost_equal(peak_snr(full_scale, adapted), 101.231, 3) + assert_almost_equal(peak_snr(full_scale, adapted), 101.2295, 3) assert_almost_equal(norm_brightness_err(full_scale, adapted), 0.041, 3) return img, adapted @@ -169,8 +170,8 @@ def test_adapthist_grayscale(): nbins=128) assert_almost_equal = np.testing.assert_almost_equal assert img.shape == adapted.shape - assert_almost_equal(peak_snr(img, adapted), 97.531, 3) - assert_almost_equal(norm_brightness_err(img, adapted), 0.0313, 3) + assert_almost_equal(peak_snr(img, adapted), 101.4070, 3) + assert_almost_equal(norm_brightness_err(img, adapted), 0.03856, 3) return data, adapted @@ -183,19 +184,36 @@ def test_adapthist_color(): hist, bin_centers = exposure.histogram(img) assert len(w) > 0 adapted = exposure.equalize_adapthist(img, clip_limit=0.01) + assert_almost_equal = np.testing.assert_almost_equal assert adapted.min() == 0 assert adapted.max() == 1.0 assert img.shape == adapted.shape full_scale = skimage.exposure.rescale_intensity(img) - assert_almost_equal(peak_snr(full_scale, adapted), 102.940, 3) + assert_almost_equal(peak_snr(full_scale, adapted), 106.865, 3) assert_almost_equal(norm_brightness_err(full_scale, adapted), - 0.0110, 3) + 0.0509, 3) return data, adapted -def test_adapthist_modes(): - '''Test adaptist `mode` parameter values. +def test_adapthist_alpha(): + '''Test an RGBA color image + ''' + img = skimage.img_as_float(data.lena()) + alpha = np.ones((img.shape[0], img.shape[1]), dtype=float) + img = np.dstack((img, alpha)) + adapted = exposure.equalize_adapthist(img) + assert adapted.shape != img.shape + img = img[:, :, :3] + full_scale = skimage.exposure.rescale_intensity(img) + assert img.shape == adapted.shape + assert_almost_equal = np.testing.assert_almost_equal + assert_almost_equal(peak_snr(full_scale, adapted), 106.862, 3) + assert_almost_equal(norm_brightness_err(full_scale, adapted), 0.0509, 3) + + +def test_adapthist_modes_scalar(): + '''Test adaptist `mode` parameter values for ndim==2 image. ''' img = skimage.img_as_float(data.moon()) full_scale = skimage.exposure.rescale_intensity(skimage.img_as_uint(img)) @@ -211,12 +229,40 @@ def test_adapthist_modes(): assert_almost_equal = np.testing.assert_almost_equal full_cropped = full_scale[:crop.shape[0], :crop.shape[1]] - assert_almost_equal(peak_snr(full_scale, ignore), 56.175, 3) - assert_almost_equal(peak_snr(full_scale, zero), 94.7988, 3) - assert_almost_equal(peak_snr(full_cropped, crop), 120.4598, 3) - assert_almost_equal(norm_brightness_err(full_scale, ignore), 0.223, 3) - assert_almost_equal(norm_brightness_err(full_scale, zero), 0.01589, 3) - assert_almost_equal(norm_brightness_err(full_cropped, crop), 0.02878, 3) + zero_cropped = zero[:crop.shape[0], :crop.shape[1]] + ignore_cropped = ignore[:crop.shape[0], :crop.shape[1]] + + for crop_img in [ignore_cropped, zero_cropped, crop]: + assert_almost_equal(peak_snr(full_cropped, crop_img), 120.456, 3) + assert_almost_equal(norm_brightness_err(full_cropped, crop_img), + 0.4398, 3) + + +def test_adapthist_modes_rgb(): + '''Test adaptist `mode` parameter values for rgb image. + ''' + img = skimage.img_as_float(data.lena()) + full_scale = skimage.exposure.rescale_intensity(skimage.img_as_uint(img)) + ignore = exposure.equalize_adapthist(img.copy(), ntiles_x=9, ntiles_y=10) + zero = exposure.equalize_adapthist(img.copy(), ntiles_x=9, ntiles_y=10, + mode='zero') + crop = exposure.equalize_adapthist(img.copy(), ntiles_x=9, ntiles_y=10, + mode='crop') + assert ignore.shape == zero.shape + assert ignore.shape != crop.shape + + assert_array_equal(zero[crop.shape[0]:, :, :], 0) + assert_array_equal(zero[:, crop.shape[1]:, :], 0) + + assert_almost_equal = np.testing.assert_almost_equal + full_cropped = full_scale[:crop.shape[0], :crop.shape[1], :] + zero_cropped = zero[:crop.shape[0], :crop.shape[1], :] + ignore_cropped = ignore[:crop.shape[0], :crop.shape[1], :] + + for crop_img in [ignore_cropped, zero_cropped, crop]: + assert np.floor(peak_snr(full_cropped, crop_img)) == 106 + assert_almost_equal(norm_brightness_err(full_cropped, crop_img), + 0.0517, 3) def peak_snr(img1, img2): @@ -260,6 +306,8 @@ def norm_brightness_err(img1, img2): nbe = ambe / dtype_range[img1.dtype.type][1] return nbe +test_adapthist_modes_rgb() + # Test Gamma Correction # ===================== @@ -431,6 +479,6 @@ def test_negative(): assert_raises(ValueError, exposure.adjust_gamma, image) -if __name__ == '__main__': - from numpy import testing - testing.run_module_suite() +#if __name__ == '__main__': +# from numpy import testing +# testing.run_module_suite() From 2bdfcc7d5e4384406632df260386d5f1f57afbd9 Mon Sep 17 00:00:00 2001 From: blink1073 Date: Mon, 17 Feb 2014 16:34:02 -0600 Subject: [PATCH 08/13] Fix failing tests and turn module level testing back ony --- skimage/exposure/tests/test_exposure.py | 19 +++++++++---------- 1 file changed, 9 insertions(+), 10 deletions(-) diff --git a/skimage/exposure/tests/test_exposure.py b/skimage/exposure/tests/test_exposure.py index 01d10219..3ee2a492 100644 --- a/skimage/exposure/tests/test_exposure.py +++ b/skimage/exposure/tests/test_exposure.py @@ -170,8 +170,8 @@ def test_adapthist_grayscale(): nbins=128) assert_almost_equal = np.testing.assert_almost_equal assert img.shape == adapted.shape - assert_almost_equal(peak_snr(img, adapted), 101.4070, 3) - assert_almost_equal(norm_brightness_err(img, adapted), 0.03856, 3) + assert_almost_equal(peak_snr(img, adapted), 105.669, 3) + assert_almost_equal(norm_brightness_err(img, adapted), 0.02470, 3) return data, adapted @@ -190,9 +190,9 @@ def test_adapthist_color(): assert adapted.max() == 1.0 assert img.shape == adapted.shape full_scale = skimage.exposure.rescale_intensity(img) - assert_almost_equal(peak_snr(full_scale, adapted), 106.865, 3) + assert_almost_equal(peak_snr(full_scale, adapted), 105.50517, 3) assert_almost_equal(norm_brightness_err(full_scale, adapted), - 0.0509, 3) + 0.0544, 3) return data, adapted @@ -208,8 +208,8 @@ def test_adapthist_alpha(): full_scale = skimage.exposure.rescale_intensity(img) assert img.shape == adapted.shape assert_almost_equal = np.testing.assert_almost_equal - assert_almost_equal(peak_snr(full_scale, adapted), 106.862, 3) - assert_almost_equal(norm_brightness_err(full_scale, adapted), 0.0509, 3) + assert_almost_equal(peak_snr(full_scale, adapted), 105.50198, 3) + assert_almost_equal(norm_brightness_err(full_scale, adapted), 0.0544, 3) def test_adapthist_modes_scalar(): @@ -306,7 +306,6 @@ def norm_brightness_err(img1, img2): nbe = ambe / dtype_range[img1.dtype.type][1] return nbe -test_adapthist_modes_rgb() # Test Gamma Correction @@ -479,6 +478,6 @@ def test_negative(): assert_raises(ValueError, exposure.adjust_gamma, image) -#if __name__ == '__main__': -# from numpy import testing -# testing.run_module_suite() +if __name__ == '__main__': + from numpy import testing + testing.run_module_suite() From cc68b0c392defab55f6a6252bf9bc290eff73554 Mon Sep 17 00:00:00 2001 From: blink1073 Date: Fri, 11 Jul 2014 13:28:00 -0500 Subject: [PATCH 09/13] Start removing optional edge pixel handling --- skimage/exposure/_adapthist.py | 24 +++--------------------- 1 file changed, 3 insertions(+), 21 deletions(-) diff --git a/skimage/exposure/_adapthist.py b/skimage/exposure/_adapthist.py index 0887224f..d583424d 100644 --- a/skimage/exposure/_adapthist.py +++ b/skimage/exposure/_adapthist.py @@ -42,8 +42,6 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, contrast). nbins : int, optional Number of gray bins for histogram ("dynamic range"). - mode : {'unchanged', 'zero', 'crop'}, optional - How to treat any pixels falling outside of the tiles. See the notes. Returns ------- @@ -52,18 +50,14 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, Notes ----- - * The algorithm relies on an image whose rows and columns are even - multiples of the number of tiles, so the extra rows and columns are not - affected by the algorithm. The handling of those outlier pixels is - determined by the `mode` parameter. If mode is 'unchanged', the values - the same as the input values. If mode is 'zero', they are set to zero. - If mode is 'trim', only the portion of the image that was equalized will - be returned. * For color images, the following steps are performed: - The image is converted to HSV color space - The CLAHE algorithm is run on the V (Value) channel - The image is converted back to RGB space and returned * For RGBA images, the original alpha channel is removed. + * The CLAHE algorithm relies on image blocks of equal size. This results + in extra border pixels that are not handled. Extra blocks are created + around the border to handle these pixels. References ---------- @@ -81,10 +75,6 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, image = skimage.img_as_uint(image) image = rescale_intensity(image, out_range=(0, NR_OF_GREY - 1)) out = _clahe(image, ntiles_x, ntiles_y, clip_limit * nbins, nbins) - if mode == 'crop': - image = image[:out.shape[0], :out.shape[1]] - if ndim == 3: - hsv_img = hsv_img[:out.shape[0], :out.shape[1], :] image[:out.shape[0], :out.shape[1]] = out image = skimage.img_as_float(image) if ndim == 3: @@ -92,14 +82,6 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, image = color.hsv2rgb(hsv_img) else: image = rescale_intensity(image) - if mode == 'zero': - mask = np.zeros(image.shape) - if ndim == 3: - for i in range(3): - mask[:out.shape[0], :out.shape[1], i] = 1 - else: - mask[:out.shape[0], :out.shape[1]] = 1 - image[mask == 0] = 0 return image From 5cafbd95e681bbb555ec4d6b1d7097af4eb3400f Mon Sep 17 00:00:00 2001 From: Tony S Yu Date: Fri, 11 Jul 2014 16:30:18 -0500 Subject: [PATCH 10/13] Add tests for intensity_range and rename parameter --- skimage/exposure/exposure.py | 1 + 1 file changed, 1 insertion(+) diff --git a/skimage/exposure/exposure.py b/skimage/exposure/exposure.py index 3a2accbc..1b42ecf8 100644 --- a/skimage/exposure/exposure.py +++ b/skimage/exposure/exposure.py @@ -3,6 +3,7 @@ import numpy as np from skimage import img_as_float from skimage.util.dtype import dtype_range, dtype_limits +from skimage._shared.utils import deprecation_warning __all__ = ['histogram', 'cumulative_distribution', 'equalize', From a924e55b5276dfd37d30bb505471f8fa468dc725 Mon Sep 17 00:00:00 2001 From: blink1073 Date: Sat, 12 Jul 2014 09:40:30 -0500 Subject: [PATCH 11/13] Fix handling of border pixels and update tests --- skimage/exposure/_adapthist.py | 42 ++++++++++++------ skimage/exposure/tests/test_exposure.py | 59 ++----------------------- 2 files changed, 33 insertions(+), 68 deletions(-) diff --git a/skimage/exposure/_adapthist.py b/skimage/exposure/_adapthist.py index d583424d..c89ba35b 100644 --- a/skimage/exposure/_adapthist.py +++ b/skimage/exposure/_adapthist.py @@ -55,9 +55,10 @@ def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, - The CLAHE algorithm is run on the V (Value) channel - The image is converted back to RGB space and returned * For RGBA images, the original alpha channel is removed. - * The CLAHE algorithm relies on image blocks of equal size. This results - in extra border pixels that are not handled. Extra blocks are created - around the border to handle these pixels. + * The CLAHE algorithm relies on image blocks of equal size. This may + result in extra border pixels that would not be handled. In that case, + we pad the image with a repeat of the border pixels, apply the + algorithm, and then trim the image to original size. References ---------- @@ -120,19 +121,36 @@ def _clahe(image, ntiles_x, ntiles_y, clip_limit, nbins=128): if clip_limit == 1.0: return image # is OK, immediately returns original image. - map_array = np.zeros((ntiles_y, ntiles_x, nbins), dtype=int) - y_res = image.shape[0] - image.shape[0] % ntiles_y x_res = image.shape[1] - image.shape[1] % ntiles_x + # make the tile size divisible by 2 while y_res % (2 * ntiles_y): y_res -= 1 while x_res % (2 * ntiles_x): x_res -= 1 - image = image[: y_res, : x_res] - x_size = image.shape[1] // ntiles_x # Actual size of contextual regions - y_size = image.shape[0] // ntiles_y + orig_shape = image.shape + x_size = x_res // ntiles_x # Actual size of contextual regions + y_size = y_res // ntiles_y + + if y_res != image.shape[0]: + ntiles_y += 1 + if x_res != image.shape[1]: + ntiles_x += 1 + if y_res != image.shape[1] or x_res != image.shape[0]: + hgt = y_size * ntiles_y - image.shape[0] + wid = x_size * ntiles_x - image.shape[1] + image = np.vstack((image, image[-hgt:, :])) + image = np.hstack((image, image[:, -wid:])) + y_res, x_res = image.shape + + bin_size = 1 + NR_OF_GREY / nbins + aLUT = np.arange(NR_OF_GREY) + aLUT //= bin_size + img_blocks = view_as_blocks(image, (y_size, x_size)) + + map_array = np.zeros((ntiles_y, ntiles_x, nbins), dtype=int) n_pixels = x_size * y_size if clip_limit > 0.0: # Calculate actual cliplimit @@ -142,11 +160,6 @@ def _clahe(image, ntiles_x, ntiles_y, clip_limit, nbins=128): else: clip_limit = NR_OF_GREY # Large value, do not clip (AHE) - bin_size = 1 + NR_OF_GREY / nbins - aLUT = np.arange(NR_OF_GREY) - aLUT //= bin_size - img_blocks = view_as_blocks(image, (y_size, x_size)) - # Calculate greylevel mappings for each contextual region for y in range(ntiles_y): for x in range(ntiles_x): @@ -203,6 +216,9 @@ def _clahe(image, ntiles_x, ntiles_y, clip_limit, nbins=128): ystart += ystep + if image.shape != orig_shape: + image = image[:orig_shape[0], :orig_shape[1]] + return image diff --git a/skimage/exposure/tests/test_exposure.py b/skimage/exposure/tests/test_exposure.py index 3ee2a492..bb584f1a 100644 --- a/skimage/exposure/tests/test_exposure.py +++ b/skimage/exposure/tests/test_exposure.py @@ -16,6 +16,8 @@ import matplotlib.pyplot as plt # Test histogram equalization # =========================== +np.random.seed(0) + # squeeze image intensities to lower image contrast test_img = skimage.img_as_float(data.camera()) test_img = exposure.rescale_intensity(test_img / 5. + 100) @@ -170,8 +172,8 @@ def test_adapthist_grayscale(): nbins=128) assert_almost_equal = np.testing.assert_almost_equal assert img.shape == adapted.shape - assert_almost_equal(peak_snr(img, adapted), 105.669, 3) - assert_almost_equal(norm_brightness_err(img, adapted), 0.02470, 3) + assert_almost_equal(peak_snr(img, adapted), 104.3168, 3) + assert_almost_equal(norm_brightness_err(img, adapted), 0.0265, 3) return data, adapted @@ -212,59 +214,6 @@ def test_adapthist_alpha(): assert_almost_equal(norm_brightness_err(full_scale, adapted), 0.0544, 3) -def test_adapthist_modes_scalar(): - '''Test adaptist `mode` parameter values for ndim==2 image. - ''' - img = skimage.img_as_float(data.moon()) - full_scale = skimage.exposure.rescale_intensity(skimage.img_as_uint(img)) - ignore = exposure.equalize_adapthist(img.copy(), ntiles_x=9, ntiles_y=9) - zero = exposure.equalize_adapthist(img.copy(), ntiles_x=9, ntiles_y=9, - mode='zero') - crop = exposure.equalize_adapthist(img.copy(), ntiles_x=9, ntiles_y=9, - mode='crop') - assert ignore.shape == zero.shape - assert ignore.shape != crop.shape - assert_array_equal(zero[crop.shape[0]:, :], 0) - assert_array_equal(zero[:, crop.shape[1]:], 0) - - assert_almost_equal = np.testing.assert_almost_equal - full_cropped = full_scale[:crop.shape[0], :crop.shape[1]] - zero_cropped = zero[:crop.shape[0], :crop.shape[1]] - ignore_cropped = ignore[:crop.shape[0], :crop.shape[1]] - - for crop_img in [ignore_cropped, zero_cropped, crop]: - assert_almost_equal(peak_snr(full_cropped, crop_img), 120.456, 3) - assert_almost_equal(norm_brightness_err(full_cropped, crop_img), - 0.4398, 3) - - -def test_adapthist_modes_rgb(): - '''Test adaptist `mode` parameter values for rgb image. - ''' - img = skimage.img_as_float(data.lena()) - full_scale = skimage.exposure.rescale_intensity(skimage.img_as_uint(img)) - ignore = exposure.equalize_adapthist(img.copy(), ntiles_x=9, ntiles_y=10) - zero = exposure.equalize_adapthist(img.copy(), ntiles_x=9, ntiles_y=10, - mode='zero') - crop = exposure.equalize_adapthist(img.copy(), ntiles_x=9, ntiles_y=10, - mode='crop') - assert ignore.shape == zero.shape - assert ignore.shape != crop.shape - - assert_array_equal(zero[crop.shape[0]:, :, :], 0) - assert_array_equal(zero[:, crop.shape[1]:, :], 0) - - assert_almost_equal = np.testing.assert_almost_equal - full_cropped = full_scale[:crop.shape[0], :crop.shape[1], :] - zero_cropped = zero[:crop.shape[0], :crop.shape[1], :] - ignore_cropped = ignore[:crop.shape[0], :crop.shape[1], :] - - for crop_img in [ignore_cropped, zero_cropped, crop]: - assert np.floor(peak_snr(full_cropped, crop_img)) == 106 - assert_almost_equal(norm_brightness_err(full_cropped, crop_img), - 0.0517, 3) - - def peak_snr(img1, img2): '''Peak signal to noise ratio of two images From 0392299fae08d05f2d77b10d6d12db748477bde5 Mon Sep 17 00:00:00 2001 From: blink1073 Date: Sat, 12 Jul 2014 18:54:50 -0500 Subject: [PATCH 12/13] Refactor to use np.pad, naming changes, update test. --- skimage/exposure/_adapthist.py | 65 ++++++++++++------------- skimage/exposure/tests/test_exposure.py | 2 +- 2 files changed, 33 insertions(+), 34 deletions(-) diff --git a/skimage/exposure/_adapthist.py b/skimage/exposure/_adapthist.py index c89ba35b..22b36542 100644 --- a/skimage/exposure/_adapthist.py +++ b/skimage/exposure/_adapthist.py @@ -26,7 +26,7 @@ NR_OF_GREY = 16384 # number of grayscale levels to use in CLAHE algorithm def equalize_adapthist(image, ntiles_x=8, ntiles_y=8, clip_limit=0.01, - nbins=256, mode='ignore'): + nbins=256): """Contrast Limited Adaptive Histogram Equalization. Parameters @@ -121,40 +121,39 @@ def _clahe(image, ntiles_x, ntiles_y, clip_limit, nbins=128): if clip_limit == 1.0: return image # is OK, immediately returns original image. - y_res = image.shape[0] - image.shape[0] % ntiles_y - x_res = image.shape[1] - image.shape[1] % ntiles_x + h_inner = image.shape[0] - image.shape[0] % ntiles_y + w_inner = image.shape[1] - image.shape[1] % ntiles_x # make the tile size divisible by 2 - while y_res % (2 * ntiles_y): - y_res -= 1 - while x_res % (2 * ntiles_x): - x_res -= 1 + while h_inner % (2 * ntiles_y): + h_inner -= 1 + while w_inner % (2 * ntiles_x): + w_inner -= 1 orig_shape = image.shape - x_size = x_res // ntiles_x # Actual size of contextual regions - y_size = y_res // ntiles_y + width = w_inner // ntiles_x # Actual size of contextual regions + height = h_inner // ntiles_y - if y_res != image.shape[0]: + if h_inner != image.shape[0]: ntiles_y += 1 - if x_res != image.shape[1]: + if w_inner != image.shape[1]: ntiles_x += 1 - if y_res != image.shape[1] or x_res != image.shape[0]: - hgt = y_size * ntiles_y - image.shape[0] - wid = x_size * ntiles_x - image.shape[1] - image = np.vstack((image, image[-hgt:, :])) - image = np.hstack((image, image[:, -wid:])) - y_res, x_res = image.shape + if h_inner != image.shape[1] or w_inner != image.shape[0]: + h_pad = height * ntiles_y - image.shape[0] + w_pad = width * ntiles_x - image.shape[1] + image = np.pad(image, ((0, h_pad), (0, w_pad)), 'reflect') + h_inner, w_inner = image.shape bin_size = 1 + NR_OF_GREY / nbins - aLUT = np.arange(NR_OF_GREY) - aLUT //= bin_size - img_blocks = view_as_blocks(image, (y_size, x_size)) + lut = np.arange(NR_OF_GREY) + lut //= bin_size + img_blocks = view_as_blocks(image, (height, width)) map_array = np.zeros((ntiles_y, ntiles_x, nbins), dtype=int) - n_pixels = x_size * y_size + n_pixels = width * height if clip_limit > 0.0: # Calculate actual cliplimit - clip_limit = int(clip_limit * (x_size * y_size) / nbins) + clip_limit = int(clip_limit * (width * height) / nbins) if clip_limit < 1: clip_limit = 1 else: @@ -164,7 +163,7 @@ def _clahe(image, ntiles_x, ntiles_y, clip_limit, nbins=128): for y in range(ntiles_y): for x in range(ntiles_x): sub_img = img_blocks[y, x] - hist = aLUT[sub_img.ravel()] + hist = lut[sub_img.ravel()] hist = np.bincount(hist) hist = np.append(hist, np.zeros(nbins - hist.size, dtype=int)) hist = clip_histogram(hist, clip_limit) @@ -176,29 +175,29 @@ def _clahe(image, ntiles_x, ntiles_y, clip_limit, nbins=128): for y in range(ntiles_y + 1): xstart = 0 if y == 0: # special case: top row - ystep = y_size / 2.0 + ystep = height / 2.0 yU = 0 yB = 0 elif y == ntiles_y: # special case: bottom row - ystep = y_size / 2.0 + ystep = height / 2.0 yU = ntiles_y - 1 yB = yU else: # default values - ystep = y_size + ystep = height yU = y - 1 yB = yB + 1 for x in range(ntiles_x + 1): if x == 0: # special case: left column - xstep = x_size / 2.0 + xstep = width / 2.0 xL = 0 xR = 0 elif x == ntiles_x: # special case: right column - xstep = x_size / 2.0 + xstep = width / 2.0 xL = ntiles_x - 1 xR = xL else: # default values - xstep = x_size + xstep = width xL = x - 1 xR = xL + 1 @@ -210,7 +209,7 @@ def _clahe(image, ntiles_x, ntiles_y, clip_limit, nbins=128): xslice = np.arange(xstart, xstart + xstep) yslice = np.arange(ystart, ystart + ystep) interpolate(image, xslice, yslice, - mapLU, mapRU, mapLB, mapRB, aLUT) + mapLU, mapRU, mapLB, mapRB, lut) xstart += xstep # set pointer on next matrix */ @@ -305,7 +304,7 @@ def map_histogram(hist, min_val, max_val, n_pixels): def interpolate(image, xslice, yslice, - mapLU, mapRU, mapLB, mapRB, aLUT): + mapLU, mapRU, mapLB, mapRB, lut): """Find the new grayscale level for a region using bilinear interpolation. Parameters @@ -316,7 +315,7 @@ def interpolate(image, xslice, yslice, Indices of the region. map* : ndarray Mappings of greylevels from histograms. - aLUT : ndarray + lut : ndarray Maps grayscale levels in image to histogram levels. Returns @@ -338,7 +337,7 @@ def interpolate(image, xslice, yslice, view = image[int(yslice[0]):int(yslice[-1] + 1), int(xslice[0]):int(xslice[-1] + 1)] - im_slice = aLUT[view] + im_slice = lut[view] new = ((y_inv_coef * (x_inv_coef * mapLU[im_slice] + x_coef * mapRU[im_slice]) + y_coef * (x_inv_coef * mapLB[im_slice] diff --git a/skimage/exposure/tests/test_exposure.py b/skimage/exposure/tests/test_exposure.py index bb584f1a..c793c2c8 100644 --- a/skimage/exposure/tests/test_exposure.py +++ b/skimage/exposure/tests/test_exposure.py @@ -172,7 +172,7 @@ def test_adapthist_grayscale(): nbins=128) assert_almost_equal = np.testing.assert_almost_equal assert img.shape == adapted.shape - assert_almost_equal(peak_snr(img, adapted), 104.3168, 3) + assert_almost_equal(peak_snr(img, adapted), 104.307, 3) assert_almost_equal(norm_brightness_err(img, adapted), 0.0265, 3) return data, adapted From 70a3760bca2f1e8f724f90f3e3703bb5733b9763 Mon Sep 17 00:00:00 2001 From: blink1073 Date: Sat, 12 Jul 2014 19:43:57 -0500 Subject: [PATCH 13/13] Remove unused deprecation import. --- skimage/exposure/exposure.py | 1 - 1 file changed, 1 deletion(-) diff --git a/skimage/exposure/exposure.py b/skimage/exposure/exposure.py index 1b42ecf8..3a2accbc 100644 --- a/skimage/exposure/exposure.py +++ b/skimage/exposure/exposure.py @@ -3,7 +3,6 @@ import numpy as np from skimage import img_as_float from skimage.util.dtype import dtype_range, dtype_limits -from skimage._shared.utils import deprecation_warning __all__ = ['histogram', 'cumulative_distribution', 'equalize',