Files
Bioxel/bioxelnodes/bioxel/scipy/_interpolation.py
Ma Nan ce397b35aa feat: more time sequence choice
fix: contrast formulation cause leak of high value
fix: empty series_id cause fail to read
chore: restructure folder
chore!: change layer vars
BREAKING CHANGE: not support 4.1 any more
2024-08-01 15:50:41 +08:00

314 lines
11 KiB
Python

import numpy as np
import warnings
from . import _ni_support
from . import _nd_image
from ._utils import normalize_axis_index
def _prepad_for_spline_filter(input, mode, cval):
if mode in ['nearest', 'grid-constant']:
npad = 12
if mode == 'grid-constant':
padded = np.pad(input, npad, mode='constant',
constant_values=cval)
elif mode == 'nearest':
padded = np.pad(input, npad, mode='edge')
else:
# other modes have exact boundary conditions implemented so
# no prepadding is needed
npad = 0
padded = input
return padded, npad
def spline_filter1d(input, order=3, axis=-1, output=np.float64,
mode='mirror'):
"""
Calculate a 1-D spline filter along the given axis.
The lines of the array along the given axis are filtered by a
spline filter. The order of the spline must be >= 2 and <= 5.
Parameters
----------
%(input)s
order : int, optional
The order of the spline, default is 3.
axis : int, optional
The axis along which the spline filter is applied. Default is the last
axis.
output : ndarray or dtype, optional
The array in which to place the output, or the dtype of the returned
array. Default is ``numpy.float64``.
%(mode_interp_mirror)s
Returns
-------
spline_filter1d : ndarray
The filtered input.
See Also
--------
spline_filter : Multidimensional spline filter.
Notes
-----
All of the interpolation functions in `ndimage` do spline interpolation of
the input image. If using B-splines of `order > 1`, the input image
values have to be converted to B-spline coefficients first, which is
done by applying this 1-D filter sequentially along all
axes of the input. All functions that require B-spline coefficients
will automatically filter their inputs, a behavior controllable with
the `prefilter` keyword argument. For functions that accept a `mode`
parameter, the result will only be correct if it matches the `mode`
used when filtering.
For complex-valued `input`, this function processes the real and imaginary
components independently.
.. versionadded:: 1.6.0
Complex-valued support added.
Examples
--------
We can filter an image using 1-D spline along the given axis:
>>> from scipy.ndimage import spline_filter1d
>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> orig_img = np.eye(20) # create an image
>>> orig_img[10, :] = 1.0
>>> sp_filter_axis_0 = spline_filter1d(orig_img, axis=0)
>>> sp_filter_axis_1 = spline_filter1d(orig_img, axis=1)
>>> f, ax = plt.subplots(1, 3, sharex=True)
>>> for ind, data in enumerate([[orig_img, "original image"],
... [sp_filter_axis_0, "spline filter (axis=0)"],
... [sp_filter_axis_1, "spline filter (axis=1)"]]):
... ax[ind].imshow(data[0], cmap='gray_r')
... ax[ind].set_title(data[1])
>>> plt.tight_layout()
>>> plt.show()
"""
if order < 0 or order > 5:
raise RuntimeError('spline order not supported')
input = np.asarray(input)
complex_output = np.iscomplexobj(input)
output = _ni_support._get_output(output, input,
complex_output=complex_output)
if complex_output:
spline_filter1d(input.real, order, axis, output.real, mode)
spline_filter1d(input.imag, order, axis, output.imag, mode)
return output
if order in [0, 1]:
output[...] = np.array(input)
else:
mode = _ni_support._extend_mode_to_code(mode)
axis = normalize_axis_index(axis, input.ndim)
_nd_image.spline_filter1d(input, order, axis, output, mode)
return output
def spline_filter(input, order=3, output=np.float64, mode='mirror'):
"""
Multidimensional spline filter.
Parameters
----------
%(input)s
order : int, optional
The order of the spline, default is 3.
output : ndarray or dtype, optional
The array in which to place the output, or the dtype of the returned
array. Default is ``numpy.float64``.
%(mode_interp_mirror)s
Returns
-------
spline_filter : ndarray
Filtered array. Has the same shape as `input`.
See Also
--------
spline_filter1d : Calculate a 1-D spline filter along the given axis.
Notes
-----
The multidimensional filter is implemented as a sequence of
1-D spline filters. The intermediate arrays are stored
in the same data type as the output. Therefore, for output types
with a limited precision, the results may be imprecise because
intermediate results may be stored with insufficient precision.
For complex-valued `input`, this function processes the real and imaginary
components independently.
.. versionadded:: 1.6.0
Complex-valued support added.
Examples
--------
We can filter an image using multidimentional splines:
>>> from scipy.ndimage import spline_filter
>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> orig_img = np.eye(20) # create an image
>>> orig_img[10, :] = 1.0
>>> sp_filter = spline_filter(orig_img, order=3)
>>> f, ax = plt.subplots(1, 2, sharex=True)
>>> for ind, data in enumerate([[orig_img, "original image"],
... [sp_filter, "spline filter"]]):
... ax[ind].imshow(data[0], cmap='gray_r')
... ax[ind].set_title(data[1])
>>> plt.tight_layout()
>>> plt.show()
"""
if order < 2 or order > 5:
raise RuntimeError('spline order not supported')
input = np.asarray(input)
complex_output = np.iscomplexobj(input)
output = _ni_support._get_output(output, input,
complex_output=complex_output)
if complex_output:
spline_filter(input.real, order, output.real, mode)
spline_filter(input.imag, order, output.imag, mode)
return output
if order not in [0, 1] and input.ndim > 0:
for axis in range(input.ndim):
spline_filter1d(input, order, axis, output=output, mode=mode)
input = output
else:
output[...] = input[...]
return output
def zoom(input, zoom, output=None, order=3, mode='constant', cval=0.0,
prefilter=True, *, grid_mode=False):
"""
Zoom an array.
The array is zoomed using spline interpolation of the requested order.
Parameters
----------
%(input)s
zoom : float or sequence
The zoom factor along the axes. If a float, `zoom` is the same for each
axis. If a sequence, `zoom` should contain one value for each axis.
%(output)s
order : int, optional
The order of the spline interpolation, default is 3.
The order has to be in the range 0-5.
%(mode_interp_constant)s
%(cval)s
%(prefilter)s
grid_mode : bool, optional
If False, the distance from the pixel centers is zoomed. Otherwise, the
distance including the full pixel extent is used. For example, a 1d
signal of length 5 is considered to have length 4 when `grid_mode` is
False, but length 5 when `grid_mode` is True. See the following
visual illustration:
.. code-block:: text
| pixel 1 | pixel 2 | pixel 3 | pixel 4 | pixel 5 |
|<-------------------------------------->|
vs.
|<----------------------------------------------->|
The starting point of the arrow in the diagram above corresponds to
coordinate location 0 in each mode.
Returns
-------
zoom : ndarray
The zoomed input.
Notes
-----
For complex-valued `input`, this function zooms the real and imaginary
components independently.
.. versionadded:: 1.6.0
Complex-valued support added.
Examples
--------
>>> from scipy import ndimage, datasets
>>> import matplotlib.pyplot as plt
>>> fig = plt.figure()
>>> ax1 = fig.add_subplot(121) # left side
>>> ax2 = fig.add_subplot(122) # right side
>>> ascent = datasets.ascent()
>>> result = ndimage.zoom(ascent, 3.0)
>>> ax1.imshow(ascent, vmin=0, vmax=255)
>>> ax2.imshow(result, vmin=0, vmax=255)
>>> plt.show()
>>> print(ascent.shape)
(512, 512)
>>> print(result.shape)
(1536, 1536)
"""
if order < 0 or order > 5:
raise RuntimeError('spline order not supported')
input = np.asarray(input)
if input.ndim < 1:
raise RuntimeError('input and output rank must be > 0')
zoom = _ni_support._normalize_sequence(zoom, input.ndim)
output_shape = tuple(
[int(round(ii * jj)) for ii, jj in zip(input.shape, zoom)])
complex_output = np.iscomplexobj(input)
output = _ni_support._get_output(output, input, shape=output_shape,
complex_output=complex_output)
if complex_output:
# import under different name to avoid confusion with zoom parameter
from scipy.ndimage._interpolation import zoom as _zoom
kwargs = dict(order=order, mode=mode, prefilter=prefilter)
_zoom(input.real, zoom, output=output.real,
cval=np.real(cval), **kwargs)
_zoom(input.imag, zoom, output=output.imag,
cval=np.imag(cval), **kwargs)
return output
if prefilter and order > 1:
padded, npad = _prepad_for_spline_filter(input, mode, cval)
filtered = spline_filter(padded, order, output=np.float64, mode=mode)
else:
npad = 0
filtered = input
if grid_mode:
# warn about modes that may have surprising behavior
suggest_mode = None
if mode == 'constant':
suggest_mode = 'grid-constant'
elif mode == 'wrap':
suggest_mode = 'grid-wrap'
if suggest_mode is not None:
warnings.warn(
(f"It is recommended to use mode = {suggest_mode} instead of {mode} "
f"when grid_mode is True."),
stacklevel=2
)
mode = _ni_support._extend_mode_to_code(mode)
zoom_div = np.array(output_shape)
zoom_nominator = np.array(input.shape)
if not grid_mode:
zoom_div -= 1
zoom_nominator -= 1
# Zooming to infinite values is unpredictable, so just choose
# zoom factor 1 instead
zoom = np.divide(zoom_nominator, zoom_div,
out=np.ones_like(input.shape, dtype=np.float64),
where=zoom_div != 0)
zoom = np.ascontiguousarray(zoom)
_nd_image.zoom_shift(filtered, zoom, None, output, order, mode, cval, npad,
grid_mode)
return output