Haversine formula: Difference between revisions
m
→{{header|Wren}}: Changed to Wren S/H
m (→{{header|REXX}}: changed some boxed comment wording.) |
m (→{{header|Wren}}: Changed to Wren S/H) |
||
| (199 intermediate revisions by 78 users not shown) | |||
Line 1:
{{task}}
{{Wikipedia}}
<br>
The '''haversine formula''' is an equation important in navigation, giving great-circle distances between two points on a sphere from their longitudes and latitudes.
It is a special case of a more general formula in spherical trigonometry, the '''law of haversines''', relating the sides and angles of spherical "triangles".
;Task:
Implement a great-circle distance function, or use a library function,
to show the great-circle distance between:
* Nashville International Airport (BNA) in Nashville, TN, USA, which is:
<big><big> '''N''' 36°7.2', '''W''' 86°40.2' (36.12, -86.67) </big></big> -and-
* Los Angeles International Airport (LAX) in Los Angeles, CA, USA, which is:
<big><big> '''N''' 33°56.4', '''W''' 118°24.0' (33.94, -118.40) </big></big>
<br>
<pre>
Line 31 ⟶ 37:
recommended that the latter value (r = 6372.8 km) be used (which
most of the given solutions have already adopted, anyways).
</pre>
Most of the examples below adopted Kaimbridge's recommended value of
6372.8 km for the earth radius. However, the derivation of this
[http://math.wikia.com/wiki/Ellipsoidal_quadratic_mean_radius ellipsoidal quadratic mean radius]
is wrong (the averaging over azimuth is biased). When applying these
examples in real applications, it is better to use the
[https://en.wikipedia.org/wiki/Earth_radius#Mean_radius mean earth radius],
6371 km. This value is recommended by the International Union of
Geodesy and Geophysics and it minimizes the RMS relative error between the
great circle and geodesic distance.
<br><br>
=={{header|11l}}==
{{trans|Python}}
<syntaxhighlight lang="11l">F haversine(=lat1, lon1, =lat2, lon2)
V r = 6372.8
V dLat = radians(lat2 - lat1)
V dLon = radians(lon2 - lon1)
lat1 = radians(lat1)
lat2 = radians(lat2)
V a = sin(dLat / 2) ^ 2 + cos(lat1) * cos(lat2) * sin(dLon / 2) ^ 2
V c = 2 * asin(sqrt(a))
R r * c
print(haversine(36.12, -86.67, 33.94, -118.40))</syntaxhighlight>
{{out}}
<pre>
2887.26
</pre>
=={{header|ABAP}}==
<syntaxhighlight lang="abap">
DATA: X1 TYPE F, Y1 TYPE F,
X2 TYPE F, Y2 TYPE F, YD TYPE F,
PI TYPE F,
PI_180 TYPE F,
MINUS_1 TYPE F VALUE '-1'.
PI = ACOS( MINUS_1 ).
PI_180 = PI / 180.
LATITUDE1 = 36,12 . LONGITUDE1 = -86,67 .
LATITUDE2 = 33,94 . LONGITUDE2 = -118,4 .
X1 = LATITUDE1 * PI_180.
Y1 = LONGITUDE1 * PI_180.
X2 = LATITUDE2 * PI_180.
Y2 = LONGITUDE2 * PI_180.
YD = Y2 - Y1.
DISTANCE = 20000 / PI *
ACOS( SIN( X1 ) * SIN( X2 ) + COS( X1 ) * COS( X2 ) * COS( YD ) ).
WRITE : 'Distance between given points = ' , distance , 'km .' .
</syntaxhighlight>
{{out}}
<pre>
Distance between given points = 2.884,2687 km .
</pre>
=={{header|Ada}}==
<
with Ada.Long_Float_Text_IO; use Ada.Long_Float_Text_IO;
with Ada.Numerics.Generic_Elementary_Functions;
Line 64 ⟶ 133:
DMS_To_Radians (33.0, 56.4), DMS_To_Radians (118.0, 24.0)), -- Los Angeles International Airport (LAX)
Aft=>3, Exp=>0);
end Haversine_Formula;</
=={{header|ALGOL 68}}==
Line 71 ⟶ 140:
{{works with|ALGOL 68G|Any - tested with release [http://sourceforge.net/projects/algol68/files/algol68g/algol68g-2.3.5 algol68g-2.3.5].}}
{{wont work with|ELLA ALGOL 68|Any (with appropriate job cards) - tested with release [http://sourceforge.net/projects/algol68/files/algol68toc/algol68toc-1.8.8d/algol68toc-1.8-8d.fc9.i386.rpm/download 1.8-8d] - due to extensive use of '''format'''[ted] ''transput''.}}
'''File: Haversine_formula.a68'''<
REAL r = 20 000/pi + 6.6 # km #,
Line 92 ⟶ 161:
# Americans don't know kilometers #
printf(($"dist: "g(0,1)" km ("g(0,1)" mi.)"l$, d, d / 1.609344))
)</
{{out}}
<pre>
dist: 2887.3 km (1794.1 mi.)
</pre>
=={{header|ALGOL W}}==
{{Trans|ALGOL 68}}
Using the mean radius value suggested in the task.
<syntaxhighlight lang="algolw">begin % compute the distance between places using the Haversine formula %
real procedure arcsin( real value x ) ; arctan( x / sqrt( 1 - ( x * x ) ) );
real procedure distance ( real value th1Deg, ph1Deg, th2Deg, ph2Deg ) ;
begin
real ph1, th1, th2, toRad, dz, dx, dy;
toRad := pi / 180;
ph1 := ( ph1Deg - ph2Deg ) * toRad;
th1 := th1Deg * toRad;
th2 := th2Deg * toRad;
dz := sin( th1 ) - sin( th2 );
dx := cos( ph1 ) * cos( th1 ) - cos( th2 );
dy := sin( ph1 ) * cos( th1 );
arcsin( sqrt( dx * dx + dy * dy + dz * dz ) / 2 ) * 2 * 6371
end distance ;
begin
real d;
integer mi, km;
d := distance( 36.12, -86.67, 33.94, -118.4 );
mi := round( d );
km := round( d / 1.609344 );
writeon( i_w := 4, s_w := 0, "distance: ", mi, " km (", km, " mi.)" )
end
end.</syntaxhighlight>
{{out}}
<pre>
distance: 2886 km (1794 mi.)
</pre>
=={{header|Amazing Hopper}}==
{{Trans|Python}}
<syntaxhighlight lang="c">
/* fórmula de Haversine para distancias en una
superficie esférica */
#include <basico.h>
#define MIN 60
#define SEG 3600
#define RADIO 6372.8
#define UNAMILLA 1.609344
algoritmo
números ( lat1, lon1, lat2, lon2, dlat, dlon, millas )
si ' total argumentos es (9) ' // LAT1 M LON1 M LAT2 M LON2 M
#basic{ lat1 = narg(2) + narg(3)/MIN
lon1 = narg(4) + narg(5)/MIN
lat2 = narg(6) + narg(7)/MIN
lon2 = narg(8) + narg(9)/MIN }
sino si ' total argumentos es (13) ' // LAT1 M LON1 M S LAT2 M LON2 M S
#basic{ lat1 = narg(2) + narg(3)/MIN + narg(4)/SEG
lon1 = narg(5) + narg(6)/MIN + narg(7)/SEG
lat2 = narg(8) + narg(9)/MIN + narg(10)/SEG
lon2 = narg(11) + narg(12)/MIN + narg(13)/SEG }
sino
imprimir("Modo de uso:\ndist.bas La1 M [S] Lo1 M [S] La2 M [S] Lo2 M [S]\n")
término prematuro
fin si
#basic{
dlat = sin(radian(lat2 - lat1)/2)^2
dlon = sin(radian(lon2 - lon1)/2)^2
RADIO*(2*arc sin(sqrt(dlat + cos(radian(lat1)) * cos(radian(lat2)) * dlon )))
}
---copiar en 'millas'---, " km. (",millas,dividido por (UNA MILLA)," mi.)\n"
decimales '2', imprimir
terminar
</syntaxhighlight>
{{out}}
<pre>
$ hopper3 basica/dist.bas -x -o bin/dist
Generating binary ‘bin/dist’...Ok!
Symbols: 27
Total size: 0.75 Kb
$ ./bin/dist 36 7.2 86 40.2 33 56.4 118 24
2887.26 km. (1794.06 mi.)
</pre>
=={{header|AMPL}}==
<syntaxhighlight lang="ampl">
set location;
set geo;
param coord{i in location, j in geo};
param dist{i in location, j in location};
data;
set location := BNA LAX;
set geo := LAT LON;
param coord:
LAT LON :=
BNA 36.12 -86.67
LAX 33.94 -118.4
;
let dist['BNA','LAX'] := 2 * 6372.8 * asin (sqrt(sin(atan(1)/45*(coord['LAX','LAT']-coord['BNA','LAT'])/2)^2 + cos(atan(1)/45*coord['BNA','LAT']) * cos(atan(1)/45*coord['LAX','LAT']) * sin(atan(1)/45*(coord['LAX','LON'] - coord
['BNA','LON'])/2)^2));
printf "The distance between the two points is approximately %f km.\n", dist['BNA','LAX'];
</syntaxhighlight>
{{out}}
<pre>
The distance between the two points is approximately 2887.259951 km.
</pre>
=={{header|APL}}==
<syntaxhighlight lang="apl">r←6371
hf←{(p q)←○⍺ ⍵÷180 ⋄ 2×rׯ1○(+/(2*⍨1○(p-q)÷2)×1(×/2○⊃¨p q))*÷2}
36.12 ¯86.67 hf 33.94 ¯118.40</syntaxhighlight>
{{out}}
<pre>2886.44</pre>
=={{header|AppleScript}}==
AppleScript provides no trigonometric functions.
Here we reach through a foreign function interface to a temporary instance of a JavaScript interpreter.
<syntaxhighlight lang="applescript">use AppleScript version "2.4" -- Yosemite (10.10) or later
use framework "Foundation"
use framework "JavaScriptCore"
use scripting additions
property js : missing value
-- haversine :: (Num, Num) -> (Num, Num) -> Num
on haversine(latLong, latLong2)
set {lat, lon} to latLong
set {lat2, lon2} to latLong2
set {rlat1, rlat2, rlon1, rlon2} to ¬
map(my radians, {lat, lat2, lon, lon2})
set dLat to rlat2 - rlat1
set dLon to rlon2 - rlon1
set radius to 6372.8 -- km
set asin to math("asin")
set sin to math("sin")
set cos to math("cos")
|round|((2 * radius * ¬
(asin's |λ|((sqrt(((sin's |λ|(dLat / 2)) ^ 2) + ¬
(((sin's |λ|(dLon / 2)) ^ 2) * ¬
(cos's |λ|(rlat1)) * (cos's |λ|(rlat2)))))))) * 100) / 100
end haversine
-- math :: String -> Num -> Num
on math(f)
script
on |λ|(x)
if missing value is js then ¬
set js to current application's JSContext's new()
(js's evaluateScript:("Math." & f & "(" & x & ")"))'s toDouble()
end |λ|
end script
end math
-------------------------- TEST ---------------------------
on run
set distance to haversine({36.12, -86.67}, {33.94, -118.4})
set js to missing value -- Clearing a c pointer.
return distance
end run
-------------------- GENERIC FUNCTIONS --------------------
-- map :: (a -> b) -> [a] -> [b]
on map(f, xs)
-- The list obtained by applying f
-- to each element of xs.
tell mReturn(f)
set lng to length of xs
set lst to {}
repeat with i from 1 to lng
set end of lst to |λ|(item i of xs, i, xs)
end repeat
return lst
end tell
end map
-- mReturn :: First-class m => (a -> b) -> m (a -> b)
on mReturn(f)
-- 2nd class handler function lifted into 1st class script wrapper.
if script is class of f then
f
else
script
property |λ| : f
end script
end if
end mReturn
-- radians :: Float x => Degrees x -> Radians x
on radians(x)
(pi / 180) * x
end radians
-- round :: a -> Int
on |round|(n)
round n
end |round|
-- sqrt :: Num -> Num
on sqrt(n)
if n ≥ 0 then
n ^ (1 / 2)
else
missing value
end if
end sqrt</syntaxhighlight>
{{Out}}
<pre>2887.26</pre>
=={{header|Arturo}}==
<syntaxhighlight lang="rebol">radians: function [x]-> x * pi // 180
haversine: function [src,tgt][
dLat: radians tgt\0 - src\0
dLon: radians tgt\1 - src\1
lat1: radians src\0
lat2: radians tgt\0
a: add product @[cos lat1, cos lat2, sin dLon/2, sin dLon/2] (sin dLat/2) ^ 2
c: 2 * asin sqrt a
return 6372.8 * c
]
print haversine @[36.12 neg 86.67] @[33.94, neg 118.40]</syntaxhighlight>
{{out}}
<pre>2887.259950607111</pre>
=={{header|ATS}}==
<syntaxhighlight lang="ats">
#include
"share/atspre_staload.hats"
Line 136 ⟶ 459:
$extfcall(void, "printf", "dist: %.1f km (%.1f mi.)\n", d, d / 1.609344)
end // end of [main0]
</syntaxhighlight>
{{out}}
Line 144 ⟶ 467:
=={{header|AutoHotkey}}==
<
GreatCircleDist(La1, La2, Lo1, Lo2, R, U) {
Line 156 ⟶ 479:
Rad(Deg) {
return, Deg * 4 * ATan(1) / 180
}</
{{out}}
<pre>2887.259951 km</pre>
=={{header|AWK}}==
<syntaxhighlight lang="awk">
# syntax: GAWK -f HAVERSINE_FORMULA.AWK
# converted from Python
Line 180 ⟶ 503:
return degree * (3.1415926 / 180.)
}
</syntaxhighlight>
{{out}}
<pre>distance: 2887.2599 km</pre>
=={{header|BASIC}}==
==={{header|Applesoft BASIC}}===
{{works with|Chipmunk Basic}}
{{works with|GW-BASIC}}
{{works with|Minimal BASIC}}
{{works with|MSX BASIC}}
{{trans|Commodore BASIC}}
<syntaxhighlight lang="qbasic">100 HOME : rem 100 CLS for GW-BASIC and MSX BASIC : DELETE for Minimal BASIC
110 LET P = ATN(1)*4
120 LET D = P/180
130 LET M = 36.12
140 LET K = -86.67
150 LET N = 33.94
160 LET L = -118.4
170 LET R = 6372.8
180 PRINT " DISTANCIA DE HAVERSINE ENTRE BNA Y LAX = ";
190 LET A = SIN((L-K)*D/2)
200 LET A = A*A
210 LET B = COS(M*D)*COS(N*D)
220 LET C = SIN((N-M)*D/2)
230 LET C = C*C
240 LET D = SQR(C+B*A)
250 LET E = D/SQR(1-D*D)
260 LET F = ATN(E)
270 PRINT 2*R*F;"KM"
280 END</syntaxhighlight>
==={{header|BASIC256}}===
<syntaxhighlight lang="basic256">
global radioTierra # radio de la tierra en km
radioTierra = 6372.8
function Haversine(lat1, long1, lat2, long2 , radio)
d_long = radians(long1 - long2)
theta1 = radians(lat1)
theta2 = radians(lat2)
dx = cos(d_long) * cos(theta1) - cos(theta2)
dy = sin(d_long) * cos(theta1)
dz = sin(theta1) - sin(theta2)
return asin(sqr(dx*dx + dy*dy + dz*dz) / 2) * radio * 2
end function
print
print " Distancia de Haversine entre BNA y LAX = ";
print Haversine(36.12, -86.67, 33.94, -118.4, radioTierra); " km"
end
</syntaxhighlight>
{{out}}
<pre> Distancia de Haversine entre BNA y LAX = 2887.25994877 km.</pre>
==={{header|Chipmunk Basic}}===
{{works with|Chipmunk Basic|3.6.4}}
<syntaxhighlight lang="qbasic">100 cls
110 pi = arctan(1)*4 : rem define pi = 3.1415...
120 deg2rad = pi/180 : rem define grados a radianes 0.01745..
130 lat1 = 36.12
140 long1 = -86.67
150 lat2 = 33.94
160 long2 = -118.4
170 radio = 6372.8
180 print " Distancia de Haversine entre BNA y LAX = ";
190 d_long = deg2rad*(long1-long2)
200 theta1 = deg2rad*(lat1)
210 theta2 = deg2rad*(lat2)
220 dx = cos(d_long)*cos(theta1)-cos(theta2)
230 dy = sin(d_long)*cos(theta1)
240 dz = sin(theta1)-sin(theta2)
250 print (asin(sqr(dx*dx+dy*dy+dz*dz)/2)*radio*2);"km"
260 end</syntaxhighlight>
==={{header|Gambas}}===
<syntaxhighlight lang="vbnet">Public deg2rad As Float = Pi / 180 ' define grados a radianes 0.01745..
Public radioTierra As Float = 6372.8 ' radio de la tierra en km
Public Sub Main()
Print "\n Distancia de Haversine entre BNA y LAX = "; Haversine(36.12, -86.67, 33.94, -118.4, radioTierra); " km"
End
Function Haversine(lat1 As Float, long1 As Float, lat2 As Float, long2 As Float, radio As Float) As Float
Dim d_long As Float = deg2rad * (long1 - long2)
Dim theta1 As Float = deg2rad * lat1
Dim theta2 As Float = deg2rad * lat2
Dim dx As Float = Cos(d_long) * Cos(theta1) - Cos(theta2)
Dim dy As Float = Sin(d_long) * Cos(theta1)
Dim dz As Float = Sin(theta1) - Sin(theta2)
Return ASin(Sqr(dx * dx + dy * dy + dz * dz) / 2) * radio * 2
End Function</syntaxhighlight>
{{out}}
<pre>Same as FreeBASIC entry.</pre>
==={{header|GW-BASIC}}===
{{works with|Applesoft BASIC}}
{{works with|BASICA}}
{{works with|Chipmunk Basic}}
{{works with|MSX BASIC}}
{{works with|Minimal BASIC}}
{{works with|PC-BASIC|any}}
{{trans|Commodore BASIC}}
<syntaxhighlight lang="qbasic">100 CLS : rem 100 HOME for Applesoft BASIC : DELETE for Minimal BASIC
110 LET P = ATN(1)*4
120 LET D = P/180
130 LET M = 36.12
140 LET K = -86.67
150 LET N = 33.94
160 LET L = -118.4
170 LET R = 6372.8
180 PRINT " DISTANCIA DE HAVERSINE ENTRE BNA Y LAX = ";
190 LET A = SIN((L-K)*D/2)
200 LET A = A*A
210 LET B = COS(M*D)*COS(N*D)
220 LET C = SIN((N-M)*D/2)
230 LET C = C*C
240 LET D = SQR(C+B*A)
250 LET E = D/SQR(1-D*D)
260 LET F = ATN(E)
270 PRINT 2*R*F;"KM"
280 END</syntaxhighlight>
==={{header|Minimal BASIC}}===
{{works with|Applesoft BASIC}}
{{works with|Chipmunk Basic}}
{{works with|GW-BASIC}}
{{works with|MSX BASIC}}
{{trans|Commodore BASIC}}
<syntaxhighlight lang="qbasic">110 LET P = ATN(1)*4
120 LET D = P/180
130 LET M = 36.12
140 LET K = -86.67
150 LET N = 33.94
160 LET L = -118.4
170 LET R = 6372.8
180 PRINT " DISTANCIA DE HAVERSINE ENTRE BNA Y LAX = ";
190 LET A = SIN((L-K)*D/2)
200 LET A = A*A
210 LET B = COS(M*D)*COS(N*D)
220 LET C = SIN((N-M)*D/2)
230 LET C = C*C
240 LET D = SQR(C+B*A)
250 LET E = D/SQR(1-D*D)
260 LET F = ATN(E)
270 PRINT 2*R*F;"KM"
280 END</syntaxhighlight>
==={{header|MSX Basic}}===
The [[#GW-BASIC|GW-BASIC]] solution works without any changes.
==={{header|QBasic}}===
{{works with|QBasic|1.1}}
{{works with|QuickBasic|4.5}}
<syntaxhighlight lang="basic">
CONST pi = 3.141593 ' define pi
CONST radio = 6372.8 ' radio de la tierra en km
PRINT : PRINT " Distancia de Haversine:";
PRINT Haversine!(36.12, -86.67, 33.94, -118.4); "km"
END
FUNCTION Haversine! (lat1!, long1!, lat2!, long2!)
deg2rad! = pi / 180 ' define grados a radianes 0.01745..
dLong! = deg2rad! * (long1! - long2!)
theta1! = deg2rad! * lat1!
theta2! = deg2rad! * lat2!
dx! = COS(dLong!) * COS(theta1!) - COS(theta2!)
dy! = SIN(dLong!) * COS(theta1!)
dz! = SIN(theta1!) - SIN(theta2!)
Haversine! = (SQR(dx! * dx! + dy! * dy! + dz! * dz!) / 2) * radio * 2
END FUNCTION</syntaxhighlight>
{{out}}
<pre> Distancia de Haversine: 2862.63 km</pre>
==={{header|Quite BASIC}}===
<syntaxhighlight lang="qbasic">100 CLS
110 LET p = atan(1)*4
120 LET d = p/180
130 LET k = 36.12
140 LET m = -86.67
150 LET l = 33.94
160 LET n = -118.4
170 LET r = 6372.8
180 PRINT " Distancia de Haversine entre BNA y LAX = ";
190 LET g = d*(m-n)
200 LET t = d*(k)
210 LET s = d*(l)
220 LET x = COS(g)*COS(t)-COS(s)
230 LET y = SIN(g)*COS(t)
240 LET z = SIN(t)-SIN(s)
250 PRINT (ASIN(SQR(x*x+y*y+z*z)/2)*r*2);"km"
260 END</syntaxhighlight>
==={{header|True BASIC}}===
<syntaxhighlight lang="basic">DEF Haversine (lat1, long1, lat2, long2)
OPTION ANGLE RADIANS
LET R = 6372.8 !radio terrestre en km.
LET dLat = RAD(lat2-lat1)
LET dLong = RAD(long2-long1)
LET lat1 = RAD(lat1)
LET lat2 = RAD(lat2)
LET Haversine = R *2 * ASIN(SQR(SIN(dLat/2)^2 + SIN(dLong/2)^2 *COS(lat1) * COS(lat2)))
END DEF
PRINT
PRINT "Distancia de Haversine:"; Haversine(36.12, -86.67, 33.94, -118.4); "km"
END</syntaxhighlight>
{{out}}
<pre>Distancia de Haversine: 2887.26 km</pre>
==={{header|Yabasic}}===
{{trans|FreeBASIC}}
<syntaxhighlight lang="yabasic">//pi está predefinido en Yabasic
deg2rad = pi / 180 // define grados a radianes 0.01745..
radioTierra = 6372.8 // radio de la tierra en km
sub Haversine(lat1, long1, lat2, long2 , radio)
d_long = deg2rad * (long1 - long2)
theta1 = deg2rad * lat1
theta2 = deg2rad * lat2
dx = cos(d_long) * cos(theta1) - cos(theta2)
dy = sin(d_long) * cos(theta1)
dz = sin(theta1) - sin(theta2)
return asin(sqr(dx*dx + dy*dy + dz*dz) / 2) * radio * 2
end sub
print " Distancia de Haversine entre BNA y LAX = ", Haversine(36.12, -86.67, 33.94, -118.4, radioTierra), " km"
end</syntaxhighlight>
{{out}}
<pre> Distancia de Haversine entre BNA y LAX = 259.478 km</pre>
=={{header|BBC BASIC}}==
Line 190 ⟶ 745:
Uses BBC BASIC's '''MOD(array())''' function which calculates
the square-root of the sum of the squares of the elements of an array.
<
END
Line 198 ⟶ 753:
\ SINRAD(e1-e2) * COSRAD(n1), \
\ SINRAD(n1) - SINRAD(n2)
= ASN(MOD(d()) / 2) * 6372.8 * 2</
{{out}}
<pre>
Distance = 2887.25995 km
</pre>
=={{header|bc}}==
{{works with|GNU_bc|GNU bc(1)|1.07.1-2}}
{{works with|OpenBSD_bc|dc(1)-based MirBSD bc(1)|as of 2021-06-11}}
{{works with|Bc|bc(1)|POSIX}}
… supplied with a small POSIX shell wrapper to feed the arguments to <code>bc</code>. ([https://edugit.org/-/snippets/27 see also])
<syntaxhighlight lang="sh">
#!/bin/sh
#-
# © 2021 mirabilos Ⓕ CC0; implementation of Haversine GCD from public sources
#
# now developed online:
# https://evolvis.org/plugins/scmgit/cgi-bin/gitweb.cgi?p=useful-scripts/mirkarte.git;a=blob;f=geo.sh;hb=HEAD
if test "$#" -ne 4; then
echo >&2 "E: syntax: $0 lat1 lon1 lat2 lon2"
exit 1
fi
set -e
# make GNU bc use POSIX mode and shut up
BC_ENV_ARGS=-qs
export BC_ENV_ARGS
# assignment of constants, variables and functions
# p: multiply with to convert from degrees to radians (π/180)
# r: earth radius in metres
# d: distance
# h: haversine intermediate
# i,j: (lat,lon) point 1
# x,y: (lat,lon) point 2
# k: delta lat
# l: delta lon
# m: sin(k/2) (square root of hav(k))
# n: sin(l/2) ( partial haversine )
# n(x): arcsin(x)
# r(x,n): round x to n decimal digits
# v(x): sign (Vorzeichen)
# w(x): min(1, sqrt(x)) (Wurzel)
bc -l <<-EOF
scale=64
define n(x) {
if (x == -1) return (-2 * a(1))
if (x == 1) return (2 * a(1))
return (a(x / sqrt(1 - x*x)))
}
define v(x) {
if (x < 0) return (-1)
if (x > 0) return (1)
return (0)
}
define r(x, n) {
auto o
o = scale
if (scale < (n + 1)) scale = (n + 1)
x += v(x) * 0.5 * A^-n
scale = n
x /= 1
scale = o
return (x)
}
define w(x) {
if (x >= 1) return (1)
return (sqrt(x))
}
/* WGS84 reference ellipsoid: große Halbachse (metres), Abplattung */
i = 6378137.000
x = 1/298.257223563
/* other axis */
j = i * (1 - x)
/* mean radius resulting */
r = (2 * i + j) / 3
/* coordinates */
p = (4 * a(1) / 180)
i = (p * $1)
j = (p * $2)
x = (p * $3)
y = (p * $4)
/* calculation */
k = (x - i)
l = (y - j)
m = s(k / 2)
n = s(l / 2)
h = ((m * m) + (c(i) * c(x) * n * n))
d = 2 * r * n(w(h))
r(d, 3)
EOF
# output is in metres, rounded to millimetres, error < ¼% in WGS84
</syntaxhighlight>
{{out}}
<pre>$ sh dist.sh 36.12 -86.67 33.94 -118.4
2886448.430</pre>
Note I used a more precise earth radius; this matches the other implementations when choosing the same.
=={{header|C}}==
<
#include <stdlib.h>
#include <math.h>
Line 230 ⟶ 883:
return 0;
}</
=={{header|c sharp|C#}}==
{{trans|Groovy}}
<
public static double calculate(double lat1, double lon1, double lat2, double lon2) {
var R = 6372.8; // In kilometers
Line 257 ⟶ 910:
// Returns: The distance between coordinates 36.12,-86.67 and 33.94,-118.4 is: 2887.25995060711
</syntaxhighlight>
=={{header|C++}}==
<syntaxhighlight lang="cpp">
#define _USE_MATH_DEFINES
#include <math.h>
#include <iostream>
const static double EarthRadiusKm = 6372.8;
inline double DegreeToRadian(double angle)
{
return M_PI * angle / 180.0;
}
class Coordinate
{
public:
Coordinate(double latitude ,double longitude):myLatitude(latitude), myLongitude(longitude)
{}
double Latitude() const
{
return myLatitude;
}
double Longitude() const
{
return myLongitude;
}
private:
double myLatitude;
double myLongitude;
};
double HaversineDistance(const Coordinate& p1, const Coordinate& p2)
{
double latRad1 = DegreeToRadian(p1.Latitude());
double latRad2 = DegreeToRadian(p2.Latitude());
double lonRad1 = DegreeToRadian(p1.Longitude());
double lonRad2 = DegreeToRadian(p2.Longitude());
double diffLa = latRad2 - latRad1;
double doffLo = lonRad2 - lonRad1;
double computation = asin(sqrt(sin(diffLa / 2) * sin(diffLa / 2) + cos(latRad1) * cos(latRad2) * sin(doffLo / 2) * sin(doffLo / 2)));
return 2 * EarthRadiusKm * computation;
}
int main()
{
Coordinate c1(36.12, -86.67);
Coordinate c2(33.94, -118.4);
std::cout << "Distance = " << HaversineDistance(c1, c2) << std::endl;
return 0;
}
</syntaxhighlight>
=={{header|clojure|Clojure}}==
{{trans|Java}}
<
(defn haversine
[{lon1 :longitude lat1 :latitude} {lon2 :longitude lat2 :latitude}]
Line 274 ⟶ 987:
(haversine {:latitude 36.12 :longitude -86.67} {:latitude 33.94 :longitude -118.40})
;=> 2887.2599506071106
</syntaxhighlight>
=={{header|CoffeeScript}}==
{{trans|JavaScript}}
<
R = 6372.8; # km
radians = args.map (deg) -> deg/180.0 * Math.PI
Line 287 ⟶ 1,000:
R * 2 * Math.asin(Math.sqrt(a))
console.log haversine(36.12, -86.67, 33.94, -118.40)</
{{out}}
<pre>2887.2599506071124</pre>
=={{header|Commodore BASIC}}==
PETSCII has the pi symbol <tt>π</tt> in place of the ASCII tilde <tt>~</tt>; Commodore BASIC interprets this symbol as the mathematical constant.
<syntaxhighlight lang="basic">10 REM================================
15 REM HAVERSINE FORMULA
20 REM
25 REM 2021-09-24
30 REM EN.WIKIPEDIA.ORG/WIKI/HAVERSINE_FORMULA
35 REM
40 REM C64 HAS PI CONSTANT
45 REM X1 LONGITUDE 1
50 REM Y1 LATITUDE 1
55 REM X2 LONGITUDE 2
60 REM Y2 LATITUDE 2
65 REM
70 REM V1, 2021-10-02, ALVALONGO
75 REM ===============================
100 REM MAIN
105 DR=π/180:REM DEGREES TO RADIANS
110 PRINT CHR$(147);CHR$(5);"HAVERSINE FORMULA"
120 PRINT "GREAT-CIRCLE DISTANCE"
130 R=6372.8:REM AVERAGE EARTH RADIUS IN KILOMETERS
200 REM GET DATA
210 PRINT
220 INPUT "LONGITUDE 1=";X1
230 INPUT "LATITUDE 1=";Y1
240 PRINT
250 INPUT "LONGITUDE 2=";X2
260 INPUT "LATITUDE 2=";Y2
270 GOSUB 500
280 PRINT
290 PRINT "DISTANCE=";D;"KM"
300 GET K$:IF K$="" THEN 300
310 GOTO 210
490 END
500 REM HAVERSINE FORMULA ------------
520 A=SIN((X2-X1)*DR/2)
530 A=A*A
540 B=COS(Y1*DR)*COS(Y2*DR)
550 C=SIN((Y2-Y1)*DR/2)
560 C=C*C
570 D=SQR(C+B*A)
580 E=D/SQR(1-D*D)
590 F=ATN(E)
600 D=2*R*F
610 RETURN</syntaxhighlight>
{{Out}}
<pre>HAVERSINE FORMULA
GREAT-CIRCLE DISTANCE
LONGITUDE 1=? -86.67
LATITUDE 1=? 36.12
LONGITUDE 2=? -118.40
LATITUDE 2=? 33.94
DISTANCE= 2887.25995 KM</pre>
=={{header|Common Lisp}}==
<
(defparameter *rad-conv* (/ pi 180))
Line 312 ⟶ 1,084:
(deg->rad lng1)
(deg->rad lat2)
(deg->rad lng2)))</
{{out}}
<pre>CL-USER> (format t "~%The distance between BNA and LAX is about ~$ km.~%"
Line 318 ⟶ 1,090:
The distance between BNA and LAX is about 2887.26 km.</pre>
=={{header|Crystal}}==
{{trans|Python}}
<syntaxhighlight lang="ruby">include Math
def haversine(lat1, lon1, lat2, lon2)
r = 6372.8 # Earth radius in kilometers
deg2rad = PI/180 # convert degress to radians
dLat = (lat2 - lat1) * deg2rad
dLon = (lon2 - lon1) * deg2rad
lat1 = lat1 * deg2rad
lat2 = lat2 * deg2rad
a = sin(dLat / 2)**2 + cos(lat1) * cos(lat2) * sin(dLon / 2)**2
c = 2 * asin(sqrt(a))
r * c
end
puts "distance is #{haversine(36.12, -86.67, 33.94, -118.40)} km "
</syntaxhighlight>
{{out}}
<pre>
distance is 2887.2599506071106 km
</pre>
=={{header|D}}==
<
real haversineDistance(in real dth1, in real dph1,
Line 343 ⟶ 1,140:
writefln("Haversine distance: %.1f km",
haversineDistance(36.12, -86.67, 33.94, -118.4));
}</
{{out}}
<pre>Haversine distance: 2887.3 km</pre>
Line 350 ⟶ 1,147:
An alternate direct implementation of the haversine formula as shown at [[wp:Haversine formula|wikipedia]]. The same length, but perhaps a little more clear about what is being done.
<
real toRad(in real degrees) pure nothrow @safe @nogc {
Line 376 ⟶ 1,173:
greatCircleDistance(36.12, -86.67, 33.94, -118.4,
earthRadius));
}</
{{out}}
<pre>Great circle distance: 2887.3 km
</pre>
=={{header|Dart}}==
{{trans|Java}}
<syntaxhighlight lang="dart">import 'dart:math';
class Haversine {
static final R = 6372.8; // In kilometers
static double haversine(double lat1, lon1, lat2, lon2) {
double dLat = _toRadians(lat2 - lat1);
double dLon = _toRadians(lon2 - lon1);
lat1 = _toRadians(lat1);
lat2 = _toRadians(lat2);
double a = pow(sin(dLat / 2), 2) + pow(sin(dLon / 2), 2) * cos(lat1) * cos(lat2);
double c = 2 * asin(sqrt(a));
return R * c;
}
static double _toRadians(double degree) {
return degree * pi / 180;
}
static void main() {
print(haversine(36.12, -86.67, 33.94, -118.40));
}
}
</syntaxhighlight>
{{out}}
<pre>2887.2599506071106</pre>
=={{header|Delphi}}==
<
uses Math;
Line 402 ⟶ 1,227:
begin
Writeln('Haversine distance: ', HaversineDist(36.12, -86.67, 33.94, -118.4):7:2, ' km.');
end.</
{{out}}
<pre>Haversine distance: 2887.26 km.</pre>
=={{header|EasyLang}}==
{{trans|C}}
<syntaxhighlight>
func dist th1 ph1 th2 ph2 .
r = 6371
ph1 -= ph2
dz = sin th1 - sin th2
dx = cos ph1 * cos th1 - cos th2
dy = sin ph1 * cos th1
return 2 * r * pi / 180 * asin (sqrt (dx * dx + dy * dy + dz * dz) / 2)
.
print dist 36.12 -86.67 33.94 -118.4
</syntaxhighlight>
=={{header|Elena}}==
ELENA 4.x:
<syntaxhighlight lang="elena">import extensions;
import system'math;
Haversine(lat1,lon1,lat2,lon2)
{
var R := 6372.8r;
var dLat := (lat2 - lat1).Radian;
var dLon := (lon2 - lon1).Radian;
var dLat1 := lat1.Radian;
var dLat2 := lat2.Radian;
var a := (dLat / 2).sin() * (dLat / 2).sin() + (dLon / 2).sin() * (dLon / 2).sin() * dLat1.cos() * dLat2.cos();
^ R * 2 * a.sqrt().arcsin()
}
public program()
{
console.printLineFormatted("The distance between coordinates {0},{1} and {2},{3} is: {4}", 36.12r, -86.67r, 33.94r, -118.40r,
Haversine(36.12r, -86.67r, 33.94r, -118.40r))
}</syntaxhighlight>
{{out}}
<pre>
The distance between coordinates 36.12,-86.67 and 33.94,-118.4 is: 2887.259950607
</pre>
=={{header|Elixir}}==
<syntaxhighlight lang="elixir">defmodule Haversine do
@v :math.pi / 180
@r 6372.8 # km for the earth radius
def distance({lat1, long1}, {lat2, long2}) do
dlat = :math.sin((lat2 - lat1) * @v / 2)
dlong = :math.sin((long2 - long1) * @v / 2)
a = dlat * dlat + dlong * dlong * :math.cos(lat1 * @v) * :math.cos(lat2 * @v)
@r * 2 * :math.asin(:math.sqrt(a))
end
end
bna = {36.12, -86.67}
lax = {33.94, -118.40}
IO.puts Haversine.distance(bna, lax)</syntaxhighlight>
{{out}}
<pre>
2887.2599506071106
</pre>
=={{header|Elm}}==
<syntaxhighlight lang="elm">haversine : ( Float, Float ) -> ( Float, Float ) -> Float
haversine ( lat1, lon1 ) ( lat2, lon2 ) =
let
r =
6372.8
dLat =
degrees (lat2 - lat1)
dLon =
degrees (lon2 - lon1)
a =
(sin (dLat / 2))
^ 2
+ (sin (dLon / 2))
^ 2
* cos (degrees lat1)
* cos (degrees lat2)
in
r * 2 * asin (sqrt a)
view =
Html.div []
[ Html.text (toString (haversine ( 36.12, -86.67 ) ( 33.94, -118.4 )))
]
</syntaxhighlight>
{{out}}
<pre>2887.2599506071106
</pre>
=={{header|Erlang}}==
<
-module(haversine).
-export([main/0]).
Line 424 ⟶ 1,347:
C = 2 * math:asin(math:sqrt(A)),
R*C.
</syntaxhighlight>
{{out}}
<pre>2887.2599506071106
Line 430 ⟶ 1,353:
=={{header|ERRE}}==
<
PROGRAM HAVERSINE_DEMO
Line 461 ⟶ 1,384:
PRINT("HAVERSINE DISTANCE: ";RES;" KM.")
END PROGRAM
</syntaxhighlight>
Using double-precision variables output is 2887.260209071741 km, while using single-precision variable output is 2887.261 Km.
Line 489 ⟶ 1,412:
</pre>
=={{header|Excel}}==
===LAMBDA===
Binding the name HAVERSINE to the following lambda expression in the Name Manager of the Excel workbook:
(See [https://www.microsoft.com/en-us/research/blog/lambda-the-ultimatae-excel-worksheet-function/ LAMBDA: The ultimate Excel worksheet function])
{{Works with|Office 365 betas 2021}}
<syntaxhighlight lang="lisp">HAVERSINE
=LAMBDA(lla,
LAMBDA(llb,
LET(
REM, "Approximate radius of Earth in km.",
earthRadius, 6372.8,
sinHalfDeltaSquared, LAMBDA(x, SIN(x / 2) ^ 2)(
RADIANS(llb - lla)
),
2 * earthRadius * ASIN(
SQRT(
INDEX(sinHalfDeltaSquared, 1) + (
PRODUCT(COS(RADIANS(
CHOOSE({1,2},
INDEX(lla, 1),
INDEX(llb, 1)
)
)))
) * INDEX(sinHalfDeltaSquared, 2)
)
)
)
)
)</syntaxhighlight>
Each of the two arguments in the example below is an Excel dynamic array of two adjacent values. The # character yields a reference to the array with the given top-left grid address.
Cell B2 is formatted to display only two decimal places.
{{Out}}
{| class="wikitable"
|-
|||style="text-align:right; font-family:serif; font-style:italic; font-size:120%;"|fx
! colspan="9" style="text-align:left; vertical-align: bottom; font-family:Arial, Helvetica, sans-serif !important;"|=HAVERSINE(E2#)(H2#)
|- style="text-align:center; font-family:Arial, Helvetica, sans-serif !important; background-color:#000000; color:#ffffff;"
|
| A
| B
| C
| D
| E
| F
| G
| H
| I
|-
| style="text-align:center; font-family:Arial, Helvetica, sans-serif !important; background-color:#000000; color:#ffffff" | 1
|
| style="text-align:left; font-weight:bold" | Distance
|
|
| colspan="2" style="font-weight:bold" | BNA
|
| colspan="2" style="font-weight:bold" | LAX
|-
| style="text-align:center; font-family:Arial, Helvetica, sans-serif !important; background-color:#000000; color:#ffffff" | 2
|
| style="text-align:left; background-color:#cbcefb" | 2887.26
| km
|
| style="text-align:left" | 36.12
| style="text-align:left" | -86.67
|
| style="text-align:left" | 33.94
| style="text-align:left" | -118.4
|}
=={{header|F_Sharp|F#}}==
{{trans|Go}} using units of measure
<
[<Measure>] type deg
Line 516 ⟶ 1,513:
let main argv =
printfn "%A" (hsDist (pos(36.12<deg>, -86.67<deg>)) (pos(33.94<deg>, -118.40<deg>)))
0</
{{out}}
<pre>2887.259951</pre>
=={{header|Factor}}==
{{trans|J}}
<
: haversin ( x -- y ) cos 1 swap - 2 / ;
Line 533 ⟶ 1,529:
2bi
v.
haversininv R_earth * ;</
<
2887.259950607113</
=={{header|FBSL}}==
Based on the Fortran and Groovy versions.
<
PRINT "Distance = ", Haversine(36.12, -86.67, 33.94, -118.4), " km"
Line 554 ⟶ 1,550:
RETURN radius * c
END FUNCTION
</syntaxhighlight>
{{out}}
Distance = 2887.25995060711 km
Press any key to continue...
=={{header|FOCAL}}==
{{trans|ZX Spectrum Basic}}
<syntaxhighlight lang="focal">1.01 S BA = 36.12; S LA = -86.67
1.02 S BB = 33.94; S LB = -118.4
1.03 S DR = 3.1415926536 / 180; S D = 2 * 6372.8
1.04 S TA = (LB - LA) * DR
1.05 S TB = DR * BA
1.06 S TC = DR * BB
1.07 S DZ = FSIN(TB) - FSIN(TC)
1.08 S DX = FCOS(TA) * FCOS(TB) - FCOS(TC)
1.09 S DY = FSIN(TA) * FCOS(TB)
1.10 S AS = DX * DX + DY * DY + DZ * DZ
1.11 S AS = FSQT(AS) / 2
1.12 S HDIST = D * FATN(AS / FSQT(1 - AS^2))
1.13 T %6.2,"Haversine distance ",HDIST,!</syntaxhighlight>
{{output}}
<pre>Haversine distance = 2887.26</pre>
Note that FOCAL lacks a built-in arcsine function, but appendix D of the FOCAL manual shows how to compute it using arctangent and square root instead.
=={{header|Forth}}==
<
: deg>rad 174532925199433e-16 f* ;
: difference f- deg>rad 2 s>f f/ fsin fdup f* ;
Line 574 ⟶ 1,589:
;
36.12e -86.67e 33.94e -118.40e haversine cr f.</
{{out}}
<pre>
Line 581 ⟶ 1,596:
=={{header|Fortran}}==
<syntaxhighlight lang="fortran">
program example
implicit none
Line 616 ⟶ 1,631:
end program example
</syntaxhighlight>
=={{header|Free Pascal}}==
Here is a Free Pascal version, works in most Pascal dialects, but also note the Delphi entry that also works in Free Pascal.
<syntaxhighlight lang="pascal">program HaversineDemo;
uses
Math;
function HaversineDistance(const lat1, lon1, lat2, lon2:double):double;inline;
const
rads = pi / 180;
dia = 2 * 6372.8;
begin
HaversineDistance := dia * arcsin(sqrt(sqr(cos(rads * (lon1 - lon2)) * cos(rads * lat1)
- cos(rads * lat2)) + sqr(sin(rads * (lon1 - lon2))
* cos(rads * lat1)) + sqr(sin(rads * lat1) - sin(rads * lat2))) / 2);
end;
begin
Writeln('Haversine distance between BNA and LAX: ', HaversineDistance(36.12, -86.67, 33.94, -118.4):7:2, ' km.');
end.</syntaxhighlight>
=={{header|FreeBASIC}}==
<syntaxhighlight lang="freebasic">' version 09-10-2016
' compile with: fbc -s console
' Nashville International Airport (BNA) in Nashville, TN, USA,
' N 36°07.2', W 86°40.2' (36.12, -86.67)
' Los Angeles International Airport (LAX) in Los Angeles, CA, USA,
' N 33°56.4', W 118°24.0' (33.94, -118.40).
' 6372.8 km is an approximation of the radius of the average circumference
#Define Pi Atn(1) * 4 ' define Pi = 3.1415..
#Define deg2rad Pi / 180 ' define deg to rad 0.01745..
#Define earth_radius 6372.8 ' earth radius in km.
Function Haversine(lat1 As Double, long1 As Double, lat2 As Double, _
long2 As Double , radius As Double) As Double
Dim As Double d_long = deg2rad * (long1 - long2)
Dim As Double theta1 = deg2rad * lat1
Dim As Double theta2 = deg2rad * lat2
Dim As Double dx = Cos(d_long) * Cos(theta1) - Cos(theta2)
Dim As Double dy = Sin(d_long) * Cos(theta1)
Dim As Double dz = Sin(theta1) - Sin(theta2)
Return Asin(Sqr(dx*dx + dy*dy + dz*dz) / 2) * radius * 2
End Function
Print
Print " Haversine distance between BNA and LAX = "; _
Haversine(36.12, -86.67, 33.94, -118.4, earth_radius); " km."
' empty keyboard buffer
While Inkey <> "" : Wend
Print : Print "hit any key to end program"
Sleep
End</syntaxhighlight>
{{out}}
<pre> Haversine distance between BNA and LAX = 2887.259950607111 km.</pre>
=={{header|Frink}}==
Frink has built-in constants for the radius of the earth, whether it is the mean radius <CODE>earthradius</CODE>, the equatorial radius <CODE>earthradius_equatorial</CODE>, or the polar radius <CODE>earthradius_polar</CODE>. Below calculates the distance between the points using the haversine formula on a sphere using the mean radius, but we can do better:
<syntaxhighlight lang="frink">haversine[theta] := (1-cos[theta])/2
dist[lat1, long1, lat2, long2] := 2 earthradius arcsin[sqrt[haversine[lat2-lat1] + cos[lat1] cos[lat2] haversine[long2-long1]]]
d = dist[36.12 deg, -86.67 deg, 33.94 deg, -118.40 deg]
println[d-> "km"]</syntaxhighlight>
{{out}}
<pre>
2886.4489734366999158 km
</pre>
Note that physical constants like degrees, kilometers, and the average radius of the earth (as well as the polar and equatorial radii) are already known to Frink. Also note that units of measure are tracked throughout all calculations, and results can be displayed in a huge number of units of distance (miles, km, furlongs, chains, feet, statutemiles, etc.) by changing the final <code>"km"</code> to something like <code>"miles"</code>.
However, Frink's library/sample program [http://futureboy.us/fsp/colorize.fsp?fileName=navigation.frink navigation.frink] (included in larger distributions) contains a much higher-precision calculation that uses ellipsoidal (not spherical) calculations to determine the distance on earth's geoid with far greater accuracy
The calculations are due to:
"Direct and Inverse Solutions of Geodesics on the Ellipsoid with Application
of Nested Equations", T.Vincenty, ''Survey Review XXII'', 176, April 1975.
http://www.ngs.noaa.gov/PUBS_LIB/inverse.pdf
There is also a slightly higher-accuracy algorithm (closer to nanometers instead of fractions of millimeters LOL):
"Algorithms for geodesics", Charles F. F. Karney, ''Journal of Geodesy'', January 2013, Volume 87, Issue 1, pp 43-55
http://link.springer.com/article/10.1007%2Fs00190-012-0578-z
<syntaxhighlight lang="frink">use navigation.frink
d = earthDistance[36.12 deg North, 86.67 deg West, 33.94 deg North, 118.40 deg West]
println[d-> "km"]</syntaxhighlight>
{{out}}
<pre>
2892.7769573807044975 km
</pre>
Which should be treated as the closest-to-right answer for the actual distance on the earth's geoid, based on the WGS84 geoid datum.
=={{header|FunL}}==
<
def haversin( theta ) = (1 - cos( theta ))/2
Line 651 ⟶ 1,746:
2R asin( sqrt(h) )
println( haversine((36.12, -86.67), (33.94, -118.40)) )</
{{out}}
Line 657 ⟶ 1,752:
<pre>
2887.259950607111
</pre>
=={{header|FutureBasic}}==
Note: The Haversine function returns an approximate theoretical value of the Great Circle Distance between two points because it does not factor the ellipsoidal shape of Earth -- fat in the middle from centrifugal force, and squashed at the ends. Navigators once relied on trigonometric functions like versine (versed sine) where angle A is 1-cos(A), and haversine (half versine) or ( 1-cos(A) ) / 2.
Also, the radius of the Earth varies, at least depending on who you talk to. Here's NASA's take on it: http://nssdc.gsfc.nasa.gov/planetary/factsheet/earthfact.html
Since it was trivial, this functions returns the distance in miles and kilometers.
<syntaxhighlight lang="futurebasic">window 1
local fn Haversine( lat1 as double, lon1 as double, lat2 as double, lon2 as double, miles as ^double, kilometers as ^double )
double deg2rad, dLat, dLon, a, c, earth_radius_miles, earth_radius_kilometers
earth_radius_miles = 3959.0 // Radius of the Earth in miles
earth_radius_kilometers = 6372.8 // Radius of the Earth in kilometers
deg2rad = Pi / 180 // Pi is predefined in FutureBasic
dLat = deg2rad * ( lat2 - lat1 )
dLon = deg2rad * ( lon2 - lon1 )
a = sin( dLat / 2 ) * sin( dLat / 2 ) + cos( deg2rad * lat1 ) * cos( deg2rad * lat2 ) * sin( dLon / 2 ) * sin( dLon / 2 )
c = 2 * asin( sqr(a) )
miles.nil# = earth_radius_miles * c
kilometers.nil# = earth_radius_kilometers * c
end fn
double miles, kilometers
fn Haversine( 36.12, -86.67, 33.94, -118.4, @miles, @kilometers )
print "Distance in miles between BNA and LAX: "; using "####.####"; miles; " miles."
print "Distance in kilometers between BNA LAX: "; using "####.####"; kilometers; " km."
HandleEvents</syntaxhighlight>
Output:
<pre>
Distance in miles between BNA and LAX: 1793.6640 miles.
Distance in kilometers between BNA LAX: 2887.2600 km.
</pre>
=={{header|Go}}==
<
import (
Line 689 ⟶ 1,821:
func main() {
fmt.Println(hsDist(degPos(36.12, -86.67), degPos(33.94, -118.40)))
}</
{{out}}
<pre>
Line 696 ⟶ 1,828:
=={{header|Groovy}}==
<
def R = 6372.8
// In kilometers
Line 711 ⟶ 1,843:
haversine(36.12, -86.67, 33.94, -118.40)
> 2887.25995060711</
=={{header|Haskell}}==
<syntaxhighlight lang="haskell">import Control.Monad (join)
import Data.Bifunctor (bimap)
import Text.Printf (printf)
-------------------- HAVERSINE FORMULA -------------------
-- The haversine of an angle.
haversine :: Float -> Float
haversine = (^ 2) . sin . (/ 2)
-- The approximate distance
--
-- The latitude and longtitude are assumed to be in degrees.
greatCircleDistance ::
(Float, Float) ->
(Float, Float) ->
Float
greatCircleDistance = distDeg 6371
where
distDeg radius p1 p2 =
distRad
radius
(deg2rad p1)
(deg2rad p2)
distRad radius (lat1, lng1) (lat2, lng2) =
(2 * radius)
* asin
( min
1.0
( sqrt $
haversine (lat2 - lat1)
+ ( (cos lat1 * cos lat2)
* haversine (lng2 - lng1)
)
)
)
deg2rad = join bimap ((/ 180) . (pi *))
--------------------------- TEST -------------------------
main :: IO ()
main =
printf
(greatCircleDistance bna lax)
where
bna = (36.12, -86.67)
lax = (33.94, -118.40)</syntaxhighlight>
{{Out}}
<pre>The distance between BNA and LAX is about 2886 km.</pre>
=={{header|Icon}} and {{header|Unicon}}==
{{trans|C}}
<
procedure main() #: Haversine formula
Line 762 ⟶ 1,912:
dy := sin(a[2]) * cos(a[1])
return asin(sqrt(dx * dx + dy * dy + dz * dz) / 2) * 2 * 6371
end</
{{libheader|Icon Programming Library}}
Line 769 ⟶ 1,919:
{{out}}
<pre>BNA to LAX is 2886 km (1793 miles)</pre>
=={{header|Idris}}==
{{trans|Haskell}}
<syntaxhighlight lang="idris">module Main
-- The haversine of an angle.
hsin : Double -> Double
hsin t = let u = sin (t/2) in u*u
-- The distance between two points, given by latitude and longtitude, on a
-- circle. The points are specified in radians.
distRad : Double -> (Double, Double) -> (Double, Double) -> Double
distRad radius (lat1, lng1) (lat2, lng2) =
let hlat = hsin (lat2 - lat1)
hlng = hsin (lng2 - lng1)
root = sqrt (hlat + cos lat1 * cos lat2 * hlng)
in 2 * radius * asin (min 1.0 root)
-- The distance between two points, given by latitude and longtitude, on a
-- circle. The points are specified in degrees.
distDeg : Double -> (Double, Double) -> (Double, Double) -> Double
distDeg radius p1 p2 = distRad radius (deg2rad p1) (deg2rad p2)
where
d2r : Double -> Double
d2r t = t * pi / 180
deg2rad (t, u) = (d2r t, d2r u)
-- The approximate distance, in kilometers, between two points on Earth.
-- The latitude and longtitude are assumed to be in degrees.
earthDist : (Double, Double) -> (Double, Double) -> Double
earthDist = distDeg 6372.8
main : IO ()
main = putStrLn $ "The distance between BNA and LAX is about " ++ show (floor dst) ++ " km."
where
bna : (Double, Double)
bna = (36.12, -86.67)
lax : (Double, Double)
lax = (33.94, -118.40)
dst : Double
dst = earthDist bna lax
</syntaxhighlight>
{{out}}
<pre>The distance between BNA and LAX is about 2887 km.</pre>
=={{header|IS-BASIC}}==
<syntaxhighlight lang="is-basic">100 PROGRAM "Haversine.bas"
110 PRINT "Haversine distance:";HAVERSINE(36.12,-86.67,33.94,-118.4);"km"
120 DEF HAVERSINE(LAT1,LON1,LAT2,LON2)
130 OPTION ANGLE RADIANS
140 LET R=6372.8
150 LET DLAT=RAD(LAT2-LAT1):LET DLON=RAD(LON2-LON1)
160 LET LAT1=RAD(LAT1):LET LAT2=RAD(LAT2)
170 LET HAVERSINE=R*2*ASIN(SQR(SIN(DLAT/2)^2+SIN(DLON/2)^2*COS(LAT1)*COS(LAT2)))
190 END DEF</syntaxhighlight>
=={{header|J}}==
'''Solution:'''
<
haversin=: 0.5 * 1 - cos
Rearth=: 6372.8
haversineDist=: Rearth * haversin^:_1@((1 , *&(cos@{.)) +/ .* [: haversin -)&rfd
</syntaxhighlight>
Note: J derives the inverse haversin ( <code>haversin^:_1</code> )
from the definition of haversin.
'''Example Use:'''
<
2887.26</
=={{header|Java}}==
{{trans|Groovy}}
<
public static final double R = 6372.8; // In kilometers
public static double haversine(double lat1, double lon1, double lat2, double lon2) {
lat1 = Math.toRadians(lat1);
lat2 = Math.toRadians(lat2);
double dLat = lat2 - lat1;
double dLon = Math.toRadians(lon2 - lon1);
double a = Math.
double c = 2 * Math.asin(Math.sqrt(a));
return R * c;
}
public static void main(String[] args) {
System.out.println(haversine(36.12, -86.67, 33.94, -118.40));
}
}</
{{out}}
<pre>2887.2599506071106</pre>
=={{header|JavaScript}}==
===ES5===
{{trans|Java}}
<
var radians = Array.prototype.map.call(arguments, function(deg) { return deg/180.0 * Math.PI; });
var lat1 = radians[0], lon1 = radians[1], lat2 = radians[2], lon2 = radians[3];
Line 816 ⟶ 2,027:
return R * c;
}
console.log(haversine(36.12, -86.67, 33.94, -118.40));</
{{out}}
<pre>2887.2599506071124</pre>
===ES6===
<syntaxhighlight lang="javascript">((x, y) => {
'use strict';
// haversine :: (Num, Num) -> (Num, Num) -> Num
const haversine = ([lat1, lon1], [lat2, lon2]) => {
// Math lib function names
const [pi, asin, sin, cos, sqrt, pow, round] = [
'PI', 'asin', 'sin', 'cos', 'sqrt', 'pow', 'round'
]
.map(k => Math[k]),
// degrees as radians
[rlat1, rlat2, rlon1, rlon2] = [lat1, lat2, lon1, lon2]
.map(x => x / 180 * pi),
dLat = rlat2 - rlat1,
dLon = rlon2 - rlon1,
radius = 6372.8; // km
// km
return round(
radius * 2 * asin(
sqrt(
pow(sin(dLat / 2), 2) +
pow(sin(dLon / 2), 2) *
cos(rlat1) * cos(rlat2)
)
) * 100
) / 100;
};
// TEST
return haversine(x, y);
// --> 2887.26
})([36.12, -86.67], [33.94, -118.40]);</syntaxhighlight>
{{Out}}
<pre>2887.26</pre>
=={{header|jq}}==
<
def radians: . * (1|atan)/45;
def sind: radians|sin;
Line 829 ⟶ 2,081:
(((lat2 - lat1)/2) | sind | sq) as $dlat
| (((lon2 - lon1)/2) | sind | sq) as $dlon
| 2 * 6372.8 * (( $dlat + (lat1|cosd) * (lat2|cosd) * $dlon ) | sqrt | asin) ;</
'''Example''':
haversine(36.12; -86.67; 33.94; -118.4)
# 2887.2599506071106
=={{header|Jsish}}==
From Javascript, ES5, except the ''arguments'' value is an Array in jsish, not an Object.
<syntaxhighlight lang="javascript">/* Haversine formula, in Jsish */
function haversine() {
var radians = arguments.map(function(deg) { return deg/180.0 * Math.PI; });
var lat1 = radians[0], lon1 = radians[1], lat2 = radians[2], lon2 = radians[3];
var R = 6372.8; // km
var dLat = lat2 - lat1;
var dLon = lon2 - lon1;
var a = Math.sin(dLat / 2) * Math.sin(dLat /2) + Math.sin(dLon / 2) * Math.sin(dLon /2) * Math.cos(lat1) * Math.cos(lat2);
var c = 2 * Math.asin(Math.sqrt(a));
return R * c;
}
;haversine(36.12, -86.67, 33.94, -118.40);
/*
=!EXPECTSTART!=
haversine(36.12, -86.67, 33.94, -118.40) ==> 2887.259950607112
=!EXPECTEND!=
*/</syntaxhighlight>
{{out}}
<pre>prompt$ jsish -u haversineFormula.jsi
[PASS] haversineFormula.jsi</pre>
=={{header|Julia}}==
{{works with|Julia|0.6}}
<syntaxhighlight lang="julia">haversine(lat1, lon1, lat2, lon2) =
2 * 6372.8 * asin(sqrt(sind((lat2 - lat1) / 2) ^ 2 +
cosd(lat1) * cosd(lat2) * sind((lon2 - lon1) / 2) ^ 2))
@show haversine(36.12, -86.67, 33.94, -118.4)</syntaxhighlight>
{{out}}
<pre>haversine(36.12, -86.67, 33.94, -118.4) = 2887.2599506071106</pre>
=={{header|Kotlin}}==
{{trans|Groovy}}
Use Unicode characters.
<syntaxhighlight lang="scala">import java.lang.Math.*
const val R = 6372.8 // in kilometers
fun haversine(lat1: Double, lon1: Double, lat2: Double, lon2: Double): Double {
val λ1 = toRadians(lat1)
val λ2 = toRadians(lat2)
val Δλ = toRadians(lat2 - lat1)
val Δφ = toRadians(lon2 - lon1)
return 2 * R * asin(sqrt(pow(sin(Δλ / 2), 2.0) + pow(sin(Δφ / 2), 2.0) * cos(λ1) * cos(λ2)))
}
fun main(args: Array<String>) = println("result: " + haversine(36.12, -86.67, 33.94, -118.40))</syntaxhighlight>
=={{header|Lambdatalk}}==
{{trans|Python}}
<syntaxhighlight lang="scheme">
{def haversine
{def diameter {* 6372.8 2}}
{def radians {lambda {:a} {* {/ {PI} 180} :a}}}
{lambda {:lat1 :lon1 :lat2 :lon2}
{let { {:dLat {radians {- :lat2 :lat1}}}
{:dLon {radians {- :lon2 :lon1}}}
{:lat1 {radians :lat1}}
{:lat2 {radians :lat2}}
} {* {diameter}
{asin {sqrt {+ {pow {sin {/ :dLat 2}} 2}
{* {cos :lat1}
{cos :lat2}
{pow {sin {/ :dLon 2}} 2} }}}}}}}}
-> haversine
{haversine 36.12 -86.67 33.94 -118.40}
-> 2887.2599506071106
or, using
{def deg2dec
{lambda {:s :w}
{let { {:s {if {or {W.equal? :s W}
{W.equal? :s S}} then - else +}}
{:dm {S.replace ° by space in
{S.replace ' by in :w}}}
} :s{S.get 0 :dm}.{round {* {/ 100 60} {S.get 1 :dm}}}}}}
-> deg2dec
we can just write
{haversine
{deg2dec N 36°7.2'}
{deg2dec W 86°40.2'}
{deg2dec N 33°56.4'}
{deg2dec W 118°24.0'}}
-> 2887.2599506071106
</syntaxhighlight>
=={{header|Liberty BASIC}}==
<
end
function havDist( th1, ph1, th2, ph2)
Line 853 ⟶ 2,197:
dy = sin( LgD) * cos( th1)
havDist = asn( ( dx^2 +dy^2 +dz^2)^0.5 /2) *diameter
end function</
<pre>Haversine distance: 2887.25995060711 km.</pre>
=={{header|LiveCode}}==
<
return n * (3.1415926 / 180)
end radians
Line 883 ⟶ 2,227:
return (radiusEarth * (2.0 * asin(sqrt(haver))))
end haversine</
Test
<
2887.259923</
=={{header|Lua}}==
<syntaxhighlight lang="lua">local function haversine(x1, y1, x2, y2)
r=0.017453292519943295769236907684886127;
x1= x1*r; x2= x2*r; y1= y1*r; y2= y2*r; dy = y2-y1; dx = x2-x1;
a = math.pow(math.sin(dx/2),2) + math.cos(x1) * math.cos(x2) * math.pow(math.sin(dy/2),2); c = 2 * math.asin(math.sqrt(a)); d = 6372.8 * c;
return d;
end</syntaxhighlight>
Usage:
<syntaxhighlight lang="lua">print(haversine(36.12, -86.67, 33.94, -118.4));</syntaxhighlight>
Output:
<pre>2887.2599506071</pre>
=={{header|Maple}}==
Inputs assumed to be in radians.
<syntaxhighlight lang="maple">distance := (theta1, phi1, theta2, phi2)->2*6378.14*arcsin( sqrt((1-cos(theta2-theta1))/2 + cos(theta1)*cos(theta2)*(1-cos(phi2-phi1))/2) );</syntaxhighlight>
If you prefer, you can define a haversine function to clarify the definition:<syntaxhighlight lang="maple">haversin := theta->(1-cos(theta))/2;
distance := (theta1, phi1, theta2, phi2)->2*6378.14*arcsin( sqrt(haversin(theta2-theta1) + cos(theta1)*cos(theta2)*haversin(phi2-phi1)) );</syntaxhighlight>
Usage:
<pre>distance(0.6304129261, -1.512676863, 0.5923647483, -2.066469834)</pre>
{{out}}
<pre>2889.679287</pre>
=={{header|Mathematica}} / {{header|Wolfram Language}}==
Inputs assumed in degrees. Sin and Haversine expect arguments in radians; the built-in variable 'Degree' converts from degrees to radians.
<syntaxhighlight lang="mathematica">
distance[{theta1_, phi1_}, {theta2_, phi2_}] :=
2*6378.14 ArcSin@
Sqrt[Haversine[(theta2 - theta1) Degree] +
Cos[theta1*Degree] Cos[theta2*Degree] Haversine[(phi2 - phi1) Degree]]
</syntaxhighlight>
Usage:
<pre>distance[{36.12, -86.67}, {33.94, -118.4}]</pre>
Line 903 ⟶ 2,270:
=={{header|MATLAB}} / {{header|Octave}}==
<
% degrees to radians
rad = degree .* pi / 180;
Line 919 ⟶ 2,286:
end;
[a,c,dlat,dlon] = haversine(36.12,-86.67,33.94,-118.40); % BNA to LAX</
{{out}}
<pre>distance: 2887.2600 km</pre>
=={{header|Maxima}}==
<
great_circle_distance(lat1, long1, lat2, long2) :=
Line 935 ⟶ 2,302:
great_circle_distance(dms( 36, 7, 28.10), -dms( 86, 40, 41.50),
dms( 33, 56, 32.98), -dms(118, 24, 29.05)), numer;
/* 2886.326609413624 */</
=={{header|МК-61/52}}==
<syntaxhighlight lang="text">П3 -> П2 -> П1 -> П0
пи 1 8 0 / П4
ИП1 МГ ИП3 МГ - ИП4 * П1 ИП0 МГ ИП4 * П0 ИП2 МГ ИП4 * П2
Line 944 ⟶ 2,311:
ИП1 cos ИП0 cos * ИП2 cos - П6
ИП1 sin ИП0 cos * П7
ИП6 x^2 ИП7 x^2 ИП8 x^2 + + КвКор 2 / arcsin 2 * ИП5 * С/П</
''Input'': 6371,1 as a radius of the Earth, taken as the ball, or 6367,554 as an average radius of the Earth, or 6367,562 as an approximation of the radius of the average circumference (by Krasovsky's ellipsoid) to Р5; В/О ''lat<sub>1</sub>'' С/П ''long<sub>1</sub>'' С/П ''lat<sub>2</sub>'' С/П ''long<sub>2</sub>'' С/П; the coordinates must be entered as ''degrees,minutes'' (example: 46°50' as 46,5).
Line 958 ⟶ 2,325:
: ''Output'': 7357,4526.
=={{header|
<syntaxhighlight lang="mysql">DELIMITER $$
CREATE FUNCTION haversine (
lat1 FLOAT, lon1 FLOAT,
lat2 FLOAT, lon2 FLOAT
) RETURNS FLOAT
NO SQL DETERMINISTIC
BEGIN
DECLARE r FLOAT unsigned DEFAULT 6372.8;
DECLARE dLat FLOAT unsigned;
DECLARE dLon FLOAT unsigned;
DECLARE a FLOAT unsigned;
DECLARE c FLOAT unsigned;
SET dLat = ABS(RADIANS(lat2 - lat1));
SET dLon = ABS(RADIANS(lon2 - lon1));
SET lat1 = RADIANS(lat1);
SET lat2 = RADIANS(lat2);
SET a = POW(SIN(dLat / 2), 2) + COS(lat1) * COS(lat2) * POW(SIN(dLon / 2), 2);
SET c = 2 * ASIN(SQRT(a));
RETURN (r * c);
END$$
DELIMITER ;</syntaxhighlight>
Usage:
<pre>SELECT haversine(36.12, -86.67, 33.94, -118.4);</pre>
{{out}}
<pre>2887.260009765625</pre>
=={{header|Nim}}==
<syntaxhighlight lang="nim">import std/math
proc haversine(lat1, lon1, lat2, lon2: float): float =
const r = 6372.8 # Earth radius in kilometers
let
dLat =
dLon =
lat1 =
lat2 =
a = sin(dLat / 2) * sin(dLat / 2) + cos(lat1) * cos(lat2) * sin(dLon / 2) * sin(dLon / 2)
c = 2 * arcsin(sqrt(a))
result = r * c
echo haversine(36.12, -86.67, 33.94, -118.40)</
{{out}}
<pre>
=={{header|Oberon-2}}==
Works with oo2c version2
<
MODULE Haversines;
IMPORT
Line 1,010 ⟶ 2,408:
Out.LongRealFix(Distance(36.12,-86.67,33.94,-118.4),6,10);Out.Ln
END Haversines.
</syntaxhighlight>
Output:
<pre>
Line 1,017 ⟶ 2,415:
=={{header|Objeck}}==
<
bundle Default {
class Haversine {
Line 1,038 ⟶ 2,436:
}
}
</syntaxhighlight>
{{out}}
<pre>
Line 1,045 ⟶ 2,443:
=={{header|Objective-C}}==
<
lat2:(double)lat2 lon2:(double)lon2 {
//degrees to radians
Line 1,061 ⟶ 2,459:
double R = 6372.8;
return R * c;
}</
=={{header|OCaml}}==
Line 1,068 ⟶ 2,466:
but with an eye toward generality and reuse,
this is how I might start:
<
let pi = 4. *. atan 1.
let radians_of_degrees = ( *. ) (pi /. 180.)
Line 1,099 ⟶ 2,497:
and lax = LatLong.of_angles (Deg 33.94) (Deg (-118.4))
in
earth_dist bna lax;;</
If the above is fed to the REPL, the last line will produce this:
Line 1,109 ⟶ 2,507:
=={{header|Oforth}}==
<syntaxhighlight lang="oforth">import: math
: haversine(lat1, lon1, lat2, lon2)
| lat lon |
Line 1,117 ⟶ 2,516:
lon 2 / sin sq lat1 asRadian cos * lat2 asRadian cos *
lat 2 / sin sq + sqrt asin 2 * 6372.8 * ;
haversine(36.12, -86.67, 33.94, -118.40) println</
{{out}}
Line 1,130 ⟶ 2,528:
{{trans|REXX}}
The rxmath library provides the required functions.
<
say " Nashville: north 36º 7.2', west 86º 40.2' = 36.12º, -86.67º"
say "Los Angles: north 33º 56.4', west 118º 24.0' = 33.94º, -118.40º"
Line 1,155 ⟶ 2,553:
asin: Return RxCalcArcSin(arg(1),,'R')
sqrt: Return RxCalcSqrt(arg(1))
::requires rxMath library</
{{out}}
<pre> Nashville: north 36º 7.2', west 86º 40.2' = 36.12º, -86.67º
Line 1,163 ⟶ 2,561:
or 1794.06 statute miles,
or 1559.00 nautical or air miles.</pre>
=={{header|PARI/GP}}==
<
my(v=[cos(ph)*cos(th1)-cos(th2),sin(ph)*cos(th1),sin(th1)-sin(th2)]);
asin(sqrt(norml2(v))/2)
Line 1,173 ⟶ 2,571:
d*dist(th1*deg, th2*deg, (ph1-ph2)*deg)
};
distEarth(36.12, -86.67, 33.94, -118.4)</
{{out}}
<pre>%1 = 2886.44444</pre>
Line 1,179 ⟶ 2,577:
=={{header|Pascal}}==
{{works with|Free_Pascal}} {{libheader|Math}}
<
uses
Line 1,202 ⟶ 2,600:
begin
writeln ('Haversine distance: ', haversineDist(36.12, -86.67, 33.94, -118.4):7:2, ' km.');
end.</
{{out}}
<pre>Haversine distance: 2887.26 km.
Line 1,208 ⟶ 2,606:
=={{header|Perl}}==
===Low-Level===
{{libheader|ntheory}}
<syntaxhighlight lang="perl">use ntheory qw/Pi/;
sub asin { my $x = shift; atan2($x, sqrt(1-$x*$x)); }
Line 1,224 ⟶ 2,626:
return $radius * $c;
}
my @BNA = (36.12, -86.67);
printf "Distance: %.3f km\n", surfacedist(@BNA, @LAX);</syntaxhighlight>
{{out}}
<pre>Distance: 2887.260 km</pre>
===Idiomatic===
Contrary to ntheory, Math::Trig is part of the Perl core distribution.
It comes with a great circle distance built-in.
<syntaxhighlight lang="perl">use Math::Trig qw(great_circle_distance deg2rad);
# Notice the 90 - latitude: phi zero is at the North Pole.
# Parameter order is: LON, LAT
my @BNA = (deg2rad(-86.67), deg2rad(90 - 36.12));
my @LAX = (deg2rad(-118.4), deg2rad(90 - 33.94));
print "Distance: ", great_circle_distance(@BNA, @LAX, 6372.8), " km\n";</syntaxhighlight>
{{out}}
<pre>Distance: 2887.25995060711 km</pre>
=={{header|Phix}}==
<!--<syntaxhighlight lang="phix">(phixonline)-->
<span style="color: #008080;">function</span> <span style="color: #000000;">haversine</span><span style="color: #0000FF;">(</span><span style="color: #004080;">atom</span> <span style="color: #000000;">lat1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">long1</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">lat2</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">long2</span><span style="color: #0000FF;">)</span>
<span style="color: #008080;">constant</span> <span style="color: #000000;">MER</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">6371</span><span style="color: #0000FF;">,</span> <span style="color: #000080;font-style:italic;">-- mean earth radius(km)</span>
<span style="color: #000000;">DEG_TO_RAD</span> <span style="color: #0000FF;">=</span> <span style="color: #004600;">PI</span><span style="color: #0000FF;">/</span><span style="color: #000000;">180</span>
<span style="color: #000000;">lat1</span> <span style="color: #0000FF;">*=</span> <span style="color: #000000;">DEG_TO_RAD</span>
<span style="color: #000000;">lat2</span> <span style="color: #0000FF;">*=</span> <span style="color: #000000;">DEG_TO_RAD</span>
<span style="color: #000000;">long1</span> <span style="color: #0000FF;">*=</span> <span style="color: #000000;">DEG_TO_RAD</span>
<span style="color: #000000;">long2</span> <span style="color: #0000FF;">*=</span> <span style="color: #000000;">DEG_TO_RAD</span>
<span style="color: #008080;">return</span> <span style="color: #000000;">MER</span><span style="color: #0000FF;">*</span><span style="color: #7060A8;">arccos</span><span style="color: #0000FF;">(</span><span style="color: #7060A8;">sin</span><span style="color: #0000FF;">(</span><span style="color: #000000;">lat1</span><span style="color: #0000FF;">)*</span><span style="color: #7060A8;">sin</span><span style="color: #0000FF;">(</span><span style="color: #000000;">lat2</span><span style="color: #0000FF;">)+</span><span style="color: #7060A8;">cos</span><span style="color: #0000FF;">(</span><span style="color: #000000;">lat1</span><span style="color: #0000FF;">)*</span><span style="color: #7060A8;">cos</span><span style="color: #0000FF;">(</span><span style="color: #000000;">lat2</span><span style="color: #0000FF;">)*</span><span style="color: #7060A8;">cos</span><span style="color: #0000FF;">(</span><span style="color: #000000;">long2</span><span style="color: #0000FF;">-</span><span style="color: #000000;">long1</span><span style="color: #0000FF;">))</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">function</span>
<span style="color: #004080;">atom</span> <span style="color: #000000;">d</span> <span style="color: #0000FF;">=</span> <span style="color: #000000;">haversine</span><span style="color: #0000FF;">(</span><span style="color: #000000;">36.12</span><span style="color: #0000FF;">,-</span><span style="color: #000000;">86.67</span><span style="color: #0000FF;">,</span><span style="color: #000000;">33.94</span><span style="color: #0000FF;">,-</span><span style="color: #000000;">118.4</span><span style="color: #0000FF;">)</span>
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">"Distance is %f km (%f miles)\n"</span><span style="color: #0000FF;">,{</span><span style="color: #000000;">d</span><span style="color: #0000FF;">,</span><span style="color: #000000;">d</span><span style="color: #0000FF;">/</span><span style="color: #000000;">1.609344</span><span style="color: #0000FF;">})</span>
<!--</syntaxhighlight>-->
{{out}}
<pre>
Distance is 2886.444443 km (1793.553425 miles)
</pre>
=={{header|PHP}}==
<
private $latitude;
private $longitude;
public function __construct($latitude, $longitude) {
$this->latitude = deg2rad($latitude);
$this->longitude = deg2rad($longitude);
}
public function
return $this->latitude;
}
public function getLongitude() {
return $this->longitude;
}
public function getDistanceInMetersTo(POI $other) {
$radiusOfEarth =
$diffLatitude = $other->getLatitude() - $this->latitude;
$diffLongitude = $other->getLongitude() - $this->longitude;
$a = sin($diffLatitude / 2) ** 2 +
cos($other->getLatitude()) *
sin($diffLongitude / 2) ** 2;
$c = 2 * asin(sqrt($a));
$distance = $radiusOfEarth * $c;
return $distance;
}
}</
'''Test:'''
<syntaxhighlight lang="php">$bna = new POI(36.12, -86.67); // Nashville International Airport
$
{{out}}
<pre>2886.44 km</pre>
=={{header|PicoLisp}}==
<
(load "@lib/math.l")
Line 1,304 ⟶ 2,727:
(/
(sqrt (+ (* DX DX) (* DY DY) (* DZ DZ)))
2 ) ) ) ) )</
Test:
<
"Haversine distance: "
(round (haversine 36.12 -86.67 33.94 -118.4))
" km" )</
{{out}}
<pre>Haversine distance: 2,886.444 km</pre>
=={{header|PL/I}}==
<
declare d float;
Line 1,342 ⟶ 2,765:
end haversine;
end test;</
{{out}}
<pre>
distance: 2887.260 km
</pre>
=={{header|PowerShell}}==
{{works with|PowerShell|3}}
<syntaxhighlight lang="powershell">
Add-Type -AssemblyName System.Device
$BNA = New-Object System.Device.Location.GeoCoordinate 36.12, -86.67
$LAX = New-Object System.Device.Location.GeoCoordinate 33.94, -118.40
$BNA.GetDistanceTo( $LAX ) / 1000
</syntaxhighlight>
{{out}}
<pre>
2888.93627213254
</pre>
{{works with|PowerShell|2}}
<syntaxhighlight lang="powershell">
function Get-GreatCircleDistance ( $Coord1, $Coord2 )
{
# Convert decimal degrees to radians
$Lat1 = $Coord1[0] / 180 * [math]::Pi
$Long1 = $Coord1[1] / 180 * [math]::Pi
$Lat2 = $Coord2[0] / 180 * [math]::Pi
$Long2 = $Coord2[1] / 180 * [math]::Pi
# Mean Earth radius (km)
$R = 6371
# Haversine formula
$ArcLength = 2 * $R *
[math]::Asin(
[math]::Sqrt(
[math]::Sin( ( $Lat1 - $Lat2 ) / 2 ) *
[math]::Sin( ( $Lat1 - $Lat2 ) / 2 ) +
[math]::Cos( $Lat1 ) *
[math]::Cos( $Lat2 ) *
[math]::Sin( ( $Long1 - $Long2 ) / 2 ) *
[math]::Sin( ( $Long1 - $Long2 ) / 2 ) ) )
return $ArcLength
}
$BNA = 36.12, -86.67
$LAX = 33.94, -118.40
Get-GreatCircleDistance $BNA $LAX
</syntaxhighlight>
{{out}}
<pre>
2886.44444283799
</pre>
=={{header|Pure Data}}==
Up until now there is no 64bit float in Pure Data, so the result of the calculation might not be completely accurate.
<pre>
#N canvas 527 1078 450 686 10;
#X obj 28 427 atan2;
#X obj 28 406 sqrt;
#X obj 62 405 sqrt;
#X obj 28 447 * 2;
#X obj 62 384 -;
#X msg 62 362 1 \$1;
#X obj 28 339 t f f;
#X obj 28 210 sin;
#X obj 83 207 sin;
#X obj 138 206 cos;
#X obj 193 206 cos;
#X obj 28 179 / 2;
#X obj 83 182 / 2;
#X obj 28 74 unpack f f;
#X obj 28 98 t f f;
#X obj 28 301 expr $f1 + ($f2 * $f3 * $f4);
#X obj 28 148 deg2rad;
#X obj 83 149 deg2rad;
#X obj 138 148 deg2rad;
#X obj 193 149 deg2rad;
#X obj 28 232 t f f;
#X obj 28 257 *;
#X obj 83 232 t f f;
#X obj 83 257 *;
#X obj 83 98 t f b;
#X obj 28 542 * 6372.8;
#X obj 193 120 f 33.94;
#X obj 28 125 - 33.94;
#X msg 28 45 36.12 -86.67;
#X obj 83 123 - -118.4;
#X floatatom 28 577 8 0 0 0 - - -, f 8;
#X connect 0 0 3 0;
#X connect 1 0 0 0;
#X connect 2 0 0 1;
#X connect 3 0 25 0;
#X connect 4 0 2 0;
#X connect 5 0 4 0;
#X connect 6 0 1 0;
#X connect 6 1 5 0;
#X connect 7 0 20 0;
#X connect 8 0 22 0;
#X connect 9 0 15 2;
#X connect 10 0 15 3;
#X connect 11 0 7 0;
#X connect 12 0 8 0;
#X connect 13 0 14 0;
#X connect 13 1 24 0;
#X connect 14 0 27 0;
#X connect 14 1 18 0;
#X connect 15 0 6 0;
#X connect 16 0 11 0;
#X connect 17 0 12 0;
#X connect 18 0 9 0;
#X connect 19 0 10 0;
#X connect 20 0 21 0;
#X connect 20 1 21 1;
#X connect 21 0 15 0;
#X connect 22 0 23 0;
#X connect 22 1 23 1;
#X connect 23 0 15 1;
#X connect 24 0 29 0;
#X connect 24 1 26 0;
#X connect 25 0 30 0;
#X connect 26 0 19 0;
#X connect 27 0 16 0;
#X connect 28 0 13 0;
#X connect 29 0 17 0;
</pre>
=={{header|PureBasic}}==
{{trans|Pascal}}
<syntaxhighlight lang="purebasic">#DIA=2*6372.8
Procedure.d Haversine(th1.d,ph1.d,th2.d,ph2.d)
Define dx.d,
dy.d,
dz.d
ph1=Radian(ph1-ph2)
th1=Radian(th1)
th2=Radian(th2)
dz=Sin(th1)-Sin(th2)
dx=Cos(ph1)*Cos(th1)-Cos(th2)
dy=Sin(ph1)*Cos(th1)
ProcedureReturn ASin(Sqr(Pow(dx,2)+Pow(dy,2)+Pow(dz,2))/2)*#DIA
EndProcedure
OpenConsole("Haversine distance")
Print("Haversine distance: ")
Print(StrD(Haversine(36.12,-86.67,33.94,-118.4),7)+" km.")
Input()</syntaxhighlight>
{{out}}
<pre>Haversine distance: 2887.2599506 km.</pre>
=={{header|Python}}==
<
def haversine(lat1, lon1, lat2, lon2):
R = 6372.8 # Earth radius in kilometers
dLat = radians(lat2 - lat1)
dLon = radians(lon2 - lon1)
lat1 = radians(lat1)
lat2 = radians(lat2)
a = sin(dLat / 2)**2 + cos(lat1) * cos(lat2) * sin(dLon / 2)**2
return R * c
>>> haversine(36.12, -86.67, 33.94, -118.40)
2887.2599506071106
>>> </
=={{header|QB64}}==
{{trans|BASIC}}
<syntaxhighlight lang="qb64">
SCREEN _NEWIMAGE(800, 100, 32)
'*** Units: K=kilometers M=miles N=nautical miles
DIM UNIT AS STRING
DIM Distance AS STRING
DIM Result AS DOUBLE
DIM ANSWER AS DOUBLE
'*** Change the To/From Latittude/Logitudes for your run
'*** LAT/LON for Nashville International Airport (BNA)
lat1 = 36.12
Lon1 = -86.67
'*** LAT/LONG for Los Angeles International Airport (LAX)
Lat2 = 33.94
Lon2 = -118.40
'*** Initialize Values
UNIT = "K"
Distance = ""
'Radius = 6378.137
Radius = 6372.8
'*** Calculate distance using Haversine Function
lat1 = (lat1 * _PI / 180)
Lon1 = (Lon1 * _PI / 180)
Lat2 = (Lat2 * _PI / 180)
Lon2 = (Lon2 * _PI / 180)
DLon = Lon1 - Lon2
ANSWER = _ACOS(SIN(lat1) * SIN(Lat2) + COS(lat1) * COS(Lat2) * COS(DLon)) * Radius
'*** Adjust Answer based on Distance Unit (kilometers, miles, nautical miles)
SELECT CASE UNIT
CASE "M"
Result = ANSWER * 0.621371192
Distance = "miles"
CASE "N"
Result = ANSWER * 0.539956803
Distance = "nautical miles"
CASE ELSE
Result = ANSWER
Distance = "kilometers"
END SELECT
'*** Change PRINT statement with your labels for FROM/TO locations
PRINT "The distance from Nashville International to Los Angeles International in "; Distance;
PRINT USING " is: ##,###.##"; Result;
PRINT "."
END
</syntaxhighlight>
=={{header|R}}==
<
# Volumetric mean radius is 6371 km, see http://nssdc.gsfc.nasa.gov/planetary/factsheet/earthfact.html
Line 1,389 ⟶ 3,019:
dms_to_rad(33, 56, 32.98), dms_to_rad(118, 24, 29.05)) # Los Angeles International Airport (LAX)
# Output: 2886.327</
=={{header|Racket}}==
Almost the same as the Scheme version.
<
#lang racket
(require math)
Line 1,409 ⟶ 3,039:
(distance (deg-to-rad 36 7.2 0) (deg-to-rad 86 40.2 0)
(deg-to-rad 33 56.4 0) (deg-to-rad 118 24.0 0))
</syntaxhighlight>
{{out}}
<pre>
2886.444442837984
</pre>
=={{header|Raku}}==
(formerly Perl 6)
<syntaxhighlight lang="raku" line>class EarthPoint {
has $.lat; # latitude
has $.lon; # longitude
has $earth_radius = 6371; # mean earth radius
has $radian_ratio = pi / 180;
# accessors for radians
method latR { $.lat * $radian_ratio }
method lonR { $.lon * $radian_ratio }
method haversine-dist(EarthPoint $p) {
my EarthPoint $arc .= new(
lat => $!lat - $p.lat,
lon => $!lon - $p.lon );
my $a = sin($arc.latR/2) ** 2 + sin($arc.lonR/2) ** 2
* cos($.latR) * cos($p.latR);
my $c = 2 * asin( sqrt($a) );
return $earth_radius * $c;
}
}
my EarthPoint $BNA .= new(lat => 36.12, lon => -86.67);
my EarthPoint $LAX .= new(lat => 33.94, lon => -118.4);
say $BNA.haversine-dist($LAX); # 2886.44444099822</syntaxhighlight>
=={{header|Raven}}==
{{trans|Groovy}}
<
-1 acos
Line 1,441 ⟶ 3,103:
}
-118.40 33.94 -86.67 36.12 haversine "haversine: %.15g\n" print</
{{out}}
<pre>haversine: 2887.25995060711</pre>
=={{header|REXX}}==
The use of normalization for angles isn't required for the Haversine formula, but those normalization functions were included
<br>herein anyway (to support normalization of input arguments to the trigonometric functions for the general case).
<syntaxhighlight lang="rexx">/*REXX program calculates the distance between Nashville and Los Angles airports.*/
call pi; numeric digits length(pi) % 2 /*use half of PI dec. digits for output*/
say " Nashville: north 36º 7.2', west 86º 40.2' = 36.12º, -86.67º"
say " Los Angles: north 33º 56.4', west 118º 24.0' = 33.94º, -118.40º"
@using_radius= 'using the mean radius of the earth as ' /*a literal for SAY.*/
do radius=1 while radii.radius\==. /*calc. distance using specific radii. */
d= surfaceDist( 36.12, -86.67, 33.94, -118.4, radii.radius); say
say ' or ' format(d/m*5280/6076.1,,2) " nautical (or air miles)."
end /*radius*/ /*show──┘ 2 dec. digs past dec. point*/
/*──────────────────────────────────────────────────────────────────────────────────────*/
surfaceDist: parse arg th1,ph1
th1 = d2r(th1) /* " " " " " */
th2 = d2r(th2)
cosTH1= cos(th1)
y = sin(ph1) * cosTH1 /* " Y " */
z = sin(th1) - sin(th2) /* " Z " */
return Asin(sqrt(x*x + y*y + z*z)*.5) *r*2 /*compute the arcsin and return value. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
r2d: return d2d( (arg(1) * 180 / pi())) /*normalize and convert rad ──► degrees*/
r2r: return arg(1) // (pi() * 2) /*normalize radians to a unit circle. */
pi: pi= 3.141592653589793238462643383279502884197169399375105820975; return pi
/*──────────────────────────────────────────────────────────────────────────────────────*/
Asin: procedure; parse arg x 1 z 1 o 1 p; a= abs(x); aa= a * a
if a >= sqrt(2) * .5 then return sign(x) * Acos( sqrt(1 - aa) )
do j=2 by 2 until p=z; p= z; o= o * aa * (j-1) / j; z= z + o / (j+1)
end /*j*/; return z /* [↑] compute until no more noise. */
/*──────────────────────────────────────────────────────────────────────────────────────*/
cos: procedure; parse arg x; x= r2r(x); a= abs(x); Hpi= pi * .5
numeric fuzz min(6, digits() - 3) ; if a=pi then return -1
if a=Hpi | a=Hpi*3 then return 0 ; if a=pi/3 then return .5
if a=pi* 2/3 then return -.5; q= x*x; p= 1; z= 1; _= 1
do k=2 by 2; _= -_*q / (k*(k-1)); z= z+_; if z=p then leave; p=z; end; return z
/*──────────────────────────────────────────────────────────────────────────────────────*/
sin: procedure; parse arg x; x= r2r(x); numeric fuzz min(5, digits() - 3)
if abs(x)=pi then return 0; q= x*x; p= x; z= x; _= x
do k=2 by 2; _= -_*q / (k*(k+1)); z= z+_; if z=p then leave; p=z; end; return z
/*──────────────────────────────────────────────────────────────────────────────────────*/
sqrt: procedure; parse arg x; if x=0 then return 0; d=digits(); m.=9; numeric form; h=d+6
numeric digits; parse value format(x,2,1,,0) 'E0' with g "E" _ .; g=g * .5'e'_ % 2
do j=0 while h>9; m.j= h; h= h%2 + 1; end /*j*/
do k=j+5 to 0 by -1; numeric digits m.k; g= (g+x/g)*.5; end /*k*/; return g</syntaxhighlight>
REXX doesn't have most of the higher math functions, so they are included here (above) as subroutines (functions).
╔════════════════════════════════════════════════════════════════════════╗
║ A note on built─in functions: REXX doesn't have a lot of mathematical ║
║ or (particularly) trigonometric functions, so REXX programmers have ║
║ to write their own. Usually, this is done once, or most likely, one ║
║ is borrowed from another program. Knowing this, the one that is used ║
║ has a lot of boilerplate in it. ║
║ ║
║ Programming note: the "general 1─liner" subroutines are taken from ║
║ other programs that I wrote, but I broke up their one line of source ║
║ so it can be viewed without shifting the viewing window. ║
║ ║
║ The pi constant (as used here) is actually a much more robust ║
║ function and will return up to one million digits in the real version. ║
║ ║
║ One bad side effect is that, like a automobile without a hood, you see ║
║ all the dirty stuff going on. Also, don't visit a sausage factory. ║
╚════════════════════════════════════════════════════════════════════════╝
{{out|output|text= when using the in-line defaults:}}
<pre>
Nashville: north 36º 7.2', west 86º 40.2' = 36.12º, -86.67º
Los Angles: north 33º 56.4', west 118º 24.0' = 33.94º, -118.40º
─────────using the mean radius of the earth as 6372.8 kilometers─────────
Distance between: 2887.26 kilometers,
──────────using the mean radius of the earth as 6371 kilometers──────────
Distance between: 2886.44 kilometers,
or 1558.56 nautical (or air miles).
</pre>
=={{header|Ring}}==
<syntaxhighlight lang="ring">
decimals(8)
see haversine(36.12, -86.67, 33.94, -118.4) + nl
func haversine x1, y1, x2, y2
r=0.01745
x1= x1*r
x2= x2*r
y1= y1*r
y2= y2*r
dy = y2-y1
dx = x2-x1
a = pow(sin(dx/2),2) + cos(x1) * cos(x2) * pow(sin(dy/2),2)
c = 2 * asin(sqrt(a))
d = 6372.8 * c
return d
</syntaxhighlight>
=={{header|RPL}}==
{{works with|Halcyon Calc|4.2.7}}
{| class="wikitable"
! Code
! Comments
|-
|
≪
ROT - 2 / DEG SIN SQ OVER COS * 3 PICK COS *
ROT ROT - 2 / SIN SQ + √ RAD ASIN 6372.8 * 2 *
≫ 'AHAV' STO
|
''( lat1 lon1 lat2 lon2 -- distance )''
Start by the end of the formula, in degree mode
Switch to radian mode to compute Arcsin
|}
The following line of command delivers what is required:
36.12 -86.67 33.94 -118.4 AHAV
Due to the uncertainty in values of Earth radius and airports coordinates, the result shall be announced as 2887 ± 1 km even if the calculation provides many digits after the decimal point
{{out}}
<pre>
1: 2887.25995061
</pre>
=={{header|Ruby}}==
<
Radius =
def spherical_distance(start_coords, end_coords)
Line 1,566 ⟶ 3,260:
lax = [33.94, -118.4]
puts "%.1f" % spherical_distance(bna, lax)</
{{out}}
<pre>
Alternatively:
{{trans|Python}}
<syntaxhighlight lang="ruby">include Math
def haversine(lat1, lon1, lat2, lon2)
r = 6372.8 # Earth radius in kilometers
deg2rad = PI/180 # convert degress to radians
dLat = (lat2 - lat1) * deg2rad
dLon = (lon2 - lon1) * deg2rad
lat1 = lat1 * deg2rad
lat2 = lat2 * deg2rad
a = sin(dLat / 2)**2 + cos(lat1) * cos(lat2) * sin(dLon / 2)**2
c = 2 * asin(sqrt(a))
r * c
end
puts "distance is #{haversine(36.12, -86.67, 33.94, -118.40)} km "
</syntaxhighlight>
{{out}}
<pre>
distance is 2887.2599506071106 km
</pre>
=={{header|Run BASIC}}==
<
diam = 2 * 6372.8
Lg1m2 = ((-86.67)-(-118.4)) * D2R
Line 1,591 ⟶ 3,311:
' Directions (destination).
' 36.12.-86.66999
' Distance is 35.37 inches.</
=={{header|Rust}}==
<syntaxhighlight lang="rust">
struct Point {
lat: f64,
lon: f64,
}
fn haversine(origin: Point, destination: Point) -> f64 {
const R: f64 = 6372.8;
let lat1 = origin.lat.to_radians();
let lat2 = destination.lat.to_radians();
let d_lat = lat2 - lat1;
let d_lon = (destination.lon - origin.lon).to_radians();
let a = (d_lat / 2.0).sin().powi(2) + (d_lon / 2.0).sin().powi(2) * lat1.cos() * lat2.cos();
let c = 2.0 * a.sqrt().asin();
R * c
}
#[cfg(test)]
mod test {
use super::{Point, haversine};
#[test]
fn test_haversine() {
let origin: Point = Point {
lat: 36.12,
lon: -86.67
};
let destination: Point = Point {
lat: 33.94,
lon: -118.4
};
let d: f64 = haversine(origin, destination);
println!("Distance: {} km ({} mi)", d, d / 1.609344);
assert_eq!(d, 2887.2599506071106);
}
}
</syntaxhighlight>Output <pre>Distance: 2887.2599506071106 km (1794.060157807846 mi)</pre>
=={{header|SAS}}==
<syntaxhighlight lang="sas">
options minoperator;
Line 1,646 ⟶ 3,408:
%haver(36.12, -86.67, 33.94, -118.40);
</syntaxhighlight>
{{out}}
<pre>Distance is 2887.2599506 K</pre>
=={{header|Scala}}==
<
object Haversine {
Line 1,668 ⟶ 3,430:
println(haversine(36.12, -86.67, 33.94, -118.40))
}
}</
{{out}}
<pre>2887.2599506071106</pre>
=={{header|Scheme}}==
<
(define pi (acos -1))
Line 1,684 ⟶ 3,446:
(distance (deg-to-rad 36 7.2 0) (deg-to-rad 86 40.2 0)
(deg-to-rad 33 56.4 0) (deg-to-rad 118 24.0 0))
; 2886.444442837984</
=={{header|Seed7}}==
<
include "float.s7i";
include "math.s7i";
Line 1,712 ⟶ 3,474:
degToRad(33.94), degToRad(-118.4)) # Los Angeles International Airport (LAX)
digits 2);
end func;</
{{out}}
Line 1,720 ⟶ 3,482:
=={{header|Sidef}}==
{{trans|
<
# accessors for radians
method latR { self.lat * radian_ratio }
method lonR { self.lon * radian_ratio }
method haversine_dist(EarthPoint p
var arc =
self.lat - p.lat,
self.lon - p.lon,
)
var a =
(arc.lonR / 2).sin*
)
earth_radius *
}
}
var BNA = EarthPoint.new(lat: 36.12, lon: -86.67)
var LAX = EarthPoint.new(lat: 33.94, lon: -118.4)
say BNA.haversine_dist(LAX) #=> 2886.444442837983299747157823945746716...</syntaxhighlight>
=={{header|smart BASIC}}==
{{trans|BASIC}}
<syntaxhighlight lang="smart basic">
'*** LAT/LONG for Nashville International Airport (BNA)
lat1=36.12
Lon1=-86.67
'*** LAT/LONG for Los Angeles International Airport (LAX)
Lat2=33.94
Lon2=-118.40
'*** Units: K=kilometers M=miles N=nautical miles
Unit$ = "K"
Result=HAVERSINE(Lat1,Lon1,Lat2,Lon2,Unit$)
R$=STR$(Result,"#,###.##")
PRINT "The distance between Nashville International Airport and Los Angeles International Airport in kilometers is: "&R$
STOP
DEF HAVERSINE(Lat1,Lon1,Lat2,Lon2,Unit$)
'---------------------------------------------------------------
'*** Haversine Formula - Calculate distances by LAT/LONG
'
'*** Pass to it the LAT/LONG of the two locations, and then unit of measure
'*** Usage: X=HAVERSINE(Lat1,Lon1,Lat2,Lon2,Unit$)
PI=3.14159265358979323846
Radius=6372.8
Lat1=(Lat1*PI/180)
Lon1=(Lon1*PI/180)
Lat2=(Lat2*PI/180)
Lon2=(Lon2*PI/180)
DLon=Lon1-Lon2
Answer=ACOS(SIN(Lat1)*SIN(Lat2)+COS(Lat1)*COS(Lat2)*COS(DLon))*Radius
IF UNIT$="M" THEN Answer=Answer*0.621371192
IF UNIT$="N" THEN Answer=Answer*0.539956803
RETURN Answer
ENDDEF
</syntaxhighlight>
{{out}}
<pre>
The distance between Nashville International Airport and Los Angeles International Airport in kilometers is: 2,887.26
</pre>
=={{header|Stata}}==
First, a program to add a distance variable to a dataset, given variables for LAT/LON of two points.
<syntaxhighlight lang="stata">program spheredist
version 15.0
syntax varlist(min=4 max=4 numeric), GENerate(namelist max=1) ///
[Radius(real 6371) ALTitude(real 0) LABel(string)]
confirm new variable `generate'
local lat1 : word 1 of `varlist'
local lon1 : word 2 of `varlist'
local lat2 : word 3 of `varlist'
local lon2 : word 4 of `varlist'
local r=2*(`radius'+`altitude'/1000)
local k=_pi/180
gen `generate'=`r'*asin(sqrt(sin((`lat2'-`lat1')*`k'/2)^2+ ///
cos(`lat1'*`k')*cos(`lat2'*`k')*sin((`lon2'-`lon1')*`k'/2)^2))
if `"`label'"' != "" {
label variable `generate' `"`label'"'
}
end</syntaxhighlight>
Illustration with a sample dataset.
<syntaxhighlight lang="stata">import delimited airports.csv, clear
format %9.4f l*
list
+----------------------------------------------------------------------------------------------------+
| iata airport city country lat lon |
|----------------------------------------------------------------------------------------------------|
1. | AMS Amsterdam Airport Schiphol Amsterdam Netherlands 52.3086 4.7639 |
2. | BNA Nashville International Airport Nashville United States 36.1245 -86.6782 |
3. | CDG Charles de Gaulle International Airport Paris France 49.0128 2.5500 |
4. | CGN Cologne Bonn Airport Cologne Germany 50.8659 7.1427 |
5. | LAX Los Angeles International Airport Los Angeles United States 33.9425 -118.4080 |
|----------------------------------------------------------------------------------------------------|
6. | MEM Memphis International Airport Memphis United States 35.0424 -89.9767 |
+----------------------------------------------------------------------------------------------------+</syntaxhighlight>
MEM/CGN joins two Fedex Express hubs. The line AMS/LAX is operated by KLM Royal Dutch Airlines.
We will compute the distance between each pair of airports, both at sea level and at typical cruising flight level (35000 ft).
Bear in mind that the actual route of an airliner is usually not a piece of great circle, so this will only give an idea. For instance, according to [http://flightaware.com/ FlightAware], the route of a Fedex flight from Memphis to Paris is 7852 km long, at FL300 altitude (9150 m). The program given here would yield 7328.33 km instead.
<syntaxhighlight lang="stata">keep iata lat lon
rename (iata lat lon) =2
gen k=0
tempfile tmp
save "`tmp'"
rename *2 *1
joinby k using `tmp'
drop if iata1>=iata2
drop k
list
+-----------------------------------------------------------+
| iata1 lat1 lon1 iata2 lat2 lon2 |
|-----------------------------------------------------------|
1. | AMS 52.3086 4.7639 BNA 36.1245 -86.6782 |
2. | AMS 52.3086 4.7639 CGN 50.8659 7.1427 |
3. | AMS 52.3086 4.7639 LAX 33.9425 -118.4080 |
4. | AMS 52.3086 4.7639 CDG 49.0128 2.5500 |
5. | AMS 52.3086 4.7639 MEM 35.0424 -89.9767 |
|-----------------------------------------------------------|
6. | BNA 36.1245 -86.6782 CGN 50.8659 7.1427 |
7. | BNA 36.1245 -86.6782 CDG 49.0128 2.5500 |
8. | BNA 36.1245 -86.6782 LAX 33.9425 -118.4080 |
9. | BNA 36.1245 -86.6782 MEM 35.0424 -89.9767 |
10. | CDG 49.0128 2.5500 LAX 33.9425 -118.4080 |
|-----------------------------------------------------------|
11. | CDG 49.0128 2.5500 MEM 35.0424 -89.9767 |
12. | CDG 49.0128 2.5500 CGN 50.8659 7.1427 |
13. | CGN 50.8659 7.1427 LAX 33.9425 -118.4080 |
14. | CGN 50.8659 7.1427 MEM 35.0424 -89.9767 |
15. | LAX 33.9425 -118.4080 MEM 35.0424 -89.9767 |
+-----------------------------------------------------------+</syntaxhighlight>
Now compute the distances and print the result.
<syntaxhighlight lang="stata">spheredist lat1 lon1 lat2 lon2, gen(dist) lab(Distance at sea level)
spheredist lat1 lon1 lat2 lon2, gen(fl350) alt(10680) lab(Distance at FL350 altitude)
format %9.2f dist fl350
list iata* dist fl350
+-----------------------------------+
| iata1 iata2 dist fl350 |
|-----------------------------------|
1. | AMS CGN 229.64 230.03 |
2. | AMS CDG 398.27 398.94 |
3. | AMS MEM 7295.19 7307.56 |
4. | AMS BNA 7004.61 7016.48 |
5. | AMS LAX 8955.95 8971.13 |
|-----------------------------------|
6. | BNA LAX 2886.32 2891.21 |
7. | BNA CGN 7222.75 7234.99 |
8. | BNA CDG 7018.39 7030.29 |
9. | BNA MEM 321.62 322.16 |
10. | CDG LAX 9102.51 9117.94 |
|-----------------------------------|
11. | CDG CGN 387.82 388.48 |
12. | CDG MEM 7317.82 7330.23 |
13. | CGN LAX 9185.47 9201.04 |
14. | CGN MEM 7514.96 7527.70 |
15. | LAX MEM 2599.71 2604.12 |
+-----------------------------------+</syntaxhighlight>
Notice that the distance from Nashville to Los Angeles is given as 2886.32 km, which is slightly different from the task description. The coordinates come from [https://openflights.org/html/apsearch OpenFlights] and are supposably more accurate. Using the data in the task description, one gets 2886.44 as expected.
=={{header|Swift}}==
{{trans|Objective-C}}
<
func haversine(lat1:Double, lon1:Double, lat2:Double, lon2:Double) -> Double {
let lat1rad = lat1 *
let lon1rad = lon1 *
let lat2rad = lat2 *
let lon2rad = lon2 *
let dLat = lat2rad - lat1rad
Line 1,769 ⟶ 3,692:
}
{{out}}
<pre>
2887.25995060711
</pre>
=={{header|Symsyn}}==
<syntaxhighlight lang="symsyn">
lat1 : 36.12
lon1 : -86.67
lat2 : 33.94
lon2 : -118.4
dx : 0.
dy : 0.
dz : 0.
kms : 0.
{degtorad(lon2 - lon1)} lon1
{degtorad lat1} lat1
{degtorad lat2} lat2
{sin lat1 - sin lat2} dz
{cos lon1 * cos lat1 - cos lat2} dx
{sin lon1 * cos lat1} dy
{arcsin(sqrt(dx^2 + dy^2 + dz^2)/2) * 12745.6} kms
"'Haversine distance: ' kms ' kms'" []
</syntaxhighlight>
{{out}}
<pre>
Haversine distance: 2887.259951 kms
</pre>
=={{header|tbas}}==
<syntaxhighlight lang="qbasic">
option angle radians ' the default
sub haversine(lat1, lon1, lat2, lon2)
dim EarthRadiusKm = 6372.8 ' Earth radius in kilometers
dim latRad1 = RAD(lat1)
dim latRad2 = RAD(lat2)
dim lonRad1 = RAD(lon1)
dim lonRad2 = RAD(lon2)
dim _diffLa = latRad2 - latRad1
dim _doffLo = lonRad2 - lonRad1
dim sinLaSqrd = sin(_diffLa / 2) ^ 2
dim sinLoSqrd = sin(_doffLo / 2) ^ 2
dim computation = asin(sqrt(sinLaSqrd + cos(latRad1) * cos(latRad2) * sinLoSqrd))
return 2 * EarthRadiusKm * computation
end sub
print using "Nashville International Airport to Los Angeles International Airport ####.########### km", haversine(36.12, -86.67, 33.94, -118.40)
print using "Perth, WA Australia to Baja California, Mexico #####.########### km", haversine(-31.95, 115.86, 31.95, -115.86)
</syntaxhighlight>
<pre>
Nashville International Airport to Los Angeles International Airport 2887.25995060712 km
Perth, WA Australia to Baja California, Mexico 15188.70229560390 km
</pre>
=={{header|Tcl}}==
{{trans|Groovy}}
<
proc haversineFormula {lat1 lon1 lat2 lon2} {
set rads [expr atan2(0,-1)/180]
Line 1,793 ⟶ 3,772:
# Don't bother with too much inappropriate accuracy!
puts [format "distance=%.1f km" [haversineFormula 36.12 -86.67 33.94 -118.40]]</
{{out}}
<pre>distance=2887.3 km</pre>
=={{header|TechBASIC}}==
<syntaxhighlight lang="techbasic">
{{trans|BASIC}}
FUNCTION HAVERSINE
!---------------------------------------------------------------
!*** Haversine Formula - Calculate distances by LAT/LONG
!
!*** LAT/LON of the two locations and Unit of measure are GLOBAL
!*** as they are defined in the main logic of the program, so they
!*** available for use in the Function.
!*** Usage: X=HAVERSINE
Radius=6378.137
Lat1=(Lat1*MATH.PI/180)
Lon1=(Lon1*MATH.PI/180)
Lat2=(Lat2*MATH.PI/180)
Lon2=(Lon2*MATH.PI/180)
DLon=Lon1-Lon2
ANSWER=ACOS(SIN(Lat1)*SIN(Lat2)+COS(Lat1)*COS(Lat2)*COS(DLon))*Radius
DISTANCE="kilometers"
SELECT CASE UNIT
CASE "M"
HAVERSINE=ANSWER*0.621371192
Distance="miles"
CASE "N"
HAVERSINE=ANSWER*0.539956803
Distance="nautical miles"
END SELECT
END FUNCTION
</syntaxhighlight>
The following is the main code that invokes the function. It takes your location and determines how far away you are from Tampa, Florida. You can change UNIT to either M=Miles, N=Nautical Miles, or K (or leave blank) as default is in Kilometers:
<pre>
!*** In techBASIC, all variables defined in the main program act as GLOBAL
!*** variables and are available to all SUBROUTINES and FUNCTIONS. So in the
!*** HAVERSINE Function being used, no paramaters need to be passed to it, so
!*** it acts as a variable when I use it as Result=HAVERSINE. The way that
!*** the Function is setup, it returns its value back as HAVERSINE.
BASE 1
!*** Get the GPS LAT/LONG of current location
location = sensors.location(30)
Lat1=location(1)
Lon1=location(2)
!*** LAT/LONG For Tampa, FL
Lat2=27.9506
Lon2=-82.4572
!*** Units: K=kilometers M=miles N=nautical miles
DIM UNIT AS STRING
DIM Distance AS STRING
DIM Result AS SINGLE
UNIT = "M"
!*** Calculate distance using Haversine Function
Result=HAVERSINE
PRINT "The distance from your current location to Tampa, FL in ";Distance;" is: ";
PRINT USING "#,###.##";Result;"."
STOP
</pre>
<B>OUTPUT:</B> *** NOTE: When I run this, I am in my house in Venice, Florida, and that distance is correct (as the crow flies). ***
<pre>
The distance from your current location to Tampa, FL in miles is: 57.94
</pre>
=={{header|Teradata Stored Procedure}}==
<syntaxhighlight lang="sql">
# syntax: call SP_HAVERSINE(36.12,33.94,-86.67,-118.40,x);
CREATE PROCEDURE SP_HAVERSINE
(
IN lat1 FLOAT,
IN lat2 FLOAT,
IN lon1 FLOAT,
IN lon2 FLOAT,
OUT distance FLOAT)
BEGIN
DECLARE dLat FLOAT;
DECLARE dLon FLOAT;
DECLARE c FLOAT;
DECLARE a FLOAT;
DECLARE km FLOAT;
SET dLat = RADIANS(lat2-lat1);
SET dLon = RADIANS(lon2-lon1);
SET a = SIN(dLat / 2) * SIN(dLat / 2) + SIN(dLon / 2) * SIN(dLon / 2) * COS(RADIANS(lat1)) * COS(RADIANS(lat2));
SET c = 2 * ASIN(SQRT(a));
SET km = 6372.8 * c;
select km into distance;
END;
</syntaxhighlight>
{{out}}
<pre>
distance: 2887.2599 km
</pre>
=={{header|Transact-SQL}}==
{{trans|C#}}
<syntaxhighlight lang="sql">CREATE FUNCTION [dbo].[Haversine](@Lat1 AS DECIMAL(9,7), @Lon1 AS DECIMAL(10,7), @Lat2 AS DECIMAL(9,7), @Lon2 AS DECIMAL(10,7))
RETURNS DECIMAL(12,7)
AS
BEGIN
DECLARE @R DECIMAL(11,7);
DECLARE @dLat DECIMAL(9,7);
DECLARE @dLon DECIMAL(10,7);
DECLARE @a DECIMAL(10,7);
DECLARE @c DECIMAL(10,7);
SET @R = 6372.8;
SET @dLat = RADIANS(@Lat2 - @Lat1);
SET @dLon = RADIANS(@Lon2 - @Lon1);
SET @Lat1 = RADIANS(@Lat1);
SET @Lat2 = RADIANS(@Lat2);
SET @a = SIN(@dLat / 2) * SIN(@dLat / 2) + SIN(@dLon / 2) * SIN(@dLon / 2) * COS(@Lat1) * COS(@Lat2);
SET @c = 2 * ASIN(SQRT(@a));
RETURN @R * @c;
END
GO
SELECT dbo.Haversine(36.12,-86.67,33.94,-118.4)
</syntaxhighlight>
{{out}}
<pre>
2887.2594934
</pre>
=={{header|TypeScript}}==
{{trans|Matlab}}
<syntaxhighlight lang="javascript">
let radians = function (degree: number) {
// degrees to radians
let rad: number = degree * Math.PI / 180;
return rad;
}
export const haversine = (lat1: number, lon1: number, lat2: number, lon2: number) => {
// var dlat: number, dlon: number, a: number, c: number, R: number;
let dlat, dlon, a, c, R: number;
R = 6372.8; // km
dlat = radians(lat2 - lat1);
dlon = radians(lon2 - lon1);
lat1 = radians(lat1);
lat2 = radians(lat2);
a = Math.sin(dlat / 2) * Math.sin(dlat / 2) + Math.sin(dlon / 2) * Math.sin(dlon / 2) * Math.cos(lat1) * Math.cos(lat2)
c = 2 * Math.asin(Math.sqrt(a));
return R * c;
}
console.log("Distance:" + haversine(36.12, -86.67, 33.94, -118.40));
</syntaxhighlight>
{{out}}
<pre>
Distance: 2887.2599506071106
</pre>
=={{header|UBASIC}}==
<
10 Point 7 'Sets decimal display to 32 places (0+.1^56)
20 Rf=#pi/180 'Degree -> Radian Conversion
Line 1,812 ⟶ 3,966:
2887.2599506 km
OK
</syntaxhighlight>
=={{header|VBA}}==
{{trans|Phix}}<syntaxhighlight lang="vb">Const MER = 6371 '-- mean earth radius(km)
Public DEG_TO_RAD As Double
Function haversine(lat1 As Double, long1 As Double, lat2 As Double, long2 As Double) As Double
lat1 = lat1 * DEG_TO_RAD
lat2 = lat2 * DEG_TO_RAD
long1 = long1 * DEG_TO_RAD
long2 = long2 * DEG_TO_RAD
haversine = MER * WorksheetFunction.Acos(Sin(lat1) * Sin(lat2) + Cos(lat1) * Cos(lat2) * Cos(long2 - long1))
End Function
Public Sub main()
DEG_TO_RAD = WorksheetFunction.Pi / 180
d = haversine(36.12, -86.67, 33.94, -118.4)
Debug.Print "Distance is "; Format(d, "#.######"); " km ("; Format(d / 1.609344, "#.######"); " miles)."
End Sub</syntaxhighlight>{{out}}
<pre>Distance is 2886,444443 km (1793,553425 miles).</pre>
=={{header|Visual Basic .NET}}==
{{trans|C#}}If you read the fine print in the Wikipedia article, you will find that the Haversine method of finding distances may have an error of up to 0.5%. This could lead one to believe that discussion about whether to use 6371.0 km or 6372.8 km for an approximation of the Earth's radius is moot.
<syntaxhighlight lang="vbnet">Imports System.Math
Module Module1
Const deg2rad As Double = PI / 180
Structure AP_Loc
Public IATA_Code As String, Lat As Double, Lon As Double
Public Sub New(ByVal iata_code As String, ByVal lat As Double, ByVal lon As Double)
Me.IATA_Code = iata_code : Me.Lat = lat * deg2rad : Me.Lon = lon * deg2rad
End Sub
Public Overrides Function ToString() As String
Return String.Format("{0}: ({1}, {2})", IATA_Code, Lat / deg2rad, Lon / deg2rad)
End Function
End Structure
Function Sin2(ByVal x As Double) As Double
Return Pow(Sin(x / 2), 2)
End Function
Function calculate(ByVal one As AP_Loc, ByVal two As AP_Loc) As Double
Dim R As Double = 6371, ' In kilometers, (as recommended by the International Union of Geodesy and Geophysics)
a As Double = Sin2(two.Lat - one.Lat) + Sin2(two.Lon - one.Lon) * Cos(one.Lat) * Cos(two.Lat)
Return R * 2 * Asin(Sqrt(a))
End Function
Sub ShowOne(pntA As AP_Loc, pntB as AP_Loc)
Dim adst As Double = calculate(pntA, pntB), sfx As String = "km"
If adst < 1000 Then adst *= 1000 : sfx = "m"
Console.WriteLine("The approximate distance between airports {0} and {1} is {2:n2} {3}.", pntA, pntB, adst, sfx)
Console.WriteLine("The uncertainty is under 0.5%, or {0:n1} {1}." & vbLf, adst / 200, sfx)
End Sub
' Airport coordinate data excerpted from the data base at http://www.partow.net/miscellaneous/airportdatabase/
' The four additional airports are the furthest and closest pairs, according to the "Fun Facts..." section.
' KBNA, BNA, NASHVILLE INTERNATIONAL, NASHVILLE, USA, 036, 007, 028, N, 086, 040, 041, W, 00183, 36.124, -86.678
' KLAX, LAX, LOS ANGELES INTERNATIONAL, LOS ANGELES, USA, 033, 056, 033, N, 118, 024, 029, W, 00039, 33.942, -118.408
' SKNV, NVA, BENITO SALAS, NEIVA, COLOMBIA, 002, 057, 000, N, 075, 017, 038, W, 00439, 2.950, -75.294
' WIPP, PLM, SULTAN MAHMUD BADARUDDIN II, PALEMBANG, INDONESIA, 002, 053, 052, S, 104, 042, 004, E, 00012, -2.898, 104.701
' LOWL, LNZ, HORSCHING INTERNATIONAL AIRPORT (AUS - AFB), LINZ, AUSTRIA, 048, 014, 000, N, 014, 011, 000, E, 00096, 48.233, 14.183
' LOXL, N/A, LINZ, LINZ, AUSTRIA, 048, 013, 059, N, 014, 011, 015, E, 00299, 48.233, 14.188
Sub Main()
ShowOne(New AP_Loc("BNA", 36.124, -86.678), New AP_Loc("LAX", 33.942, -118.408))
ShowOne(New AP_Loc("NVA", 2.95, -75.294), New AP_Loc("PLM", -2.898, 104.701))
ShowOne(New AP_Loc("LNZ", 48.233, 14.183), New AP_Loc("N/A", 48.233, 14.188))
End Sub
End Module</syntaxhighlight>
{{out}}
<pre>The approximate distance between airports BNA: (36.124, -86.678) and LAX: (33.942, -118.408) is 2,886.36 km.
The uncertainty is under 0.5%, or 14.4 km.
The approximate distance between airports NVA: (2.95, -75.294) and PLM: (-2.898, 104.701) is 20,009.28 km.
The uncertainty is under 0.5%, or 100.0 km.
The approximate distance between airports LNZ: (48.233, 14.183) and N/A: (48.233, 14.188) is 370.34 m.
The uncertainty is under 0.5%, or 1.9 m.</pre>Looking at the altitude difference between the last two airports, (299 - 96 = 203), the reported distance of 370 meters ought to be around 422 meters if you actually went there and saw it for yourself.
=={{header|V (Vlang)}}==
{{trans|go}}
<syntaxhighlight lang="v (vlang)">import math
fn haversine(h f64) f64 {
return .5 * (1 - math.cos(h))
}
struct Pos {
lat f64 // latitude, radians
long f64 // longitude, radians
}
fn deg_pos(lat f64, lon f64) Pos {
return Pos{lat * math.pi / 180, lon * math.pi / 180}
}
const r_earth = 6372.8 // km
fn hs_dist(p1 Pos, p2 Pos) f64 {
return 2 * r_earth * math.asin(math.sqrt(haversine(p2.lat-p1.lat)+
math.cos(p1.lat)*math.cos(p2.lat)*haversine(p2.long-p1.long)))
}
fn main() {
println(hs_dist(deg_pos(36.12, -86.67), deg_pos(33.94, -118.40)))
}</syntaxhighlight>
{{out}}
<pre>2887.2599506071</pre>
=={{header|Wren}}==
{{trans|Julia}}
<syntaxhighlight lang="wren">var R = 6372.8 // Earth's approximate radius in kilometers.
/* Class containing trig methods which work with degrees rather than radians. */
class D {
static deg2Rad(deg) { (deg*Num.pi/180 + 2*Num.pi) % (2*Num.pi) }
static sin(d) { deg2Rad(d).sin }
static cos(d) { deg2Rad(d).cos }
}
var haversine = Fn.new { |lat1, lon1, lat2, lon2|
var dlat = lat2 - lat1
var dlon = lon2 - lon1
return 2 * R * (D.sin(dlat/2).pow(2) + D.cos(lat1) * D.cos(lat2) * D.sin(dlon/2).pow(2)).sqrt.asin
}
System.print(haversine.call(36.12, -86.67, 33.94, -118.4))</syntaxhighlight>
{{out}}
<pre>
2887.2599506071
</pre>
=={{header|X86 Assembly}}==
Assemble with tasm /m /l; tlink /t
<
0000 .code
.486
Line 1,901 ⟶ 4,192:
;(TASM isn't smart enough to do floating point constant calculations)
end start
</syntaxhighlight>
{{out}}
<pre>
Line 1,908 ⟶ 4,199:
=={{header|XPL0}}==
<
func real Haversine(Ang);
Line 1,921 ⟶ 4,212:
def D2R = 3.141592654/180.0; \degrees to radians
RlOut(0, Dist(36.12*D2R, 33.94*D2R, -86.67*D2R, -118.40*D2R ));</
{{out}}
Line 1,927 ⟶ 4,218:
2887.25995
</pre>
=={{header|XQuery}}==
<syntaxhighlight lang="xquery">declare namespace xsd = "http://www.w3.org/2001/XMLSchema";
declare namespace math = "http://www.w3.org/2005/xpath-functions/math";
declare function local:haversine($lat1 as xsd:float, $lon1 as xsd:float, $lat2 as xsd:float, $lon2 as xsd:float)
as xsd:float
{
let $dlat := ($lat2 - $lat1) * math:pi() div 180
let $dlon := ($lon2 - $lon1) * math:pi() div 180
let $rlat1 := $lat1 * math:pi() div 180
let $rlat2 := $lat2 * math:pi() div 180
let $a := math:sin($dlat div 2) * math:sin($dlat div 2) + math:sin($dlon div 2) * math:sin($dlon div 2) * math:cos($rlat1) * math:cos($rlat2)
let $c := 2 * math:atan2(math:sqrt($a), math:sqrt(1-$a))
return xsd:float($c * 6371.0)
};
local:haversine(36.12, -86.67, 33.94, -118.4)</syntaxhighlight>
{{out}}
<pre>
2886.444
</pre>
=={{header|Zig}}==
{{trans|R}}
When a Zig <em>struct</em> type can be inferred then anonymous structs .{} can be used for initialisation.
This can be seen on the lines where the constants <em>bna</em> and <em>lax</em> are instantiated.
A Zig <em>struct</em> can have methods, the same as an <em>enum</em> and or a <em>union</em>.
They are only namespaced functions that can be called with dot syntax.
<syntaxhighlight lang="zig">
const std = @import("std");
const math = std.math; // Save some typing, reduce clutter. Otherwise math.sin() would be std.math.sin() etc.
pub fn main() !void {
// Coordinates are found here:
// http://www.airport-data.com/airport/BNA/
// http://www.airport-data.com/airport/LAX/
const bna = LatLong{
.lat = .{ .d = 36, .m = 7, .s = 28.10 },
.long = .{ .d = 86, .m = 40, .s = 41.50 },
};
const lax = LatLong{
.lat = .{ .d = 33, .m = 56, .s = 32.98 },
.long = .{ .d = 118, .m = 24, .s = 29.05 },
};
const distance = calcGreatCircleDistance(bna, lax);
std.debug.print("Output: {d:.6} km\n", .{distance});
// Output: 2886.326609 km
}
const LatLong = struct { lat: DMS, long: DMS };
/// degrees, minutes, decimal seconds
const DMS = struct {
d: f64,
m: f64,
s: f64,
fn toRadians(self: DMS) f64 {
return (self.d + self.m / 60 + self.s / 3600) * math.pi / 180;
}
};
// Volumetric mean radius is 6371 km, see http://nssdc.gsfc.nasa.gov/planetary/factsheet/earthfact.html
// The diameter is thus 12742 km
fn calcGreatCircleDistance(lat_long1: LatLong, lat_long2: LatLong) f64 {
const lat1 = lat_long1.lat.toRadians();
const lat2 = lat_long2.lat.toRadians();
const long1 = lat_long1.long.toRadians();
const long2 = lat_long2.long.toRadians();
const a = math.sin(0.5 * (lat2 - lat1));
const b = math.sin(0.5 * (long2 - long1));
return 12742 * math.asin(math.sqrt(a * a + math.cos(lat1) * math.cos(lat2) * b * b));
}
</syntaxhighlight>
=={{header|zkl}}==
{{trans|Erlang}}
<
fcn haversine(Lat1, Long1, Lat2, Long2){
Line 1,943 ⟶ 4,320:
C := 2.0 * A.sqrt().asin();
R*C;
}</
{{out}}
<pre>
2887.26
</pre>
=={{header|ZX Spectrum Basic}}==
{{trans|Run_BASIC}}
<syntaxhighlight lang="zxbasic">10 LET diam=2*6372.8
20 LET Lg1m2=FN r((-86.67)-(-118.4))
30 LET Lt1=FN r(36.12)
40 LET Lt2=FN r(33.94)
50 LET dz=SIN (Lt1)-SIN (Lt2)
60 LET dx=COS (Lg1m2)*COS (Lt1)-COS (Lt2)
70 LET dy=SIN (Lg1m2)*COS (Lt1)
80 LET hDist=ASN ((dx*dx+dy*dy+dz*dz)^0.5/2)*diam
90 PRINT "Haversine distance: ";hDist;" km."
100 STOP
1000 DEF FN r(a)=a*0.017453293: REM convert degree to radians</syntaxhighlight>
[[Category:Geometry]]
| |||