use strict;

sub distance {
    my ( $lat1, $lon1, $lat2, $lon2, $unit ) = @_;
    my $theta = $lon1 - $lon2;
    my $dist
        = sin( deg2rad($lat1) ) * sin( deg2rad($lat2) )
        + cos( deg2rad($lat1) ) * cos( deg2rad($lat2) ) *
        cos( deg2rad($theta) );
    $dist = acos($dist);
    $dist = rad2deg($dist);
    my $miles = $dist * 60 * 1.1515;
    $unit = uc $unit;
    if ( $unit eq "K" ) {
        return ( $miles * 1.609344 );
    } ## end if ( $unit eq "K" )
    elsif ( $unit eq "N" ) {
        return ( $miles * 0.8684 );
    } ## end elsif ( $unit eq "N" )
    else {
        return $miles;
    } ## end else [ if ( $unit eq "K" )
} ## end sub distance

sub acos { atan2( sqrt( 1 - $_[0] * $_[0] ), $_[0] ) }

sub deg2rad { $_[0] * 0.0174532925 }

sub rad2deg { $_[0] * 57.2957795 }


##</code><code>##

print distance(60.1843160, -128.5069940, 29.46786, -98.53506, "M"), " miles\n";
print distance(32.9697, -96.80322, 29.46786, -98.53506, "K"), " kilometers\n";
print distance(32.9697, -96.80322, 29.46786, -98.53506, "N"), " nautical miles\n";
##</code><code>##

2535.33852270174 miles
422.738930838725 kilometers
228.109395841006 nautical miles
##</code><code>##

use GIS::Distance;
my $gis = GIS::Distance->new();
$gis->formula('Cosine');
print $gis->distance( 60.1843160, -128.5069940 => 29.46786, -98.53506 )
    ->value('mile'), " miles\n";
print $gis->distance( 32.9697, -96.80322 => 29.46786, -98.53506 )
    ->value('kilometer'), " kilometers\n";
print $gis->distance( 32.9697, -96.80322 => 29.46786, -98.53506 )
    ->value('nautical_mile'), " nautical miles\n";