Applying affine transformation (rotate, scale and translate) to align 2 shapes

s(:,1)=[20 40 40 20]';

s(:,2)=[20 20 40 40]';

plot(s(:,1),s(:,2),'.g');title('Reference Shape');

refv=s(:);

d(:,1)=[30 40 30 20];

d(:,2)=[20 30 40 30];

plot(d(:,1),d(:,2),'.g');title('Shape to be aligned');

v=d(:);

[x, y] = COG( v ); %[x,y]=[30,30]

v = Translate( v, x, y); %v-(x,y)

[x, y] = COG( refv ); % here COG is same for both

refv = Translate( refv, x, y);

vSize = Norm2( v ); %square all of them and then sum and find square root.

refvSize = Norm2( refv );

v = Scale( v, refvSize / vSize ); % refvSize / vSize gives how much reference shape is more in scale %than v and so multiply v with that it is less of .

t=v+30;

plot(t(1:4),t(5:8),'.g');

r=GetRotation( refv, v )

re = Rotate( v, r )

re=re+30;

plot(re(1:4),re(5:8),'.g');

the below figure shows that the 2 shapes are aligned.

===================================================

function [x, y] = COG( v )

n = size(v, 1) / 2;

x = sum(v(1:n, 1)) / n;

y = sum(v(n+1:2*n, 1)) / n;

function r = Translate( v, x, y)

[h, w] = size(v);

n = h / 2;

r = zeros(h, w);

r(1:n, 1) = v(1:n, 1) - x;

r(n+1:2*n, 1) = v(n+1:2*n, 1) - y;

function r = Norm2( v )

r = sqrt(sum(v.^2'));

function r = Scale( v, scale )

r = v.* scale;

====================================================

To filter out rotational effects, the following Singular Value Decomposition is used as suggested by Bookstein [1]:

1. Arrange the size and position of aligned shapes x1 and x2 as n x k matrices.

2. Calculate the SVD, UDVT, of  VUT

3. Then the rotation matrix needed to optimally superimpose x1upon x2 is VUT. In the planar case:

function r = GetRotation( refv, v )

                          

n = size( refv, 1 ) / 2;

refs = reshape( refv, n, 2 );

s = reshape( v, n, 2 );

res = refs' * s; %(reference shape ‘ *shape) we get a 2x2 matrix

[U,S,V] = svd( res );

res = V * U';

cos_theta = res(1,1);

sin_theta = res(2,1);

epsilon = 1e-12;

if 1 - cos_theta < epsilon

   r = 0;

elseif abs(cos_theta) < epsilon

   r = pi / 2;

elseif 1 + cos_theta < epsilon

   r = pi;

else

   a_cos = acos( cos_theta );

   a_sin = asin( sin_theta );

   if a_sin < 0

       r = -a_cos;

   else

       r = a_cos;

   end

end
========================================================

 function r = Rotate( v, theta)
rm = [ cos( theta ) sin( theta ); -sin( theta ) cos( theta ) ];
n = size( v, 1) / 2;
s = reshape( v, n, 2 );
for i = 1 : n
    s(i,:) = s(i,:) * rm;
end
r = reshape( s, n*2, 1 );

[1] M. B. Stegmann, Active Appearance Models: Theory, Extensions and Cases, pp. 262, Informatics and Mathematical Modelling, Technical University of Denmark, DTU, 2000.

All the above functions are written by WANG Lei, CG&CV Lab, Hunan University, Changsha

% $Id: ModelAAM.m, v 1.2 2004-4-16 15:17 Lei$