I have 2D arrays representing surfaces. I have created them by scanning some theoretically flat objects and I ended up with what can be called "imperfections", or else, deviation from the perfectly flat surface.
I hope to get information on the dominant wavelengths of imperfections and their amplitude, so I employed numpy's fft2. Since I am a first timer, I did some digging and tried to play around with some basic examples to get a grip.
Here is one:
import numpy as np import matplotlib.pyplot as plt """ fft2 playground. """ # Initialise an empty array field = np.empty([29, 481]) # Set the wave amplitude amp = 2.5 # Create the synthetic sinusoidal field for n, i in enumerate(np.linspace(0, 2 * np.pi, field.shape)): for m, j in enumerate(np.linspace(0, 2 * np.pi, field.shape)): field[n, m] = amp * np.sin(j) # Perform 2D fourier and shift the result to centre f = np.fft.fft2(field) fshift = np.fft.fftshift(f) # Calculate the magnitude and phase spectra magnitude_spectrum = 20*np.log(np.abs(fshift)) phase_spectrum = np.angle(fshift) # Reconstruct the initial field f_ishift = np.fft.ifftshift(fshift) re_field = np.abs(np.fft.ifft2(f_ishift)) # Plot fig = plt.figure() fig.add_subplot(411) plt.imshow(field, cmap='gray') plt.title('Field'), plt.xticks(), plt.yticks() plt.colorbar() fig.add_subplot(412) plt.imshow(magnitude_spectrum, cmap='gray') plt.title('Magnitude spectrum'), plt.xticks(), plt.yticks() fig.add_subplot(413) plt.imshow(phase_spectrum, cmap='gray') plt.title('Phase spectrum'), plt.xticks(), plt.yticks() fig.add_subplot(414) plt.imshow(re_field, cmap='gray') plt.title('Reconstructed field'), plt.xticks(), plt.yticks() plt.colorbar() plt.show()
Obviously the returned field is not the same as the input. The amplitude is correct but it is mirrored and 90 degrees off-phase. I feel I am missing something but I cannot nail it.
Also, I would expect that the magnitude spectrum would be something like 2 dots, one on either side of the central pixel and nothing else. What I understand is that the central pixel is the zero-frequency and moving to the edges are the higher frequencies, up to the Nyquist.
Worth noting, if I change the dimensions of the input field I get different results on the magnitude and phase spectra.
Could someone point out what am I missing please? How would I get back the initial field from the fft and why the magnitude spectrum is the way it is?