Add image enhancement measure (EME)

doc/examples/filters/plot_enhancement_measure.py

@@ -0,0 +1,61 @@
+"""
+=========================================
+Quantitative Measure of Image Enhancement
+=========================================
+In this example, we will quantify the effect of applying contrast enhancement
+to an image using a measure called EME [1]_.
+
+It is defined as the log of the ratio of local maximum to minimum within a
+sliding window of given size. This measure aims to provide a perceptual
+estimate of image quality, useful when comparing image enhancement algorithms.
+
+"""
+import numpy as np
+import matplotlib.pyplot as plt
+from skimage import data
+from skimage.filters import gaussian
+from skimage.exposure import equalize_adapthist
+from skimage.measure import simple_metrics
+
+img_gray = data.camera()
+img_rgb = data.astronaut()
+
+
+def compare_side_by_side(before: np.ndarray, after: np.ndarray):
+    fig, ax = plt.subplots(nrows=1, ncols=2, figsize=(16, 14))
+    for idx, image in enumerate([before, after]):
+        ax[idx].imshow(image)
+        ax[idx].axis('off')
+        ax[idx].set_title("EME = {:.2f} ".format(
+            simple_metrics.enhancement_measure(image, size=3)))
+    fig.show()
+
+
+#####################################################################
+# We will create a function that takes an original image as `before`
+# and applies a `transform` to it. We will use a gaussian filter to
+# strongly blur the image and an adaptive histogram equalization as
+# transforms.
+
+
+compare_side_by_side(before=img_gray, after=gaussian(image=img_gray,
+                                                     sigma=10))
+compare_side_by_side(before=img_rgb, after=equalize_adapthist(image=img_rgb,
+                                                              nbins=256 * 3))
+
+############################################################################
+# You can see that the greyscale camera image is very blurred due to
+# strong gaussian filtering. This resulted in smoothing a lot of pixel
+# differences and, as a result, to lower value of EME.
+# On the other hand, the `astronaut` image was enhanced by adaptive
+# histogram equalization with `256*3` bins. Histogram equalization
+# brought out details and the EME value is higher than before the transform.
+#
+# References
+# ----------
+#
+# .. [1] Sos S. Agaian, Karen Panetta, and Artyom M. Grigoryan.
+#        "A new measure of image enhancement.",
+#        IASTED International Conference on Signal Processing
+#        & Communication, Citeseer, 2000,
+#        :DOI:`10.1.1.35.4021`

skimage/measure/__init__.py

@@ -1,9 +1,14 @@
 from ._find_contours import find_contours
 from ._marching_cubes_lewiner import marching_cubes_lewiner
 from ._marching_cubes_classic import (marching_cubes_classic,
-                                      mesh_surface_area)
+                                      mesh_surface_area,
+                                      _correct_mesh_orientation)
+from ._regionprops import regionprops, perimeter
+from .simple_metrics import (compare_mse,
+                             compare_nrmse,
+                             compare_psnr,
+                             enhancement_measure)
 from ._regionprops import regionprops, perimeter, regionprops_table
-from .simple_metrics import compare_mse, compare_nrmse, compare_psnr
 from ._structural_similarity import compare_ssim
 from ._polygon import approximate_polygon, subdivide_polygon
 from .pnpoly import points_in_poly, grid_points_in_poly
@@ -47,5 +52,6 @@
            'compare_mse',
            'compare_nrmse',
            'compare_psnr',
+           'enhancement_measure',
            'shannon_entropy',
-]
+           ]

skimage/measure/simple_metrics.py

@@ -1,12 +1,19 @@
+import numpy as np
+from ..util.dtype import dtype_range
+from .._shared.utils import skimage_deprecation, warn
+from ..util import img_as_float
+from ..exposure import rescale_intensity
+from scipy.ndimage import maximum_filter, minimum_filter
 from warnings import warn
 import numpy as np
 from ..metrics.simple_metrics import (mean_squared_error,
-                                       peak_signal_noise_ratio,
-                                       normalized_root_mse)
+                                      peak_signal_noise_ratio,
+                                      normalized_root_mse)
 
 __all__ = ['compare_mse',
            'compare_nrmse',
            'compare_psnr',
+           'enhancement_measure'
            ]
 
 
@@ -23,10 +30,10 @@ def compare_mse(im1, im2):
     -----
     Deprecated:
         .. versionadded:: 0.16
-
-        This function is deprecated and will be removed in scikit-image 0.18.
-        Please use the function named ``mean_squared_error`` from the
-        ``metrics`` module instead.
+          This function is deprecated and will be removed
+          in scikit-image 0.18.
+          Please use the function named ``mean_squared_error`` from the
+          ``metrics`` module instead.
 
     See also
     --------
