fitAffineTransform3d

 Compute the affine transform that best register two 3D point sets.

   TRANS = fitAffineTransform3d(REF, SRC)
   Returns the affine transform matrix that minimizes the distance between
   the reference point set REF and the point set SRC after transformation.
   Both REF and SRC must by N-by-3 arrays with the same number of rows,
   and the points must be in correspondence.
   The function minimizes the sum of the squared distances:
   CRIT = sum(distancePoints3d(REF, transformPoint3d(PTS, TRANSFO)).^2);

   Example
     N = 50;
     pts = rand(N, 3)*10;
     trans = createRotationOx([5 4 3], pi/4);
     pts2 = transformPoint3d(pts, trans);
     pts3 = pts2 + randn(N, 3)*2;
     fitted = fitAffineTransform3d(pts, pts2);
   

   See also
     transforms3d, transformPoint3d, transformVector3d,
     fitAffineTransform2d, registerPoints3dAffine

0001 function trans = fitAffineTransform3d(ref, src)
0002 % Compute the affine transform that best register two 3D point sets.
0003 %
0004 %   TRANS = fitAffineTransform3d(REF, SRC)
0005 %   Returns the affine transform matrix that minimizes the distance between
0006 %   the reference point set REF and the point set SRC after transformation.
0007 %   Both REF and SRC must by N-by-3 arrays with the same number of rows,
0008 %   and the points must be in correspondence.
0009 %   The function minimizes the sum of the squared distances:
0010 %   CRIT = sum(distancePoints3d(REF, transformPoint3d(PTS, TRANSFO)).^2);
0011 %
0012 %   Example
0013 %     N = 50;
0014 %     pts = rand(N, 3)*10;
0015 %     trans = createRotationOx([5 4 3], pi/4);
0016 %     pts2 = transformPoint3d(pts, trans);
0017 %     pts3 = pts2 + randn(N, 3)*2;
0018 %     fitted = fitAffineTransform3d(pts, pts2);
0019 %
0020 %
0021 %   See also
0022 %     transforms3d, transformPoint3d, transformVector3d,
0023 %     fitAffineTransform2d, registerPoints3dAffine
0024 %
0025 
0026 % ------
0027 % Author: David Legland
0028 % e-mail: david.legland@inrae.fr
0029 % Created: 2009-07-31,    using Matlab 7.7.0.471 (R2008b)
0030 % Copyright 2009 INRAE - Cepia Software Platform.
0031 
0032 
0033 % number of points
0034 N = size(src, 1);
0035 if size(ref, 1) ~= N
0036     error('Requires the same number of points for both arrays');
0037 end
0038 
0039 % main matrix of the problem
0040 tmp = [src(:,1) src(:,2) src(:,3) ones(N,1)];
0041 A = [...
0042     tmp zeros(N, 8) ; ...
0043     zeros(N, 4) tmp zeros(N, 4) ; ...
0044     zeros(N, 8) tmp ];
0045 
0046 % conditions initialisations
0047 B = [ref(:,1) ; ref(:,2) ; ref(:,3)];
0048 
0049 % compute coefficients using least square
0050 coefs = A\B;
0051 
0052 % format to a matrix
0053 trans = [coefs(1:4)' ; coefs(5:8)'; coefs(9:12)'; 0 0 0 1];