3 people like it.

Azimuthal equidistant projection

Simple and a more optimized implementation of the azimuthal equidistant projection. Input is expected in degrees.

 1: 
 2: 
 3: 
 4: 
 5: 
 6: 
 7: 
 8: 
 9: 
10: 
11: 
12: 
13: 
14: 
15: 
16: 
17: 
18: 
19: 
20: 
21: 
22: 
23: 
24: 
25: 
26: 
27: 
28: 
29: 
30: 
31: 
32: 
33: 
34: 
35: 
36: 
open System
module AzimuthalEquidistantProjection =

    let inline degToRad d = 0.0174532925199433 * d; // (1.0/180.0 * Math.PI) * d

    let project centerlon centerlat lon lat =
        // http://mathworld.wolfram.com/AzimuthalEquidistantProjection.html
        // http://www.radicalcartography.net/?projectionref
        let t:float = degToRad lat
        let l:float = degToRad lon
        let t1 = degToRad centerlat // latitude center of projection
        let l0 = degToRad centerlon // longitude center of projection
        let c = Math.Acos ((sin t1) * (sin t) + (cos t1) * (cos t) * (cos (l-l0)))
        let k = c / (sin c)
        let x = k * (cos t) * (sin (l-l0))
        let y = k * (cos t1) * (sin t) - (sin t1) * (cos t) * (cos (l-l0))
        (x, y)

    let project_optimized l0 t1 cost1 lon lat =
        let t:float = degToRad lat
        let l:float = degToRad lon

        let costcosll0 = (cos t) * (cos (l-l0))
        let sint = sin t
        let sint1 = sin t1
        let c = Math.Acos ((sint1) * (sint) + (cost1) * costcosll0)
        let k = c / (sin c)
        let x = k * (cos t) * (sin (l-l0))
        let y = k * (cost1) * (sint) - (sint1) * costcosll0
        (x, y)

    let buildProjection centerLon centerLat = 
        let t1 = degToRad centerLat
        let l0 = degToRad centerLon
        let cost1 = cos t1
        project_optimized l0 t1 cost1
namespace System
module AzimuthalEquidistantProjection

from Script
val degToRad : d:float -> float

Full name: Script.AzimuthalEquidistantProjection.degToRad
val d : float
val project : centerlon:float -> centerlat:float -> lon:float -> lat:float -> float * float

Full name: Script.AzimuthalEquidistantProjection.project
val centerlon : float
val centerlat : float
val lon : float
val lat : float
val t : float
Multiple items
val float : value:'T -> float (requires member op_Explicit)

Full name: Microsoft.FSharp.Core.Operators.float

--------------------
type float = Double

Full name: Microsoft.FSharp.Core.float

--------------------
type float<'Measure> = float

Full name: Microsoft.FSharp.Core.float<_>
val l : float
val t1 : float
val l0 : float
val c : float
type Math =
  static val PI : float
  static val E : float
  static member Abs : value:sbyte -> sbyte + 6 overloads
  static member Acos : d:float -> float
  static member Asin : d:float -> float
  static member Atan : d:float -> float
  static member Atan2 : y:float * x:float -> float
  static member BigMul : a:int * b:int -> int64
  static member Ceiling : d:decimal -> decimal + 1 overload
  static member Cos : d:float -> float
  ...

Full name: System.Math
Math.Acos(d: float) : float
val sin : value:'T -> 'T (requires member Sin)

Full name: Microsoft.FSharp.Core.Operators.sin
val cos : value:'T -> 'T (requires member Cos)

Full name: Microsoft.FSharp.Core.Operators.cos
val k : float
val x : float
val y : float
val project_optimized : l0:float -> t1:float -> cost1:float -> lon:float -> lat:float -> float * float

Full name: Script.AzimuthalEquidistantProjection.project_optimized
val cost1 : float
val costcosll0 : float
val sint : float
val sint1 : float
val buildProjection : centerLon:float -> centerLat:float -> (float -> float -> float * float)

Full name: Script.AzimuthalEquidistantProjection.buildProjection
val centerLon : float
val centerLat : float
Raw view Test code New version

More information

Link:http://fssnip.net/lA
Posted:10 years ago
Author:Samuel Bosch
Tags: math , spatial