Coordinate Systems#

SpatialData provides two types of coordinate systems: Cartesian (CSCart) and georeferenced (CSGeo). The coordinate systems can be specified

CSCart#

CSCart provides 1D, 2D, and 3D Cartesian coordinate systems.

Pyre User Interface

See CSCart component.

Examples#

We specify a 2D coordinate system with the coordinates given in meters.

Listing 1 Specifying a CSCart in a spatial database file.#
cs-data = cartesian {
    to-meters = 1000.0
    space-dim = 2
}
Listing 2 Create a CSCart object in C++.#
#include "spatialdata/geocoords/CSCart.hh"

spatialdata::geocoords:CSCart cs;
cs.setToMeters(1000.0);
cs.setSpaceDim(2);
Listing 3 Create a CSCart object in Python.#
from spatialdata.geocoords.CSCart import CSCart

cs = CSCart()
cs.units = "km"
cs.spaceDim = 2
cs._configure()

CSGeo#

CSGeo provides 2D and 3D georeferenced coordinate systems. We use Proj to perform coordinate system operations, so any coordinate system supported by Proj can be used. The coordinate systems can be specified using EPSG codes, Well-Known Text (WKT), or Proj parameters. See the Proj documentation for more information.

Table 3 Frequently used EPSG codes.#

EPSG Code

Description

EPSG:4326

WGS84 (latitude, longitude)

EPSG:4269

NAD83 (latitude, longitude)

EPSG:4267

NAD27 (latitude, longitude)

EPSG:32610[1]

UTM Zone 10N WGS84 datum (easting, northing)

EPSG:26910[1]

UTM Zone 10N NAD83 datum (easting, northing)

EPSG:26710[1]

UTM Zone 10N NAD27 datum (easting, northing)

Important

Georeferencing has become much more accurate since the availability of the Global Positional System and its successors. As a result, geographic coordinate systems or Coordinate Reference Systems (CRS) have become more complicated. Traditional Proj parameters (used exclusively in Proj versions 5.0 and earlier) do not provide complete specification of a CRS; specifying CRS using EPSG codes or WKT are preferable to using Proj parameters.

Pyre User Interface

See CSGeo component.

Examples#

We create a 3D georeferenced coordinate system corresponding to UTM Zone 10N with the WGS84 datum (EPSG:32610).

Listing 4 Specifying a CSGeo in a spatial database file.#
cs-data = geographic {
    crs-string = EPSG:32610
    space-dim = 3
}
Listing 5 Create a CSGeo object in C++.#
#include "spatialdata/geocoords/CSGeo.hh"

spatialdata::geocoords:CSGeo cs;
cs.setString("EPSG:32610");
cs.setSpaceDim(3);
Listing 6 Create a CSGeo object in Python.#
from spatialdata.geocoords.CSGeo import CSGeo

cs = CSGeo()
cs.crsString = "EPSG:32610"
cs.spaceDim = 3
cs._configure()

CSGeoLocal#

Note

New in v3.1.0.

CSGeoLocal provides the flexibility of adding a local origin and rotation to a georeferenced coordinate system. As with CSGeo, we use Proj to perform georeferenced coordinate system operations, so any geographic coordinate system supported by Proj can be used. See the Proj documentation for more information.

Fig. 1 Diagram of the relationship among the the geographic coordinate system (\(xy\)), the local unrotated coordinate system (\(x'y'\)), and the local rotated coordinate system (\(x''y''\)). CSGeoLocal corresponds to the local rotated coordinate system with the geographic coordinate system specified by the Proj string.#

Transformation from geographic coordinates to local coordinates:

(2)#\[\begin{align} x' &= x - x_0 \\ y' &= y - y_0 \\ x'' &= \cos\theta \, x' - \sin\theta \, y' \\ y'' &= sin\theta \, x' + cos\theta \, y' \\ \end{align}\]

Transformation from local coordinates to geographic coordinates:

(3)#\[\begin{align} x' &= \cos\theta \, x'' + \sin\theta \, y'' \\ y' &= -sin\theta \, x'' + cos\theta \, y'' \\ x &= x' + x_0 \\ y &= y' + y_0 \\ \end{align}\]

Pyre User Interface

See CSGeoLocal component.

Examples#

We create a 2D georeferenced coordinate system in UTM Zone 11N with the WGS84 datum (EPSG:32611) that has a local origin at 500000 easting and 3750000 northing and rotated 30 degrees clockwise from north.

Listing 7 Specifying a CSGeoLocal in a spatial database file.#
cs-data = local-geographic {
    crs-string = EPSG:32611
    space-dim = 2
    origin-x = 500000.0
    origin-y = 3750000.0
    y-azimuth = 30.0
}
Listing 8 Create a CSGeoLocal object in C++.#
#include "spatialdata/geocoords/CSGeo.hh"

spatialdata::geocoords:CSGeo cs;
cs.setString("EPSG:32611");
cs.setSpaceDim(2);
cs.setLocal(500000.0, 3750000.0, 30.0);
Listing 9 Create a CSGeoLocal object in Python.#
from spatialdata.geocoords.CSGeoLocal import CSGeoLocal

cs = CSGeoLocal()
cs.crsString = "EPSG:32611"
cs.spaceDim = 2
cs.originX = 500000.0
cs.originY = 3750000.0
cs.yAzimuth = 30.0
cs._configure()

Converting between coordinate systems#

We convert from 34.5N 118.0W in the WGS84 datum (EPSG:4326) to UTM Zone 10N in the NAD27 datum.

Examples#

Listing 10 Using the Proj cs2cs command line tool to convert between coordinate systems.#
$ echo "37.5 -121.0" | cs2cs EPSG:4326 EPSG:26710 -f "%.1f"
676885.9        4152025.1 0.0

The points are converted in place (array values are changed).

#include "spatialdata/geocoords/CSGeo.hh"
#include "spatialdata/geocoords/Converter.hh"

const size_t numPoints = 1;
const size_t spaceDim = 2;
const double points[numPoints*spaceDim] = { 37.5, -121.0 };

spatialdata::geocoords::CSGeo csSrc;
csSrc.setString("EPSG:4326");
csSrc.setSpaceDim(spaceDim);

spatialdata::geocoords::CSGeo csDest;
csSrc.setString("EPSG:26710");
csSrc.setSpaceDim(spaceDim);

spatialdata::geocoords::Converter converter;
converter.convert(points, numPoints, spaceDim, &csDest, &csSrc);

The points are converted in place (array values are changed).

import numpy

from spatialdata.geocoords.CSGeo import CSGeo
from spatialdata.geocoords.Converter import convert

csSrc = CSGeo()
csSrc.crsString = "EPSG:4326"
csSrc.spaceDim = 2
csSrc._configure()

csDest = CSGeo()
csDest.crsString = "EPSG:26710"
csDest.spaceDim = 2
csDest._configure()

points = numpy.array([[37.5, -121.0]], dtype=numpy.float64)
convert(points, csDest, csSrc)