Simpe example if using NFFT3 library

Question

I am trying to use the NFFT3 library to perform a NFFT from non-equidistant spatial data to the corresponding equidistant Fourier coefficients. I'm implementing my project in C, in a Linux environment.

To test the use of the library I am implementing a transform from nodes [-0.1, 0.1] with values [1.0, 1.0] to N=10 Fourier coefficients. The expected result should be 2.0*cos(pi*nu/5), where nu is the frequency.

I have some experience with FFTW, but I can't seem to figure out how to get NFFT3 to work. Any advice on how to improve my C-code to get the expected results using NFFT3 would be very much appreciated.

The C-code I have so far is as follows:

#include <stdio.h>
#include <math.h>
#include <string.h>
#include <stdlib.h>
#include <complex.h>
#include "nfft3.h"

static void test_nfft_1d(void)
{
  /** problem size */
  int N=10; // # of Fourier coefficient sample points - equispaced
  int M=2;  // # of non-equispaced nodes

  /** init one dimensional plan */
  nfft_plan p;
  nfft_init_1d(&p, N, M);

  /** init non-equidistant nodes */
  p.x[0] = - 0.1;
  p.x[1] = 0.1;

  /** precompute psi */
  if(p.flags & PRE_ONE_PSI)
      nfft_precompute_one_psi(&p);

  /** init spacial values */
  p.f[0] = 1.0;
  p.f[1] = 1.0;

  /** perform adjoint nfft */
  nfft_adjoint(&p);

  /** show the result */
  for (int idx_k = 0; idx_k < N; idx_k++) {
    printf("%#Zg\n",p.f_hat[idx_k]);
  }

  /** finalise the one dimensional plan */
  nfft_finalize(&p);
}

int main(void)
{
  // computing a one dimensional nfft
  test_nfft_1d();

  return 1;
}

As output I get:

Where I expect to get:

The following Matlab code produces the expected results:

%% Define Signal Parameters
N = 10; % # of Fourier coefficient sample points - equispaced
nu_range = 40; 
nu = linspace(-nu_range,nu_range,N);

%% Define Non-Equidistant Spacial Signal
x_nu = [-0.1,0.1]; % non-equidistant nodes
sig_nu = [1.0,1.0]; % signal at given nodes
 
%% Compute the NFFT
NFFT_mat = nufft(double(sig_nu),double(x_nu),kx);
NFFT_an = ndft(double(sig_nu),double(x_nu),kx);
FFT_an = 2.0*cos(pi*kx/5);

% NFFT results
figure(1)
hold on;
plot(kx,real(NFFT_mat),"--");
plot(kx,real(NFFT_an),":");
plot(kx,real(FFT_an),".-.");
legend("Matlab NFFT","Analytic NFFT","Rect function FT");
xlabel("frequency");
ylabel("real part of Fourier coefficient");
title("NFFT Results [Real Part]")
hold off;

%% Analytical expression for NDFT
function [out] = ndft(X,t,f)
  out = zeros(size(f));
  for j = 1:length(X)
    out(:) = out(:) + X(j) * exp(-1j*2*pi*t(j)*f(:));
  end
end

Simpe example if using NFFT3 library

0 Answers0