Skip to content

Instantly share code, notes, and snippets.

@josip
Created September 30, 2012 23:19
Show Gist options
  • Select an option

  • Save josip/3808709 to your computer and use it in GitHub Desktop.

Select an option

Save josip/3808709 to your computer and use it in GitHub Desktop.
Renders a map of GPS-coordinates from a file (LAT,LNG\n)
import Data.List
import Data.List.Split
import System.Environment
import System.Directory
import System.FilePath
import Control.Monad
import Graphics.Rendering.Cairo as C
type Coord = (Double, Double)
type Pixel = (Double, Double)
type Extent = (Double, Double, Double, Double)
type Dimensions = (Double, Double)
type RGB = (Double, Double, Double)
earthCirc = 40075.04 -- [KM]
earthRadius = 6378.14 -- [KM]
main = do
paths <- getDirectoryContents "txt"
--mapM drawPath (filter (\x -> length x > 3) paths)
drawPath "txt/hills.txt"
putStrLn "Fin"
drawPath p = renderMap (takeBaseName p) 600.0 600.0
crgb (r,g,b) = (r/255, g/255, b/255)
renderMap filename w h = do
putStrLn $ filename ++ "..."
coordsFile <- readFile $ "txt/" ++ filename ++ ".txt"
let (coords, alts) = unzip $ map parseLine $ lines coordsFile
centre = mapCentre coords
rStart = calculateRadius centre coords
(r, extent, pixels) = fitMap centre (w, h) coords rStart
centrePx = coordToPixel extent (w, h) centre
C.withImageSurface C.FormatARGB32 (round w) (round h) $ \surface -> do
C.renderWith surface $ drawRegularMap w h pixels
C.surfaceWriteToPNG surface $ "png/basic_" ++ filename ++ ".png"
C.renderWith surface $ drawAltMap w h pixels alts
C.surfaceWriteToPNG surface $ "png/alt_" ++ filename ++ ".png"
return ()
parseLine :: String -> (Coord, Double)
parseLine xs = (coord, alt)
where (lat:lng:alts:[]) = splitOn "," xs
coord = (read lat :: Double, read lng :: Double)
alt = read alts :: Double
drawRegularMap w h pixels = do
drawBackground w h
C.setSourceRGB (15/255) (18/255) (38/255) -- #732412
C.setLineWidth 1
C.moveTo (fst . head $ pixels) (snd . head $ pixels)
mapM (\(x, y) -> C.lineTo x y) (tail pixels)
C.stroke
showPx 3 (head pixels) (217, 59, 59)
showPx 5 (last pixels) (217, 59, 59)
drawAltMap w h pixels alts = do
let maxAlt = maximum alts
maxAltI = case elemIndex maxAlt alts of
(Just x) -> x
otherwise -> 0
minAlt = minimum alts
minAltI = case elemIndex minAlt alts of
(Just x) -> x
otherwise -> 0
nAlts = map (\x -> (x - minAlt)/maxAlt) alts
avgAlt = (sum nAlts) / (fromIntegral (length nAlts) :: Double)
liftIO $ print (minimum nAlts, avgAlt, maximum nAlts)
drawBackground w h
C.setLineWidth 1
C.moveTo (fst . head $ pixels) (snd . head $ pixels)
mapM (drawPointWithAlt avgAlt) $ zip (tail pixels) nAlts
showPx 3 (pixels !! minAltI) (0, 0, 255)
showPx 3 (pixels !! maxAltI) (255, 0, 0)
drawPointWithAlt avg ((x, y), alt) = do
let d = abs $ avg - alt
--(r, g, b) = gradient (0, 0, 255) (255, 0, 0) alt
setColor r g b = C.setSourceRGB (r/255) (g/255) (b/255)
case alt of
x | x >= 0.9 -> setColor 255 0 0
x | x >= 0.5 -> setColor 0 255 0
otherwise -> setColor 0 0 255
--setColor 0 0 0 $ (0.5 + d)
C.lineTo x y
C.stroke
C.moveTo x y
splitLine :: String -> Coord
splitLine xs = (read lat :: Double, read (drop 1 lng) :: Double)
where (lat, lng) = break (== ',') xs
{-# INLINE validPixel #-}
validPixel :: Dimensions -> Pixel -> Bool
validPixel (w, h) (x, y) = (x >= 0 && x <= w) && (y >= 0 && y <= h)
showPx rad (x, y) (r, g, b) = do
C.save
C.setSourceRGB (r/255) (g/255) (b/255)
C.arc x y rad 0 (2 * pi)
C.fill
C.restore
drawBackground w h = do
C.save
C.setSourceRGB (235/255) (230/255) (220/255) -- #EBE6DC
C.rectangle 0 0 w h
C.fill
drawGrid 10 w h
C.restore
drawTriangle d (x, y) (r, g, b) = do
C.save
C.moveTo (x - d) (y + d)
C.lineTo x (y - d)
C.lineTo (x + d) (y + d)
C.lineTo (x - d) (y + d)
C.setSourceRGB (r/255) (g/255) (b/255)
C.closePath
C.fill
C.restore
drawRect d (x, y) (r, g, b) = do
C.save
C.setSourceRGB r g b
C.rectangle x y (d - 0.5) (d - 0.5)
C.fill
C.restore
drawGrid d w h = do
C.save
C.setSourceRGBA (69/255) (192/255) (146/255) 0.2 -- #45C092
C.setLineWidth 1
mapM (drawHorizontalLine d w) [1..w/d]
mapM (drawVerticalLine d h) [1..h/d]
C.restore
drawHorizontalLine d w n = do
C.moveTo 0 (n*d - 0.5)
C.lineTo w (n*d - 0.5)
C.stroke
drawVerticalLine d h n = do
C.moveTo (n*d - 0.5) 0
C.lineTo (n*d - 0.5) h
C.stroke
gradient :: RGB -> RGB -> Double -> RGB
gradient (r1,g1,b1) (r2,g2,b2) n = (c r1 r2, c g1 g2, c b1 b2)
where c x y = (x + n * (y - x) / 255)
{-
gradient :: RGB -> RGB -> Double -> RGB
gradient (r1, g1, b1) (r2, g2, b2) step = (r/255, g/255, b/255)
where combos = [(r, g, b) | r <- [r1..r2], g <- [g1..g2], b <- [b1..b2]]
len = fromIntegral $ length combos
(r,g,b) = combos !! (round $ len * step)
-}
fitMap :: Coord -> Dimensions -> [Coord] -> Double -> (Double, Extent, [Pixel])
fitMap centre dim cs r = if allFit then (r, ex, pxs)
else fitMap centre dim cs (r + r/2)
where ex = calculateExtent r centre
pxs = map (coordToPixel ex dim) cs
allFit = all (validPixel dim) pxs
mapCentre :: [Coord] -> Coord
mapCentre coords = (lats/n, lngs/n)
where n = fromIntegral (length coords) :: Double
sumCoords (x, y) (xs, ys) = (x + xs, y + ys)
(lats, lngs) = foldr1 sumCoords coords
calculateExtent :: Double -> Coord -> Extent
calculateExtent r (lat, lng) = (lat + sLat, lng + sLng, lat - sLat, lng - sLng)
where sLat = (r / (earthCirc / pi)) / d2r
sLng = sLat / (cos $ lat * d2r)
d2r = 0.0174532925
calculateRadius :: Coord -> [Coord] -> Double
calculateRadius centre cs = (ds !! maxD) / 80
where ds = map (coordDistance centre) cs
maxD = case elemIndex (maximum ds) ds of
(Just d) -> d
otherwise -> 0
coordToPixel :: Extent -> Dimensions -> Coord -> Pixel
coordToPixel (n, e, s, w) (width, height) (lat, lng) = (x, y)
where dLat = n - s
dLng = e - w
x = width * (lng - w) / dLng
y = height * (1 - (lat - s) / dLat)
-- Haversine formula
coordDistance :: Coord -> Coord -> Double
coordDistance (lat1, lng1) (lat2, lng2) = earthRadius * c
where dLat = (lat2 - lat1) / 2
dLng = (lng2 - lng1) / 2
a = (sin dLat)^2 + (cos lat1) * (cos lat2) * ((sin dLng)^2)
c = 2.0 * (atan2 (sqrt a) (sqrt (1.0 - a)))
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment