#
# Trigonometric functions, mostly inherited from Math::Complex.
-# -- Jarkko Hietaniemi, April 1997
+# -- Jarkko Hietaniemi, since April 1997
# -- Raphael Manfredi, September 1996 (indirectly: because of Math::Complex)
#
use vars qw($VERSION $PACKAGE
@ISA
- @EXPORT);
+ @EXPORT @EXPORT_OK %EXPORT_TAGS);
@ISA = qw(Exporter);
@EXPORT = (@{$Math::Complex::EXPORT_TAGS{'trig'}},
@angcnv);
-use constant pi2 => 2 * pi;
-use constant DR => pi2/360;
-use constant RD => 360/pi2;
-use constant DG => 400/360;
-use constant GD => 360/400;
-use constant RG => 400/pi2;
-use constant GR => pi2/400;
+my @rdlcnv = qw(cartesian_to_cylindrical
+ cartesian_to_spherical
+ cylindrical_to_cartesian
+ cylindrical_to_spherical
+ spherical_to_cartesian
+ spherical_to_cylindrical);
+
+@EXPORT_OK = (@rdlcnv, 'great_circle_distance');
+
+%EXPORT_TAGS = ('radial' => [ @rdlcnv ]);
+
+use constant pi2 => 2 * pi;
+use constant pip2 => pi / 2;
+use constant DR => pi2/360;
+use constant RD => 360/pi2;
+use constant DG => 400/360;
+use constant GD => 360/400;
+use constant RG => 400/pi2;
+use constant GR => pi2/400;
#
# Truncating remainder.
sub grad2rad ($) { remt(GR * $_[0], pi2) }
+sub cartesian_to_spherical {
+ my ( $x, $y, $z ) = @_;
+
+ my $rho = sqrt( $x * $x + $y * $y + $z * $z );
+
+ return ( $rho,
+ atan2( $y, $x ),
+ $rho ? acos( $z / $rho ) : 0 );
+}
+
+sub spherical_to_cartesian {
+ my ( $rho, $theta, $phi ) = @_;
+
+ return ( $rho * cos( $theta ) * sin( $phi ),
+ $rho * sin( $theta ) * sin( $phi ),
+ $rho * cos( $phi ) );
+}
+
+sub spherical_to_cylindrical {
+ my ( $x, $y, $z ) = spherical_to_cartesian( @_ );
+
+ return ( sqrt( $x * $x + $y * $y ), $_[1], $z );
+}
+
+sub cartesian_to_cylindrical {
+ my ( $x, $y, $z ) = @_;
+
+ return ( sqrt( $x * $x + $y * $y ), atan2( $y, $x ), $z );
+}
+
+sub cylindrical_to_cartesian {
+ my ( $rho, $theta, $z ) = @_;
+
+ return ( $rho * cos( $theta ), $rho * sin( $theta ), $z );
+}
+
+sub cylindrical_to_spherical {
+ return ( cartesian_to_spherical( cylindrical_to_cartesian( @_ ) ) );
+}
+
+sub great_circle_distance {
+ my ( $theta0, $phi0, $theta1, $phi1, $rho ) = @_;
+
+ $rho = 1 unless defined $rho; # Default to the unit sphere.
+
+ my $lat0 = pip2 - $phi0;
+ my $lat1 = pip2 - $phi1;
+
+ return $rho *
+ acos(cos( $lat0 ) * cos( $lat1 ) * cos( $theta0 - $theta1 ) +
+ sin( $lat0 ) * sin( $lat1 ) );
+}
+
+=pod
+
=head1 NAME
Math::Trig - trigonometric functions
The tangent
- tan
+=over 4
+
+=item B<tan>
+
+=back
The cofunctions of the sine, cosine, and tangent (cosec/csc and cotan/cot
are aliases)
- csc cosec sec cot cotan
+B<csc>, B<cosec>, B<sec>, B<sec>, B<cot>, B<cotan>
The arcus (also known as the inverse) functions of the sine, cosine,
and tangent
- asin acos atan
+B<asin>, B<acos>, B<atan>
The principal value of the arc tangent of y/x
- atan2(y, x)
+B<atan2>(y, x)
The arcus cofunctions of the sine, cosine, and tangent (acosec/acsc
and acotan/acot are aliases)
- acsc acosec asec acot acotan
+B<acsc>, B<acosec>, B<asec>, B<acot>, B<acotan>
The hyperbolic sine, cosine, and tangent
- sinh cosh tanh
+B<sinh>, B<cosh>, B<tanh>
The cofunctions of the hyperbolic sine, cosine, and tangent (cosech/csch
and cotanh/coth are aliases)
- csch cosech sech coth cotanh
+B<csch>, B<cosech>, B<sech>, B<coth>, B<cotanh>
The arcus (also known as the inverse) functions of the hyperbolic
sine, cosine, and tangent
- asinh acosh atanh
+B<asinh>, B<acosh>, B<atanh>
The arcus cofunctions of the hyperbolic sine, cosine, and tangent
(acsch/acosech and acoth/acotanh are aliases)
- acsch acosech asech acoth acotanh
+B<acsch>, B<acosech>, B<asech>, B<acoth>, B<acotanh>
The trigonometric constant B<pi> is also defined.
- $pi2 = 2 * pi;
+$pi2 = 2 * B<pi>;
=head2 ERRORS DUE TO DIVISION BY ZERO
The following functions
- tan
- sec
- csc
- cot
- asec
+ acoth
acsc
- tanh
- sech
- csch
- coth
- atanh
- asech
acsch
- acoth
+ asec
+ asech
+ atanh
+ cot
+ coth
+ csc
+ csch
+ sec
+ sech
+ tan
+ tanh
cannot be computed for all arguments because that would mean dividing
by zero or taking logarithm of zero. These situations cause fatal
That is, a complex number with the real part of approximately C<1.571>
and the imaginary part of approximately C<-1.317>.
-=head1 ANGLE CONVERSIONS
+=head1 PLANE ANGLE CONVERSIONS
(Plane, 2-dimensional) angles may be converted with the following functions.
The full circle is 2 I<pi> radians or I<360> degrees or I<400> gradians.
+=head1 RADIAL COORDINATE CONVERSIONS
+
+B<Radial coordinate systems> are the B<spherical> and the B<cylindrical>
+systems, explained shortly in more detail.
+
+You can import radial coordinate conversion functions by using the
+C<:radial> tag:
+
+ use Math::Trig ':radial';
+
+ ($rho, $theta, $z) = cartesian_to_cylindrical($x, $y, $z);
+ ($rho, $theta, $phi) = cartesian_to_spherical($x, $y, $z);
+ ($x, $y, $z) = cylindrical_to_cartesian($rho, $theta, $z);
+ ($rho_s, $theta, $phi) = cylindrical_to_spherical($rho_c, $theta, $z);
+ ($x, $y, $z) = spherical_to_cartesian($rho, $theta, $phi);
+ ($rho_c, $theta, $z) = spherical_to_cylindrical($rho_s, $theta, $phi);
+
+B<All angles are in radians>.
+
+=head2 COORDINATE SYSTEMS
+
+B<Cartesian> coordinates are the usual rectangular I<(x, y,
+z)>-coordinates.
+
+Spherical coordinates, I<(rho, theta, pi)>, are three-dimensional
+coordinates which define a point in three-dimensional space. They are
+based on a sphere surface. The radius of the sphere is B<rho>, also
+known as the I<radial> coordinate. The angle in the I<xy>-plane
+(around the I<z>-axis) is B<theta>, also known as the I<azimuthal>
+coordinate. The angle from the I<z>-axis is B<phi>, also known as the
+I<polar> coordinate. The `North Pole' is therefore I<0, 0, rho>, and
+the `Bay of Guinea' (think of the missing big chunk of Africa) I<0,
+pi/2, rho>.
+
+B<Beware>: some texts define I<theta> and I<phi> the other way round,
+some texts define the I<phi> to start from the horizontal plane, some
+texts use I<r> in place of I<rho>.
+
+Cylindrical coordinates, I<(rho, theta, z)>, are three-dimensional
+coordinates which define a point in three-dimensional space. They are
+based on a cylinder surface. The radius of the cylinder is B<rho>,
+also known as the I<radial> coordinate. The angle in the I<xy>-plane
+(around the I<z>-axis) is B<theta>, also known as the I<azimuthal>
+coordinate. The third coordinate is the I<z>, pointing up from the
+B<theta>-plane.
+
+=head2 3-D ANGLE CONVERSIONS
+
+Conversions to and from spherical and cylindrical coordinates are
+available. Please notice that the conversions are not necessarily
+reversible because of the equalities like I<pi> angles being equal to
+I<-pi> angles.
+
+=over 4
+
+=item cartesian_to_cylindrical
+
+ ($rho, $theta, $z) = cartesian_to_cylindrical($x, $y, $z);
+
+=item cartesian_to_spherical
+
+ ($rho, $theta, $phi) = cartesian_to_spherical($x, $y, $z);
+
+=item cylindrical_to_cartesian
+
+ ($x, $y, $z) = cylindrical_to_cartesian($rho, $theta, $z);
+
+=item cylindrical_to_spherical
+
+ ($rho_s, $theta, $phi) = cylindrical_to_spherical($rho_c, $theta, $z);
+
+Notice that when C<$z> is not 0 C<$rho_s> is not equal to C<$rho_c>.
+
+=item spherical_to_cartesian
+
+ ($x, $y, $z) = spherical_to_cartesian($rho, $theta, $phi);
+
+=item spherical_to_cylindrical
+
+ ($rho_c, $theta, $z) = spherical_to_cylindrical($rho_s, $theta, $phi);
+
+Notice that when C<$z> is not 0 C<$rho_c> is not equal to C<$rho_s>.
+
+=back
+
+=head1 GREAT CIRCLE DISTANCES
+
+You can compute spherical distances, called B<great circle distances>,
+by importing the C<great_circle_distance> function:
+
+ use Math::Trig 'great_circle_distance'
+
+ $distance = great_circle_distance($theta0, $phi0, $theta1, $phi, [, $rho]);
+
+The I<great circle distance> is the shortest distance between two
+points on a sphere. The distance is in C<$rho> units. The C<$rho> is
+optional, it defaults to 1 (the unit sphere), therefore the distance
+defaults to radians.
+
+=head EXAMPLES
+
+To calculate the distance between London (51.3N 0.5W) and Tokyo (35.7N
+139.8E) in kilometers:
+
+ use Math::Trig qw(great_circle_distance deg2rad);
+
+ # Notice the 90 - latitude: phi zero is at the North Pole.
+ @L = (deg2rad(-0.5), deg2rad(90 - 51.3));
+ @T = (deg2rad(139.8),deg2rad(90 - 35.7));
+
+ $km = great_circle_distance(@L, @T, 6378);
+
+The answer may be off by up to 0.3% because of the irregular (slightly
+aspherical) form of the Earth.
+
=head1 BUGS
Saying C<use Math::Trig;> exports many mathematical routines in the
+++ /dev/null
-package Math::Trig::Radial;
-
-use strict;
-use vars qw(@ISA @EXPORT);
-@ISA = qw(Exporter);
-
-@EXPORT =
- qw(
- cartesian_to_cylindrical
- cartesian_to_spherical
- cylindrical_to_cartesian
- cylindrical_to_spherical
- spherical_to_cartesian
- spherical_to_cylindrical
- great_circle_distance
- );
-
-use Math::Trig;
-
-sub pip2 { pi/2 }
-
-=pod
-
-=head1 NAME
-
-Math::Trig::Radial - spherical and cylindrical trigonometry
-
-=head1 SYNOPSIS
-
- use Math::Trig::Radial;
-
- ($rho, $theta, $z) = cartesian_to_cylindrical($x, $y, $z);
- ($rho, $theta, $phi) = cartesian_to_spherical($x, $y, $z);
- ($x, $y, $z) = cylindrical_to_cartesian($rho, $theta, $z);
- ($rho_s, $theta, $phi) = cylindrical_to_spherical($rho_c, $theta, $z);
- ($x, $y, $z) = spherical_to_cartesian($rho, $theta, $phi);
- ($rho_c, $theta, $z) = spherical_to_cylindrical($rho_s, $theta, $phi);
-
-=head1 DESCRIPTION
-
-This module contains a few basic spherical and cylindrical
-trigonometric formulas. B<All angles are in radians>, if needed
-use C<Math::Trig> angle unit conversions.
-
-=head2 COORDINATE SYSTEMS
-
-B<Cartesian> coordinates are the usual rectangular I<xyz>-coordinates.
-
-Spherical coordinates are three-dimensional coordinates which define a
-point in three-dimensional space. They are based on a sphere surface.
-The radius of the sphere is B<rho>, also known as the I<radial>
-coordinate. The angle in the I<xy>-plane (around the I<z>-axis) is
-B<theta>, also known as the I<azimuthal> coordinate. The angle from
-the I<z>-axis is B<phi>, also known as the I<polar> coordinate. The
-`North Pole' is therefore I<0, 0, rho>, and the `Bay of Guinea' (think
-Africa) I<0, pi/2, rho>.
-
-Cylindrical coordinates are three-dimensional coordinates which define
-a point in three-dimensional space. They are based on a cylinder
-surface. The radius of the cylinder is B<rho>, also known as the
-I<radial> coordinate. The angle in the I<xy>-plane (around the
-I<z>-axis) is B<theta>, also known as the I<azimuthal> coordinate.
-The third coordinate is the I<z>.
-
-=head2 CONVERSIONS
-
-Conversions to and from spherical and cylindrical coordinates are
-available. Please notice that the conversions are not necessarily
-reversible because of the equalities like I<pi> angles equals I<-pi>
-angles.
-
-=over 4
-
-=item cartesian_to_cylindrical
-
- ($rho, $theta, $z) = cartesian_to_cylindrical($x, $y, $z);
-
-=item cartesian_to_spherical
-
- ($rho, $theta, $phi) = cartesian_to_spherical($x, $y, $z);
-
-=item cylindrical_to_cartesian
-
- ($x, $y, $z) = cylindrical_to_cartesian($rho, $theta, $z);
-
-=item cylindrical_to_spherical
-
- ($rho_s, $theta, $phi) = cylindrical_to_spherical($rho_c, $theta, $z);
-
-Notice that when C<$z> is not 0 C<$rho_s> is not equal to C<$rho_c>.
-
-=item spherical_to_cartesian
-
- ($x, $y, $z) = spherical_to_cartesian($rho, $theta, $phi);
-
-=item spherical_to_cylindrical
-
- ($rho_c, $theta, $z) = spherical_to_cylindrical($rho_s, $theta, $phi);
-
-Notice that when C<$z> is not 0 C<$rho_c> is not equal to C<$rho_s>.
-
-=back
-
-=head2 GREAT CIRCLE DISTANCE
-
- $distance = great_circle_distance($x0, $y0, $z0, $x1, $y1, $z1 [, $rho]);
-
-The I<great circle distance> is the shortest distance between two
-points on a sphere. The distance is in C<$rho> units. The C<$rho> is
-optional, it defaults to 1 (the unit sphere), therefore the distance
-defaults to radians. The coordinates C<$x0> ... C<$z1> are in
-cartesian coordinates.
-
-=head EXAMPLES
-
-To calculate the distance between London (51.3N 0.5W) and Tokyo (35.7N
-139.8E) in kilometers:
-
- use Math::Trig::Radial;
- use Math::Trig;
-
- my @L = spherical_to_cartesian(1, map { deg2rad $_ } qw(51.3 -0.5));
- my @T = spherical_to_cartesian(1, map { deg2rad $_ } qw(35.7 139.8));
-
- $km = great_circle_distance(@L, @T, 6378);
-
-The answer may be off by up to 0.3% because of the irregular (slightly
-aspherical) form of the Earth.
-
-=head2 AUTHOR
-
-Jarkko Hietaniemi F<E<lt>jhi@iki.fiE<gt>>
-
-=cut
-
-sub cartesian_to_spherical {
- my ( $x, $y, $z ) = @_;
-
- my $rho = sqrt( $x * $x + $y * $y + $z * $z );
-
- return ( $rho,
- atan2( $y, $x ),
- $rho ? acos( $z / $rho ) : 0 );
-}
-
-sub spherical_to_cartesian {
- my ( $rho, $theta, $phi ) = @_;
-
- return ( $rho * cos( $theta ) * sin( $phi ),
- $rho * sin( $theta ) * sin( $phi ),
- $rho * cos( $phi ) );
-}
-
-sub spherical_to_cylindrical {
- my ( $x, $y, $z ) = spherical_to_cartesian( @_ );
-
- return ( sqrt( $x * $x + $y * $y ), $_[1], $z );
-}
-
-sub cartesian_to_cylindrical {
- my ( $x, $y, $z ) = @_;
-
- return ( sqrt( $x * $x + $y * $y ), atan2( $y, $x ), $z );
-}
-
-sub cylindrical_to_cartesian {
- my ( $rho, $theta, $z ) = @_;
-
- return ( $rho * cos( $theta ), $rho * sin( $theta ), $z );
-}
-
-sub cylindrical_to_spherical {
- return ( cartesian_to_spherical( cylindrical_to_cartesian( @_ ) ) );
-}
-
-sub great_circle_distance {
- my ( $x0, $y0, $z0, $x1, $y1, $z1, $rho ) = @_;
-
- $rho = 1 unless defined $rho; # Default to the unit sphere.
-
- my ( $r0, $theta0, $phi0 ) = cartesian_to_spherical( $x0, $y0, $z0 );
- my ( $r1, $theta1, $phi1 ) = cartesian_to_spherical( $x1, $y1, $z1 );
-
- my $lat0 = pip2 - $phi0;
- my $lat1 = pip2 - $phi1;
-
- return $rho *
- acos(cos( $lat0 ) * cos( $lat1 ) * cos( $theta0 - $theta1 ) +
- sin( $lat0 ) * sin( $lat1 ) );
-}
-
-1;
-
abs($_[0] - $_[1]) < (defined $_[2] ? $_[2] : $eps);
}
-print "1..7\n";
+print "1..20\n";
$x = 0.9;
print 'not ' unless (near(tan($x), sin($x) / cos($x)));
print 'not ' unless (near(rad2deg(pi), 180));
print "ok 7\n";
+use Math::Trig ':radial';
+
+{
+ my ($r,$t,$z) = cartesian_to_cylindrical(1,1,1);
+
+ print 'not ' unless (near($r, sqrt(2))) and
+ (near($t, deg2rad(45))) and
+ (near($z, 1));
+ print "ok 8\n";
+
+ ($x,$y,$z) = cylindrical_to_cartesian($r, $t, $z);
+
+ print 'not ' unless (near($x, 1)) and
+ (near($y, 1)) and
+ (near($z, 1));
+ print "ok 9\n";
+
+ ($r,$t,$z) = cartesian_to_cylindrical(1,1,0);
+
+ print 'not ' unless (near($r, sqrt(2))) and
+ (near($t, deg2rad(45))) and
+ (near($z, 0));
+ print "ok 10\n";
+
+ ($x,$y,$z) = cylindrical_to_cartesian($r, $t, $z);
+
+ print 'not ' unless (near($x, 1)) and
+ (near($y, 1)) and
+ (near($z, 0));
+ print "ok 11\n";
+}
+
+{
+ my ($r,$t,$f) = cartesian_to_spherical(1,1,1);
+
+ print 'not ' unless (near($r, sqrt(3))) and
+ (near($t, deg2rad(45))) and
+ (near($f, atan2(sqrt(2), 1)));
+ print "ok 12\n";
+
+ ($x,$y,$z) = spherical_to_cartesian($r, $t, $f);
+
+ print 'not ' unless (near($x, 1)) and
+ (near($y, 1)) and
+ (near($z, 1));
+ print "ok 13\n";
+
+ ($r,$t,$f) = cartesian_to_spherical(1,1,0);
+
+ print 'not ' unless (near($r, sqrt(2))) and
+ (near($t, deg2rad(45))) and
+ (near($f, deg2rad(90)));
+ print "ok 14\n";
+
+ ($x,$y,$z) = spherical_to_cartesian($r, $t, $f);
+
+ print 'not ' unless (near($x, 1)) and
+ (near($y, 1)) and
+ (near($z, 0));
+ print "ok 15\n";
+}
+
+{
+ my ($r,$t,$z) = cylindrical_to_spherical(spherical_to_cylindrical(1,1,1));
+
+ print 'not ' unless (near($r, 1)) and
+ (near($t, 1)) and
+ (near($z, 1));
+ print "ok 16\n";
+
+ ($r,$t,$z) = spherical_to_cylindrical(cylindrical_to_spherical(1,1,1));
+
+ print 'not ' unless (near($r, 1)) and
+ (near($t, 1)) and
+ (near($z, 1));
+ print "ok 17\n";
+}
+
+{
+ use Math::Trig 'great_circle_distance';
+
+ print 'not '
+ unless (near(great_circle_distance(0, 0, 0, pi/2), pi/2));
+ print "ok 18\n";
+
+ print 'not '
+ unless (near(great_circle_distance(0, 0, pi, pi), pi));
+ print "ok 19\n";
+
+ # London to Tokyo.
+ my @L = (deg2rad(-0.5), deg2rad(90 - 51.3));
+ my @T = (deg2rad(139.8),deg2rad(90 - 35.7));
+
+ my $km = great_circle_distance(@L, @T, 6378);
+
+ print 'not ' unless (near($km, 9605.26637021388));
+ print "ok 20\n";
+}
+
# eof
projects / p5sagit/p5-mst-13.2.git / commitdiff