Re: Help with geometry problem required.
- To: mathgroup at smc.vnet.net
- Subject: [mg20830] Re: Help with geometry problem required.
- From: Michael Ellis <michael at digsci.demon.co.uk>
- Date: Wed, 17 Nov 1999 03:40:44 -0500 (EST)
- Organization: Digital Scientific
- References: <80nfvr$24m@smc.vnet.net> <382FDFCD.540F711E@informatik.uni-stuttgart.de>
- Sender: owner-wri-mathgroup at wolfram.com
Thanks for your help and I have indeed found your solution to work. However,
when I substitute unbound variables for the points the solution is very long
winded and contains a great many common subexpressions. See:
In[105]:=
GetTransformation[p1_, p2_, p3_, q1_, q2_, q3_] :=
Module[{pcom, qcom, sp1, sp2, sp3, sq1, sq2, sq3, sp1Xsp2, sq1Xsq2, tp, tq,
m, m4},
pcom = (p1 + p2 + p3)/3;
qcom = (q1 + q2 + q3)/3;
sp1 = p1 - pcom;
sp2 = p2 - pcom;
sp3 = p3 - pcom;
sp1Xsp2 = Simplify[Cross[sp1, sp2]];
tp = {{1, 0, 0, -pcom[[1]]}, {0, 1, 0, -pcom[[2]]}, {0, 0,
1, -pcom[[3]]}, {0, 0, 0, 1}};
sq1 = q1 - qcom;
sq2 = q2 - qcom;
sq3 = q3 - qcom;
sq1Xsq2 = Simplify[Cross[sq1, sq2]];
tq = {{1, 0, 0, qcom[[1]]}, {0, 1, 0, qcom[[2]]}, {0, 0, 1,
qcom[[3]]}, {0, 0, 0, 1}};
m = Transpose[{sq1, sq2, sq1Xsq2}].Simplify[
Inverse[Transpose[{sp1, sp2, sp1Xsp2}]]];
m4 = tq.{{m[[1, 1]], m[[1, 2]], m[[1, 3]], 0}, {m[[2, 1]], m[[2, 2]],
m[[2, 3]], 0}, {m[[3, 1]], m[[3, 2]], m[[3, 3]], 0}, {0, 0, 0,
1}}.tp;
m4
]
Calling GetTransformation with hard coded values for the points quickly
yields the correct solution:
In[112]:=
p1 = {3, 1, 5}; p2 = {-4, 1, 2}; p3 = {1, 3, -5};
q1 = {-1, 3, 5}; q2 = {-1, -4, 2}; q3 = {-3, 1, -5};
M = GetTransformation[p1, p2, p3, q1, q2, q3]
Out[114]=
{{0, -1, 0, 0}, {1, 0, 0, 0}, {0, 0, 1, 0}, {0, 0, 0, 1}}
Calling GetTransformation with the following unbound points takes a little
longer but does indeed yield a very lenghthy solution.
p1 = {x1, y1, z1}; p2 = {x2, y2, z2}; p3 = {x3, y3, z3};
q1 = {X1, Y1, Z1}; q2 = {X2, Y2, Z2}; q3 = {X3, Y3, Z3};
M = GetTransformation[p1, p2, p3, q1, q2, q3]
<< Solution ommitted; it's about 400 lines long >>
I want to implement GeTtransformation in my C application by coding the
result above. I observe that there are many common subexpressions in the
solution, can I get Mathematica to eliminate common sub-expressions?
Note that calling Simplify doesnt remove all common sub expressions and
Simplify on my 400MHz PowerMac with the MatherKernal given 80MB of memory
evntually runs out of memory when asked to simplify the result.
Many thnaks -- Michael Ellis
Martin Kraus wrote:
> Michael Ellis wrote:
> >
> > I am new to this news group so please forgive me if this is an
> > inappropriate posting.
>
> I guess it is appropriate as we had a similar problem some
> weeks ago. :)
>
> > My Problem:
> > I have three points marked on a piece of rigid card at positions p1, p2
> > and p3. The card is moved, by translation and or rotation, but not
> > otherwise distorted to a new location. The three points p1, p2 and p3
> > are now at new positions say P1, P2 and P3. My question: Is there a 4 by
> > 4 Transform M that uniquely describes this relocation and if so how can
> > I derive M given p1, p2, p3, P1, P2 and P3.
> [code deleted]
>
> I assume that the points are given in a three-dimensional space
> (in spite of this card thing).
>
> I am not sure how to do this the simplest way, but one solutions
> goes as follows: Build up the matrix by
> 1) moving the center of mass (com) of the three points p1, p2, and p3
> to the origin.
> 2) rotate such that two points, say p1 - com and p2 - com, are mapped
> to the corresponding points P1 - Com and P2 - Com.
> 3) move the origin to the center of mass of P1, P2, and P3 (= Com)
>
> The idea is to rotate only around the origin. And only care about
> 2 points, as the third has to come out correctly if it was calculated
> correctly. How do we calculate this rotation matrix? In order to avoid
> sine and cosine, we use the following interpretation of column vectors
> of a matrix: The column vectors are the images of the unit vectors.
> How does this help? Well, if we want to map p1 - com and p2 - com
> to P1 - com and P2 - com we can do so with one additional step:
> mapping to the unit vectors. Thus, we map p1 - com and p2 - com to
> the unit vectors. But that's just the inverse of mapping the unit
> vectors
> to p1 - com and p2 - com! Thus, we have to compute the inverse of
> a matrix with column vectors p1 - com and p2 - com. But wait! We need
> a third column. What about using p3 - com as a third vector?
> Actually, that's a bad idea, because it is not linear independent.
> (Because we have moved the center of mass to the origin: Think of it
> as moving the center of mass of a triangle into the origin. Then the
> vertices of the triangle are in one plane containing the origin, i.e.
> not
> linear independent.)
> Thus, we just construct a linear independent vector with the cross
> product.
> So we have our matrix to map p1 - com and p2 - com to the unit vectors.
> In a second step the unit vectors are mapped to P1 - com and P2 - com
> (and the cross product of p1 - com and p2 - com to the cross product
> of P1 - com and P2 - com). That's really easy as all we have to do is
> to write P1 - com and company into the columns of a matrix.
> We multiply the matrix and multiply with translation matrices as
> mentioned
> above and we are done.
>
> Code with test:
>
> p1 = {3, 1, 5}; p2 = {-4, 1, 2}; p3 = {1,
> 3, -5}; alpha = 0.2; beta = 0.3; gamma = 0.4; rot = {{Cos[alpha],
> Sin[alpha], 0}, {-Sin[alpha], Cos[alpha], 0}, {0, 0, 1}}; rot =
> rot.{{Cos[beta], 0, Sin[beta]}, {0, 1, 0}, {-Sin[beta], 0,
> Cos[beta]}}; t = {3, 1, 2}; rot =
> rot.{{1, 0, 0}, {0, Cos[gamma], Sin[gamma]}, {0, -Sin[gamma],
> Cos[gamma]}}; P1 = rot.(p1 + t); P2 = rot.(p2 + t); P3 = rot.(p3
> + t);
>
> com = (p1 + p2 + p3)/3;
> Com = (P1 + P2 + P3)/3;
>
> pp1 = p1 - com;
> pp2 = p2 - com;
> pp3 = p3 - com;
> np = Cross[pp1, pp2];
> Pp1 = P1 - Com;
> Pp2 = P2 - Com;
> Pp3 = P3 - Com;
> Np = Cross[Pp1, Pp2];
>
> m = Transpose[{Pp1, Pp2, Np}].Inverse[Transpose[{pp1, pp2, np}]];
>
> m4 = {{1, 0, 0, Com[[1]]}, {0, 1, 0, Com[[2]]}, {0, 0, 1, Com[[3]]}, {0,
> 0, 0,
> 1}}.{{m[[1, 1]], m[[1, 2]], m[[1, 3]], 0}, {m[[2, 1]], m[[2,
> 2]],
> m[[2, 3]], 0}, {m[[3, 1]], m[[3, 2]], m[[3, 3]], 0}, {0, 0, 0,
> 1}}.{{1, 0, 0, -com[[1]]}, {0, 1, 0, -com[[2]]}, {0, 0,
> 1, -com[[3]]}, {0, 0, 0, 1}};
>
> m4.Append[p1, 1]
>
> P1
>
> Note that the points may not be on one straight line for this solution.
> (That's a special case which needs extra code...)
>
> Hope that helps
>
> Martin Kraus