nufft package

Core implementation of the non-uniform Fourier transform (NUFFT), based on [Fessler2003].

nufft.compute_interp_coeffs(om, Nd, Jd, Kd, n_shift, varargin)

Compute the interpolation coefficients involved in the direct NUFFT operator (adapted from original code associated with [Fessler2003]).

Parameters
  • om (double[:, d]) – Non-uniform coordinates of the target of the data points in the \(d\)-dimensional Fourier domain.

  • Nd (int[1, d]) – Input tensor size along each dimension.

  • Jd (int[1, d]) – Size of the interpolation kernel along each dimension.

  • Kd (int[1, d]) – Size of the Fourier domain along each dimension.

  • n_shift (int[1, d]) – Phase-shift along each dimension in Fourier space (expressed in number of samples in the Fourier domain).

  • varargin ({string} or {string, cell, cell}) – User-defined type of kernel, to be selected among the options 'minmax:kb' (used by default), 'minmax:tuned' and 'kaiser'. Possible options are documented below.

  • ktype ('minmax:kb') – Minmax interpolator with excellent KB scaling ('minmax:kb' is recommended, and used by default).

  • ktype ('minmax:tuned') – Minmax interpolator, somewhat numerically tuned.

  • ktype ('kaiser') – Kaiser-bessel (KB) interpolator (using default minmax best alpha, m).

  • ktype, alpha, m ('kaiser', double[:], double[:]) – Kaiser-bessel (KB) interpolator using the specified values for the parameters alpha and m).

Returns

st – Structure containing the NUFFT scaling coefficients (st.sn) and Auxiliary parameters to build the NUFFT interpolation matrix (st.alpha, st.beta, st.uu).

Return type

struct

nufft.createG(om, J, N, K, nshift, varargin)

Create the gridding matrix G to perform the 2D NUFFT.

Parameters
  • om (double[M, 2]) – Normalized coordinates of the target of the data points in the \(d\)-dimensional Fourier domain (uv frequencies) (M data points).

  • J (int[1, 2]) – Size of the interpolation kernel along each dimension.

  • N (int[1, 2]) – Image size along each dimension.

  • K (int[1, 2]) – Size of the Fourier domain along each dimension.

  • nshift (int[1, 2]) – Phase-shift along each dimension in Fourier space (expressed in number of samples in the Fourier domain).

  • ktype (string (varargin)) – Selected interpolator (among 'minmax:kb', 'minmax:tuned' and 'kaiser', see nufft.compute_interp_coeffs_extended() for associated kernel documentation). Recommended default is 'minmax:kb'.

Returns

  • G (sparse complex[M, prod(K)]) – Sparse de-gridding interpolation matrix involed in the 2D NUFFT.

  • scale (double[:, :] [N]) – Scaling weights involed in the 2D NUFFT, of size N.

nufft.createG_calib(D, om, N, K, J, S, nshift, varargin)

Create the gridding matrix G to perform the 2D NUFFT including DDE kernels.

Parameters
  • D (double[:, :]) – DDEs kernels for a single time instant [S2, na].

  • om (double[M, 2]) – Normalized coordinates of the target of the data points in the \(d\)-dimensional Fourier domain (uv frequencies).

  • N (int[1, 2]) – Image size along each dimension.

  • K (int[1, 2]) – Size of the Fourier domain along each dimension.

  • J (int) – Size of the square interpolation kernel along a single dimension.

  • S (int) – Size of the square DDE kernels (in the spatial Fourier domain) along a single dimension. In the following, S2 = S^2.

  • nshift (int[1, 2]) – Phase-shift along each dimension in Fourier space (expressed in number of samples in the Fourier domain).

  • ktype (string (varargin)) – Selected interpolator (among 'minmax:kb', 'minmax:tuned' and 'kaiser', see nufft.compute_interp_coeffs_extended() for associated kernel documentation). Recommended default is 'minmax:kb'.

Returns

  • G (sparse complex[M, prod(K)]) – Sparse de-gridding interpolation matrix involed in the 2D NUFFT.

  • scale (double[:, :] [N]) – Scaling weights involed in the 2D NUFFT, of size N.

  • ll (int[:]) – Position of the nonzero values in each row of G.

  • v (double[Q2, M]) – convolutions between the DDE kernels (Q2 = (2 * S - 1)^2).

  • V (double[Q2, J2, M]) – values contained in G (J2 = J^2).

nufft.kaiser_bessel(x, J, alpha, kb_m, K_N)

Generalized Kaiser-Bessel function for x in support [-J/2,J/2] shape parameter “alpha” (default 2.34 J) order parameter “kb_m” (default 0) see (A1) in lewitt:90:mdi, JOSA-A, Oct. 1990

Parameters
  • x (double[M, 1]) – Arguments.

  • J (int[:]) – Number of neighbors.

  • alpha – Parameters of the interpolation kernel.

  • kb_m – Parameters of the interpolation kernel.

  • K_N – Parameters of the interpolation kernel.

Returns

  • kb ([M, 1]) – KB function values, if x is numbers, or string for kernel(k,J), if x is a string, or inline function, if x is inline.

  • alpha – Parameters of the interpolation kernel.

  • kb_m – Parameters of the interpolation kernel.

nufft.kaiser_bessel_ft(u, J, alpha, kb_m, d)

Fourier transform of generalized Kaiser-Bessel function, in dimension d.

Parameters
  • u – frequency arguments or ‘string’ to return function string or ‘handle’ to return function handle

  • J (int) – Default 6.

  • alpha (double) – Shape parameter. Default 2.34 * J.

  • kb_m (int) – Order parameter. Default 0.

  • d (int) – Number of dimensions considered. Default 1.

