Path: news.mathworks.com!not-for-mail
From: <HIDDEN>
Newsgroups: comp.soft-sys.matlab
Subject: Re: Directional Cosine Matrix
Date: Fri, 20 Jun 2008 06:54:02 +0000 (UTC)
Organization: The MathWorks, Inc.
Lines: 71
Message-ID: <g3fk6a$co7$1@fred.mathworks.com>
References: <g3do4d$667$1@fred.mathworks.com>
Reply-To: <HIDDEN>
NNTP-Posting-Host: webapp-02-blr.mathworks.com
Content-Type: text/plain; charset="ISO-8859-1"
Content-Transfer-Encoding: 8bit
X-Trace: fred.mathworks.com 1213944842 13063 172.30.248.37 (20 Jun 2008 06:54:02 GMT)
X-Complaints-To: news@mathworks.com
NNTP-Posting-Date: Fri, 20 Jun 2008 06:54:02 +0000 (UTC)
X-Newsreader: MATLAB Central Newsreader 1187260
Xref: news.mathworks.com comp.soft-sys.matlab:474901



"Chia Khai Low" <kaiserwulf@gmail.com> wrote in message <g3do4d$667
$1@fred.mathworks.com>...
> I have two vectors in different excel files. I have ........

  It is not entirely evident from your description what you asking, Chia.  For 
example, the statement, "These 3d vectors are in cartesian coordinates" is not 
completely clear to me.  I am guessing what you mean is that you have two n 
by 3 arrays, A and B, in which the corresponding rows in them are x,y,z 
(Cartesian) coordinates in two different coordinate systems, one presumably 
rotated with respect to the other about a common origin.  It is this 3 x 3 
rotation matrix that you are apparently seeking as your "directional cosine 
matrix".  Have I guessed correctly?

  You didn't mention a possible translation in aligning the two sets of points.  
If that were necessary, you could subtract the mean values of each so that the 
sum of the coordinates in A and B would then be zero.  The necessary 
translation would then be the vector difference between these two mean value 
vectors.  I have therefore assumed below that the sums of each coordinate 
column in A and B are zero: sum(A,1) = sum(B,1) = [0,0,0].

  I do not regard the method I am about to give as entirely satisfactory.  It has 
two fundamental weaknesses.  The first weakness is that if A and B are not 
accurate rotations of one another, the method may not give a perfect least 
squares optimum rotation matrix for them.  In fact, if they vary too far from 
being rotations of one another, the method could conceivably fail altogether.

  The second weakness has to do with the eigenvalues (singular values in this 
case) that are computed below.  If it were to happen that two or three of them 
were essentially equal, the method again falls apart.

  If we (temporarily) ignore these difficulties, do this:

 [U1,S1,V1] = svd(A'*A,0); % Find eigenvalues and eigenvectors
 [U2,S2,V2] = svd(B'*B,0); % for A'*A and B'*B
 mn = inf;
 II = [-1,-1,-1;-1,-1,1;-1,1,-1;-1,1,1;1,-1,-1;1,-1,1;1,1,-1;1,1,1];
 for k = 1:8
  Rt = V1*diag(II(k,:))*V2'; % One of these eight
  X = B-A*Rt; %  possibilities should be the right one
  tt = norm(X(:)); % Choose the closest fit to find it
  if tt < mn, mn = tt; R = Rt; end
 end

Then R is the desired rotation matrix from the n points in A to those in B.  We 
should have essential equality, B = A*R, provided the points in B are actually 
rotated versions of those in A.

  The reasoning here is that the eigenvalues of A'*A and B'*B should 
(hopefully) be the same, and then the matrix V1*V2', with one of eight 
possible combinations of sign changes in the columns of V1, should be the 
required transformation from A to B.  We have only to search through the 
eight possibilities to find it.

  As stated above, if two or three of the eigenvalues in S1 or S2 should 
happen to be equal, the logic in the previous paragraph fails.  There would 
have to be a search through infinitely many possible matrices between V1 and 
V2', which of course is impossible.  The method would break down in this 
case.

  I should mention here that if some kind of reflection has occurred between 
the points of A and those of B, this method will come up with an R matrix 
whose determinant is -1, that is, it is a reflection, not a rotation, matrix.

  I am sure there exists a superior method of solving your problem that would 
avoid the two aforementioned weaknesses, but I have not been able, up to 
this time, to come up with a practical one.  If I manage to think of a better 
method, I will inform you, or perhaps some other kind soul will present one 
here.

Roger Stafford