Geo-CheapRuler
view release on metacpan or search on metacpan
lib/Geo/CheapRuler.pm view on Meta::CPAN
# ABSTRACT: Geo::CheapRuler Fast Geodesic (Lat,Lon) approximations for city scale / not near the poles. Port of mapbox/cheap-ruler.
package Geo::CheapRuler;
#
# how to update README.md
# pod2markdown CheapRuler.pm > ../../README.md
#
# $VERSION appears twice here, and once on GITHUB
#
our $VERSION = 'v0.1.1';
=head1 NAME
Geo::CheapRuler - Fast GPS Geodesic functions, port of mapbox/cheap-ruler
=head1 VERSION
v0.1.1
=head1 SYNOPSIS
A collection of very fast approximations to common geodesic measurements. Useful for performance-sensitive code that measures things on a city scale (less than 500km, not near the poles).
Can be an order of magnitude faster than Haversine based methods.
A Perl port of Mapbox's cheap-ruler v4.0.0 https://github.com/mapbox/cheap-ruler
Very fast as they use just 1 trig function per call.
=head1 MATHS MODEL
The Maths model is based upon an approximation to Vincenty's formulae, which uses the Earth's actual shape, an oblate ellipsoid (squashed sphere). For 'city' scale work, it is still more accurate than
the Haversine formulae (which uses several trig calls based upon a spherical Earth). For an explanation, see
https://blog.mapbox.com/fast-geodesic-approximations-with-cheap-ruler-106f229ad016
=head1 EXPORT
Nothing
=head1 WEBSITE
https://github.com/aavmurphy/CheapRuler
=head1 USAGE
This module uses "geojson" style GPS geometery. Points are [lon, lat]. Polygons are a series of rings. The first ring is exterior and clockwise. Subsequent rings are interior (holes) and anticlockwise.
The latitude (lat) parameter passed to the constructor should be the 'median' of the lat's used, i.e. (max-lat + min-lat)/2.
Some methods have units, e.g. "expand a bounding box by 10 meters/miles/kilometers". The default 'units' are 'kilometers', you may which to use 'meters'.
Data is passed / retured as arrayrefs, e.g. $p = [ 0.1, 54.1 ];
=head1 EXAMPLES
In the examples below, $p is a point, $a and $b are a line segment.
$p = [ -1, 57 ];
$a = [0.1, 54.1];
$b = [0.2, 54.2];
$ruler = Cheap::Ruler->new( ( 54.1 + 54.2 )/2, 'meters' );
# so the 'units' below are meters
$distance = $ruler->distance( $a, $b );
# return meters
$bearing = $ruler->bearing( $a, $b );
# returns degrees
$point = $ruler->destination( $a, 1000, 90);
# returns a new point, 1000 units away at 90 degrees
$point = $ruler->offset( $p, 100, 200 );
# returns a point 100 units east, 200 units north
$distance = $ruler->lineDistance( [ $p, $a, $b ] );
# length of the line
$area = $ruler->area( [
[-67.031, 50.458], [-67.031, 50.534], [-66.929, 50.534], [-66.929, 50.458], [-67.031, 50.458]
] ); # area of a polygon
$point = $ruler->along( [ [-67.031, 50.458], [-67.031, 50.534], [-66.929, 50.534] ], 2.5);
# returns a point 2.5 units along the line
$distance = $ruler->pointToSegmentDistance( $p, $a, $b );
# distance from point to a 2 point line segment
=cut
use strict;
lib/Geo/CheapRuler.pm view on Meta::CPAN
=head1 API
=head2 CheapRuler::fromTile( $y, $z, $units='kilometers' )
Creates a ruler object from Google web mercator tile coordinates (y and z). That's correct, y and z, not x.
See 'new' below for available units.
Example
$ruler = CheapRuler::fromTile( 11041, 15, 'meters');
=cut
sub fromTile( $y, $z, $units='kilometers') {
my $n = pi * (1 - 2 * ($y + 0.5) / ( 2**$z ));
my $lat = atan(0.5 * ( exp($n) - exp(-$n) )) / $RAD;
return CheapRuler->new($lat, $units);
}
=head2 units()
Multipliers for converting between units.
See 'new' below for available units.
Example : convert 50 meters to yards
$units = CheapRuler::units();
$yards = 50 * $units->{yards} / $units->{meters};
=cut
sub CheapRuler::units() {
return { %FACTORS };
}
=head2 CheapRuler->new( $lat, $units='kilometers' )
Create a ruler instance for very fast approximations to common geodesic measurements around a certain latitude.
param latitude, e.g. 54.31
param units (optional), one of: kilometers miles nauticalmiles meters metres yards feet inches
Example
$ruler = CheapRuler->new(35.05, 'miles');
=cut
sub new ( $class, $lat, $units='kilometers' ) {
croak( 'No latitude given.') if ! defined $lat;
croak( "Unknown unit $units. Use one of: " + join( ', ', keys(%FACTORS) ) ) if $units && ! $FACTORS{ $units };
my $self = bless {};
# Curvature formulas from https://en.wikipedia.org/wiki/Earth_radius#Meridional
my $m = $RAD * $RE * ( $units ? $FACTORS{ $units } : 1 );
my $coslat = cos( $lat * $RAD );
my $w2 = 1 / (1 - $E2 * (1 - $coslat**2) );
my $w = sqrt($w2);
# multipliers for converting longitude and latitude degrees into distance
$self->{kx} = $m * $w * $coslat; # based on normal radius of curvature
$self->{ky} = $m * $w * $w2 * (1 - $E2); # based on meridonal radius of curvature
return $self;
}
=head2 distance( $a, $b )
Given two points of the form [longitude, latitude], returns the distance in 'ruler' units.
param $a, point [longitude, latitude]
param $b, point [longitude, latitude]
returns distance (in chosen units)
Example
$distance = $ruler->distance([30.5, 50.5], [30.51, 50.49]);
=cut
sub distance( $self, $a, $b) {
my $dx = &wrap( $a->[0] - $b->[0]) * $self->{kx};
my $dy = ( $a->[1] - $b->[1]) * $self->{ky};
return sqrt( $dx * $dx + $dy * $dy);
}
=head2 bearing( $a, $b )
Returns the bearing between two points in degrees
param $a, point [longitude, latitude]
param $b, point [longitude, latitude]
returns $bearing (degrees)
Example
$bearing = $ruler->bearing([30.5, 50.5], [30.51, 50.49]);
=cut
sub bearing($self, $a, $b) {
my $dx = &wrap($b->[0] - $a->[0]) * $self->{kx};
my $dy = ($b->[1] - $a->[1]) * $self->{ky};
return atan2($dx, $dy) / $RAD;
}
=head2 destination( $point, $distance, $bearing)
Returns a new point given distance and bearing from the starting point.
lib/Geo/CheapRuler.pm view on Meta::CPAN
param $b, point, [ lon, lat ]
=cut
sub equals($a, $b) {
return ( $a->[0] == $b->[0] && $a->[1] == $b->[1] ) ? 1 : 0;
}
=head2 CheapRuler::interpolate( $a, $b, $t )
Returns a point along a line segment from $a to $b
a function not a method!
param $a, point, [lon, lat]
param $b, point, [lon, lat]
param $t, ratio (0 <= $t < 1 ), of the way along the line segment
returns $p, point [ lon, lat]
=cut
sub interpolate($a, $b, $t) {
my $dx = &wrap($b->[0] - $a->[0]);
my $dy = $b->[1] - $a->[1];
return [
$a->[0] + $dx * $t,
$a->[1] + $dy * $t
];
}
=head2 CheapRuler::normalize( $degrees )
Normalize a lon degree value into [-180..180] range
a function not a method!
param $degrees
return $degrees
=cut
sub wrap( $deg) {
while ( $deg < -180) { $deg += 360; }
while ( $deg > 180) { $deg -= 360; }
return $deg;
}
=head1 SEE ALSO
GIS::Distance
https://metacpan.org/pod/GIS::Distance
Group of modules with XS versions which implement Haversine and Vincenty distance formulas
Consider for longer than 'city scale' distances, or near the Poles.
GIS::Distance::Vincenty
The more accurate formulae which this module approximates. There is an XS version of it.
https://metacpan.org/pod/GIS::Distance::Vincenty
=head1 SUPPORT
You can find documentation for this module with the perldoc command.
perldoc Geo::CheapRuler
Or at
https://github.com/aavmurphy/CheapRuler
=head1 BUGS
Please report any bugs or feature requests of this port to
https://github.com/aavmurphy/CheapRuler
For the original, please see
https://github.com/mapbox/cheap-ruler
For testing, see
https://github.com/aavmurphy/CheapRuler/blob/main/test/test.pl
=head1 AUTHOR
Andrew Murphy, C<< <aavm at perl.org> >>
=head1 LICENSE AND COPYRIGHT
The original, mapbox/cheap-ruler, is (c) Mapbox.
This port is Copyright (c) 2025 by Andrew Murphy.
This is free software, licensed under:
The Artistic License 2.0 (GPL Compatible)
=head1 ACKNOWLEDGEMENTS
This module is a direct port of mapbox/cheap-ruler
=head1 DEVELOPER NOTES
README.md is auto-generated from Perl POD
The original's test were also ported. See the Github ./tests directory.
=cut
1; # End of Geo::CheapRuler
( run in 0.448 second using v1.01-cache-2.11-cpan-39bf76dae61 )