module Data.Geo.Jord.GreatCircle
(
GreatCircle
, through
, headingOn
, MinorArc
, minorArc
, minorArcStart
, minorArcEnd
, Side(..)
, alongTrackDistance
, alongTrackDistance'
, angularDistance
, crossTrackDistance
, crossTrackDistance'
, destination
, distance
, enclosedBy
, finalBearing
, initialBearing
, interpolated
, intersection
, intersections
, mean
, projection
, side
, turn
) where
import Data.Geo.Jord.Angle (Angle)
import qualified Data.Geo.Jord.Angle as Angle
import Data.Geo.Jord.Ellipsoid (equatorialRadius)
import Data.Geo.Jord.Geodetic (HorizontalPosition)
import qualified Data.Geo.Jord.Geodetic as Geodetic
import Data.Geo.Jord.Length (Length)
import qualified Data.Geo.Jord.Length as Length
import Data.Geo.Jord.Math3d (V3)
import qualified Data.Geo.Jord.Math3d as Math3d
import Data.Geo.Jord.Model (Spherical, surface)
data GreatCircle a =
GreatCircle !V3 !a String
deriving (GreatCircle a -> GreatCircle a -> Bool
(GreatCircle a -> GreatCircle a -> Bool)
-> (GreatCircle a -> GreatCircle a -> Bool) -> Eq (GreatCircle a)
forall a. Eq a => GreatCircle a -> GreatCircle a -> Bool
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: GreatCircle a -> GreatCircle a -> Bool
$c/= :: forall a. Eq a => GreatCircle a -> GreatCircle a -> Bool
== :: GreatCircle a -> GreatCircle a -> Bool
$c== :: forall a. Eq a => GreatCircle a -> GreatCircle a -> Bool
Eq)
instance (Spherical a) => Show (GreatCircle a) where
show :: GreatCircle a -> String
show (GreatCircle V3
_ a
_ String
s) = String
s
through :: (Spherical a) => HorizontalPosition a -> HorizontalPosition a -> Maybe (GreatCircle a)
through :: HorizontalPosition a
-> HorizontalPosition a -> Maybe (GreatCircle a)
through HorizontalPosition a
p1 HorizontalPosition a
p2 = (V3 -> GreatCircle a) -> Maybe V3 -> Maybe (GreatCircle a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap (\V3
n -> V3 -> a -> String -> GreatCircle a
forall a. V3 -> a -> String -> GreatCircle a
GreatCircle V3
n (HorizontalPosition a -> a
forall a. Model a => HorizontalPosition a -> a
Geodetic.model HorizontalPosition a
p1) String
dscr) (HorizontalPosition a -> HorizontalPosition a -> Maybe V3
forall a.
Spherical a =>
HorizontalPosition a -> HorizontalPosition a -> Maybe V3
arcNormal HorizontalPosition a
p1 HorizontalPosition a
p2)
where
dscr :: String
dscr = String
"Great Circle { through " String -> ShowS
forall a. [a] -> [a] -> [a]
++ HorizontalPosition a -> String
forall a. Show a => a -> String
show HorizontalPosition a
p1 String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" & " String -> ShowS
forall a. [a] -> [a] -> [a]
++ HorizontalPosition a -> String
forall a. Show a => a -> String
show HorizontalPosition a
p2 String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" }"
headingOn :: (Spherical a) => HorizontalPosition a -> Angle -> GreatCircle a
headingOn :: HorizontalPosition a -> Angle -> GreatCircle a
headingOn HorizontalPosition a
p Angle
b = V3 -> a -> String -> GreatCircle a
forall a. V3 -> a -> String -> GreatCircle a
GreatCircle (V3 -> V3 -> V3
Math3d.subtract V3
n' V3
e') a
m String
dscr
where
m :: a
m = HorizontalPosition a -> a
forall a. Model a => HorizontalPosition a -> a
Geodetic.model HorizontalPosition a
p
v :: V3
v = HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p
e :: V3
e = V3 -> V3 -> V3
Math3d.cross (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector (HorizontalPosition a -> V3)
-> (a -> HorizontalPosition a) -> a -> V3
forall b c a. (b -> c) -> (a -> b) -> a -> c
. a -> HorizontalPosition a
forall a. Model a => a -> HorizontalPosition a
Geodetic.northPole (a -> V3) -> a -> V3
forall a b. (a -> b) -> a -> b
$ a
m) V3
v
n :: V3
n = V3 -> V3 -> V3
Math3d.cross V3
v V3
e
e' :: V3
e' = V3 -> Double -> V3
Math3d.scale V3
e (Angle -> Double
Angle.cos Angle
b Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ V3 -> Double
Math3d.norm V3
e)
n' :: V3
n' = V3 -> Double -> V3
Math3d.scale V3
n (Angle -> Double
Angle.sin Angle
b Double -> Double -> Double
forall a. Fractional a => a -> a -> a
/ V3 -> Double
Math3d.norm V3
n)
dscr :: String
dscr = String
"Great Circle { by " String -> ShowS
forall a. [a] -> [a] -> [a]
++ HorizontalPosition a -> String
forall a. Show a => a -> String
show HorizontalPosition a
p String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" & heading on " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Angle -> String
forall a. Show a => a -> String
show Angle
b String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" }"
data MinorArc a =
MinorArc !V3 (HorizontalPosition a) (HorizontalPosition a)
deriving (MinorArc a -> MinorArc a -> Bool
(MinorArc a -> MinorArc a -> Bool)
-> (MinorArc a -> MinorArc a -> Bool) -> Eq (MinorArc a)
forall a. Model a => MinorArc a -> MinorArc a -> Bool
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: MinorArc a -> MinorArc a -> Bool
$c/= :: forall a. Model a => MinorArc a -> MinorArc a -> Bool
== :: MinorArc a -> MinorArc a -> Bool
$c== :: forall a. Model a => MinorArc a -> MinorArc a -> Bool
Eq)
instance (Spherical a) => Show (MinorArc a) where
show :: MinorArc a -> String
show (MinorArc V3
_ HorizontalPosition a
s HorizontalPosition a
e) = String
"Minor Arc { from: " String -> ShowS
forall a. [a] -> [a] -> [a]
++ HorizontalPosition a -> String
forall a. Show a => a -> String
show HorizontalPosition a
s String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
", to: " String -> ShowS
forall a. [a] -> [a] -> [a]
++ HorizontalPosition a -> String
forall a. Show a => a -> String
show HorizontalPosition a
e String -> ShowS
forall a. [a] -> [a] -> [a]
++ String
" }"
minorArc :: (Spherical a) => HorizontalPosition a -> HorizontalPosition a -> Maybe (MinorArc a)
minorArc :: HorizontalPosition a -> HorizontalPosition a -> Maybe (MinorArc a)
minorArc HorizontalPosition a
p1 HorizontalPosition a
p2 = (V3 -> MinorArc a) -> Maybe V3 -> Maybe (MinorArc a)
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap (\V3
n -> V3 -> HorizontalPosition a -> HorizontalPosition a -> MinorArc a
forall a.
V3 -> HorizontalPosition a -> HorizontalPosition a -> MinorArc a
MinorArc V3
n HorizontalPosition a
p1 HorizontalPosition a
p2) (HorizontalPosition a -> HorizontalPosition a -> Maybe V3
forall a.
Spherical a =>
HorizontalPosition a -> HorizontalPosition a -> Maybe V3
arcNormal HorizontalPosition a
p1 HorizontalPosition a
p2)
minorArcStart :: (Spherical a) => MinorArc a -> HorizontalPosition a
minorArcStart :: MinorArc a -> HorizontalPosition a
minorArcStart (MinorArc V3
_ HorizontalPosition a
s HorizontalPosition a
_) = HorizontalPosition a
s
minorArcEnd :: (Spherical a) => MinorArc a -> HorizontalPosition a
minorArcEnd :: MinorArc a -> HorizontalPosition a
minorArcEnd (MinorArc V3
_ HorizontalPosition a
_ HorizontalPosition a
e) = HorizontalPosition a
e
data Side
= LeftOf
| RightOf
| None
deriving (Side -> Side -> Bool
(Side -> Side -> Bool) -> (Side -> Side -> Bool) -> Eq Side
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: Side -> Side -> Bool
$c/= :: Side -> Side -> Bool
== :: Side -> Side -> Bool
$c== :: Side -> Side -> Bool
Eq, Int -> Side -> ShowS
[Side] -> ShowS
Side -> String
(Int -> Side -> ShowS)
-> (Side -> String) -> ([Side] -> ShowS) -> Show Side
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [Side] -> ShowS
$cshowList :: [Side] -> ShowS
show :: Side -> String
$cshow :: Side -> String
showsPrec :: Int -> Side -> ShowS
$cshowsPrec :: Int -> Side -> ShowS
Show)
alongTrackDistance :: (Spherical a) => HorizontalPosition a -> MinorArc a -> Length
alongTrackDistance :: HorizontalPosition a -> MinorArc a -> Length
alongTrackDistance HorizontalPosition a
p (MinorArc V3
n HorizontalPosition a
s HorizontalPosition a
_) = HorizontalPosition a -> HorizontalPosition a -> V3 -> Length
forall a.
Spherical a =>
HorizontalPosition a -> HorizontalPosition a -> V3 -> Length
alongTrackDistance'' HorizontalPosition a
p HorizontalPosition a
s V3
n
alongTrackDistance' ::
(Spherical a) => HorizontalPosition a -> HorizontalPosition a -> Angle -> Length
alongTrackDistance' :: HorizontalPosition a -> HorizontalPosition a -> Angle -> Length
alongTrackDistance' HorizontalPosition a
p HorizontalPosition a
s Angle
b = HorizontalPosition a -> HorizontalPosition a -> V3 -> Length
forall a.
Spherical a =>
HorizontalPosition a -> HorizontalPosition a -> V3 -> Length
alongTrackDistance'' HorizontalPosition a
p HorizontalPosition a
s V3
n
where
(GreatCircle V3
n a
_ String
_) = HorizontalPosition a -> Angle -> GreatCircle a
forall a.
Spherical a =>
HorizontalPosition a -> Angle -> GreatCircle a
headingOn HorizontalPosition a
s Angle
b
angularDistance ::
(Spherical a)
=> HorizontalPosition a
-> HorizontalPosition a
-> Maybe (HorizontalPosition a)
-> Angle
angularDistance :: HorizontalPosition a
-> HorizontalPosition a -> Maybe (HorizontalPosition a) -> Angle
angularDistance HorizontalPosition a
p1 HorizontalPosition a
p2 Maybe (HorizontalPosition a)
n = V3 -> V3 -> Maybe V3 -> Angle
signedAngleBetween V3
v1 V3
v2 Maybe V3
vn
where
v1 :: V3
v1 = HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p1
v2 :: V3
v2 = HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p2
vn :: Maybe V3
vn = (HorizontalPosition a -> V3)
-> Maybe (HorizontalPosition a) -> Maybe V3
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector Maybe (HorizontalPosition a)
n
crossTrackDistance :: (Spherical a) => HorizontalPosition a -> GreatCircle a -> Length
crossTrackDistance :: HorizontalPosition a -> GreatCircle a -> Length
crossTrackDistance HorizontalPosition a
p (GreatCircle V3
n a
_ String
_) =
Angle -> Length -> Length
Angle.arcLength (Angle -> Angle -> Angle
Angle.subtract Angle
a (Double -> Angle
Angle.decimalDegrees Double
90)) (HorizontalPosition a -> Length
forall a. Spherical a => HorizontalPosition a -> Length
radius HorizontalPosition a
p)
where
a :: Angle
a = V3 -> V3 -> Angle
angleBetween V3
n (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p)
crossTrackDistance' ::
(Spherical a) => HorizontalPosition a -> HorizontalPosition a -> Angle -> Length
crossTrackDistance' :: HorizontalPosition a -> HorizontalPosition a -> Angle -> Length
crossTrackDistance' HorizontalPosition a
p HorizontalPosition a
s Angle
b = HorizontalPosition a -> GreatCircle a -> Length
forall a.
Spherical a =>
HorizontalPosition a -> GreatCircle a -> Length
crossTrackDistance HorizontalPosition a
p (HorizontalPosition a -> Angle -> GreatCircle a
forall a.
Spherical a =>
HorizontalPosition a -> Angle -> GreatCircle a
headingOn HorizontalPosition a
s Angle
b)
destination :: (Spherical a) => HorizontalPosition a -> Angle -> Length -> HorizontalPosition a
destination :: HorizontalPosition a -> Angle -> Length -> HorizontalPosition a
destination HorizontalPosition a
p Angle
b Length
d
| Length
d Length -> Length -> Bool
forall a. Eq a => a -> a -> Bool
== Length
Length.zero = HorizontalPosition a
p
| Bool
otherwise = V3 -> a -> HorizontalPosition a
forall a. Model a => V3 -> a -> HorizontalPosition a
Geodetic.nvectorPos' V3
nvd a
m
where
m :: a
m = HorizontalPosition a -> a
forall a. Model a => HorizontalPosition a -> a
Geodetic.model HorizontalPosition a
p
nv :: V3
nv = HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p
ed :: V3
ed = V3 -> V3
Math3d.unit (V3 -> V3 -> V3
Math3d.cross (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector (HorizontalPosition a -> V3)
-> (a -> HorizontalPosition a) -> a -> V3
forall b c a. (b -> c) -> (a -> b) -> a -> c
. a -> HorizontalPosition a
forall a. Model a => a -> HorizontalPosition a
Geodetic.northPole (a -> V3) -> a -> V3
forall a b. (a -> b) -> a -> b
$ a
m) V3
nv)
nd :: V3
nd = V3 -> V3 -> V3
Math3d.cross V3
nv V3
ed
r :: Length
r = HorizontalPosition a -> Length
forall a. Spherical a => HorizontalPosition a -> Length
radius HorizontalPosition a
p
ta :: Angle
ta = Length -> Length -> Angle
Angle.central Length
d Length
r
de :: V3
de = V3 -> V3 -> V3
Math3d.add (V3 -> Double -> V3
Math3d.scale V3
nd (Angle -> Double
Angle.cos Angle
b)) (V3 -> Double -> V3
Math3d.scale V3
ed (Angle -> Double
Angle.sin Angle
b))
nvd :: V3
nvd = V3 -> V3 -> V3
Math3d.add (V3 -> Double -> V3
Math3d.scale V3
nv (Angle -> Double
Angle.cos Angle
ta)) (V3 -> Double -> V3
Math3d.scale V3
de (Angle -> Double
Angle.sin Angle
ta))
distance :: (Spherical a) => HorizontalPosition a -> HorizontalPosition a -> Length
distance :: HorizontalPosition a -> HorizontalPosition a -> Length
distance HorizontalPosition a
p1 HorizontalPosition a
p2 = Angle -> Length -> Length
Angle.arcLength Angle
a (HorizontalPosition a -> Length
forall a. Spherical a => HorizontalPosition a -> Length
radius HorizontalPosition a
p1)
where
a :: Angle
a = V3 -> V3 -> Angle
angleBetween (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p1) (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p2)
enclosedBy :: (Spherical a) => HorizontalPosition a -> [HorizontalPosition a] -> Bool
enclosedBy :: HorizontalPosition a -> [HorizontalPosition a] -> Bool
enclosedBy HorizontalPosition a
p [HorizontalPosition a]
ps
| [HorizontalPosition a] -> Bool
forall (t :: * -> *) a. Foldable t => t a -> Bool
null [HorizontalPosition a]
ps = Bool
False
| HorizontalPosition a
p HorizontalPosition a -> [HorizontalPosition a] -> Bool
forall (t :: * -> *) a. (Foldable t, Eq a) => a -> t a -> Bool
`elem` [HorizontalPosition a]
ps = Bool
False
| [HorizontalPosition a] -> HorizontalPosition a
forall a. [a] -> a
head [HorizontalPosition a]
ps HorizontalPosition a -> HorizontalPosition a -> Bool
forall a. Eq a => a -> a -> Bool
== [HorizontalPosition a] -> HorizontalPosition a
forall a. [a] -> a
last [HorizontalPosition a]
ps = HorizontalPosition a -> [HorizontalPosition a] -> Bool
forall a.
Spherical a =>
HorizontalPosition a -> [HorizontalPosition a] -> Bool
enclosedBy HorizontalPosition a
p ([HorizontalPosition a] -> [HorizontalPosition a]
forall a. [a] -> [a]
init [HorizontalPosition a]
ps)
| [HorizontalPosition a] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [HorizontalPosition a]
ps Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
< Int
3 = Bool
False
| Bool
otherwise =
let aSum :: Angle
aSum =
(Angle -> (V3, V3) -> Angle) -> Angle -> [(V3, V3)] -> Angle
forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl
(\Angle
a (V3, V3)
v' -> Angle -> Angle -> Angle
Angle.add Angle
a ((V3 -> V3 -> Maybe V3 -> Angle) -> (V3, V3) -> Maybe V3 -> Angle
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry V3 -> V3 -> Maybe V3 -> Angle
signedAngleBetween (V3, V3)
v' (V3 -> Maybe V3
forall a. a -> Maybe a
Just V3
v)))
Angle
Angle.zero
([V3] -> [(V3, V3)]
egdes ((V3 -> V3) -> [V3] -> [V3]
forall a b. (a -> b) -> [a] -> [b]
map (V3 -> V3 -> V3
Math3d.subtract V3
v) [V3]
vs))
in Double -> Double
forall a. Num a => a -> a
abs (Angle -> Double
Angle.toDecimalDegrees Angle
aSum) Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
> Double
180.0
where
v :: V3
v = HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p
vs :: [V3]
vs = (HorizontalPosition a -> V3) -> [HorizontalPosition a] -> [V3]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector [HorizontalPosition a]
ps
finalBearing :: (Spherical a) => HorizontalPosition a -> HorizontalPosition a -> Maybe Angle
finalBearing :: HorizontalPosition a -> HorizontalPosition a -> Maybe Angle
finalBearing HorizontalPosition a
p1 HorizontalPosition a
p2
| HorizontalPosition a
p1 HorizontalPosition a -> HorizontalPosition a -> Bool
forall a. Eq a => a -> a -> Bool
== HorizontalPosition a
p2 = Maybe Angle
forall a. Maybe a
Nothing
| Bool
otherwise = Angle -> Maybe Angle
forall a. a -> Maybe a
Just (Angle -> Angle -> Angle
Angle.normalise (HorizontalPosition a -> HorizontalPosition a -> Angle
forall a.
Spherical a =>
HorizontalPosition a -> HorizontalPosition a -> Angle
initialBearing' HorizontalPosition a
p2 HorizontalPosition a
p1) (Double -> Angle
Angle.decimalDegrees Double
180))
initialBearing :: (Spherical a) => HorizontalPosition a -> HorizontalPosition a -> Maybe Angle
initialBearing :: HorizontalPosition a -> HorizontalPosition a -> Maybe Angle
initialBearing HorizontalPosition a
p1 HorizontalPosition a
p2
| HorizontalPosition a
p1 HorizontalPosition a -> HorizontalPosition a -> Bool
forall a. Eq a => a -> a -> Bool
== HorizontalPosition a
p2 = Maybe Angle
forall a. Maybe a
Nothing
| Bool
otherwise = Angle -> Maybe Angle
forall a. a -> Maybe a
Just (HorizontalPosition a -> HorizontalPosition a -> Angle
forall a.
Spherical a =>
HorizontalPosition a -> HorizontalPosition a -> Angle
initialBearing' HorizontalPosition a
p1 HorizontalPosition a
p2)
interpolated ::
(Spherical a)
=> HorizontalPosition a
-> HorizontalPosition a
-> Double
-> HorizontalPosition a
interpolated :: HorizontalPosition a
-> HorizontalPosition a -> Double -> HorizontalPosition a
interpolated HorizontalPosition a
p0 HorizontalPosition a
p1 Double
f
| Double
f Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
< Double
0 Bool -> Bool -> Bool
|| Double
f Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
> Double
1 = String -> HorizontalPosition a
forall a. HasCallStack => String -> a
error (String
"fraction must be in range [0..1], was " String -> ShowS
forall a. [a] -> [a] -> [a]
++ Double -> String
forall a. Show a => a -> String
show Double
f)
| Double
f Double -> Double -> Bool
forall a. Eq a => a -> a -> Bool
== Double
0 = HorizontalPosition a
p0
| Double
f Double -> Double -> Bool
forall a. Eq a => a -> a -> Bool
== Double
1 = HorizontalPosition a
p1
| Bool
otherwise = V3 -> a -> HorizontalPosition a
forall a. Model a => V3 -> a -> HorizontalPosition a
Geodetic.nvectorPos' V3
iv (HorizontalPosition a -> a
forall a. Model a => HorizontalPosition a -> a
Geodetic.model HorizontalPosition a
p0)
where
nv0 :: V3
nv0 = HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p0
nv1 :: V3
nv1 = HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p1
iv :: V3
iv = V3 -> V3
Math3d.unit (V3 -> V3 -> V3
Math3d.add V3
nv0 (V3 -> Double -> V3
Math3d.scale (V3 -> V3 -> V3
Math3d.subtract V3
nv1 V3
nv0) Double
f))
intersection :: (Spherical a) => MinorArc a -> MinorArc a -> Maybe (HorizontalPosition a)
intersection :: MinorArc a -> MinorArc a -> Maybe (HorizontalPosition a)
intersection a1 :: MinorArc a
a1@(MinorArc V3
n1 HorizontalPosition a
s1 HorizontalPosition a
e1) a2 :: MinorArc a
a2@(MinorArc V3
n2 HorizontalPosition a
s2 HorizontalPosition a
e2) =
case V3 -> V3 -> a -> Maybe (HorizontalPosition a, HorizontalPosition a)
forall a.
Spherical a =>
V3 -> V3 -> a -> Maybe (HorizontalPosition a, HorizontalPosition a)
intersections' V3
n1 V3
n2 (HorizontalPosition a -> a
forall a. Model a => HorizontalPosition a -> a
Geodetic.model HorizontalPosition a
s1) of
Maybe (HorizontalPosition a, HorizontalPosition a)
Nothing -> Maybe (HorizontalPosition a)
forall a. Maybe a
Nothing
(Just (HorizontalPosition a
i1, HorizontalPosition a
i2))
| HorizontalPosition a -> MinorArc a -> Bool
forall a. Spherical a => HorizontalPosition a -> MinorArc a -> Bool
onMinorArc HorizontalPosition a
pot MinorArc a
a1 Bool -> Bool -> Bool
&& HorizontalPosition a -> MinorArc a -> Bool
forall a. Spherical a => HorizontalPosition a -> MinorArc a -> Bool
onMinorArc HorizontalPosition a
pot MinorArc a
a2 -> HorizontalPosition a -> Maybe (HorizontalPosition a)
forall a. a -> Maybe a
Just HorizontalPosition a
pot
| Bool
otherwise -> Maybe (HorizontalPosition a)
forall a. Maybe a
Nothing
where mid :: V3
mid =
[V3] -> V3
meanV
[ HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
s1
, HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
e1
, HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
s2
, HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
e2
]
pot :: HorizontalPosition a
pot =
if V3 -> V3 -> Double
Math3d.dot (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
i1) V3
mid Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
> Double
0
then HorizontalPosition a
i1
else HorizontalPosition a
i2
intersections ::
(Spherical a)
=> GreatCircle a
-> GreatCircle a
-> Maybe (HorizontalPosition a, HorizontalPosition a)
intersections :: GreatCircle a
-> GreatCircle a
-> Maybe (HorizontalPosition a, HorizontalPosition a)
intersections (GreatCircle V3
n1 a
m String
_) (GreatCircle V3
n2 a
_ String
_) = V3 -> V3 -> a -> Maybe (HorizontalPosition a, HorizontalPosition a)
forall a.
Spherical a =>
V3 -> V3 -> a -> Maybe (HorizontalPosition a, HorizontalPosition a)
intersections' V3
n1 V3
n2 a
m
mean :: (Spherical a) => [HorizontalPosition a] -> Maybe (HorizontalPosition a)
mean :: [HorizontalPosition a] -> Maybe (HorizontalPosition a)
mean [] = Maybe (HorizontalPosition a)
forall a. Maybe a
Nothing
mean [HorizontalPosition a
p] = HorizontalPosition a -> Maybe (HorizontalPosition a)
forall a. a -> Maybe a
Just HorizontalPosition a
p
mean [HorizontalPosition a]
ps =
if (HorizontalPosition a -> Bool) -> [HorizontalPosition a] -> Bool
forall (t :: * -> *) a. Foldable t => (a -> Bool) -> t a -> Bool
any (HorizontalPosition a -> [HorizontalPosition a] -> Bool
forall (t :: * -> *) a. (Foldable t, Eq a) => a -> t a -> Bool
`elem` [HorizontalPosition a]
ps) [HorizontalPosition a]
antipodes
then Maybe (HorizontalPosition a)
forall a. Maybe a
Nothing
else HorizontalPosition a -> Maybe (HorizontalPosition a)
forall a. a -> Maybe a
Just
(V3 -> a -> HorizontalPosition a
forall a. Model a => V3 -> a -> HorizontalPosition a
Geodetic.nvectorPos'
([V3] -> V3
meanV ((HorizontalPosition a -> V3) -> [HorizontalPosition a] -> [V3]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector [HorizontalPosition a]
ps))
(HorizontalPosition a -> a
forall a. Model a => HorizontalPosition a -> a
Geodetic.model (HorizontalPosition a -> a)
-> ([HorizontalPosition a] -> HorizontalPosition a)
-> [HorizontalPosition a]
-> a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. [HorizontalPosition a] -> HorizontalPosition a
forall a. [a] -> a
head ([HorizontalPosition a] -> a) -> [HorizontalPosition a] -> a
forall a b. (a -> b) -> a -> b
$ [HorizontalPosition a]
ps))
where
antipodes :: [HorizontalPosition a]
antipodes = (HorizontalPosition a -> HorizontalPosition a)
-> [HorizontalPosition a] -> [HorizontalPosition a]
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap HorizontalPosition a -> HorizontalPosition a
forall a. Model a => HorizontalPosition a -> HorizontalPosition a
Geodetic.antipode [HorizontalPosition a]
ps
projection :: (Spherical a) => HorizontalPosition a -> MinorArc a -> Maybe (HorizontalPosition a)
projection :: HorizontalPosition a -> MinorArc a -> Maybe (HorizontalPosition a)
projection HorizontalPosition a
p ma :: MinorArc a
ma@(MinorArc V3
na HorizontalPosition a
mas HorizontalPosition a
mae) =
case Maybe V3
mnb of
Maybe V3
Nothing -> Maybe (HorizontalPosition a)
forall a. Maybe a
Nothing
(Just V3
nb) ->
case Maybe (HorizontalPosition a, HorizontalPosition a)
is of
Maybe (HorizontalPosition a, HorizontalPosition a)
Nothing -> Maybe (HorizontalPosition a)
forall a. Maybe a
Nothing
(Just (HorizontalPosition a
i1, HorizontalPosition a
i2)) ->
if HorizontalPosition a -> MinorArc a -> Bool
forall a. Spherical a => HorizontalPosition a -> MinorArc a -> Bool
onMinorArc HorizontalPosition a
pot MinorArc a
ma
then HorizontalPosition a -> Maybe (HorizontalPosition a)
forall a. a -> Maybe a
Just HorizontalPosition a
pot
else Maybe (HorizontalPosition a)
forall a. Maybe a
Nothing
where mid :: V3
mid =
[V3] -> V3
meanV
[ HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
mas
, HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
mae
, HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
nap
, HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
nbp
]
pot :: HorizontalPosition a
pot =
if V3 -> V3 -> Double
Math3d.dot V3
mid (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
i1) Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
> Double
0
then HorizontalPosition a
i1
else HorizontalPosition a
i2
where nbp :: HorizontalPosition a
nbp = V3 -> a -> HorizontalPosition a
forall a. Model a => V3 -> a -> HorizontalPosition a
Geodetic.nvectorPos' V3
nb a
m
is :: Maybe (HorizontalPosition a, HorizontalPosition a)
is = V3 -> V3 -> a -> Maybe (HorizontalPosition a, HorizontalPosition a)
forall a.
Spherical a =>
V3 -> V3 -> a -> Maybe (HorizontalPosition a, HorizontalPosition a)
intersections' (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
nap) (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
nbp) a
m
where
m :: a
m = HorizontalPosition a -> a
forall a. Model a => HorizontalPosition a -> a
Geodetic.model HorizontalPosition a
p
nap :: HorizontalPosition a
nap = V3 -> a -> HorizontalPosition a
forall a. Model a => V3 -> a -> HorizontalPosition a
Geodetic.nvectorPos' V3
na a
m
mnb :: Maybe V3
mnb = HorizontalPosition a -> HorizontalPosition a -> Maybe V3
forall a.
Spherical a =>
HorizontalPosition a -> HorizontalPosition a -> Maybe V3
arcNormal HorizontalPosition a
nap HorizontalPosition a
p
side ::
(Spherical a) => HorizontalPosition a -> HorizontalPosition a -> HorizontalPosition a -> Side
side :: HorizontalPosition a
-> HorizontalPosition a -> HorizontalPosition a -> Side
side HorizontalPosition a
p0 HorizontalPosition a
p1 HorizontalPosition a
p2
| HorizontalPosition a
p0 HorizontalPosition a -> HorizontalPosition a -> Bool
forall a. Eq a => a -> a -> Bool
== HorizontalPosition a
p1 Bool -> Bool -> Bool
|| HorizontalPosition a
p0 HorizontalPosition a -> HorizontalPosition a -> Bool
forall a. Eq a => a -> a -> Bool
== HorizontalPosition a
p2 = Side
None
| Bool
otherwise = Side -> (Double -> Side) -> Maybe Double -> Side
forall b a. b -> (a -> b) -> Maybe a -> b
maybe Side
None Double -> Side
toSide Maybe Double
ms
where
ms :: Maybe Double
ms = (V3 -> Double) -> Maybe V3 -> Maybe Double
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap (V3 -> V3 -> Double
Math3d.dot (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p0)) (HorizontalPosition a -> HorizontalPosition a -> Maybe V3
forall a.
Spherical a =>
HorizontalPosition a -> HorizontalPosition a -> Maybe V3
arcNormal HorizontalPosition a
p1 HorizontalPosition a
p2)
turn ::
(Spherical a)
=> HorizontalPosition a
-> HorizontalPosition a
-> HorizontalPosition a
-> Angle
turn :: HorizontalPosition a
-> HorizontalPosition a -> HorizontalPosition a -> Angle
turn HorizontalPosition a
a HorizontalPosition a
b HorizontalPosition a
c =
case (Maybe V3, Maybe V3)
ns of
(Just V3
n1, Just V3
n2) -> V3 -> V3 -> Maybe V3 -> Angle
signedAngleBetween V3
n1 V3
n2 (V3 -> Maybe V3
forall a. a -> Maybe a
Just (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
b))
(Maybe V3
_, Maybe V3
_) -> Angle
Angle.zero
where
ns :: (Maybe V3, Maybe V3)
ns = ((V3 -> V3) -> Maybe V3 -> Maybe V3
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap V3 -> V3
Math3d.unit (HorizontalPosition a -> HorizontalPosition a -> Maybe V3
forall a.
Spherical a =>
HorizontalPosition a -> HorizontalPosition a -> Maybe V3
arcNormal HorizontalPosition a
a HorizontalPosition a
b), (V3 -> V3) -> Maybe V3 -> Maybe V3
forall (f :: * -> *) a b. Functor f => (a -> b) -> f a -> f b
fmap V3 -> V3
Math3d.unit (HorizontalPosition a -> HorizontalPosition a -> Maybe V3
forall a.
Spherical a =>
HorizontalPosition a -> HorizontalPosition a -> Maybe V3
arcNormal HorizontalPosition a
b HorizontalPosition a
c))
alongTrackDistance'' ::
(Spherical a) => HorizontalPosition a -> HorizontalPosition a -> V3 -> Length
alongTrackDistance'' :: HorizontalPosition a -> HorizontalPosition a -> V3 -> Length
alongTrackDistance'' HorizontalPosition a
p HorizontalPosition a
s V3
n = Angle -> Length -> Length
Angle.arcLength Angle
a (HorizontalPosition a -> Length
forall a. Spherical a => HorizontalPosition a -> Length
radius HorizontalPosition a
s)
where
a :: Angle
a =
V3 -> V3 -> Maybe V3 -> Angle
signedAngleBetween
(HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
s)
(V3 -> V3 -> V3
Math3d.cross (V3 -> V3 -> V3
Math3d.cross V3
n (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p)) V3
n)
(V3 -> Maybe V3
forall a. a -> Maybe a
Just V3
n)
egdes :: [V3] -> [(V3, V3)]
egdes :: [V3] -> [(V3, V3)]
egdes [V3]
ps = [V3] -> [V3] -> [(V3, V3)]
forall a b. [a] -> [b] -> [(a, b)]
zip [V3]
ps ([V3] -> [V3]
forall a. [a] -> [a]
tail [V3]
ps [V3] -> [V3] -> [V3]
forall a. [a] -> [a] -> [a]
++ [[V3] -> V3
forall a. [a] -> a
head [V3]
ps])
intersections' ::
(Spherical a) => V3 -> V3 -> a -> Maybe (HorizontalPosition a, HorizontalPosition a)
intersections' :: V3 -> V3 -> a -> Maybe (HorizontalPosition a, HorizontalPosition a)
intersections' V3
n1 V3
n2 a
s
| (V3 -> Double
Math3d.norm V3
i :: Double) Double -> Double -> Bool
forall a. Eq a => a -> a -> Bool
== Double
0.0 = Maybe (HorizontalPosition a, HorizontalPosition a)
forall a. Maybe a
Nothing
| Bool
otherwise
, let ni :: HorizontalPosition a
ni = V3 -> a -> HorizontalPosition a
forall a. Model a => V3 -> a -> HorizontalPosition a
Geodetic.nvectorPos' (V3 -> V3
Math3d.unit V3
i) a
s = (HorizontalPosition a, HorizontalPosition a)
-> Maybe (HorizontalPosition a, HorizontalPosition a)
forall a. a -> Maybe a
Just (HorizontalPosition a
ni, HorizontalPosition a -> HorizontalPosition a
forall a. Model a => HorizontalPosition a -> HorizontalPosition a
Geodetic.antipode HorizontalPosition a
ni)
where
i :: V3
i = V3 -> V3 -> V3
Math3d.cross V3
n1 V3
n2
initialBearing' :: (Spherical a) => HorizontalPosition a -> HorizontalPosition a -> Angle
initialBearing' :: HorizontalPosition a -> HorizontalPosition a -> Angle
initialBearing' HorizontalPosition a
p1 HorizontalPosition a
p2 = Angle -> Angle -> Angle
Angle.normalise Angle
a (Double -> Angle
Angle.decimalDegrees Double
360)
where
m :: a
m = HorizontalPosition a -> a
forall a. Model a => HorizontalPosition a -> a
Geodetic.model HorizontalPosition a
p1
v1 :: V3
v1 = HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p1
v2 :: V3
v2 = HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p2
gc1 :: V3
gc1 = V3 -> V3 -> V3
Math3d.cross V3
v1 V3
v2
gc2 :: V3
gc2 = V3 -> V3 -> V3
Math3d.cross V3
v1 (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector (HorizontalPosition a -> V3)
-> (a -> HorizontalPosition a) -> a -> V3
forall b c a. (b -> c) -> (a -> b) -> a -> c
. a -> HorizontalPosition a
forall a. Model a => a -> HorizontalPosition a
Geodetic.northPole (a -> V3) -> a -> V3
forall a b. (a -> b) -> a -> b
$ a
m)
a :: Angle
a = V3 -> V3 -> Maybe V3 -> Angle
signedAngleBetween V3
gc1 V3
gc2 (V3 -> Maybe V3
forall a. a -> Maybe a
Just V3
v1)
arcNormal :: (Spherical a) => HorizontalPosition a -> HorizontalPosition a -> Maybe V3
arcNormal :: HorizontalPosition a -> HorizontalPosition a -> Maybe V3
arcNormal HorizontalPosition a
p1 HorizontalPosition a
p2
| HorizontalPosition a
p1 HorizontalPosition a -> HorizontalPosition a -> Bool
forall a. Eq a => a -> a -> Bool
== HorizontalPosition a
p2 = Maybe V3
forall a. Maybe a
Nothing
| HorizontalPosition a -> HorizontalPosition a
forall a. Model a => HorizontalPosition a -> HorizontalPosition a
Geodetic.antipode HorizontalPosition a
p1 HorizontalPosition a -> HorizontalPosition a -> Bool
forall a. Eq a => a -> a -> Bool
== HorizontalPosition a
p2 = Maybe V3
forall a. Maybe a
Nothing
| Bool
otherwise = V3 -> Maybe V3
forall a. a -> Maybe a
Just (V3 -> V3 -> V3
Math3d.cross (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p1) (HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p2))
onMinorArc :: (Spherical a) => HorizontalPosition a -> MinorArc a -> Bool
onMinorArc :: HorizontalPosition a -> MinorArc a -> Bool
onMinorArc HorizontalPosition a
p (MinorArc V3
_ HorizontalPosition a
s HorizontalPosition a
e) =
V3 -> V3 -> Double
Math3d.squaredDistance V3
v0 V3
v1 Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
<= Double
l Bool -> Bool -> Bool
&& V3 -> V3 -> Double
Math3d.squaredDistance V3
v0 V3
v2 Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
<= Double
l
where
v0 :: V3
v0 = HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
p
v1 :: V3
v1 = HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
s
v2 :: V3
v2 = HorizontalPosition a -> V3
forall a. HasCoordinates a => a -> V3
Geodetic.nvector HorizontalPosition a
e
l :: Double
l = V3 -> V3 -> Double
Math3d.squaredDistance V3
v1 V3
v2
radius :: (Spherical a) => HorizontalPosition a -> Length
radius :: HorizontalPosition a -> Length
radius = Ellipsoid -> Length
equatorialRadius (Ellipsoid -> Length)
-> (HorizontalPosition a -> Ellipsoid)
-> HorizontalPosition a
-> Length
forall b c a. (b -> c) -> (a -> b) -> a -> c
. a -> Ellipsoid
forall a. Model a => a -> Ellipsoid
surface (a -> Ellipsoid)
-> (HorizontalPosition a -> a) -> HorizontalPosition a -> Ellipsoid
forall b c a. (b -> c) -> (a -> b) -> a -> c
. HorizontalPosition a -> a
forall a. Model a => HorizontalPosition a -> a
Geodetic.model
angleBetween :: V3 -> V3 -> Angle
angleBetween :: V3 -> V3 -> Angle
angleBetween V3
v1 V3
v2 = V3 -> V3 -> Maybe V3 -> Angle
signedAngleBetween V3
v1 V3
v2 Maybe V3
forall a. Maybe a
Nothing
signedAngleBetween :: V3 -> V3 -> Maybe V3 -> Angle
signedAngleBetween :: V3 -> V3 -> Maybe V3 -> Angle
signedAngleBetween V3
v1 V3
v2 Maybe V3
n = Double -> Double -> Angle
Angle.atan2 Double
sinO Double
cosO
where
sign :: Double
sign = Double -> (V3 -> Double) -> Maybe V3 -> Double
forall b a. b -> (a -> b) -> Maybe a -> b
maybe Double
1 (Double -> Double
forall a. Num a => a -> a
signum (Double -> Double) -> (V3 -> Double) -> V3 -> Double
forall b c a. (b -> c) -> (a -> b) -> a -> c
. V3 -> V3 -> Double
Math3d.dot (V3 -> V3 -> V3
Math3d.cross V3
v1 V3
v2)) Maybe V3
n
sinO :: Double
sinO = Double
sign Double -> Double -> Double
forall a. Num a => a -> a -> a
* V3 -> Double
Math3d.norm (V3 -> V3 -> V3
Math3d.cross V3
v1 V3
v2)
cosO :: Double
cosO = V3 -> V3 -> Double
Math3d.dot V3
v1 V3
v2
meanV :: [Math3d.V3] -> V3
meanV :: [V3] -> V3
meanV [V3]
vs = V3 -> V3
Math3d.unit (V3 -> V3) -> V3 -> V3
forall a b. (a -> b) -> a -> b
$ (V3 -> V3 -> V3) -> V3 -> [V3] -> V3
forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl V3 -> V3 -> V3
Math3d.add V3
Math3d.zero [V3]
vs
toSide :: Double -> Side
toSide :: Double -> Side
toSide Double
s
| Double
s Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
< Double
0 = Side
RightOf
| Double
s Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
> Double
0 = Side
LeftOf
| Bool
otherwise = Side
None