25 lines
694 B
FSharp
25 lines
694 B
FSharp
open System
|
|
|
|
[<Measure>] type deg
|
|
[<Measure>] type rad
|
|
[<Measure>] type km
|
|
|
|
let haversine (θ: float<rad>) = 0.5 * (1.0 - Math.Cos(θ/1.0<rad>))
|
|
|
|
let radPerDeg = (Math.PI / 180.0) * 1.0<rad/deg>
|
|
|
|
type pos(latitude: float<deg>, longitude: float<deg>) =
|
|
member this.φ = latitude * radPerDeg
|
|
member this.ψ = longitude * radPerDeg
|
|
|
|
let rEarth = 6372.8<km>
|
|
|
|
let hsDist (p1: pos) (p2: pos) =
|
|
2.0 * rEarth *
|
|
Math.Asin(Math.Sqrt(haversine(p2.φ - p1.φ)+
|
|
Math.Cos(p1.φ/1.0<rad>)*Math.Cos(p2.φ/1.0<rad>)*haversine(p2.ψ - p1.ψ)))
|
|
|
|
[<EntryPoint>]
|
|
let main argv =
|
|
printfn "%A" (hsDist (pos(36.12<deg>, -86.67<deg>)) (pos(33.94<deg>, -118.40<deg>)))
|
|
0
|