Code covered by the BSD License

### Highlights from Absolute Orientation

5.0
5.0 | 8 ratings Rate this file 13 Downloads (last 30 days) File Size: 2.38 KB File ID: #22422 Version: 1.3

# Absolute Orientation

### Christian Wengert (view profile)

12 Dec 2008 (Updated )

Computes the transformation to register two corresponding 3D point sets.

File Information
Description

[s R T error] = absoluteOrientationQuaternion( A, B, doScale)

Computes the orientation and position (and optionally the uniform scale factor) for the transformation between two corresponding 3D point sets Ai and Bi such as they are related by:

Bi = sR*Ai+T

Implementation is based on the paper by Berthold K.P. Horn:
"Closed-from solution of absolute orientation using unit quaternions"
http://people.csail.mit.edu/bkph/papers/Absolute_Orientation.pdf

Authors:
Dr. Christian Wengert, Dr. Gerald Bianchi

ETH Zurich, Computer Vision Laboratory, Switzerland

Parameters:
A 3xN matrix representing the N 3D points
B 3xN matrix representing the N 3D points
doScale Flag indicating whether to estimate the uniform scale factor as well [default=0]

Return:
s The scale factor
R The 3x3 rotation matrix
T The 3x1 translation vector
err Residual error (optional)

Notes: Minimum 3D point number is N > 4

The residual error is being computed as the sum of the residuals:

for i=1:Npts
d = (B(:,i) - (s*R*A(:,i) + T));
err = err + norm(d);
end

Example:

s=0.7;
R = [0.36 0.48 -0.8 ; -0.8 0.6 0 ; 0.48 0.64 0.6];
T= [45 -78 98]';
X = [ 0.272132 0.538001 0.755920 0.582317;
0.728957 0.089360 0.507490 0.100513;
0.578818 0.779569 0.136677 0.785203];
Y = s*R*X+repmat(T,1,4);

%Compute
[s2 R2 T2 error] = absoluteOrientationQuaternion( X, Y, 1);

error = 0;

Noise = [
-0.23 -0.01 0.03 -0.06;
0.07 -0.09 -0.037 -0.08;
0.009 0.09 -0.056 0.012];

Y = Y+Noise;
%Compute
[s2 R2 T2 error] = absoluteOrientationQuaternion( X, Y, 1);

error = 0.33

Acknowledgements

This file inspired Absolute Orientation Horn's Method.

MATLAB release MATLAB 6.5 (R13)
05 Feb 2016 Cezary Sieluzycki

### Cezary Sieluzycki (view profile)

Hi Christian

Thanks a lot for this very useful script. I have the following question, if I may: When I apply your script on two random \$3 \ times 4\$ matrices, I get larger error when scale s is taken into account vs when it is not. How is it possible? Of course, I run both scenarios on exactly the same matrices.

M1 =

-0.79816 -0.71453 -0.58903 -1.1201
1.0187 1.3514 -0.29375 2.526
-0.13322 -0.22477 -0.84793 1.6555

>> M2 = randn(3, 4)

M2 =

0.30754 -0.17653 -2.3299 0.39135
-1.2571 0.79142 -1.4491 0.45168
-0.86547 -1.332 0.33351 -0.13028

>> [s R T e] = absoluteOrientationQuaternion(M1, M2, 0)

s =

1

R =

0.46424 0.54524 0.69799
-0.83965 0.52176 0.15087
-0.28192 -0.65611 0.70003

T =

-0.78375
-1.6594
-0.049414

e =

4.3475

>> [s R T e] = absoluteOrientationQuaternion(M1, M2, 1)

s =

1.5861

R =

0.46424 0.54524 0.69799
-0.83965 0.52176 0.15087
-0.28192 -0.65611 0.70003

T =

-0.97824
-2.4175
0.21381

e =

5.6218

Comment only
22 Jul 2013 Saber

### Saber (view profile)

Great work.

06 Mar 2013 Daniell Algar

### Daniell Algar (view profile)

An extremely useful tool that works great in my current application!

Big thanks to the author

21 Feb 2013 Sina Abolhoseini

### Sina Abolhoseini (view profile)

thank you man ... it was so helpful!

13 Dec 2011 Georg Stillfried

### Georg Stillfried (view profile)

Just what I was looking for

03 Nov 2011 Peter

### Peter (view profile)

In your result, you compute the residual error as:
err =0;
for i=1:Na
d = (B(:,i) - (s*R*A(:,i) + T));
err = err + norm(d);
end

Wouldn't it be more appropriate to compute the sum of squared errors, because this is what you actually minimize? (so just use norm(d).^2 instead)

A nice addition would be to add the symmetric scale computation as mentioned later in the paper as a third option.

05 Feb 2011 Matt J

### Matt J (view profile)

@Dirk-Jan, the method that you've shown is not a constrained least squares estimation. The transformation matrix that it produces is therefore not gauranteed to be of the form
[sR,T;zeros(1,3), 1].

If you add noise to X and Y, you will see that your method does not produce the same results as the absolute orientation solver, nor will it be in the desired family of transformations.

Comment only
23 Jun 2010 Christian Wengert

### Christian Wengert (view profile)

Here is an interesting paper to the topic:

A Comparison of Four Algorithm s for Estimating tD Rigid Transformations

http://www.homepages.inf.ed.ac.uk/rbf/MY_DAI_OLD_FTP/rp765.pdf

Actually the Horn approach (unit quaternions) and the above mentioned (SVD?) are pretty much equivalent

Comment only
19 Jun 2010 Dirk-Jan Kroon

### Dirk-Jan Kroon (view profile)

What is the advantage of this absolute orientation method?

If I want the least-squares transformation matrix, with your coordinates it is simple :
s=0.7;
R = [0.36 0.48 -0.8 ; -0.8 0.6 0 ; 0.48 0.64 0.6];
T= [45 -78 98]';
X = [ 0.272132 0.538001 0.755920 0.582317;
0.728957 0.089360 0.507490 0.100513;
0.578818 0.779569 0.136677 0.785203];
Y = s*R*X+repmat(T,1,4);

>> X(4,:)=1; Y(4,:)=1;
>> Y/X

ans =

0.2520 0.3360 -0.5600 45.0000
-0.5600 0.4200 0.0000 -78.0000
0.3360 0.4480 0.4200 98.0000
0 0 0 1.0000

Comment only
29 May 2010 Bryan Murawski

### Bryan Murawski (view profile)

Great work, but I believe that there is a mistake in your error metric computation. If you want to compute the sum of the squared error like you're doing (I assume for performance reasons) you should divide by Na to compute the average squared error and then sqrt that quantity (ie. sqrt(err/Na)).

06 Mar 2009 Dave Ligthart

### Dave Ligthart (view profile)

Great implementation!

15 Dec 2008 Christian Wengert

### Christian Wengert (view profile)

I am sorry for the inconvenience caused by the missing crossprod function files. I uploaded the new version and it should be online soon.
Cheers

Comment only
15 Dec 2008 Thomas Pieper

### Thomas Pieper (view profile)

Very good implementation of Horn's paper, but the functions crossprodquaternion and crossprodquaternion2 are still missing.

15 Dec 2008 1.1

Update, included the missing function crossprodQuaternion.
Sorry for that

15 Dec 2008 1.2