Skip to content

Unitful coordinate reference systems for geographic maps in Julia

License

Notifications You must be signed in to change notification settings

JuliaEarth/CoordRefSystems.jl

Repository files navigation

Build Status Coverage

CoordRefSystems.jl provides conversions between Coordinate Reference Systems (CRS) in native Julia. It was designed to work with units fromUnitful.jl, respects projection bounds catalogued inhttps://epsg.io,and is very fast thanks to advanced parametrizations at compile-time.

This package addresses various design issues encountered in previous attempts such as Geodesy.jlandMapMaths.jl. Ourbenchmarksshow that CoordRefSystems.jl is sometimes faster thanPROJ, which is the most widely used software library for cartography in the world (written in C/C++).

Installation

Get the latest stable release with Julia's package manager:

] add CoordRefSystems

Usage

Basic usage

Consider the following conversions betweenCartesian,Spherical, CylindricalandPolarcoordinates to get started:

Cartesian <> Spherical

julia>cartesian=Cartesian(1,1,1)
Cartesian{NoDatum} coordinates
├─ x:1.0m
├─ y:1.0m
└─ z:1.0m

julia>spherical=convert(Spherical, cartesian)
Spherical{NoDatum} coordinates
├─ r:1.7320508075688772m
├─ θ:0.9553166181245093rad
└─ ϕ:0.7853981633974483rad

julia>convert(Cartesian, spherical)
Cartesian{NoDatum} coordinates
├─ x:1.0m
├─ y:0.9999999999999998m
└─ z:0.9999999999999999m

Cartesian <> Cylindrical

julia>cartesian=Cartesian(1,1,1)
Cartesian{NoDatum} coordinates
├─ x:1.0m
├─ y:1.0m
└─ z:1.0m

julia>cylindrical=convert(Cylindrical, cartesian)
Cylindrical{NoDatum} coordinates
├─ ρ:1.4142135623730951m
├─ ϕ:0.7853981633974483rad
└─ z:1.0m

julia>convert(Cartesian, cylindrical)
Cartesian{NoDatum} coordinates
├─ x:1.0000000000000002m
├─ y:1.0m
└─ z:1.0m

Cartesian <> Polar

julia>cartesian=Cartesian(1,1)
Cartesian{NoDatum} coordinates
├─ x:1.0m
└─ y:1.0m

julia>polar=convert(Polar, cartesian)
Polar{NoDatum} coordinates
├─ ρ:1.4142135623730951m
└─ ϕ:0.7853981633974483rad

julia>convert(Cartesian, polar)
Cartesian{NoDatum} coordinates
├─ x:1.0000000000000002m
└─ y:1.0m

Special syntax

Julia'sconvertmethods can be triggered with special syntax assuming that a list of coordinates is available:

julia>Mercator[LatLon(0,0),LatLon(30,30),LatLon(20,30)]
3-element Vector{Mercator}:
Mercator{WGS84Latest}(x:0.0m, y:0.0m)
Mercator{WGS84Latest}(x:3.33958e6m, y:3.48219e6m)
Mercator{WGS84Latest}(x:3.33958e6m, y:2.25842e6m)

The example above is equivalent to running convert(Mercator, latlon)for alllatlon coordinates in the list.

Advanced usage

CRS are most useful to locate objets in the physical world. Given an ellipsoid of revolution and a standardized origin (a.k.a. datum), it is possible assign coordinates to points without ambiguity.

Below is an example converting geodeticLatLoncoordinates on theWGS84Latestdatum toMercator,WebMercator,and Robinsonprojected coordinates on the same datum:

julia>latlon=LatLon(30,60)
GeodeticLatLon{WGS84Latest} coordinates
├─ lat:30.0°
└─ lon:60.0°

julia>convert(Mercator, latlon)
Mercator{WGS84Latest} coordinates
├─ x:6.679169447596414e6m
└─ y:3.482189085408618e6m

julia>convert(WebMercator, latlon)
WebMercator{WGS84Latest} coordinates
├─ x:6.679169447596414e6m
└─ y:3.5035498435043753e6m

julia>convert(Robinson, latlon)
Robinson{WGS84Latest} coordinates
├─ x:5.441866544132874e6m
└─ y:3.2085576115038935e6m

julia>latlon=LatLon(30,60)
GeodeticLatLon{WGS84Latest} coordinates
├─ lat:30.0°
└─ lon:60.0°

julia>mercator=convert(Mercator, latlon)
Mercator{WGS84Latest} coordinates
├─ x:6.679169447596414e6m
└─ y:3.482189085408618e6m

julia>convert(LatLon, mercator)
GeodeticLatLon{WGS84Latest} coordinates
├─ lat:29.999999999999996°
└─ lon:59.99999999999999°

It is also possible to convert between different datum, transparently. In the following examples, we convert coordinates between theWGS84Latest datum, currently an alias toWGS84{1762},and theITRF{2008}datum:

julia>latlon=LatLon{WGS84Latest}(30,60)
GeodeticLatLon{WGS84Latest} coordinates
├─ lat:30.0°
└─ lon:60.0°

julia>convert(LatLon{ITRF{2008}}, latlon)
GeodeticLatLon{ITRF{2008}} coordinates
├─ lat:30.00000000081754°
└─ lon:59.99999999999999°

julia>mercator=convert(Mercator{WGS84Latest}, latlon)
Mercator{WGS84Latest} coordinates
├─ x:6.679169447596414e6m
└─ y:3.482189085408618e6m

julia>convert(WebMercator{WGS84Latest}, mercator)
WebMercator{WGS84Latest} coordinates
├─ x:6.679169447596414e6m
└─ y:3.5035498435043753e6m

julia>convert(WebMercator{ITRF{2008}}, mercator)
WebMercator{ITRF{2008}} coordinates
├─ x:6.679169447596414e6m
└─ y:3.5035498436094625e6m

EPSG/ERSI codes

CRS are catalogued with numerical codes inhttps://epsg.io. The package providesEPSG{code}andERSI{code},and the utilityCoordRefSystems.getfunction to query the database:

julia>CRS1=CoordRefSystems.get(EPSG{4326})
GeodeticLatLon{WGS84Latest}

julia>CRS2=CoordRefSystems.get(EPSG{3395})
Mercator{WGS84Latest}

julia>CRS1(0,90)
GeodeticLatLon{WGS84Latest} coordinates
├─ lat:0.0°
└─ lon:90.0°

julia>convert(CRS2,CRS1(0,90))
Mercator{WGS84Latest} coordinates
├─ x:1.0018754171394622e7m
└─ y:0.0m

Credits

Most implementations in this package are adaptations from PROJ - Cartographic Projections and Coordinate Transformations Library and itslist of references. Our tests were designed to match their results to the last digit via theProj.jlwrapper.