FFT of shifted image in Numpy

Question

I want to prove using numpy the theory of Fourier transforms in which translation in space corresponds to a shift in the phase domain (frequency domain remains constant).

I have generated these three 2d-numpy arrays showing translated versions of the same boxes:

I then use numpy to calculate the FFT of each of these images doing:

H1 = np.fft.fft2(im1)
H2 = np.fft.fft2(im2)
H3 = np.fft.fft2(im3)

Although the frequency and phase plots of the transformed arrays look very similar I perform a comparison (MSE) between them to see the differences.

This is the MSE plot for the frequency domain:

And this is the MSE plot for the phase domain:

I'm confused as I was expecting to see identical spectra in the frequency domain but different in the phase domain.

I would appreciate if someone can provide some insight to detect if the problem comes from a theoretical misunderstanding or a problem in the code.

Here is the code to replicate the results shown here.

Answer 1

In the code the OP is comparing the real and the imaginary components of the FFT in each to determine similarity. If a phase shift occurs then the real and imaginary components will also change. What the OP should be comparing is the magnitudes of the complex values in each sample, which should be invariant under such image shifts.

Confirming this; the OP had done the following to compare the FFT results:

plt.bar(['1-1','2-1','3-1'],[np.mean(np.square(H1.imag-H1.imag)),np.mean(np.square(H1.imag-H2.imag)),np.mean(np.square(H1.imag-H3.imag))])

resulting in:

If one does the following instead:

plt.bar(['1-1','2-1','3-1'],[np.mean(np.square(np.abs(H1)-np.abs(H1))),np.mean(np.square(np.abs(H1)-np.abs(H2))),np.mean(np.square(np.abs(H1)-np.abs(H3)))])

We get a result that we expect:

The equivalent magnitude is given by the translation property of the 2D DFT where:

$$x[k-m, l-n] \mapsto X[u,v]e^{-j2\pi (mu/M + nv/N)}$$

This is similar to the time shift property of the 1D DFT (which shows fundamentally what is occurring and can be extended to a 2D DFT):

Given 1D DFT as

$$x[n] \mapsto \sum_{n=0}^{N-1}x[n]W_N^{nk} = X[k]$$

Where $W_N^n=e^{-j2\pi n/N}$

From which:

$$x[n-m] \mapsto \sum_{n=0}^{N-1}x[n-m]W_N^{nk}$$

Substitute $u= n-m$

$$\sum_{u = -m}^{N-1-m}x[u]W_N^{(u+m)k} = \sum_{u = -m}^{N-1-m}x[u]W_N^{uk}W_N^{mk}$$

$$=W_N^{mk}\sum_{u = -m}^{N-1-m}x[u]W_N^{uk}$$

for any k:

$$\sum_{u = -m}^{N-1-m}x[u]W_N^{uk} = \sum_{u=0}^{N-1}x[u]W_N^{uk}$$

Therefore

$$x[n-m] \mapsto W_N^{mk}X[k]$$

And $W_N^{mk}$ as $e^{-j2\pi mk/N}$ only modifies phase, not magnitude

Also the DFT is circularly periodic in the time domain, so any such rotation must done in circular fashion; x[mod(n-m,M)] = x [n]