Help with gridding 1D case

[schmidtijoe]

Hello,
Thanks for the nice package!
I have a bit of a different use-case. I got a sequence with various non-uniform readout line samplings but matching phase encodes. If i suppose the phase encode errors are negligible, i basically just want to grid all of my read lines in 1D. Since i need another recon method afterwards i was only looking at the KbInterpAdjoint method of the interpolation module.

• I got N batches and C coils, (data dims [N, C, k-length), trajectory dims [N, 1, k-length]).
• set the grid and image dimension to 1d tuples
• set the k-space trajectory per batch to radians / voxel.

The interpolation seems to work, however, it appears that the output pushed the k-space center towards the edges and i appear to get edge artifacts for a partial fourier region which is not acquired (i guess this is acceptable since using interpolation to fill the non acquired data with 0s is not the intended use). When i am simply using a torch.fft.ifftshift it seems to do the job, i was just wondering if i am missing something or should fftshift the input prior to application of the interpolation. Since its not mentioned in its respective doc.

Any hints are welcome !
All the best!

[mmuckley]

Hello @schmidtijoe, the fftshift is normally built into the phase of the output data from the interpolation operation:

torchkbnufft/torchkbnufft/_nufft/interp.py

Lines 679 to 686 in 20c45a2

	# phase for fftshift
	data = (
	data
	* imag_exp(
	torch.sum(omega * n_shift.unsqueeze(-1), dim=-2, keepdim=True),
	return_complex=True,
	).conj()
	)

If it helps, you can see how the interpolation operation is wrapped with the inverse FFT here:

torchkbnufft/torchkbnufft/functional/nufft.py

Lines 194 to 260 in 20c45a2

	def kb_table_nufft_adjoint(
	data: Tensor,
	scaling_coef: Tensor,
	im_size: Tensor,
	grid_size: Tensor,
	omega: Tensor,
	tables: List[Tensor],
	n_shift: Tensor,
	numpoints: Tensor,
	table_oversamp: Tensor,
	offsets: Tensor,
	norm: Optional[str] = None,
	) -> Tensor:
	"""Kaiser-Bessel NUFFT adjoint with table interpolation.
	See :py:class:`~torchkbnufft.KbNufftAdjoint` for an overall description of
	the adjoint NUFFT.
	Args:
	data: Scattered data to be iNUFFT'd to an image.
	scaling_coef: Image-domain coefficients to compensate for interpolation
	errors.
	im_size: Size of image with length being the number of dimensions.
	grid_size: Size of grid to use for interpolation, typically 1.25 to 2
	times ``im_size``.
	omega: k-space trajectory (in radians/voxel).
	tables: Interpolation tables (one table for each dimension).
	n_shift: Size for fftshift, usually ``im_size // 2``.
	numpoints: Number of neighbors to use for interpolation.
	table_oversamp: Table oversampling factor.
	offsets: A list of offsets, looping over all possible combinations of
	`numpoints`.
	norm: Whether to apply normalization with the FFT operation.
	Options are ``"ortho"`` or ``None``.
	Returns:
	``data`` transformed to an image.
	"""
	is_complex = True
	if not data.is_complex():
	if not data.shape[-1] == 2:
	raise ValueError("For real inputs, last dimension must be size 2.")
	is_complex = False
	data = torch.view_as_complex(data)
	image = ifft_and_scale(
	image=kb_table_interp_adjoint(
	data=data,
	omega=omega,
	tables=tables,
	n_shift=n_shift,
	numpoints=numpoints,
	table_oversamp=table_oversamp,
	offsets=offsets,
	grid_size=grid_size,
	),
	scaling_coef=scaling_coef,
	im_size=im_size,
	grid_size=grid_size,
	norm=norm,
	)
	if is_complex is False:
	image = torch.view_as_real(image)
	return image

So once you have your coordinates and apply interpolation, the relationship between the image and the data should be a simple FFT.

Let me know if this helps. Also might be helpful to see an image of what you mean by "pushed the k-space center towards the edges" to help any further debugging.

[schmidtijoe]

Yeah thanks, with that phrase i meant exactly that it looks like an fft shift applied along my 1 dimensional data. Which was also what i did to "mitigate" the problem. I was just wondering if i did something wrong since it wasn't obvious to me that the output should be related via the shift. Thanks for the clarification!