Coverage for install/scipp/scipy/ndimage/__init__.py: 95%
74 statements
« prev ^ index » next coverage.py v7.4.0, created at 2024-04-28 01:28 +0000
1# SPDX-License-Identifier: BSD-3-Clause
2# Copyright (c) 2023 Scipp contributors (https://github.com/scipp)
3# @author Simon Heybrock
4"""Sub-package for multidimensional image processing.
6This subpackage provides wrappers for a subset of functions from
7:py:mod:`scipy.ndimage`.
8"""
10from functools import wraps
11from typing import Callable, Dict, Optional, Union
13import scipy.ndimage
15from ...core import (
16 CoordError,
17 DataArray,
18 DimensionError,
19 Variable,
20 VariancesError,
21 empty_like,
22 islinspace,
23 ones,
24)
25from ...typing import VariableLike, VariableLikeType
28def _ndfilter(func: Callable) -> Callable:
29 @wraps(func)
30 def function(x: Union[Variable, DataArray], **kwargs) -> Union[Variable, DataArray]:
31 if 'output' in kwargs:
32 raise TypeError("The 'output' argument is not supported")
33 if x.variances is not None:
34 raise VariancesError(
35 "Filter cannot be applied to input array with variances."
36 )
37 if getattr(x, 'masks', None):
38 raise ValueError("Filter cannot be applied to input array with masks.")
39 return func(x, **kwargs)
41 return function
44def _delta_to_positional(x: Union[Variable, DataArray], dim, index, dtype):
45 if not isinstance(index, Variable):
46 return index
47 coord = x.coords[dim]
48 if not islinspace(coord, dim).value:
49 raise CoordError(
50 f"Data points not regularly spaced along {dim}. To ignore this, "
51 f"provide a plain value (convertible to {dtype.__name__}) instead of a "
52 "scalar variable. Note that this will correspond to plain positional "
53 "indices/offsets."
54 )
55 pos = (len(coord) - 1) * (index.to(unit=coord.unit) / (coord[-1] - coord[0])).value
56 return dtype(pos)
59def _require_matching_dims(index, x, name):
60 if set(index) != set(x.dims):
61 raise KeyError(
62 f"Data has dims={x.dims} but input argument '{name}' provides "
63 f"values for {tuple(index)}"
64 )
67def _positional_index(x: Union[Variable, DataArray], index, name=None, dtype=int):
68 if not isinstance(index, dict):
69 return [_delta_to_positional(x, dim, index, dtype=dtype) for dim in x.dims]
70 _require_matching_dims(index, x, name)
71 return [_delta_to_positional(x, dim, index[dim], dtype=dtype) for dim in x.dims]
74@_ndfilter
75def gaussian_filter(
76 x: VariableLikeType,
77 /,
78 *,
79 sigma: Union[float, Variable, Dict[str, Union[int, float, Variable]]],
80 order: Optional[Union[int, Dict[str, int]]] = 0,
81 **kwargs,
82) -> VariableLikeType:
83 """
84 Multidimensional Gaussian filter.
86 This is a wrapper around :py:func:`scipy.ndimage.gaussian_filter`. See there for
87 full argument description. There are two key differences:
89 - This wrapper uses explicit dimension labels in the ``sigma`` and ``order``
90 arguments. For example, instead of ``sigma=[4, 6]`` use
91 ``sigma={'time':4, 'space':6}``
92 (with appropriate dimension labels for the data).
93 - Coordinate values can be used (and should be preferred) for ``sigma``. For
94 example, instead of ``sigma=[4, 6]`` use
95 ``sigma={'time':sc.scalar(5.0, unit='ms'), 'space':sc.scalar(1.2, unit='mm')}``.
96 In this case it is required that the corresponding coordinates exist and form a
97 "linspace", i.e., are evenly spaced.
99 Warning
100 -------
101 If ``sigma`` is an integer or a mapping to integers then coordinate values are
102 ignored. That is, the filter is applied even if the data points are not evenly
103 spaced. The resulting filtered data may thus have no meaningful interpretation.
105 Parameters
106 ----------
107 x: scipp.typing.VariableLike
108 Input variable or data array.
109 sigma:
110 Standard deviation for Gaussian kernel. The standard deviations of the Gaussian
111 filter are given as a mapping from dimension labels to numbers or scalar
112 variables, or as a single number or scalar variable, in which case it is equal
113 for all axes.
114 order:
115 The order of the filter along each dimension, given as mapping from dimension
116 labels to integers, or as a single integer. An order of 0 corresponds to
117 convolution with a Gaussian kernel. A positive order corresponds to convolution
118 with that derivative of a Gaussian.
120 Returns
121 -------
122 : scipp.typing.VariableLike
123 Filtered variable or data array
125 Examples
126 --------
128 .. plot:: :context: close-figs
130 >>> from scipp.scipy.ndimage import gaussian_filter
131 >>> da = sc.data.data_xy()
132 >>> da.plot()
134 With sigma as integer:
136 .. plot:: :context: close-figs
138 >>> filtered = gaussian_filter(da, sigma=4)
139 >>> filtered.plot()
141 With sigma based on input coordinate values:
143 .. plot:: :context: close-figs
145 >>> filtered = gaussian_filter(da, sigma=sc.scalar(0.1, unit='mm'))
146 >>> filtered.plot()
148 With different sigma for different dimensions:
150 .. plot:: :context: close-figs
152 >>> filtered = gaussian_filter(da, sigma={'x':sc.scalar(0.01, unit='mm'),
153 ... 'y':sc.scalar(1.0, unit='mm')})
154 >>> filtered.plot()
155 """
156 sigma = _positional_index(x, sigma, name='sigma', dtype=float)
157 if isinstance(order, dict):
158 _require_matching_dims(order, x, 'order')
159 order = [order[dim] for dim in x.dims]
160 out = empty_like(x)
161 scipy.ndimage.gaussian_filter(
162 x.values, sigma=sigma, order=order, output=out.values, **kwargs
163 )
164 return out
167def _make_footprint(x: Union[Variable, DataArray], size, footprint) -> Variable:
168 if footprint is None:
169 size = _positional_index(x, size, name='size')
170 footprint = ones(dims=x.dims, shape=size, dtype='bool')
171 else:
172 if size is not None:
173 raise ValueError("Provide either 'size' or 'footprint', not both.")
174 if set(footprint.dims) != set(x.dims):
175 raise DimensionError(
176 f"Dimensions {footprint.dims} must match data dimensions {x.dim}"
177 )
178 return footprint
181def _make_footprint_filter(name, example=True, extra_args=''):
182 def footprint_filter(
183 x: VariableLike,
184 /,
185 *,
186 size: Optional[Union[int, Variable, Dict[str, Union[int, Variable]]]] = None,
187 footprint: Optional[Variable] = None,
188 origin: Optional[Union[int, Variable, Dict[str, Union[int, Variable]]]] = 0,
189 **kwargs,
190 ) -> VariableLike:
191 footprint = _make_footprint(x, size=size, footprint=footprint)
192 origin = _positional_index(x, origin, name='origin')
193 out = empty_like(x)
194 scipy_filter = getattr(scipy.ndimage, name)
195 scipy_filter(
196 x.values,
197 footprint=footprint.values,
198 origin=origin,
199 output=out.values,
200 **kwargs,
201 )
202 return out
204 footprint_filter.__name__ = name
205 if extra_args:
206 extra_args = f', {extra_args}'
207 doc = f"""
208 Calculate a multidimensional {name.replace('_', ' ')}.
210 This is a wrapper around :py:func:`scipy.ndimage.{name}`. See there for full
211 argument description. There are two key differences:
213 - This wrapper uses explicit dimension labels in the ``size``, ``footprint``, and
214 ``origin`` arguments. For example, instead of ``size=[4, 6]`` use
215 ``size={{'time':4, 'space':6}}`` (with appropriate dimension labels for the data).
216 - Coordinate values can be used (and should be preferred) for ``size`` and
217 ``origin``. For example, instead of ``size=[4, 6]`` use
218 ``size={{'time':sc.scalar(5.0, unit='ms'), 'space':sc.scalar(1.2, unit='mm')}}``.
219 In this case it is required that the corresponding coordinates exist and form a
220 "linspace", i.e., are evenly spaced.
222 Warning
223 -------
224 When ``size`` is an integer or a mapping to integers or when ``footprint`` is
225 given, coordinate values are ignored. That is, the filter is applied even if the
226 data points are not evenly spaced. The resulting filtered data may thus have no
227 meaningful interpretation.
229 Parameters
230 ----------
231 x: scipp.typing.VariableLike
232 Input variable or data array.
233 size:
234 Integer or scalar variable or mapping from dimension labels to integers or
235 scalar variables. Defines the footprint (see below).
236 footprint:
237 Variable with same dimension labels (but different shape) as the input data.
238 The boolean values specify (implicitly) a shape, but also which of the elements
239 within this shape will get passed to the filter function.
240 origin:
241 Integer or scalar variable or mapping from dimension labels to integers or
242 scalar variables. Controls the placement of the filter on the input array.
244 Returns
245 -------
246 : scipp.typing.VariableLike
247 Filtered variable or data array
248 """
249 if example:
250 doc += f"""
251 Examples
252 --------
254 .. plot:: :context: close-figs
256 >>> from scipp.scipy.ndimage import {name}
257 >>> da = sc.data.data_xy()
258 >>> da.plot()
260 With size as integer:
262 .. plot:: :context: close-figs
264 >>> filtered = {name}(da, size=4{extra_args})
265 >>> filtered.plot()
267 With size based on input coordinate values:
269 .. plot:: :context: close-figs
271 >>> filtered = {name}(da, size=sc.scalar(0.2, unit='mm'){extra_args})
272 >>> filtered.plot()
274 With different size for different dimensions:
276 .. plot:: :context: close-figs
278 >>> filtered = {name}(da, size={{'x':sc.scalar(0.2, unit='mm'),
279 ... {' ' * len(name)}'y':sc.scalar(1.1, unit='mm')}}{extra_args})
280 >>> filtered.plot()
281 """ # noqa: E501
282 footprint_filter.__doc__ = doc
283 return _ndfilter(footprint_filter)
286generic_filter = _make_footprint_filter('generic_filter', example=False)
287maximum_filter = _make_footprint_filter('maximum_filter')
288median_filter = _make_footprint_filter('median_filter')
289minimum_filter = _make_footprint_filter('minimum_filter')
290percentile_filter = _make_footprint_filter(
291 'percentile_filter', extra_args='percentile=80'
292)
293rank_filter = _make_footprint_filter('rank_filter', extra_args='rank=3')
295__all__ = [
296 'gaussian_filter',
297 'generic_filter',
298 'maximum_filter',
299 'median_filter',
300 'minimum_filter',
301 'percentile_filter',
302 'rank_filter',
303]