Estimate Euclidean transformation with python

Question

I want to do something similar to what in image analysis would be a standard 'image registration' using features.

I want to find the best transformation that transforms a set of 2D coordinates A in another one B. But I want to add an extra constraint being that the transformation is 'rigid/Euclidean transformation' Meaning that there is no scaling but only translation and rotation. Normally allowing scaling I would do:

 from skimage import io, transform
 destination = array([[1.0,2.0],[1.0,4.0],[3.0,3.0],[3.0,7.0]])
 source = array([[1.2,1.7],[1.1,3.8],[3.1,3.4],[2.6,7.0]])
 T = transform.estimate_transform('similarity',source,destination)

I believeestimate_transform under the hood just solves a least squares problem. But I want to add the constraint of no scaling.

Are there any function in skimage or other packages that solve this? Probably I need to write my own optimization problem with scipy, CVXOPT or cvxpy. Any help to phrase/implement this optimization problem?

EDIT: My implementation thanks to Stefan van der Walt Answer

from matplotlib.pylab import *
from scipy.optimize import *

def obj_fun(pars,x,src):
    theta, tx, ty = pars
    H = array([[cos(theta), -sin(theta), tx],\
         [sin(theta), cos(theta), ty],
         [0,0,1]])
    src1 = c_[src,ones(src.shape[0])]
    return sum( (x - src1.dot(H.T)[:,:2])**2 )

def apply_transform(pars, src):
    theta, tx, ty = pars
    H = array([[cos(theta), -sin(theta), tx],\
         [sin(theta), cos(theta), ty],
         [0,0,1]])
    src1 = c_[src,ones(src.shape[0])]
    return src1.dot(H.T)[:,:2]

res = minimize(obj_fun,[0,0,0],args=(dst,src), method='Nelder-Mead')

score 4 · Accepted Answer · answered Oct 26 '14 at 16:13

With that extra constraint you are no longer solving a linear least squares problem, so you'll have to use one of SciPy's minimization functions. The inner part of your minimization would set up a matrix H:

H = np.array([[np.cos(theta), -np.sin(theta), tx],
              [np.sin(theta), np.cos(theta), ty],
              [0, 0, 1]])

Then, you would compute the distance

|x_target - H.dot(x_source)|

for all data-points and sum the errors. Now, you have a cost function that you can send to the minimization function. You probably will also want to make use of RANSAC, which is available as skimage.measure.ransac, to reject outliers.

I left out some detail such as homogeneous coordinates, etc., so kudos for getting it to work so easily! — Stefan van der Walt, Oct 26 '14 at 22:52

score 1 · Answer 2 · answered Apr 14 '18 at 16:33

skimage now provides native support in the transform module.

http://scikit-image.org/docs/dev/api/skimage.transform.html#skimage.transform.estimate_transform

Somewhat easier than OpenCV I find. There is an extensive set of functions which covers all use cases.

Estimate Euclidean transformation with python

2 Answers2