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$
