Code covered by the BSD License  

Highlights from
Absolute Orientation - Horn's method

5.0
5.0 | 7 ratings Rate this file 109 Downloads (last 30 days) File Size: 18.7 KB File ID: #26186

Absolute Orientation - Horn's method

by

Matt J (view profile)

 

23 Dec 2009 (Updated )

Solves weighted absolute orientation problem using Horn's quaternion-based method.

| Watch this File

File Information
Description

 ABSOR - a tool for solving the absolute orientation problem using Horn's
 quaternion-based method, that is, for finding the rotation, translation, and
 optionally also the scaling, that best maps one collection of point coordinates
 to another in a least squares sense. The function works for both 2D and 3D
 coordinates, and also gives the option of weighting the coordinates non-uniformly.
 The code avoids for-loops to maximize speed.
 
 DESCRIPTION:
 
 As input data, one has
 
   A: a 2xN or 3xN matrix whos columns are the coordinates of N source points.
   B: a 2xN or 3xN matrix whos columns are the coordinates of N target points.
 
 The syntax
 
      [regParams,Bfit,ErrorStats]=absor(A,B)
 
 solves the unweighted/unscaled registration problem
 
            min. sum_i ||R*A(:,i) + t - B(:,i)||^2
 
 for unknown rotation matrix R and unknown translation vector t.
 
 This is a special case of the more general problem
 
            min. sum_i w(i)*||s*R*A(:,i) + t - B(:,i)||^2
 
 where s>=0 is an unknown global scale factor, to be estimated along with R and t,
 and w is a user-supplied length N vector of weights. One can include either
 s or w or both in the problem formulation using the syntax,
 
   [regParams,Bfit,ErrorStats]=absor(A,B,'param1',value1,'param2',value2,...)
 
 with parameter/value pair options
 
   'doScale' - Boolean flag. If TRUE, the global scale factor, s, is included.
                Default=FALSE.
 
   'weights' - the length N vector of weights, w. Default, no weighting.
 
 
 
OUTPUT:
 
 
  regParams: structure output with estimated registration parameters,
 
      regParams.R: The estimated rotation matrix, R
      regParams.t: The estimated translation vector, t
      regParams.s: The estimated scale factor (set to 1 if doScale=false).
      regParams.M: Homogenous coordinate transform matrix [s*R,t;[0 0 ... 1]].
 
      For 3D problems, the structure includes
 
         regParams.q: A unit quaternion [q0 qx qy qz] corresponding to R and
                    signed to satisfy max(q)=max(abs(q))>0
 
      For 2D problems, it includes
 
         regParams.theta: the counter-clockwise rotation angle about the
                    2D origin
 
 
   Bfit: The rotation, translation, and scaling (as applicable) of A that
         best matches B.
 
 
  ErrorStats: structure output with error statistics. In particular,
              defining err(i)=sqrt(w(i))*norm( Bfit(:,i)-B(:,i) ),
              it contains
 
       ErrorStats.errlsq = norm(err)
       ErrorStats.errmax = max(err)
 
 

Acknowledgements

Absolute Orientation inspired this file.

MATLAB release MATLAB 7.9 (R2009b)
Tags for This File   Please login to tag files.
Please login to add a comment or rating.
Comments and Ratings (19)
15 Apr 2015 Matt J

Matt J (view profile)

Hi Fox,

I'm not sure if you've been through the help documentation or not, but the input syntax is described there and in the Description section above, in particular where it says

A: a 2xN or 3xN matrix whos columns are the coordinates of N source points.
B: a 2xN or 3xN matrix whos columns are the coordinates of N target points.

So, put your points to be matched in the columns of matrices A and B. They must be in correspondng order, i.e., B(:,j) must be the rototranslation of A(:,j).

Comment only
13 Apr 2015 Fox

Fox (view profile)

I'd just like to know how to write the co-ordinates out. How would I input the 3 sets of A values and 3 sets of B values?

Thanks

Comment only
12 Apr 2015 Fox

Fox (view profile)

Thanks for the quick response. I've just sent you a message.

I've just realised this might not work for what I require.

Comment only
12 Apr 2015 Matt J

Matt J (view profile)

@Fox,
You haven't shown your call to absor(), so we don't know precisely what you're doing. However, your data example,

