P is for Practical PerlMonks

### Rotation around an axis clockwise and anti-clockwise algorithm question?

by fraizerangus (Sexton)
 on Jun 12, 2012 at 16:19 UTC Need Help??
fraizerangus has asked for the wisdom of the Perl Monks concerning the following question:

Hello Monks

I've written a program which rotates atoms (with a set rotation and translation) around an axis for a Bioinformatics program I'm working on, it works going in one direction, does anyone know of a simple tweak I could do make it rotate in the opposite (anti-clockwise direction) please, which wouldn't involve me re-writing the whole program again?

```
sub rotator {

use Math::Vector::Real;

(*xcord, *ycord, *zcord, *BMatoms, \$AxisLineUnitvector, \$XlineAxis, \$Y
+lineAxis, \$ZlineAxis) = @_;

\$transl = 0.4;

use Math::Trig;

#coefficient processing

#Breakdown axis line vector into constituent X,Y & Z co-ordinates
(\$LUVx, \$LUVy, \$LUVz) = &VectorBreakdown(\$AxisLineUnitvector);

#for loop to process all atoms
for (my \$i = 0; \$i < @BMatoms; \$i++) {

#distances of co-ordinates to axis line
\$Xco[\$i] = \$xcord[\$i] - \$XlineAxis;
\$Yco[\$i] = \$ycord[\$i] - \$YlineAxis;
\$Zco[\$i] = \$zcord[\$i] - \$ZlineAxis;

#convert diatance co-ordinates into a vector
\$Vector[\$i] = V(\$Xco[\$i], \$Yco[\$i], \$Zco[\$i]);

#dot product
\$scal[\$i] = \$LineUnitvector * \$Vector[\$i];

#scalar product
\$f[\$i] = (\$scal[\$i] * \$AxisLineUnitvector) - \$Vector[\$i];

#dot product
\$f[\$i] = \$coefl1 * \$f[\$i];

#cross product
\$s[\$i] = \$AxisLineUnitvector x \$Vector[\$i];

#dot product
\$s[\$i] = \$coefl2 * \$s[\$i];

#vector broken down into constituent parts
(\$fx[\$i], \$fy[\$i], \$fz[\$i]) = &VectorBreakdown(\$f[\$i]);
(\$sx[\$i], \$sy[\$i], \$sz[\$i]) = &VectorBreakdown(\$s[\$i]);

#new co-ordinates
\$xo[\$i] = \$xcord[\$i] + \$fx[\$i] + \$sx[\$i] + (\$transl * \$LUVx);
\$yo[\$i] = \$ycord[\$i] + \$fy[\$i] + \$sy[\$i] + (\$transl * \$LUVy);
\$zo[\$i] = \$zcord[\$i] + \$fz[\$i] + \$sz[\$i] + (\$transl * \$LUVz);

#rounding up
\$xo[\$i] = sprintf("%.3f",\$xo[\$i]);
\$yo[\$i] = sprintf("%.3f",\$yo[\$i]);
\$zo[\$i] = sprintf("%.3f",\$zo[\$i]);

}

return(\@xo, \@yo, \@zo);

• Comment on Rotation around an axis clockwise and anti-clockwise algorithm question?

Replies are listed 'Best First'.
Re: Rotation around an axis clockwise and anti-clockwise algorithm question?
by BrowserUk (Pope) on Jun 12, 2012 at 19:52 UTC

One suggestion. It appears as though you may be rotating your points multiple times; and doing so by modifying the original coordinates by each rotation. The problem with this is that as the rotations accumulate, the rounding errors inevitably introduced by each successive rotation also accumulate, but those rounding errors will be different and accumulate differently for each point on the structure being rotated. The net affect is that after many small rotations, the registration between the points making up the structure will wander.

A concrete example might clarify here. Let's say you start with a nice unit vector cube distributed equidistant around the origin:

```my \$cube = [
[-0.5, 0.5, 0.5],[ 0.5, 0.5, 0.5],[ 0.5,-0.5, 0.5],[-0.5,-0.5, 0.5
+], # front face
[-0.5, 0.5,-0.5],[ 0.5, 0.5,-0.5],[ 0.5,-0.5,-0.5],[-0.5,-0.5,-0.5
+], # back face
];

