Related to another problem I'm having, I was looking into the workings of numpy's rfft2 and irfft2. These are special versions of the FFT routine, in so far that it needs less input; because you require the real-space image to be real you only need to 'fill' half of Fourier space - due to symmetry, that's all the information you need.

Now, in the docs it is mentioned that performing the inverse transform right after the forward transform should return the original array - to within numerical precision. Time to put this to the test.

If I do the following:

import numpy as np

fourier_image = np.random.uniform(-1,1,(100,51))
fourier_recovered = np.fft.rfft2(np.fft.irfft2(fourier_image))

print 'Original power: ', sum(sum(abs(fourier_image)**2.))
print 'Recovered power: ', sum(sum(abs(fourier_recovered)**2.))

... I find the original power and the recovered power differ at the percent level. That seems a rather big difference to me. Am I doing something wrong, or is this just the best we can do?


I naively assumed any (N,N/2+1) Fourier image had the right symmetry, but this goes wrong at the 0- and Nyquist frequency because they don't have a symmetric counterpart. Here, you have to impose the symmetry explicitly.

So, if we generate an image:

N = 100
fourier_image = np.random.uniform(-1,1,(N,1+N/2))

we impose the required symmetry as follows:

fourier_image[1:N/2,0] = np.conj(fourier_image[N:N/2:-1,0])
fourier_image[1:N/2,-1] = np.conj(fourier_image[N:N/2:-1,-1])

And doing the same as in the question above, we find the input - and recovered power:

print 'Original power: ', sum(sum(abs(fourier_image)**2.))
print 'Recovered power: ', sum(sum(abs(fourier_recovered)**2.))

which are indeed the same:

Original power: 1664.3177945
Recovered power: 1664.3177945

For an odd-shaped image this will be slightly different but the idea is the same.

| improve this answer | |
  • 1
    $\begingroup$ Thanks for answering your own question! Please make this the accepted answer once you are able to. $\endgroup$ – Peter K. Oct 9 '15 at 15:06

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.