@@ -43,19 +50,19 @@ def compare_nrmse(im_true, im_test, norm_type='euclidean'):
 
 if normalized_root_mse.__doc__ is not None:
     compare_nrmse.__doc__ = normalized_root_mse.__doc__ + """
-    Warns
-    -----
-    Deprecated:
-        .. versionadded:: 0.16
+        Warns
+        -----
+        Deprecated:
+            .. versionadded:: 0.16
+              This function is deprecated and will be removed
+              in scikit-image 0.18.
+              Please use the function named ``normalized_root_mse`` from the
+              ``metrics`` module instead.
 
-        This function is deprecated and will be removed in scikit-image 0.18.
-        Please use the function named ``normalized_root_mse`` from the
-        ``metrics`` module instead.
-
-    See also
-    --------
-    skimage.metrics.normalized_root_mse
-    """
+        See also
+        --------
+        skimage.metrics.normalized_root_mse
+        """
 
 
 def compare_psnr(im_true, im_test, data_range=None):
@@ -71,12 +78,63 @@ def compare_psnr(im_true, im_test, data_range=None):
     -----
     Deprecated:
         .. versionadded:: 0.16
-
-        This function is deprecated and will be removed in scikit-image 0.18.
-        Please use the function named ``peak_signal_noise_ratio`` from the
-        ``metrics`` module instead.
+          This function is deprecated and will be removed
+          in scikit-image 0.18.
+          Please use the function named ``peak_signal_noise_ratio`` from the
+          ``metrics`` module instead.
 
     See also
     --------
     skimage.metrics.peak_signal_noise_ratio
     """
+
+
+def enhancement_measure(image: np.ndarray,
+                        size: int = 3,
+                        eps: float = 1e-6) -> float:
+    """ The image enhancement measure called EME based on [1]_.
+    It is a way of quantifying improvement of the image after enhancement.
+
+    The function uses a sliding window of user-provided size to measure
+    the mean of log of maximal and minimal intensity ratio
+    within the window.
+
+    Parameters
+    ----------
+    image : array
+        Input image of which the quality should be assessed.
+        Can be either 3-channel RGB or 1-channel grayscale.
+    size : int, optional
+        Size of the window.
+    eps : float, optional
+        Parameter to avoid division by zero.
+
+    Returns
+    -------
+    eme : float
+        The number describing image quality.
+
+    References
+    ----------
+    .. [1] Sos S. Agaian,  Karen Panetta, and Artyom M. Grigoryan.
+            "A new measure of image enhancement.",
+            IASTED International Conference on Signal Processing
+            & Communication, Citeseer, 2000,
+            :DOI:`10.1.1.35.4021`
+
+    Examples
+    --------
+    >>> from skimage.data import camera
+    >>> from skimage.exposure import equalize_hist
+    >>> img = camera()
+    >>> before = enhancement_measure(img)
+    >>> after = enhancement_measure(equalize_hist(img))
+    >>> before < after
+    True
+    """
+    image = img_as_float(image)
+    image = rescale_intensity(image, out_range=(0., 1.))
+    eme = np.divide(maximum_filter(image, size=size),
+                    minimum_filter(image, size=size) + eps)
+    eme = np.mean(20 * np.log(eme + eps))
+    return eme

skimage/measure/tests/test_simple_metrics.py

@@ -1,12 +1,20 @@
 import numpy as np
 
 import skimage.data
-from skimage.measure import compare_nrmse, compare_psnr, compare_mse
-
+from skimage.measure import (compare_psnr,
+                             compare_nrmse,
+                             compare_mse,
+                             enhancement_measure)
+from skimage.metrics import peak_signal_noise_ratio
+from skimage.metrics import normalized_root_mse
+from skimage.metrics import mean_squared_error
 from skimage._shared import testing
-from skimage._shared.testing import assert_equal, assert_almost_equal
+from skimage._shared.testing import (assert_equal,
+                                     assert_almost_equal,
+                                     assert_greater,
+                                     assert_less)
 from skimage._shared._warnings import expected_warnings
-
+from skimage.exposure import equalize_hist, adjust_gamma
 
 np.random.seed(5)
 cam = skimage.data.camera()
@@ -19,66 +27,88 @@ def test_PSNR_vs_IPOL():
     # Tests vs. imdiff result from the following IPOL article and code:
     # https://www.ipol.im/pub/art/2011/g_lmii/
     p_IPOL = 22.4497
-    with expected_warnings(['DEPRECATED']):
-        p = compare_psnr(cam, cam_noisy)
+    p = peak_signal_noise_ratio(cam, cam_noisy)
     assert_almost_equal(p, p_IPOL, decimal=4)
 
 
 def test_PSNR_float():
-    with expected_warnings(['DEPRECATED']):
-        p_uint8 = compare_psnr(cam, cam_noisy)
-        p_float64 = compare_psnr(cam / 255., cam_noisy / 255.,
-                                 data_range=1)
+    p_uint8 = peak_signal_noise_ratio(cam, cam_noisy)
+    p_float64 = peak_signal_noise_ratio(cam / 255., cam_noisy / 255.,
+                                        data_range=1)
     assert_almost_equal(p_uint8, p_float64, decimal=5)
 
     # mixed precision inputs