After rotating 360° around each axis in 1° steps -- 1080 floating point calculations, you will likely end up with values that look more like:

```[
[-0.50000000000000011, 0.49999999999999989, 0.49999999999999971],
[ 0.50000000000000011, 0.500000000023, 0.499999999987873],
[ 0.50000000000000011,-0.49999999999445, 0.5000000000123],
[-0.49999999999999576,-0.49999999999999989, 0.500000000000876], #
+front face
[-0.500000000000001, 0.49999999999999989,-0.49999999999999987],
[ 0.50000000000000011, 0.500000000000000345,-0.500000000000000234]
+,
[ 0.49999999999999921,-0.500000000000000001,-0.500000000000000022]
+,
[-0.50000000000000011,-0.4999999999999999002,-0.50000000000000001]
+, # back face
];

As you can see, the 'cube' has become irrevocably distorted. And that was starting with relatively large (compared to your atomic scale measurements) and floating point friendly (0.5 stored exactly as a FP number) values. With smaller values and values that cannot be exactly represented, the distortions are going to be even more noticeable.

An oft-used alternative is to retain the original coordinates and accumulate the rotations -- usually in the form of a transform matrix -- and only apply the rotations when fixing a position for output (eg. when displaying).

In this way, the you only have the distortions from a single calculations rounding errors at any given position, rather than the accumulated distortions of many calculations.

With the rise and rise of 'Social' network sites: 'Computers are making people easier to use everyday'
Examine what is said, not who speaks -- Silence betokens consent -- Love the truth but pardon error.
"Science is about questioning the status quo. Questioning authority".
In the absence of evidence, opinion is indistinguishable from prejudice.

The start of some sanity?

Re: Rotation around an axis clockwise and anti-clockwise algorithm question?
by dave_the_m (Prior) on Jun 12, 2012 at 16:34 UTC
What's the matter with just specifying a negative angle of rotation?

Dave.

Can that be done in the conversion to radians? I thought it couldn't? I just tried changing it and the program crashed?
Well, if your example hard-codes the rotation angle (9.4) as it appears, can't you use (360 - angle_in_degrees) for anti-clockwise?

fnord

Re: Rotation around an axis clockwise and anti-clockwise algorithm question?
by johngg (Abbot) on Jun 12, 2012 at 22:23 UTC

Not related to your problem (at least, I don't think so) but did you really mean to use the string multiplier (x) in this line of code?

```\$s[\$i] = \$AxisLineUnitvector x \$Vector[\$i];

It just seems a little out of place when all the other multipliers are the arithmetic (*) sort.

Cheers,

JohnGG

Math::Vector::Real overloads the x operator to perform the cross product of the two vectors. It is customary in maths to use x for that operation.
Re: Rotation around an axis clockwise and anti-clockwise algorithm question?
by zentara (Archbishop) on Jun 13, 2012 at 13:25 UTC
Amplifying on what others said about using matrix operations, I was wondering if you want to give your program a visualization GUI? If so, I would like to mention Tk::Zinc and Gtk2's Goo::Canvas. The Zinc canvas does matrix rotations, and is very good at maintaining accuracy of its matrix with tset and tget. I'm disappointed that its parent website has seem to have disappeared, as Zinc was great.

The Goo canvas is still alive and well and being developed. In both their rotation algorithms, you can specify negative angles. As a matter of fact, a good vector equation should accept negative angles and deal with them properly within its coordinate system.

Here is a simple Goo canvas positive and negative angle rotation example.

