diff options
Diffstat (limited to 'systems/texlive/tlnet/tlpkg/tlperl/lib/Math/Trig.pm')
-rw-r--r-- | systems/texlive/tlnet/tlpkg/tlperl/lib/Math/Trig.pm | 94 |
1 files changed, 53 insertions, 41 deletions
diff --git a/systems/texlive/tlnet/tlpkg/tlperl/lib/Math/Trig.pm b/systems/texlive/tlnet/tlpkg/tlperl/lib/Math/Trig.pm index 1d9612a41c..218ab690a5 100644 --- a/systems/texlive/tlnet/tlpkg/tlperl/lib/Math/Trig.pm +++ b/systems/texlive/tlnet/tlpkg/tlperl/lib/Math/Trig.pm @@ -15,7 +15,7 @@ require Exporter; our @ISA = qw(Exporter); -our $VERSION = 1.23; +our $VERSION = 1.62; my @angcnv = qw(rad2deg rad2grad deg2rad deg2grad @@ -47,8 +47,9 @@ my @pi = qw(pi pi2 pi4 pip2 pip4); our @EXPORT_OK = (@rdlcnv, @greatcircle, @pi, 'Inf'); # See e.g. the following pages: -# http://www.movable-type.co.uk/scripts/LatLong.html -# http://williams.best.vwh.net/avform.htm +# https://www.movable-type.co.uk/scripts/latlong.html +# https://edwilliams.org/avform.htm +# https://en.wikipedia.org/wiki/Great-circle_distance our %EXPORT_TAGS = ('radial' => [ @rdlcnv ], 'great_circle' => [ @greatcircle ], @@ -153,12 +154,19 @@ sub great_circle_distance { $rho = 1 unless defined $rho; # Default to the unit sphere. - my $lat0 = pip2 - $phi0; - my $lat1 = pip2 - $phi1; + my $dphi = $phi1 - $phi0; + my $dtheta = $theta1 - $theta0; + + # A formula that is accurate for all distances is the following special + # case of the Vincenty formula for an ellipsoid with equal major and minor + # axes. See + # https://en.wikipedia.org/wiki/Great-circle_distance#Computational_formulas - return $rho * - acos_real( cos( $lat0 ) * cos( $lat1 ) * cos( $theta0 - $theta1 ) + - sin( $lat0 ) * sin( $lat1 ) ); + my $c1 = sin($phi1) * sin($dtheta); + my $c2 = sin($phi1) * cos($dtheta); + my $c3 = sin($phi0) * cos($phi1) - cos($phi0) * $c2; + my $c4 = cos($phi0) * cos($phi1) + sin($phi0) * $c2; + return $rho * atan2(sqrt($c1 * $c1 + $c3 * $c3), $c4); } sub great_circle_direction { @@ -247,7 +255,7 @@ Math::Trig - trigonometric functions $rad = deg2rad(120); - # Import constants pi2, pip2, pip4 (2*pi, pi/2, pi/4). + # Import constants pi2, pi4, pip2, pip4 (2*pi, 4*pi, pi/2, pi/4). use Math::Trig ':pi'; # Import the conversions between cartesian/spherical/cylindrical. @@ -417,7 +425,7 @@ and the imaginary part of approximately C<-1.317>. =back The full circle is 2 I<pi> radians or I<360> degrees or I<400> gradians. -The result is by default wrapped to be inside the [0, {2pi,360,400}[ circle. +The result is by default wrapped to be inside the [0, {2pi,360,400}] circle. If you don't want this, supply a true second argument: $zillions_of_radians = deg2rad($zillions_of_degrees, 1); @@ -465,15 +473,15 @@ B<All angles are in radians>. B<Cartesian> coordinates are the usual rectangular I<(x, y, z)>-coordinates. -Spherical coordinates, I<(rho, theta, pi)>, are three-dimensional +Spherical coordinates, I<(rho, theta, phi)>, 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 Gulf of Guinea (think of the missing big chunk of Africa) I<0, -pi/2, rho>. In geographical terms I<phi> is latitude (northward +I<polar> coordinate. The North Pole is therefore I<rho, 0, 0>, and +the Gulf of Guinea (think of the missing big chunk of Africa) I<rho, +0, pi/2>. In geographical terms I<phi> is latitude (northward positive, southward negative) and I<theta> is longitude (eastward positive, westward negative). @@ -537,26 +545,19 @@ points. =head2 great_circle_distance -You can compute spherical distances, called B<great circle distances>, -by importing the great_circle_distance() function: +Returns the great circle distance between two points on a sphere. - use Math::Trig 'great_circle_distance'; + $distance = great_circle_distance($theta0, $phi0, $theta1, $phi1, [, $rho]); - $distance = great_circle_distance($theta0, $phi0, $theta1, $phi1, [, $rho]); +Where ($theta0, $phi0) and ($theta1, $phi1) are the spherical coordinates of +the two points, respectively. The distance is in C<$rho> units. The C<$rho> +is optional. It defaults to 1 (the unit sphere). -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. - -If you think geographically the I<theta> are longitudes: zero at the -Greenwhich meridian, eastward positive, westward negative -- and the -I<phi> are latitudes: zero at the North Pole, northward positive, -southward negative. B<NOTE>: this formula thinks in mathematics, not -geographically: the I<phi> zero is at the North Pole, not at the -Equator on the west coast of Africa (Bay of Guinea). You need to -subtract your geographical coordinates from I<pi/2> (also known as 90 -degrees). +If you are using geographic coordinates, latitude and longitude, you need to +adjust for the fact that latitude is zero at the equator increasing towards +the north and decreasing towards the south. Assuming ($lat0, $lon0) and +($lat1, $lon1) are the geographic coordinates in radians of the two points, +the distance can be computed with $distance = great_circle_distance($lon0, pi/2 - $lat0, $lon1, pi/2 - $lat1, $rho); @@ -617,15 +618,22 @@ The great_circle_midpoint() is just a special case of ($thetai, $phii) = great_circle_waypoint($theta0, $phi0, $theta1, $phi1, $way); -Where the $way is a value from zero ($theta0, $phi0) to one ($theta1, -$phi1). Note that antipodal points (where their distance is I<pi> -radians) do not have waypoints between them (they would have an an -"equator" between them), and therefore C<undef> is returned for -antipodal points. If the points are the same and the distance -therefore zero and all waypoints therefore identical, the first point -(either point) is returned. +Where $way indicates the position of the waypoint along the great +circle arc through the starting point ($theta0, $phi0) and the end +point ($theta1, $phi1) relative to the distance from the starting +point to the end point. So $way = 0 gives the starting point, $way = 1 +gives the end point, $way < 0 gives a point "behind" the starting +point, and $way > 1 gives a point beyond the end point. $way defaults +to 0.5 if not given. + +Note that antipodal points (where their distance is I<pi> radians) do +not have unique waypoints between them, and therefore C<undef> is +returned in such cases. If the points are the same, so the distance +between them is zero, all waypoints are identical to the starting/end +point. -The thetas, phis, direction, and distance in the above are all in radians. +The thetas, phis, direction, and distance in the above are all in +radians. You can import all the great circle formulas by @@ -661,11 +669,13 @@ straight north being zero, straight east being pi/2). The midpoint between London and Tokyo being - use Math::Trig qw(great_circle_midpoint); + use Math::Trig qw(great_circle_midpoint rad2deg); my @M = great_circle_midpoint(@L, @T); + sub SWNE { rad2deg( $_[0] ), 90 - rad2deg( $_[1] ) } + my @lonlat = SWNE(@M); -or about 69 N 89 E, in the frozen wastes of Siberia. +or about 69 N 89 E, on the Putorana Plateau of Siberia. B<NOTE>: you B<cannot> get from A to B like this: @@ -743,6 +753,8 @@ an answer instead of giving a fatal runtime error. Do not attempt navigation using these formulas. +=head1 SEE ALSO + L<Math::Complex> =head1 AUTHORS |