KbNufft
- class KbNufft(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 layer.
This object applies the FFT and interpolates a grid of Fourier data to off-grid locations using a Kaiser-Bessel kernel. Mathematically, in one dimension it estimates \(Y_m, m \in [0, ..., M-1]\) at frequency locations \(\omega_m\) from \(X_k, k \in [0, ..., K-1]\), the oversampled DFT of \(x_n, n \in [0, ..., N-1]\). To perform the estimate, this layer applies
\[X_k = \sum_{n=0}^{N-1} s_n x_n e^{-i \gamma k n} \]\[Y_m = \sum_{j=1}^J X_{\{k_m+j\}_K} u^*_j(\omega_m) \]In the first step, an image-domain signal \(x_n\) is converted to a gridded, oversampled frequency-domain signal, \(X_k\). The scaling coefficeints \(s_n\) are multiplied to precompensate for NUFFT interpolation errors. The oversampling coefficient is \(\gamma = 2\pi / K, K >= N\).
In the second step, \(u\), the Kaiser-Bessel kernel, is used to estimate \(X_k\) at off-grid frequency locations \(\omega_m\). \(k_m\) is the index of the root offset of nearest samples of \(X\) to frequency location \(\omega_m\), and \(J\) is the number of nearest neighbors to use from \(X_k\). Multiple dimensions are handled separably. For a detailed description see Nonuniform fast Fourier transforms using min-max interpolation (JA Fessler and BP Sutton).
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 forward FFT, 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 timesim_size
.numpoints
is the number of nearest neighbors to use for interpolation, i.e., \(J\).n_shift
is the FFT shift distance, typicallyim_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 timesim_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 forfftshift
. 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
>>> image = torch.randn(1, 1, 8, 8) + 1j * torch.randn(1, 1, 8, 8) >>> omega = torch.rand(2, 12) * 2 * np.pi - np.pi >>> kb_ob = tkbn.KbNufft(im_size=(8, 8)) >>> data = kb_ob(image, omega)
- forward(image, omega, interp_mats=None, smaps=None, norm=None)[source]
Apply FFT and interpolate from gridded data to scattered data.
Input tensors should be of shape
(N, C) + im_size
, whereN
is the batch size andC
is the number of sensitivity coils.omega
, the k-space trajectory, should be of size(len(grid_size), klength)
or(N, len(grid_size), klength)
, whereklength
is the length of the k-space trajectory.Note
If the batch dimension is included in
omega
, the interpolator will parallelize over the batch dimension. This is efficient for many small trajectories that might occur in dynamic imaging settings.If your tensors are real, ensure that 2 is the size of the last dimension.
- Parameters
image (
Tensor
) – Object to calculate off-grid Fourier samples from.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"
orNone
.
- Return type
- Returns
image
calculated at Fourier frequencies specified byomega
.