7
$\begingroup$

I need a function that rotates and translates a huge amount of line segments.

For example, I have a set of line segments in the form {{x0,y0,z0},{x1,y1,z1}}

data = RandomReal[{-1, 1}, {20000, 2, 3}]

and three vectors

n = {n1,n2,n3}
s = {s1,s2,s3} (* Orthogonal to n *)
p = {p1,p2,p3}

I want to perform a rigid (right handed) transformation that sends p to the origin, n to the $z$ axis and $s$ to the $y$ axis, in order to apply it to data.

I know about RotationMatrix, RotationTransform, and I have coded a (very ugly) solution but, due to inexperience and the amount of data, I'm struggling with performance.

I would be very grateful if someone can provide a fast solution and, if possible, to explain why is fast.

EDIT

Here is the ugly code, for my own embarrassment

rotZ = RotationTransform[{n, {0, 0, 1}}];
rotY = RotationTransform[{s, {0, 1, 0}}];
RT[v_]:= rotY[rotZ[-p + #]] & /@ v;

and then

dataTrans = RT /@ data

Taking

n = RandomReal[{-1, 1}, 3];
s = {-n[[2]], n[[1]], 0}; (* n and s are orthogonal *)
p = RandomReal[{-1, 1}, 3];

then

Timing[RT /@ data][[1]]
(*    3.60422    *)

In an Intel Core2 Duo CPU T8100 @2.10GHz, 4gb ram Ubuntu 12.04, Mathematica 8 distribution.

Improve this question
$\endgroup$
6
  • 3
    $\begingroup$ It might help to post the code, "ugly" or not. $\endgroup$
    – Verbeia
    Aug 11, 2012 at 21:36
  • $\begingroup$ Are n,s,p already orthogonal? $\endgroup$
    – Dr. belisarius
    Aug 11, 2012 at 22:00
  • $\begingroup$ n and s are. $\endgroup$
    – Pragabhava
    Aug 11, 2012 at 22:13
  • 1
    $\begingroup$ There's an ampersand missing in the definition of RT $\endgroup$
    – Sjoerd C. de Vries
    Aug 11, 2012 at 22:23
  • 1
    $\begingroup$ @SjoerdC.deVries Thanks for all the debugging! $\endgroup$
    – Pragabhava
    Aug 11, 2012 at 22:43

1 Answer 1

Reset to default
7
$\begingroup$

Defining n, s and p with their numeric values before defining the transformation function helps considerably, but a huge difference is made by changing the definition of RT to use the numerical matrices instead of the RotationTransform-s:

data = RandomReal[{-1, 1}, {20000, 2, 3}];
n = RandomReal[{-1, 1}, 3];
s = {-n[[2]], n[[1]], 0};(*n and s are orthogonal*)
p = RandomReal[{-1, 1}, 3];
rotZ = RotationTransform[{n, {0, 0, 1}}];
rotY = RotationTransform[{s, {0, 1, 0}}];

m = TransformationMatrix[rotY][[1 ;; 3, 1 ;; 3]].TransformationMatrix[
     rotZ][[1 ;; 3, 1 ;; 3]];

RT1[v_] := m.(-p + #) & /@ v;
RT2[v_] := rotY[rotZ[-p + #]] & /@ v;

Timing[r1 = RT1 /@ data][[1]]
Timing[r2 = RT2 /@ data][[1]]

0.328

2.824

Mathematica graphics

Improve this answer
$\endgroup$
4
  • $\begingroup$ I didn't know that one could extract the numerical matrix from RotationTransform. Is there any difference (qualitatively) than simply use RotationMatrix? $\endgroup$
    – Pragabhava
    Aug 12, 2012 at 4:08
  • 3
    $\begingroup$ Let RT3[segs_] := Transpose[m.(Transpose[segs, {3, 2, 1}] - p), {3, 2, 1}] and data = RandomReal[{-1, 1}, {200000, 2, 3}]. The timing for RT1 /@ data is 3.156 seconds; for RT3[data] (using efficient multithreaded Transpose and Dot) it is a matter of milliseconds. $\endgroup$
    – Oleksandr R.
    Aug 12, 2012 at 4:31
  • $\begingroup$ @OleksandrR. I suspected the (-p + #) & /@ approach was inefficient, but I didn't know how else to do the sum. Can you explain why the program can do the sum when dimensions are {3,2,N} but it can't when they are {N,2,3}? $\endgroup$
    – Pragabhava
    Aug 12, 2012 at 5:14
  • 1
    $\begingroup$ @Manuel that's just the way Mathematica is--which is my way of saying I don't know. Personally, I'd find exhaustive array broadcasting rules (like one finds in NumPy) very useful, but for whatever reason there hasn't been much effort in this direction. There are only so many features WRI can work on at once, I suppose, and Transpose usually gets the job done after all. $\endgroup$
    – Oleksandr R.
    Aug 12, 2012 at 5:21

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.