-- |
-- Module:      Data.Geo.Jord.Geocentric
-- Copyright:   (c) 2020 Cedric Liegeois
-- License:     BSD3
-- Maintainer:  Cedric Liegeois <ofmooseandmen@yahoo.fr>
-- Stability:   experimental
-- Portability: portable
--
-- Geocentric coordinates of points (X, Y, and Z Cartesian coordinates) in specified models.
--
-- For the Earth the coordinate system is known as ECEF (acronym for earth-centered, earth-fixed),
-- or ECR (initialism for earth-centered rotational).
--
-- In order to use this module you should start with the following imports:
--
-- @
-- import qualified Data.Geo.Jord.Geocentric as Geocentric
-- import qualified Data.Geo.Jord.Length as Length
-- import Data.Geo.Jord.Models
-- @
--
-- see "Data.Geo.Jord.Models" for supported models.
module Data.Geo.Jord.Geocentric
    ( Position(..)
    , coords
    , metresCoords
    , metresPos
    , metresPos'
    , antipode
    , northPole
    , southPole
    ) where

import Data.Geo.Jord.Ellipsoid (polarRadius)
import Data.Geo.Jord.Length (Length)
import qualified Data.Geo.Jord.Length as Length (metres, toMetres)
import qualified Data.Geo.Jord.Math3d as Math3d
import Data.Geo.Jord.Model

-- | Geocentric coordinates (cartesian X, Y, Z) of a position in a specified 'Model'.
--
-- @gx-gy@ plane is the equatorial plane, @gx@ is on the prime meridian, and @gz@ on the polar axis.
--
-- On a spherical celestial body, an /n/-vector is equivalent to a normalised version of a
-- geocentric cartesian coordinate.
data Position a =
    Position
        { Position a -> Length
gx :: Length -- ^ x-coordinate
        , Position a -> Length
gy :: Length -- ^ y-coordinate
        , Position a -> Length
gz :: Length -- ^ z-coordinate
        , Position a -> a
model :: a -- ^ model (e.g. WGS84)
        }
    deriving (Position a -> Position a -> Bool
(Position a -> Position a -> Bool)
-> (Position a -> Position a -> Bool) -> Eq (Position a)
forall a. Eq a => Position a -> Position a -> Bool
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: Position a -> Position a -> Bool
$c/= :: forall a. Eq a => Position a -> Position a -> Bool
== :: Position a -> Position a -> Bool
$c== :: forall a. Eq a => Position a -> Position a -> Bool
Eq, Int -> Position a -> ShowS
[Position a] -> ShowS
Position a -> String
(Int -> Position a -> ShowS)
-> (Position a -> String)
-> ([Position a] -> ShowS)
-> Show (Position a)
forall a. Show a => Int -> Position a -> ShowS
forall a. Show a => [Position a] -> ShowS
forall a. Show a => Position a -> String
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [Position a] -> ShowS
$cshowList :: forall a. Show a => [Position a] -> ShowS
show :: Position a -> String
$cshow :: forall a. Show a => Position a -> String
showsPrec :: Int -> Position a -> ShowS
$cshowsPrec :: forall a. Show a => Int -> Position a -> ShowS
Show)

-- | 3d vector representing the (X, Y, Z) coordinates in __metres__ of the given position.
metresCoords :: (Model a) => Position a -> Math3d.V3
metresCoords :: Position a -> V3
metresCoords Position a
p = Position a -> (Length -> Double) -> V3
forall a. Model a => Position a -> (Length -> Double) -> V3
coords Position a
p Length -> Double
Length.toMetres

