Skip to content

Commit

Permalink
Remove dimensional dependence from tests, and more updates.
Browse files Browse the repository at this point in the history
  • Loading branch information
PaulJohnson committed Dec 27, 2024
1 parent 0278cd3 commit 1d3313e
Show file tree
Hide file tree
Showing 8 changed files with 180 additions and 187 deletions.
39 changes: 28 additions & 11 deletions src/Geodetics/Ellipsoids.hs
Original file line number Diff line number Diff line change
Expand Up @@ -14,28 +14,31 @@ specific area.
-}

module Geodetics.Ellipsoids (
radiansToDegrees,
degreesToRadians,
-- ** Helmert transform between geodetic reference systems
-- * Conversion constants
degree,
arcminute,
arcsecond,
kilometer,
-- * Helmert transform between geodetic reference systems
Helmert (..),
inverseHelmert,
ECEF,
applyHelmert,
-- ** Ellipsoid models of the Geoid
-- * Ellipsoid models of the Geoid
Ellipsoid (..),
WGS84 (..),
LocalEllipsoid (..),
flattening,
minorRadius,
eccentricity2,
eccentricity'2,
-- ** Auxiliary latitudes and related Values
-- * Auxiliary latitudes and related Values
normal,
latitudeRadius,
meridianRadius,
primeVerticalRadius,
isometricLatitude,
-- ** Tiny linear algebra library for 3D vectors
-- * Tiny linear algebra library for 3D vectors
Vec3,
Matrix3,
add3,
Expand All @@ -49,15 +52,29 @@ module Geodetics.Ellipsoids (
) where


-- | All angles in this library are radians.
radiansToDegrees :: Double -> Double
radiansToDegrees = (* (180/pi))
-- | All angles in this library are in radians. This is one degree in radians.
degree :: Double
degree = pi/180

-- | One arc-minute in radians.
arcminute :: Double
arcminute = degree / 60

degreesToRadians :: Double -> Double
degreesToRadians = (* (pi/180))
-- | One arc-second in radians.
arcsecond :: Double
arcsecond = arcminute / 60


-- | All distances in this library are in meters. This is one kilometer in meters.
kilometer :: Double
kilometer = 1000


-- | Small integers, specialised to @Int@, used for raising to powers.
--
-- If you say @x^2@ then Haskell complains that the @2@ has ambiguous type, so you
-- need to say @x^(2::Int)@ to disambiguate it. This gets tedious in complex formulae.

-- | 3d vector as @(X,Y,Z)@.
type Vec3 a = (a,a,a)

Expand Down
18 changes: 9 additions & 9 deletions src/Geodetics/Geodetic.hs
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
module Geodetics.Geodetic (
-- ** Geodetic Coordinates
-- * Geodetic Coordinates
Geodetic (..),
readGroundPosition,
toLocal,
Expand All @@ -10,11 +10,11 @@ module Geodetics.Geodetic (
groundDistance,
properAngle,
showAngle,
-- ** Earth Centred Earth Fixed Coordinates
-- * Earth Centred Earth Fixed Coordinates
ECEF,
geoToEarth,
earthToGeo,
-- ** Re-exported for convenience
-- * Re-exported for convenience
WGS84 (..)
) where

Expand Down Expand Up @@ -89,7 +89,7 @@ readGroundPosition :: (Ellipsoid e) => e -> String -> Maybe (Geodetic e)
readGroundPosition e str =
case map fst $ filter (null . snd) $ readP_to_S latLong str of
[] -> Nothing
(lat,long) : _ -> Just $ groundPosition $ Geodetic (degreesToRadians lat) (degreesToRadians long) undefined e
(lat,long) : _ -> Just $ groundPosition $ Geodetic (lat * degree) (long * degree) undefined e


-- | Show an angle as degrees, minutes and seconds to two decimal places.
Expand All @@ -103,7 +103,7 @@ showAngle a
where
sgn = if a < 0 then "-" else ""
centisecs :: Integer
centisecs = abs $ round $ (radiansToDegrees a) * 360000 -- hundredths of arcsec per degree.
centisecs = abs $ round $ (a * degree * 360000) -- hundredths of arcsec per degree.
(d, m1) = centisecs `divMod` 360000
(m, s1) = m1 `divMod` 6000 -- hundredths of arcsec per arcmin
(s, ds) = s1 `divMod` 100
Expand All @@ -122,8 +122,8 @@ antipode :: (Ellipsoid e) => Geodetic e -> Geodetic e
antipode g = Geodetic lat long (geoAlt g) (ellipsoid g)
where
lat = negate $ latitude g
long' = longitude g - degreesToRadians 180
long | long' < 0 = long' + degreesToRadians 360
long' = longitude g - 180 * degree
long | long' < 0 = long' + 360 * degree
| otherwise = long'


Expand Down Expand Up @@ -230,8 +230,8 @@ groundDistance p1 p2 = do
listToMaybe $ dropWhile converging $ take 100 $ zip lambdas $ drop 1 lambdas
let
uSq = cos2Alpha * (a**2 - b**2) / b**2
bigA = 1 + uSq/16384 * 4096 + uSq * (-768) + uSq * ((320 - (175*uSq)))
bigB = uSq/1024 * 256 + uSq * (-128) + uSq * ((74 - (47* uSq)))
bigA = 1 + uSq/16384 * (4096 + uSq * ((-768) + uSq * ((320 - 175*uSq))))
bigB = uSq/1024 * (256 + uSq * ((-128) + uSq * ((74 - 47* uSq))))
deltaDelta =
bigB * sinDelta * (cos2DeltaM +
bigB/4 * (cosDelta * (2 * cos2DeltaM**2 - 1)
Expand Down
8 changes: 4 additions & 4 deletions src/Geodetics/Grid.hs
Original file line number Diff line number Diff line change
@@ -1,11 +1,11 @@
{-# LANGUAGE FunctionalDependencies #-}

module Geodetics.Grid (
-- ** Grid types
-- * Grid types
GridClass (..),
GridPoint (..),
GridOffset (..),
-- ** Grid operations
-- * Grid operations
polarOffset,
offsetScale,
offsetNegate,
Expand All @@ -14,9 +14,9 @@ module Geodetics.Grid (
offsetDistanceSq,
offsetBearing,
gridOffset,
-- ** Unsafe conversion
-- * Unsafe conversion
unsafeGridCoerce,
-- ** Utility functions for grid references
-- * Utility functions for grid references
fromGridDigits,
toGridDigits
) where
Expand Down
1 change: 0 additions & 1 deletion src/Geodetics/LatLongParser.hs
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,6 @@
-- All angles are returned in degrees.

module Geodetics.LatLongParser (

degreesMinutesSeconds,
degreesMinutesSecondsUnits,
degreesDecimalMinutes,
Expand Down
1 change: 0 additions & 1 deletion src/Geodetics/Stereographic.hs
Original file line number Diff line number Diff line change
Expand Up @@ -115,5 +115,4 @@ instance (Ellipsoid e) => GridClass (GridStereo e) e where
e2 = eccentricity2 $ gridEllipsoid grid
lats = Stream.iterate next lat1
latN = snd $ Stream.head $ Stream.dropWhile (\(v1, v2) -> abs (v1-v2) > 0.01 * arcsecond) $ Stream.zip lats $ Stream.drop 1 lats
arcsecond = degreesToRadians $ 1/3600
gridEllipsoid = ellipsoid . gridTangent
10 changes: 2 additions & 8 deletions src/Geodetics/UK.hs
Original file line number Diff line number Diff line change
Expand Up @@ -36,8 +36,6 @@ instance Ellipsoid OSGB36 where
cX = 446.448, cY = (-125.157), cZ = 542.06,
helmertScale = (-20.4894),
rX = 0.1502 * arcsecond, rY = 0.247 * arcsecond, rZ = 0.8421 * arcsecond }
where
arcsecond = degreesToRadians $ 1/3600

-- | The UK National Grid is a Transverse Mercator projection with a true origin at
-- 49 degrees North, 2 degrees West on OSGB36, and a false origin 400km West and 100 km North of
Expand All @@ -53,8 +51,8 @@ instance GridClass UkNationalGrid OSGB36 where

ukTrueOrigin :: Geodetic OSGB36
ukTrueOrigin = Geodetic {
latitude = degreesToRadians 49,
longitude = degreesToRadians (-2),
latitude = 49 * degree,
longitude = (-2) * degree,
geoAlt = 0,
ellipsoid = OSGB36
}
Expand All @@ -63,10 +61,6 @@ ukFalseOrigin :: GridOffset
ukFalseOrigin = GridOffset ((-400) * kilometer) (100 * kilometer) (0 * kilometer)


kilometer :: Double
kilometer = 1000


-- | Numerical definition of the UK national grid.
ukGrid :: GridTM OSGB36
ukGrid = mkGridTM ukTrueOrigin ukFalseOrigin
Expand Down
Loading

0 comments on commit 1d3313e

Please sign in to comment.