KbNufftAdjoint

class KbNufftAdjoint(im_size, grid_size=None, numpoints=6, n_shift=None, table_oversamp=1024, kbwidth=2.34, order=0.0, dtype=None, device=None)[source]

Non-uniform FFT adjoint layer.

This object interpolates off-grid Fourier data to on-grid locations using a Kaiser-Bessel kernel prior to inverse DFT. Mathematically, in one dimension it estimates \(x_n, n \in [0, ..., N-1]\) from a off-grid signal \(Y_m, m \in [0, ..., M-1]\) where the off-grid frequency locations are \(\omega_m\). To perform the estimate, this layer applies

\[X_k = \sum_{j=1}^J \sum_{m=0}^{M-1} Y_m u_j(\omega_m) \mathbb{1}_{\{\{k_m+j\}_K=k\}}, \]
\[x_n = s_n^* \sum_{k=0}^{K-1} X_k e^{i \gamma k n} \]

In the first step, \(u\), the Kaiser-Bessel kernel, is used to estimate \(Y\) at on-grid frequency locations from locations at \(\omega\). \(k_m\) is the index of the root offset of nearest samples of \(X\) to frequency location \(\omega_m\), \(\mathbb{1}\) is an indicator function, and \(J\) is the number of nearest neighbors to use from \(X_k, k \in [0, ..., K-1]\).

In the second step, an image-domain signal \(x_n\) is estimated from a gridded, oversampled frequency-domain signal, \(X_k\) by applying the inverse FFT, after which the complex conjugate scaling coefficients \(s_n\) are multiplied. The oversampling coefficient is \(\gamma = 2\pi / K, K >= N\). Multiple dimensions are handled separably. For a detailed description see Nonuniform fast Fourier transforms using min-max interpolation (JA Fessler and BP Sutton).

Note

This function is not the inverse of KbNufft; it is the adjoint.

When called, the parameters of this class define properties of the kernel and how the interpolation is applied.

  • im_size is the size of the base image, analagous to \(N\).

  • grid_size is the size of the grid after adjoint interpolation, analogous to \(K\). To reduce errors, NUFFT operations are done on an oversampled grid to reduce interpolation distances. This will typically be 1.25 to 2 times im_size.

  • numpoints is the number of nearest neighbors to use for interpolation, i.e., \(J\).

  • n_shift is the FFT shift distance, typically im_size // 2.

Parameters
  • im_size (Sequence[int]) – Size of image with length being the number of dimensions.

  • grid_size (Optional[Sequence[int]]) – Size of grid to use for interpolation, typically 1.25 to 2 times im_size. Default: 2 * im_size

  • numpoints (Union[int, Sequence[int]]) – Number of neighbors to use for interpolation in each dimension.

  • n_shift (Optional[Sequence[int]]) – Size for fftshift. Default: im_size // 2.

  • table_oversamp (Union[int, Sequence[int]]) – Table oversampling factor.

  • kbwidth (float) – Size of Kaiser-Bessel kernel.

  • order (Union[float, Sequence[float]]) – Order of Kaiser-Bessel kernel.

  • dtype (Optional[dtype]) – Data type for tensor buffers. Default: torch.get_default_dtype()

  • device (Optional[device]) – Which device to create tensors on. Default: torch.device('cpu')

Examples

>>> data = torch.randn(1, 1, 12) + 1j * torch.randn(1, 1, 12)
>>> omega = torch.rand(2, 12) * 2 * np.pi - np.pi
>>> adjkb_ob = tkbn.KbNufftAdjoint(im_size=(8, 8))
>>> image = adjkb_ob(data, omega)
forward(data, omega, interp_mats=None, smaps=None, norm=None)[source]

Interpolate from scattered data to gridded data and then iFFT.

Input tensors should be of shape (N, C) + klength, where N is the batch size and C is the number of sensitivity coils. omega, the k-space trajectory, should be of size (len(im_size), klength), where klength is the length of the k-space trajectory.

If your tensors are real, ensure that 2 is the size of the last dimension.

Parameters
  • data (Tensor) – Data to be gridded and then inverse FFT’d.

  • omega (Tensor) – k-space trajectory (in radians/voxel).

  • interp_mats (Optional[Tuple[Tensor, Tensor]]) – 2-tuple of real, imaginary sparse matrices to use for sparse matrix NUFFT interpolation (overrides default table interpolation).

  • smaps (Optional[Tensor]) – Sensitivity maps. If input, these will be multiplied before the forward NUFFT.

  • norm (Optional[str]) – Whether to apply normalization with the FFT operation. Options are "ortho" or None.

Return type

Tensor

Returns

data transformed to the image domain.