Thread Subject:

Defining a Coordinate System

Hey everybody:

I have a question pertaining to what should be a fairly simple linear algebra problem. Goes like this:

I am given three points in a "standard" cartesian coordinate system, A = (ax,ay,az), B = (bx,by,bz), and C = (cx,cy,cz). I would like to define a plane with these three points, and establish a new orthogonal coordinate system based on this plane, where the x- and y-axes are on the plane and the z-axis is normal to it. I want to set A to be the origin of this system (0,0,0), B to have only an x component (bx,0,0) and C to have x and y components (cx, cy, 0). What would be the best way to go about making this change of basis in Matlab?

I am new to this program, so your patience is very appreciated.

am

"Armen" wrote in message <itt2rk$bv7$1@newscl01ah.mathworks.com>...
> Hey everybody:
>
> I have a question pertaining to what should be a fairly simple linear algebra problem. Goes like this:
>
> I am given three points in a "standard" cartesian coordinate system, A = (ax,ay,az), B = (bx,by,bz), and C = (cx,cy,cz). I would like to define a plane with these three points, and establish a new orthogonal coordinate system based on this plane, where the x- and y-axes are on the plane and the z-axis is normal to it. I want to set A to be the origin of this system (0,0,0), B to have only an x component (bx,0,0) and C to have x and y components (cx, cy, 0). What would be the best way to go about making this change of basis in Matlab?
>
> I am new to this program, so your patience is very appreciated.
>
> am
- - - - - - - - -
  One way is to think of this in terms of vector analysis. First compute three unit mutually orthogonal vectors based at point A and of the kind you have described.

 v1 = B-A; v1 = v1/norm(v1);
 v3 = cross(v1,C-A); v3 = v3/norm(v3);
 v2 = cross(v3,v1); % v2 is already a unit vector

  Now for any arbitrary point P = [x,y,z] we always have:

 P-A = dot(P-A,v1)*v1 + dot(P-A,v2)*v2 + dot(P-A,v3)*v3

This means that [dot(P-A,v1),dot(P-A,v2),dot(P-A,v3)] are the coordinates of point P with respect to the v1, v2, v3 axes, and you can rewrite this as a set of linear transformation equations from the original coordinates to the new set. Finding the inverse of this transformation would allow you to express the original coordinates in terms of the new ones.

Roger Stafford

"Roger Stafford" wrote in message <itt5t5$l3n$1@newscl01ah.mathworks.com>...
> "Armen" wrote in message <itt2rk$bv7$1@newscl01ah.mathworks.com>...
> > Hey everybody:
> >
> > I have a question pertaining to what should be a fairly simple linear algebra problem. Goes like this:
> >
> > I am given three points in a "standard" cartesian coordinate system, A = (ax,ay,az), B = (bx,by,bz), and C = (cx,cy,cz). I would like to define a plane with these three points, and establish a new orthogonal coordinate system based on this plane, where the x- and y-axes are on the plane and the z-axis is normal to it. I want to set A to be the origin of this system (0,0,0), B to have only an x component (bx,0,0) and C to have x and y components (cx, cy, 0). What would be the best way to go about making this change of basis in Matlab?
> >
> > I am new to this program, so your patience is very appreciated.
> >
> > am
> - - - - - - - - -
> One way is to think of this in terms of vector analysis. First compute three unit mutually orthogonal vectors based at point A and of the kind you have described.
>
> v1 = B-A; v1 = v1/norm(v1);
> v3 = cross(v1,C-A); v3 = v3/norm(v3);
> v2 = cross(v3,v1); % v2 is already a unit vector
>
> Now for any arbitrary point P = [x,y,z] we always have:
>
> P-A = dot(P-A,v1)*v1 + dot(P-A,v2)*v2 + dot(P-A,v3)*v3
>
> This means that [dot(P-A,v1),dot(P-A,v2),dot(P-A,v3)] are the coordinates of point P with respect to the v1, v2, v3 axes, and you can rewrite this as a set of linear transformation equations from the original coordinates to the new set. Finding the inverse of this transformation would allow you to express the original coordinates in terms of the new ones.
>
> Roger Stafford

---

Hey Roger:

Thanks for your reply - I've set up the code you posted and it runs, but I'm encountering a different problem. As a check I've plugged in the vectors B and C, expecting returns of (bx,0,0) and (cx,cy,0), but instead it's just returning the original vectors I put in. At this point, my problem certainly doesn't have anything to do with Matlab but rather my linear algebra is very rusty (funny how quickly you forget things you learn in school). Any suggestions would be appreciated.

am