-- | @coords p f@ returns the 3d vector representing the (X, Y, Z) coordinates in the unit
-- of @f@. For example:
--
-- >>> Geocentric.coords (Geocentric.metresPos 3194669.145061 3194669.145061 4487701.962256 WGS84) Length.toKilometres
-- V3 {vx = 3194.669145061, vy = 3194.669145061, vz = 4487.701962256}
coords :: (Model a) => Position a -> (Length -> Double) -> Math3d.V3
coords :: Position a -> (Length -> Double) -> V3
coords (Position Length
x Length
y Length
z a
_) Length -> Double
f = Double -> Double -> Double -> V3
Math3d.vec3 (Length -> Double
f Length
x) (Length -> Double
f Length
y) (Length -> Double
f Length
z)

-- | Geocentric position from given (X, Y, Z) in __metres__ an given 'Model'.
metresPos :: (Model a) => Double -> Double -> Double -> a -> Position a
metresPos :: Double -> Double -> Double -> a -> Position a
metresPos Double
xm Double
ym Double
zm = Length -> Length -> Length -> a -> Position a
forall a. Length -> Length -> Length -> a -> Position a
Position (Double -> Length
Length.metres Double
xm) (Double -> Length
Length.metres Double
ym) (Double -> Length
Length.metres Double
zm)

-- | Geocentric position from given 3d vector (X, Y, Z) in __metres__ an given 'Model'.
metresPos' :: (Model a) => Math3d.V3 -> a -> Position a
metresPos' :: V3 -> a -> Position a
metresPos' V3
v = Double -> Double -> Double -> a -> Position a
forall a. Model a => Double -> Double -> Double -> a -> Position a
metresPos (V3 -> Double
Math3d.v3x V3
v) (V3 -> Double
Math3d.v3y V3
v) (V3 -> Double
Math3d.v3z V3
v)

-- | @antipode p@ computes the antipodal position of @p@: the position which is diametrically
-- opposite to @p@.
antipode :: (Model a) => Position a -> Position a
antipode :: Position a -> Position a
antipode Position a
p = Double -> Double -> Double -> a -> Position a
forall a. Model a => Double -> Double -> Double -> a -> Position a
metresPos (V3 -> Double
Math3d.v3x V3
avm) (V3 -> Double
Math3d.v3y V3
avm) (V3 -> Double
Math3d.v3z V3
avm) (Position a -> a
forall a. Position a -> a
model Position a
p)
  where
    c :: V3
c = Position a -> V3
forall a. Model a => Position a -> V3
metresCoords Position a
p
    avm :: V3
avm = V3 -> Double -> V3
Math3d.scale V3
c (-Double
1.0)

-- | Surface position of the North Pole in the given model.
northPole :: (Model a) => a -> Position a
northPole :: a -> Position a
northPole a
m = Double -> Double -> Double -> a -> Position a
forall a. Model a => Double -> Double -> Double -> a -> Position a
metresPos Double
0 Double
0 Double
r a
m
  where
    r :: Double
r = Length -> Double
Length.toMetres (Length -> Double) -> (a -> Length) -> a -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Ellipsoid -> Length
polarRadius (Ellipsoid -> Length) -> (a -> Ellipsoid) -> a -> Length
forall b c a. (b -> c) -> (a -> b) -> a -> c
. a -> Ellipsoid
forall a. Model a => a -> Ellipsoid
surface (a -> Double) -> a -> Double
forall a b. (a -> b) -> a -> b
$ a
m

-- | Surface position of the South Pole in the given model.
southPole :: (Model a) => a -> Position a
southPole :: a -> Position a
southPole a
m = Double -> Double -> Double -> a -> Position a
forall a. Model a => Double -> Double -> Double -> a -> Position a
metresPos Double
0 Double
0 (-Double
r) a
m
  where
    r :: Double
r = Length -> Double
Length.toMetres (Length -> Double) -> (a -> Length) -> a -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Ellipsoid -> Length
polarRadius (Ellipsoid -> Length) -> (a -> Ellipsoid) -> a -> Length
forall b c a. (b -> c) -> (a -> b) -> a -> c
. a -> Ellipsoid
forall a. Model a => a -> Ellipsoid
surface (a -> Double) -> a -> Double
forall a b. (a -> b) -> a -> b
$ a
m