Returns

y – Transform values.

Return type

complex[M, 1]

Note

Original code taken from [Fessler2003], available at https://github.com/JeffFessler/mirt. See (A3) in lewitt:90:mdi, JOSA-A, Oct. 1990.

nufft.nufft_T(N, J, K, tol, alpha, beta, use_true_diric)

Precompute the matrix T = [C' S S' C]\inv used in NUFFT. This can be precomputed, being independent of frequency location.

Parameters
  • N (int[:]) – Signal length.

  • J (int[:]) – Number of neighbors.

  • K (int[:]) – FFT length.

  • tol (double) – Tolerance for smallest eigenvalue.

  • alpha (double[L+1, 1]) – Fourier coefficient vector for scaling.

  • beta (double[:]) – Scale gamma = 2 * pi / K by this for Fourier series.

  • use_true_diric (bool) – Use true Diric function (default is to use sinc approximation).

Returns

T – Precomputed matrix.

Return type

double[J, J]

nufft.nufft_alpha_kb_fit(N, J, K, varargin)

Return the alpha and beta corresponding to LS fit of L components to optimized Kaiser-Bessel scaling factors (m=0, alpha=2.34*J). This is the best method known by J. Fessler currently for choosing alpha.

Parameters
  • N (int[:]) – Signal length.

  • J (int[:]) – Number of neighbors.

  • K (int[:]) – FFT length.

  • varargin (option, [1]) – Specify the Nmid midpoint: floor(N/2) or default (N-1)/2

Returns

  • alphas (complex[:]) – Parameters of the interpolation kernel.

  • beta (double[:]) – Parameters of the interpolation kernel.

Note

Original code taken from [Fessler2003], available at https://github.com/JeffFessler/mirt.

nufft.nufft_coef

NUFFT interpolation coefficient vectors for a given kernel function.

Parameters
  • om (double[M, 1]) – Digital frequency omega in radians.

  • J (int) – Number of neighbors used per frequency location.

  • K (_type_) – FFT size (should be >= N, the signal length).

  • kernel (function handle) – Kernel function handle.

Returns

  • coef (double[J, M]) – Coef. vector for each frequency.

  • arg (double[J, M]) – Kernel argument.

Note

Original code taken from [Fessler2003], available at https://github.com/JeffFessler/mirt.

nufft.nufft_diric(k, N, K, use_true_diric)

“Regular fourier” Dirichlet-function WITHOUT phase nufft_diric(t) = sin(pi N t / K) / ( N * sin(pi t / K) ) \approx sinc(t / (K/N))

Parameters
  • k (double[…]) – Sample locations (unitless real numbers)

  • N (int[…]) – Signal length.

  • K (int[…]) – FFT length.

  • use_true_diric (bool) – Use true Diric function (default is to use sinc approximation).

Returns

f – Corresponding function values.

Return type

double[…]

Warning

Caution: - matlab’s version is different: sin(N * x / 2) / (N * sin(x / 2)) - nufft_diric() is K-periodic for odd N but 2K-periodic for even N.

Note

Original code taken from [Fessler2003], available at https://github.com/JeffFessler/mirt.

nufft.nufft_offset(om, J, K)

Offset for NUFFT.

Parameters
  • om (double[M, 1]) – omega (radians), typically in \([-\pi, \pi)\) (not essential!)

  • J (int[1, d]) – Number of neighbors used for NUFFT interpolation.

  • K (int[1, d]) – FFT size along each domension.

Returns

k0 – Prepared for mod(k0 + [1:J], K) (+ 1 for matlab).

Return type

int[M,1]

Note

Original code taken from [Fessler2003], available at https://github.com/JeffFessler/mirt.

nufft.nufft_r(om, N, J, K, alpha, beta, use_true_diric)

Make NUFFT r vector.

Parameters
  • om (double[M, d]) – Digital frequency omega in radians.

  • N (int[1, d]) – Signal length.

  • J (int[1, d]) – Number of neighbors used per frequency location.

  • K (int[1, d]) – FFT size (should be > N)

  • alpha (complex[0:L]) – Fourier series coefficients of scaling factors.

  • beta (double[:]) – Scale gamma=2pi./K by this in Fourier series, typically is K./N [Fessler2003] or 0.5 (Liu)

  • use_true_diric (bool) – Optional flag for debugging purposes (defaults to false).

Returns

  • rr (double[J, M]) – r vector for each frequency

  • arg (double[J, M]) – Dirac argument for t=0.

Note

Original code taken from [Fessler2003], available at https://github.com/JeffFessler/mirt.

nufft.nufft_scale(Nd, Kd, alpha, beta, Nmid)

Compute scaling factors for NUFFT.

Parameters
  • Nd (int[1, d]) – Input size (along each dimension).

  • Kd (int[1, d]) – Fourier size (along each dimension).

  • alpha (cell{d, 1}) – Parameters of the interpolation kernel.

  • beta (cell{d, 1}) – Parameters of the interpolation kernel.

  • Nmid (int[1, d]) – Midpoint: floor(Nd/2) or default (Nd-1)/2.

Returns

sn – Scaling factors.

Return type

double[Nd]

Note

Original code taken from [Fessler2003], available at https://github.com/JeffFessler/mirt.

nufft.nufft_sinc(x)

J. Fessler’s version of “sinc” function, because matlab’s sinc() is in a toolbox.

Parameters

x (double[:, …]) – Input array.

Returns

y – Evaluation of sinc(x).

Return type

double[:, …]

Note

Original code taken from [Fessler2003], available at https://github.com/JeffFessler/mirt.