"Armen" wrote in message <itvnh1$323$1@newscl01ah.mathworks.com>...
> ..... I've set up the code you posted and it runs, but I'm encountering a different problem. As a check I've plugged in the vectors B and C, expecting returns of (bx,0,0) and (cx,cy,0), but instead it's just returning the original vectors I put in.
- - - - - - - - - -
  I ran the test below and the results look correct to me, Armen. The quantities TA, TB, and TC have the transformed coordinates of points A, B, and C. Note that TA is at the new origin with all zeros, TB has only one non-zero component, and TC has only two, as should be the case.

  Notice that the first component of TB will not be bx = B(1). It must be the norm of the vector B-A, that is, the distance between A and B because you wanted A to be relocated at the new origin and the axes rotated so that B-A points along the new x-axis. A similar statement holds for TC - it is not [C(1),C(2),0]. It has also undergone a rotation and translation.

  The general formula

 TP = [dot(P-A,v1),dot(P-A,v2),dot(P-A,v3)]

gives you the transform of any arbitrary point P = [x,y,z].

  You can test that the set v1,v2,v3 is orthonormal by seeing that

 [v1;v2;v3]*[v1;v2;v3]'

produces the identity matrix (except for round off errors.)

% The code
clear
format long

% Choose three random points, A, B, C
A = randn(1,3), B = randn(1,3), C = randn(1,3)

% Compute the three new axes
v1 = B-A; v1 = v1/norm(v1);
v3 = cross(v1,C-A); v3 = v3/norm(v3);
v2 = cross(v3,v1); % v2 is already a unit vector

% Get the transforms of A, B, and C
P = A;
TA = [dot(P-A,v1),dot(P-A,v2),dot(P-A,v3)]
P = B;
TB = [dot(P-A,v1),dot(P-A,v2),dot(P-A,v3)]
P = C;
TC = [dot(P-A,v1),dot(P-A,v2),dot(P-A,v3)]

% Results
A = 0.63940564208852 0.87421289486361 1.75240173032956
B = -0.32005082643214 -0.13741380814487 0.61576962808672
C = 0.97789406984520 -1.11534771220514 -0.55002144880449

TA = 0 0 0
TB = 1.79885459030876 0 0.00000000000000
TC = 2.39315386990053 1.90968397273444 -0.00000000000000

Roger Stafford

"Roger Stafford" wrote in message <iu0979$ah$1@newscl01ah.mathworks.com>...
> "Armen" wrote in message <itvnh1$323$1@newscl01ah.mathworks.com>...
> > ..... I've set up the code you posted and it runs, but I'm encountering a different problem. As a check I've plugged in the vectors B and C, expecting returns of (bx,0,0) and (cx,cy,0), but instead it's just returning the original vectors I put in.
> - - - - - - - - - -
> I ran the test below and the results look correct to me, Armen. The quantities TA, TB, and TC have the transformed coordinates of points A, B, and C. Note that TA is at the new origin with all zeros, TB has only one non-zero component, and TC has only two, as should be the case.
>
> Notice that the first component of TB will not be bx = B(1). It must be the norm of the vector B-A, that is, the distance between A and B because you wanted A to be relocated at the new origin and the axes rotated so that B-A points along the new x-axis. A similar statement holds for TC - it is not [C(1),C(2),0]. It has also undergone a rotation and translation.
>
> The general formula
>
> TP = [dot(P-A,v1),dot(P-A,v2),dot(P-A,v3)]
>
> gives you the transform of any arbitrary point P = [x,y,z].
>
> You can test that the set v1,v2,v3 is orthonormal by seeing that
>
> [v1;v2;v3]*[v1;v2;v3]'
>
> produces the identity matrix (except for round off errors.)
>
> % The code
> clear
> format long
>
> % Choose three random points, A, B, C
> A = randn(1,3), B = randn(1,3), C = randn(1,3)
>
> % Compute the three new axes
> v1 = B-A; v1 = v1/norm(v1);
> v3 = cross(v1,C-A); v3 = v3/norm(v3);
> v2 = cross(v3,v1); % v2 is already a unit vector
>
> % Get the transforms of A, B, and C
> P = A;
> TA = [dot(P-A,v1),dot(P-A,v2),dot(P-A,v3)]
> P = B;
> TB = [dot(P-A,v1),dot(P-A,v2),dot(P-A,v3)]
> P = C;
> TC = [dot(P-A,v1),dot(P-A,v2),dot(P-A,v3)]
>
> % Results
> A = 0.63940564208852 0.87421289486361 1.75240173032956
> B = -0.32005082643214 -0.13741380814487 0.61576962808672
> C = 0.97789406984520 -1.11534771220514 -0.55002144880449
>
> TA = 0 0 0
> TB = 1.79885459030876 0 0.00000000000000
> TC = 2.39315386990053 1.90968397273444 -0.00000000000000
>
> Roger Stafford

Aha! I see what I did wrong - to test I was simply calling on the original points instead of newly-defined ones, of course that would simply map it back to the original coordinate space. Thanks again for your help and patience.

am