-    with expected_warnings(['DEPRECATED']):
-        p_mixed = compare_psnr(cam / 255., np.float32(cam_noisy / 255.),
-                               data_range=1)
+    p_mixed = peak_signal_noise_ratio(cam / 255., np.float32(cam_noisy / 255.),
+                                      data_range=1)
     assert_almost_equal(p_mixed, p_float64, decimal=5)
 
     # mismatched dtype results in a warning if data_range is unspecified
-    with expected_warnings(['Inputs have mismatched dtype', 'DEPRECATED']):
-        p_mixed = compare_psnr(cam / 255., np.float32(cam_noisy / 255.))
+    with expected_warnings(['Inputs have mismatched dtype']):
+        p_mixed = peak_signal_noise_ratio(cam / 255.,
+                                          np.float32(cam_noisy / 255.))
     assert_almost_equal(p_mixed, p_float64, decimal=5)
 
 
 def test_PSNR_errors():
-    with expected_warnings(['DEPRECATED']):
-        # shape mismatch
-        with testing.raises(ValueError):
-            compare_psnr(cam, cam[:-1, :])
+    # shape mismatch
+    with testing.raises(ValueError):
+        peak_signal_noise_ratio(cam, cam[:-1, :])
 
 
 def test_NRMSE():
     x = np.ones(4)
     y = np.asarray([0., 2., 2., 2.])
-    with expected_warnings(['DEPRECATED']):
-        assert_equal(compare_nrmse(y, x, 'mean'), 1 / np.mean(y))
-        assert_equal(compare_nrmse(y, x, 'Euclidean'), 1 / np.sqrt(3))
-        assert_equal(compare_nrmse(y, x, 'min-max'), 1 / (y.max() - y.min()))
-
-        # mixed precision inputs are allowed
-        assert_almost_equal(compare_nrmse(y, np.float32(x), 'min-max'),
-                            1 / (y.max() - y.min()))
+    assert_equal(compare_nrmse(y, x, 'mean'), 1 / np.mean(y))
+    assert_equal(compare_nrmse(y, x, 'Euclidean'), 1 / np.sqrt(3))
+    assert_equal(compare_nrmse(y, x, 'min-max'), 1 / (y.max() - y.min()))
+    assert_equal(normalized_root_mse(
+        y, x, normalization='mean'), 1 / np.mean(y))
+    assert_equal(normalized_root_mse(
+        y, x, normalization='Euclidean'), 1 / np.sqrt(3))
+    assert_equal(normalized_root_mse(
+        y, x, normalization='min-max'), 1 / (y.max() - y.min()))
+
+    # mixed precision inputs are allowed
+    assert_almost_equal(normalized_root_mse(y, np.float32(x),
+                                            normalization='min-max'),
+                        1 / (y.max() - y.min()))
 
 
 def test_NRMSE_no_int_overflow():
     camf = cam.astype(np.float32)
     cam_noisyf = cam_noisy.astype(np.float32)
-    with expected_warnings(['DEPRECATED']):
-        assert_almost_equal(compare_mse(cam, cam_noisy),
-                            compare_mse(camf, cam_noisyf))
-        assert_almost_equal(compare_nrmse(cam, cam_noisy),
-                            compare_nrmse(camf, cam_noisyf))
+    assert_almost_equal(mean_squared_error(cam, cam_noisy),
+                        mean_squared_error(camf, cam_noisyf))
+    assert_almost_equal(normalized_root_mse(cam, cam_noisy),
+                        normalized_root_mse(camf, cam_noisyf))
 
 
 def test_NRMSE_errors():
     x = np.ones(4)
-    with expected_warnings(['DEPRECATED']):
-        # shape mismatch
-        with testing.raises(ValueError):
-            compare_nrmse(x[:-1], x)
-        # invalid normalization name
-        with testing.raises(ValueError):
-            compare_nrmse(x, x, norm_type='foo')
+    # shape mismatch
+    with testing.raises(ValueError):
+        normalized_root_mse(x[:-1], x)
+    # invalid normalization name
+    with testing.raises(ValueError):
+        compare_nrmse(x, x, 'foo')
+
+
+def test_EME_greyscale():
+    orig = cam
+    enhanced = equalize_hist(orig)
+    assert_less(enhancement_measure(orig),
+                enhancement_measure(enhanced))
+    assert_less(enhancement_measure(orig, size=5),
+                enhancement_measure(enhanced, size=5))
+    assert_less(enhancement_measure(orig, size=11),
+                enhancement_measure(enhanced, size=11))
+
+
+def test_EME_color():
+    orig = skimage.data.astronaut()
+    enhanced = adjust_gamma(orig, gamma=5)
+    assert_less(enhancement_measure(orig, size=3),
+                enhancement_measure(enhanced, size=3, eps=1))
+    assert_less(enhancement_measure(orig, size=7),
+                enhancement_measure(enhanced, size=7))
+    # normalized_root_mse(x, x, 'foo')