# Source code for africanus.filters.kaiser_bessel_filter

# -*- 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