# -*- coding: utf-8 -*-
from __future__ import absolute_import
from __future__ import division
from __future__ import print_function
import numpy as np
[docs]def estimate_kaiser_bessel_beta(W):
r"""
Estimate the kaiser bessel beta using the following heuristic:
.. math::
\beta = 2.34 \times W
Derived from `Nonuniform fast Fourier transforms
using min-max interpolation
<nufft-min-max-ref_>`_.
.. _nufft-min-max-ref: https://ieeexplore.ieee.org/document/1166689/
Parameters
----------
W : int
Width of the filter
Returns
-------
float
Kaiser Bessel beta shape parameter
"""
return 2.34*W
[docs]def kaiser_bessel(u, W, beta):
r"""
Compute a 1D Kaiser Bessel filter as defined
in `Selection of a Convolution Function
for Fourier Inversion Using Gridding
<kbref_>`_.
.. _kbref: https://ieeexplore.ieee.org/document/97598/
Parameters
----------
u : :class:`numpy.ndarray`
Filter positions
W : int
Width of the filter
beta : float, optional
Kaiser Bessel shape parameter
Returns
-------
:class:`numpy.ndarray`
Kaiser Bessel filter with the same shape as `u`
"""
# Sanity check
hW = W // 2
assert np.all(-hW <= u) & np.all(u <= hW)
param = 1 - (2 * u / W)**2
param[param < 0] = 0 # Zero negative numbers
return np.i0(beta * np.sqrt(param)) / np.i0(beta)
[docs]def kaiser_bessel_with_sinc(u, W, oversample, beta, normalise=True):
"""
Produces a filter composed of Kaiser Bessel multiplied by a sinc.
Accounts for the oversampling factor, as well as normalising the filter.
Parameters
----------
u : :class:`numpy.ndarray`
Filter positions
W : int
Width of the filter
oversample : int
Oversampling factor
beta : float
Kaiser Bessel shape parameter
normalise : optional, {True, False}
True if the filter should be normalised
Returns
-------
:class:`numpy.ndarray`
Filter with the same shape as `u`
"""
kb = kaiser_bessel(u, W, beta)
kb *= oversample
kb *= np.sinc(u / oversample)
if normalise:
kb /= np.trapz(kb, u)
return kb
[docs]def kaiser_bessel_fourier(x, W, beta):
r"""
Computes the Fourier Transform of a 1D Kaiser Bessel filter.
as defined in `Selection of a Convolution Function
for Fourier Inversion Using Gridding
<kbref_>`_.
.. _kbref: https://ieeexplore.ieee.org/document/97598/
Parameters
----------
x : :class:`numpy.ndarray`
Filter positions
W : int
Width of the filter.
beta : float
Kaiser bessel shape parameter
Returns
-------
:class:`numpy.ndarray`
Fourier Transform of the Kaiser Bessel,
with the same shape as `x`.
"""
term = (np.pi*W*x)**2 - beta**2
val = np.lib.scimath.sqrt(term).real
return np.sin(val)/val