Add optional pipeline noise reduction module with Gaussian and CLAHE methods

docs/_toc.yml

@@ -22,6 +22,7 @@ parts:
   chapters:
   - file: notebooks/pipeline_step_by_step
   - file: notebooks/background_correction
+  - file: notebooks/noise_reduction
   - file: notebooks/exploring_pipeline_data
   - file: notebooks/montaging
 - caption: Code docs

docs/notebooks/noise_reduction.md

@@ -0,0 +1,44 @@
+# Noise reduction
+
+For particularly noisy image data, you can add an optional noise reduction step in the pipeline.
+This can be inserted anywhere it is useful, for example between image loading and background correction,
+or just before segmentation.
+
+## Example: Gaussian noise reduction before background correction
+
+```toml
+[steps.load]
+pipeline_class = "pyopia.instrument.silcam.SilCamLoad"
+
+[steps.noisereduction]
+pipeline_class = "pyopia.noise.ReduceNoise"
+method = "gaussian"
+image_source = "imraw"
+gaussian_sigma = 1.0
+
+[steps.correctbackground]
+pipeline_class = "pyopia.background.CorrectBackgroundAccurate"
+bgshift_function = "accurate"
+average_window = 5
+image_source = "im_denoised"
+```
+
+## Example: CLAHE before segmentation
+
+```toml
+[steps.imageprep]
+pipeline_class = "pyopia.instrument.silcam.ImagePrep"
+image_level = "im_corrected"
+
+[steps.noisereduction]
+pipeline_class = "pyopia.noise.ReduceNoise"
+method = "clahe"
+image_source = "im_minimum"
+clahe_clip_limit = 0.01
+clahe_nbins = 256
+
+[steps.segmentation]
+pipeline_class = "pyopia.process.Segment"
+threshold = 0.85
+segment_source = "im_denoised"
+```

pyopia/noise.py

@@ -0,0 +1,81 @@
+'''
+Noise reduction module for optional preprocessing steps in the pipeline.
+'''
+from skimage import exposure
+from skimage import filters
+
+
+class ReduceNoise():
+    '''
+    :class:`pyopia.pipeline` compatible class for optional noise reduction.
+
+    Required keys in :class:`pyopia.pipeline.Data`:
+        - key given by ``image_source``
+
+    Parameters
+    ----------
+    method : str, optional
+        Noise reduction method. Supported values are ``'gaussian'`` and ``'clahe'``.
+        Defaults to ``'gaussian'``.
+    image_source : str, optional
+        Key in Pipeline.data containing the input image. Defaults to ``'im_corrected'``.
+    output_key : str | None, optional
+        Key where the filtered image will be stored. If ``None``, the filtered image is stored
+        in ``'im_denoised'`` (i.e. not in place). Defaults to ``None``.
+    gaussian_sigma : float, optional
+        Gaussian sigma when ``method='gaussian'``. Defaults to ``1.0``.
+    clahe_clip_limit : float, optional
+        CLAHE clip limit when ``method='clahe'``. Defaults to ``0.01``.
+    clahe_nbins : int, optional
+        Number of bins for CLAHE when ``method='clahe'``. Defaults to ``256``.
+
+    Returns
+    -------
+    data : :class:`pyopia.pipeline.Data`
+        containing filtered image in ``output_key``.
+
+    Example pipeline use
+    --------------------
+
+    .. code-block:: toml
+
+        [steps.noisereduction]
+        pipeline_class = 'pyopia.noise.ReduceNoise'
+        method = 'gaussian'
+        image_source = 'imraw'
+        output_key = 'im_denoised'
+        gaussian_sigma = 1.0
+    '''
+
+    def __init__(self,
+                 method='gaussian',
+                 image_source='im_corrected',
+                 output_key=None,
+                 gaussian_sigma=1.0,
+                 clahe_clip_limit=0.01,
+                 clahe_nbins=256):
+
+        self.method = method
+        self.image_source = image_source
+        self.output_key = 'im_denoised' if output_key is None else output_key
+        self.gaussian_sigma = gaussian_sigma
+        self.clahe_clip_limit = clahe_clip_limit
+        self.clahe_nbins = clahe_nbins
+
+    def __call__(self, data):
+        im = data[self.image_source]
+
+        if self.method == 'gaussian':
+            data[self.output_key] = filters.gaussian(im,
+                                                     sigma=self.gaussian_sigma,
+                                                     preserve_range=True,
+                                                     channel_axis=-1 if im.ndim == 3 else None)
+        elif self.method == 'clahe':
+            data[self.output_key] = exposure.equalize_adapthist(im,
+                                                                clip_limit=self.clahe_clip_limit,
+                                                                nbins=self.clahe_nbins)
+        else:
+            raise ValueError(f"Unknown noise reduction method '{self.method}'. "
+                             "Expected one of: 'gaussian', 'clahe'.")
+
+        return data

pyopia/tests/test_noise.py

@@ -0,0 +1,67 @@
+import numpy as np
+import pytest
+
+from pyopia.noise import ReduceNoise
+from pyopia.pipeline import Pipeline
+
+
+class _Data(dict):
+    pass
+
+
+def test_reduce_noise_gaussian_smooths_image():
+    data = _Data()
+    data['imraw'] = np.zeros((5, 5), dtype=float)
+    data['imraw'][2, 2] = 1.0
+
+    reducer = ReduceNoise(method='gaussian', image_source='imraw', output_key='im_denoised', gaussian_sigma=1.0)
+    out = reducer(data)
+
+    assert out['im_denoised'].shape == data['imraw'].shape
+    assert out['im_denoised'][2, 2] < 1.0
+    assert out['im_denoised'][2, 2] > 0.0
+
+
+def test_reduce_noise_clahe_returns_contrast_enhanced_float_image():
+    data = _Data()
+    data['im_corrected'] = np.full((8, 8), 0.4, dtype=float)
+    data['im_corrected'][2:6, 2:6] = 0.6
+
+    reducer = ReduceNoise(method='clahe', image_source='im_corrected')
+    out = reducer(data)
+
+    assert out['im_denoised'].shape == data['im_corrected'].shape
+    assert np.issubdtype(out['im_denoised'].dtype, np.floating)
+    assert out['im_denoised'].min() >= 0.0
+    assert out['im_denoised'].max() <= 1.0
+
+
+def test_reduce_noise_invalid_method_raises_value_error():
+    data = _Data(im_corrected=np.ones((4, 4), dtype=float))
+
+    reducer = ReduceNoise(method='invalid', image_source='im_corrected')
+
+    with pytest.raises(ValueError, match='Unknown noise reduction method'):
+        reducer(data)
+
+
+def test_reduce_noise_can_run_as_pipeline_step():
+    settings = {
+        'general': {'raw_files': None},
+        'steps': {
+            'noisereduction': {
+                'pipeline_class': 'pyopia.noise.ReduceNoise',
+                'method': 'gaussian',
+                'image_source': 'im_corrected',
+                'gaussian_sigma': 0.5,
+            }
+        }
+    }
+    pipeline = Pipeline(settings=settings, initial_steps='')
+    pipeline.data['im_corrected'] = np.zeros((6, 6), dtype=float)
+    pipeline.data['im_corrected'][3, 3] = 1.0
+
+    pipeline.run_step('noisereduction')
+
+    assert pipeline.data['im_denoised'][3, 3] < 1.0
+    assert pipeline.data['im_corrected'][3, 3] == 1.0