# 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,:img.shape] = 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

def convertfreqs(x):