```#!/usr/bin/perl -w
use strict;
use warnings;
use Goo::Canvas;
use Gtk2 '-init';
use Glib qw(TRUE FALSE);

my \$window = Gtk2::Window->new('toplevel');
\$window->signal_connect('delete_event' => sub { Gtk2->main_quit; });
\$window->set_size_request(640, 600);

my \$swin = Gtk2::ScrolledWindow->new;

my \$canvas = Goo::Canvas->new();
\$canvas->set_size_request(800, 650);
\$canvas->set_bounds(0, 0, 1000, 1000);
my \$root = \$canvas->get_root_item();

# first offset set
my \$pts_ref = [50,50,180,120,90,100,50,50];
my \$line = Goo::Canvas::Polyline->new(
\$root, TRUE,
\$pts_ref,
'stroke-color' => 'black',
'line-width' => 3,
'fill-color-rgba' => 0x3cb37180,
);

my (\$midx, \$midy) = _get_CM( @\$pts_ref );

my \$ellipse = Goo::Canvas::Ellipse->new(
\$root, \$midx-2, \$midy-2,\$midx+2, \$midy+2,
'stroke-color' => 'goldenrod',
'line-width' => 8
);

my \$ellipse1 = Goo::Canvas::Ellipse->new(
\$root, -2, -2, +2, +2,
'stroke-color' => 'black',
'line-width' => 4
);

\$ellipse1->translate(\$midx,\$midy);
# end first set

#if you have equilateral shapes it's
#possible to make at origin and translate
my \$group = Goo::Canvas::Group->new(\$root);
my \$pts_ref1 = [-60,0, 60,0, 0, 40, -60, 0];
my \$line1 = Goo::Canvas::Polyline->new(
\$group, TRUE,
\$pts_ref1,
'stroke-color' => 'black',
'line-width' => 3,
'fill-color-rgba' => 0xffb37180,
);

my (\$midx1, \$midy1) = _get_CM( @\$pts_ref1 );
print "\$midx1, \$midy1\n";

my \$ellipse2 = Goo::Canvas::Ellipse->new(
\$group, -2, -2, +2, +2,
'stroke-color' => 'black',
'line-width' => 4
);

\$ellipse2->translate(0,\$midy1 + 4);

my \$ellipse3 = Goo::Canvas::Ellipse->new(
\$group,
,-60,-60,60,60,
'stroke-color' => 'green',
'line-width' => 4
);

\$ellipse3->translate(60,60);

#move whole group
\$group->translate(400,400);

my \$id = Glib::Timeout->add (10, sub {
\$line->rotate (10, \$midx, \$midy);
\$group->rotate (-1, 0, 0 );
return 1;
});

\$window->show_all();
Gtk2->main;

#################################################################
# This sub finds the center of mass of a polygon.
# I grabbed the algorithm somewhere from the web.
# I grabbed it from Ala Qumsieh's RotCanvas :-)
sub _get_CM {
my (\$x, \$y, \$area);

my \$i = 0;

while (\$i < \$#_) {
my \$x0 = \$_[\$i];
my \$y0 = \$_[\$i+1];

my (\$x1, \$y1);
if (\$i+2 > \$#_) {
\$x1 = \$_[0];
\$y1 = \$_[1];
} else {
\$x1 = \$_[\$i+2];
\$y1 = \$_[\$i+3];
}

\$i += 2;

my \$a1 = 0.5*(\$x0 + \$x1);
my \$a2 = (\$x0**2 + \$x0*\$x1 + \$x1**2)/6;
my \$a3 = (\$x0*\$y1 + \$y0*\$x1 + 2*(\$x1*\$y1 + \$x0*\$y0))/6;
my \$b0 = \$y1 - \$y0;

\$area += \$a1 * \$b0;
\$x    += \$a2 * \$b0;
\$y    += \$a3 * \$b0;
}

return split ' ', sprintf "%.0f %0.f" => \$x/\$area, \$y/\$area;
}
####################################################################
__END__

I'm not really a human, but I play one on earth.
Old Perl Programmer Haiku ................... flash japh
Re: Rotation around an axis clockwise and anti-clockwise algorithm question?
by salva (Abbot) on Jun 13, 2012 at 14:37 UTC
I have just extended Math::Vector::Real to support a rotate_3d method:
```@r = \$axis->rotate_3d(\$rad, @vs);  # anticlockwise
@s = \$axis->rotate_3d(-\$rad, @vs); # clockwise

Create A New User
Node Status?
node history
Node Type: perlquestion [id://975842]
Front-paged by Corion
help
Chatterbox?
and all is quiet...

How do I use this? | Other CB clients
Other Users?
Others chilling in the Monastery: (5)
As of 2018-01-19 16:04 GMT
Sections?
Information?
Find Nodes?
Leftovers?
Voting Booth?
How did you see in the new year?

Results (221 votes). Check out past polls.

Notices?