A = [1,2,3]'
B = [2,1,2]'

only contains 1 point in each set. You need at least 3 points in each set to do a sensible registration. They are to be placed in the columns of A,B in corresponding order.

Comment only
11 Apr 2015 Fox

Fox (view profile)

I'm pretty new to this type of coding, could someone tell me where to put my two sets of 3D co-ordinates?

If I put for example:
A = [1,2,3]'
B = [2,1,2]'

It says that A and B are used in the function, so they get changed to - in the first line of the code.

It does run, however doesn't give me any values. Just gives like [3x1 double matrix] in writing.

Thanks.

Comment only
06 Nov 2014 ZZZZZ

ZZZZZ (view profile)

 
23 Jun 2014 Zohar Bar-Yehuda  
14 Apr 2014 Matt J

Matt J (view profile)

@Steve,
The points must be placed into the columns of a and b, not the rows.

Comment only
14 Apr 2014 Steve

Steve (view profile)

One question about the translation vector (regParams.t):
I have the following 3 points before rotation:
A1 = (400, 400, 200,)
A2 = (600, 600, 200,)
A3 = (2000, 2000, 2000)
and after rotation:
B1 = (400, 400, 200,)
B2 = (600, 600, 200,)
B3 = (2000, 2000, 2001)
so all point are the same, just B3z is 2001 instead of 2000. Then the computed translation vector using your code is: (162, -76,4001)
What doesn't make sense to me? Because just one single point changed it's z-value from 2000 to 2001 what is just a difference of 1...?
Nevertheless, the result seems to be right, because:
b(:,1) = regParams.R * a(:,1) + regParams.t
is true...

Comment only
14 Apr 2014 Steve

Steve (view profile)

Thanks for this very nice tool!

14 Nov 2013 Cong

Cong (view profile)

 
14 Nov 2013 Matt J

Matt J (view profile)

Hi Cong,

If you're asking how to apply the rotation matrix R to rotate a point P=[x;y;z], you would just use matrix multiplication Pnew=R*P.

Note, however, that ABSOR calculates a translation vector, t, as well (and sometimes also a scaling if you select the doScale option). The full transformation would be Pnew=R*P+t.

Finally, if you are working in homogeneous coordinates P=[x;y;z;1], ABSOR also returns a 4x4 total transformation matrix M. In that case, the transformation can then be done Pnew=M*P.

Comment only
14 Nov 2013 Cong

Cong (view profile)

Thank you to provide this very useful toolbox. I have a question: if I obtain the rotation matrix from your method, then how can I get the coordinates of any point in the new coordinate system (target).

Comment only
18 Apr 2013 monica

monica (view profile)

Thank you very much.

Comment only
16 Apr 2013 Matt J

Matt J (view profile)

Hi Monica. You can find the derivation of the scaling factor in the original paper by Horn, Section 2D

http://people.csail.mit.edu/bkph/papers/Absolute_Orientation.pdf

Comment only
16 Apr 2013 monica

monica (view profile)

Thank you Matt for this submission. I want to know how you got scaling factor sss=summ(right.*(R*left))/summ(left.*left). I know t+RSp-q=0 but i am not able to figure out the formula you have used for getting scaling factor. Could you please elaborate.

06 Jun 2012 Georg Wiora

Very nice work!

10 May 2012 Matt J

Matt J (view profile)

Thanks for the feedback, Georg. I'm not entirely sure why the case of 3 points has been giving you trouble, though. It is true that alternative SVD-based methods and orthogonal matrix methods had to be modified to handle coplanar point data, but I always understood that to be a fairly stable solution.

Comment only
10 May 2012 Georg Stillfried

great work! (seems to work fine and also allows matching of three points)

Updates
25 Dec 2009

1. Important typo fix in help doc. Model is s*R*A+t not s*R*A-t
2. Added version for earlier MATLAB versions lacking bsxfun

22 Jan 2010

Added absorientParams. See README.txt for details

09 Apr 2010

Added tools for purely 2D registration

29 Sep 2010

*The many previous files in the distribution have been consolidated into a single file absor.m.

*New capability:weighted least squares registration.

*Fixed minor bug in the 2D registration routines, occurring when rotation angle was 0.

Contact us