# Sub pixel shift causes weird results

I am trying to shift an image on a subpixel level, with the method described here: http://www.mathworks.com/matlabcentral/fileexchange/18401-efficient-subpixel-image-registration-by-cross-correlation/content/html/efficient_subpixel_registration.html

The results are however disappointing. Since I am actually programming this algorithm for the rotation of images, where I resample inbetween pixels, this problem keeps on bugging me: http://spin.ecn.purdue.edu/fmri/PDFLibrary/TongR_MRM_1999_41_253_256.pdf

Can somebody tell me how these strange artifacts can be resolved? According to multiple scientific articles I found, it should be possible to resample on other grids with satisfactory results. My minimal code example (in python) is below, along with the resulting images, which are not satisfactory at all.

Does someone know what causes this effect, and how it can be resolved to reproduce the results others have?

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.image as mpimg

N = 2048

f = np.zeros((N,N))
f[:img.shape[0],:img.shape[1]] = img[:,:,3]

def fft(x):
return np.fft.fft2(x)

def ifft(x):
return np.fft.ifft2(x)

exp = np.exp
pi = np.pi
j = 1j

def s(r):
plt.figure()
#plt.imshow(np.real(r), vmin=0.0, vmax=1.0, interpolation='none')
plt.imshow(np.absolute(r), interpolation='none', vmin=0.0, vmax=1.0)
plt.colorbar()

x = np.tile(np.arange(N),(N,1)).T
y = x.T
result = ifft( fft(f)* np.exp(1j*2*np.pi*(-0.5*x/N-0.5*y/N)) )
s(result)
s(f)
plt.show()

• Did you fftshift? – Mikhail Dec 15 '14 at 5:38
• Where and/or when should I use an fftshift? The formula to shift half a pixel is as following, right? result = ifft( fft(f)* np.exp(1j*2*np.pi*(-0.5*x/N-0.5*y/N)) ) – 317070 Dec 15 '14 at 14:24
• Haven't fully grokked your application, but did you mean result = ifft( fft(f) .* np.exp(1j*2*np.pi*(-0.5*x/N-0.5*y/N)) ) instead (element-wise mult)? – rickhg12hs Dec 18 '14 at 11:02
• No, the "*" is already an element-wise multiplication in python. – 317070 Dec 19 '14 at 12:22

The solution lies in the use of the x's. You need to convert the x's you use to lower frequencies, which can be either positive or negative. So only use x's between -N/2 and N/2 and the problem is solved.

def convertfreqs(x):
return np.mod((x + N/2),N) - N/2


This works because it is the higher frequencies which cause the signal to change heavily inbetween pixels. Higher frequencies means higher x's.