module Data.GPS
(
Coordinate(..)
, Location(..)
, DMSCoordinate(..)
, DMS(..)
, Heading
, Speed
, Vector
, Trail
, smoothTrails
, longestSmoothTrail
, restLocations
, closestDistance
, normalizeDMS
, KML
, trailToKML
, pointsToKML
, kmlToString
) where
import Data.Ord (comparing)
import Data.List (sort, mapAccumL, maximumBy)
import Data.Binary
import Data.Binary.Put
import Data.Binary.Get
import Text.PrettyPrint.HughesPJClass
import Text.PrettyPrint
import Numeric (showFFloat)
import Data.Time
import Data.Maybe (listToMaybe)
import Text.XML.Light
type KML = Element
kmlToString :: KML -> String
kmlToString = (++) xml_header . ppElement
kmlTop = Element (QName "kml" Nothing Nothing) [Attr (basicName "xmlns") "http://www.opengis.net/kml/2.2"] [] Nothing
addElem :: Element -> Element -> Element
addElem top child = top { elContent = Elem child : elContent top }
basicName :: String -> QName
basicName n = QName n Nothing Nothing
basicElem :: String -> Element
basicElem n = Element (basicName n) [] [] Nothing
filledElem :: String -> String -> Element
filledElem name val = Element (basicName name) [] [Text (CData CDataText val Nothing)] Nothing
trailToKML :: Coordinate a => String -> Trail a -> KML
trailToKML name ps = addElem kmlTop doc
where
doc = foldl addElem (basicElem "Document") (nameElem : trailElem : [])
nameElem = filledElem "name" name
trailName = filledElem "name" (name ++ " trail")
trailElem = foldl addElem (basicElem "Placemark") [trailName, pointsElem]
pointsElem = addElem (basicElem "LineString") (filledElem "coordinates" points)
points :: String
points = concatMap packCoordinate ps
packCoordinate :: Coordinate a => a -> String
packCoordinate x =
let (lat, lon) = getDegreeLatLon x
in concat [lon, ",", lat, ",0 "]
pointsToKML :: Coordinate a => String -> [a] -> [String] -> KML
pointsToKML name ps pointNames = addElem kmlTop doc
where
doc = foldl addElem (basicElem "Document") (nameElem : placemarkElems)
nameElem = filledElem "name" name
placemarkElems = map buildPoint (zip ps pointNames)
buildPoint (p,n) =
let (t,g) = getDegreeLatLon p
in foldl addElem (basicElem "Placemark")
[ filledElem "name" n
, foldl addElem (basicElem "LookAt")
[ filledElem "longitude" g
, filledElem "latitude" t
, filledElem "altitude" "0"
, filledElem "range" "500"
]
, addElem (basicElem "Point") (filledElem "coordinates" (g ++ "," ++ t))
]
getDegreeLatLon :: Coordinate a => a -> (String,String)
getDegreeLatLon x =
let (lat,lon) = dmsToLatLon . toDMS $ x
f = show . toDecimal
in (f lat, f lon)
radiusOfEarth :: Double
radiusOfEarth = 6378700
class Coordinate a where
toDMS :: a -> DMSCoordinate
fromDMS :: DMSCoordinate -> a
distance :: a -> a -> Distance
distance a b =
3963000 * acos( sin(lat1) * sin(lat2) + cos(lat1) * cos(lat2) * cos(lon2 lon1) )
where
(lat1, lon1) = dmsToLatLon (toDMS a)
(lat2, lon2) = dmsToLatLon (toDMS b)
heading :: a -> a -> Heading
heading a b = dmsHeading (toDMS a) (toDMS b)
getVector :: a -> a -> Vector
getVector a b = (distance a b, heading a b)
class (Coordinate coord) => Location loc coord time | loc -> coord, loc -> time where
getCoordinate :: loc -> coord
getTime :: loc -> time
getUTC :: loc -> UTCTime
speed :: loc -> loc -> Speed
speed a b = distMeter / timeSec
where
timeSec = realToFrac $ diffUTCTime (getUTC b) (getUTC a)
distMeter = distance (getDMS a) (getDMS b)
toDecimal = (*) (180 / pi)
getDMS :: (Location loc x y) => loc -> DMSCoordinate
getDMS = toDMS . getCoordinate
dmsHeading :: DMSCoordinate -> DMSCoordinate -> Heading
dmsHeading a b =
atan2 (sin (diffLon) * cos (lat2))
(cos(lat1) * sin (lat2) sin(lat1) * cos lat2 * cos (diffLon))
where
(lat1, lon1) = dmsToLatLon a
(lat2, lon2) = dmsToLatLon b
diffLon = lon1 lon2
dmsToLatLon :: DMSCoordinate -> (Double, Double)
dmsToLatLon (DMSCoord (DMS latD latM latS) (DMS lonD lonM lonS)) = (lat, lon)
where
lon = toRadians lonD lonM lonS
lat = toRadians latD latM latS
toRadians d m s = (d + m / 60 + s / 3600) * (pi / 180)
normalizeDMS :: DMSCoordinate -> DMSCoordinate
normalizeDMS (DMSCoord lat lon) =
DMSCoord lat' lon'
where
lat' = normalize lat
lon' = normalize lon
normalize (DMS d m s) =
let all = d + (m / 60) + (s / 3600)
d' = fromIntegral (floor all)
m' = (alld') * 60
in DMS d' m' 0
type Distance = Double
type Heading = Double
type Speed = Double
type Vector = (Distance, Heading)
type Trail a = [a]
data DMSCoordinate = DMSCoord
{ latitude :: DMS
, longitude :: DMS
} deriving (Eq, Ord, Show, Read)
data DMS = DMS { degrees :: Double
, minutes :: Double
, seconds :: Double }
deriving (Eq, Ord, Show, Read)
instance Coordinate DMSCoordinate where
toDMS = id
fromDMS = id
instance (Coordinate c) => Location (c, UTCTime) c UTCTime where
getCoordinate = fst
getTime = snd
getUTC = snd
instance (Coordinate c) => Coordinate (c,UTCTime) where
toDMS = toDMS . getCoordinate
fromDMS d = (fromDMS d , UTCTime (toEnum 0) 0)
data TempTrail a = T (Trail a) a
fromTemp :: TempTrail a -> Trail a
fromTemp (T ps p) = reverse (p:ps)
smoothTrails :: Location a b c => Speed -> Trail a -> [Trail a]
smoothTrails s ps =
let (trails, _) = mapAccumL go [] (linearTime ps)
in map fromTemp trails
where
go acc p =
let (i,trails) = mapAccumL (add p) 0 acc
in if i == 0
then (T [] p : trails, ())
else (trails, ())
add p i (T ts l) =
if (speed l p) <= s
then (i + 1, T (l:ts) p)
else (i, T ts l)
longestSmoothTrail :: Location a b c => Speed -> Trail a -> Trail a
longestSmoothTrail metersPerSec trail = maximumBy (comparing length) ([] : smoothTrails metersPerSec trail)
linearTime :: Location a b c => Trail a -> Trail a
linearTime [] = []
linearTime (p:ps) = go (getUTC p) ps
where
go _ [] = []
go t (p:ps) = if getUTC p < t then go t ps else p : go (getUTC p) ps
restLocations :: Location a b c => Distance -> NominalDiffTime -> Trail a -> [Trail a]
restLocations d s = filter timeSpan . span2
where
span2 :: Location a b c => [a] -> [[a]]
span2 [] = []
span2 (a:as) = let (x,y) = span ((>) d . distance (getCoordinate a) . getCoordinate) as in (a : x) : (span2 y)
timeSpan :: Location a b c => Trail a -> Bool
timeSpan [] = False
timeSpan t = s <= diffUTCTime (getUTC (last t)) (getUTC (head t))
closestDistance :: Coordinate a => Trail a -> Trail a -> Maybe Distance
closestDistance as bs = listToMaybe $ sort [distance a b | a <- as, b <- bs]
instance Pretty DMSCoordinate where
pPrint c =
let (DMSCoord (DMS d' m' _) (DMS d m _)) = normalizeDMS c
showF = showFFloat (Just 6)
showX = shows . floor
in text $ showX d' (' ' : showF m' ('\t' : showX d (' ' : showF m "")))
instance Binary DMSCoordinate where
put (DMSCoord (DMS tD tM tS) (DMS gD gM gS)) = do
put tD >> put tM >> put tS >> put gD >> put gM >> put gS
get = do
tD <- get
tM <- get
tS <- get
gD <- get
gM <- get
gS <- get
return (DMSCoord (DMS tD tM tS) (DMS gD gM gS))
putLocation :: Binary c => (c, UTCTime) -> Put
putLocation (c,t) = put c >> put t
getLocation :: Binary c => Get (c, UTCTime)
getLocation = do
c <- get
t <- get
return (c,t)
instance Binary UTCTime where
get = do
day <- get
diff <- get
return (UTCTime (toEnum day) (fromRational diff))
put (UTCTime day diff) = do
put (fromEnum day)
put (toRational diff)