def gaussian2D(window):
"""Real 2D Gaussian interpolation for sub pixel shift"""
#ref on paper
w = np.ones((3, 3))*(1./9)
rhs = np.zeros(6)
M = np.zeros((6,6))
for i in [-1, 0, 1]:
for j in [-1, 0, 1]:
rhs = rhs + np.array([i*w[i+1, j+1]*np.log(np.abs(window[i+1, j+1])),
j*w[i+1, j+1]*np.log(np.abs(window[i+1, j+1])),
i*j*w[i+1, j+1]*np.log(np.abs(window[i+1, j+1])),
i*i*w[i+1, j+1]*np.log(np.abs(window[i+1, j+1])),
j*j*w[i+1, j+1]*np.log(np.abs(window[i+1, j+1])),
w[i+1, j+1]*np.log(np.abs(window[i+1, j+1]))], dtype='float')
M = M + w[i+1, j+1]*np.array([[ i*i, i*j, i*i*j, i*i*i, i*j*j, i],
[ i*j, j*j, i*j*j, i*i*j, j*j*j, j],
[i*i*j, i*j*j, i*i*j*j, i*i*i*j, i*j*j*j, i*j],
[i*i*i, i*i*j, i*i*i*j, i*i*i*i, i*i*j*j, i*i],
[i*j*j, j*j*j, i*j*j*j, i*i*j*j, j*j*j*j, j*j],
[ i, j, i*j, i*i, j*j, 1]], dtype='float')
solution = nl.solve(M, rhs)
dx = ( solution[2]*solution[1] - 2.0*solution[0]*solution[4])/ \
(4.0*solution[3]*solution[4] - solution[2]*solution[2])
dy = ( solution[2]*solution[0] - 2.0*solution[1]*solution[3])/ \
(4.0*solution[3]*solution[4] - solution[2]*solution[2])
return dx, dy
评论列表
文章目录