This code illustrates how to find the different momentum components in the output of fftrn. The initial array of size (,) is generated with only two sinusoidal modes excited, the real part of and the imaginary part of . The loop at the end shows how to find the different momenta. The excitation only shows up once but the mode and its complex conjugate are redundantly included in the output. Note that no matter what you put in as input the output components , , and will be real.
#define N1 8 // Number of real points in the data - first dimension #define N2 4 // Second dimension void main() { int a,b,ka,kb; float f[N1][N2]; // An array of real numbers float fnyquist[2*N1]; // For storing k_2 = N2/2 frequencies int size[2]={N1,N2}; // Size of data in all dimensions float pi=2.*asin(1.),x1=2.*pi/(float)N1,x2=2.*pi/(float)N2; for(a=0;a<N1;a++) for(b=0;b<N2;b++) f[a][b]=5.*cos(3.*x1*(float)a + x2*(float)b)+sin(2.*x1*(float)a); fftrn((float *)f,(float *)fnyquist,2,size,1); for(a=0;a<N1;a++) { ka = (a<=N1/2 ? a : a-N1); // ka is the frequency stored at position a for(b=0;b<N2/2;b++) { kb = b; // kb is the frequency stored at position b printf("Frequency={%d,%d}, Transform value={%f + i %f}\n",ka,kb,f[a][2*b],f[a][2*b+1]); } printf("Frequency={%d,%d}, Transform value={%f + i %f}\n",ka,N2/2,fnyquist[2*a],fnyquist[2*a+1]); } } OUTPUT: Frequency={0,0}, Transform value={-0.000002 + i 0.000000} Frequency={0,1}, Transform value={0.000000 + i 0.000003} ... Frequency={2,0}, Transform value={0.000002 + i 15.999999} ... Frequency={3,1}, Transform value={79.999992 + i -0.000023} ... Frequency={-2,0}, Transform value={0.000002 + i -15.999999} ... Frequency={-1,1}, Transform value={0.000004 + i 0.000004} Frequency={-1,2}, Transform value={-0.000000 + i 0.000001}