Map transformations¶
Most of the georeferencing machinery for gridded datasets is
handled by the Grid
class: its capacity to handle
gridded datasets in a painless manner was one of the primary
motivations to develop Salem.
Grids¶
A point on earth can be defined unambiguously in two ways:
- DATUM (lon, lat, datum)
- longitudes and latitudes are angular coordinates of a point on an ellipsoid (often called “datum”)
- PROJ (x, y, projection)
- x (eastings) and y (northings) are cartesian coordinates of a point in a map projection (the unit of x, y is usually meter)
Salem adds a third coordinate reference system (crs) to this list:
- GRID (i, j, Grid)
- on a structured grid, the (x, y) coordinates are distant of a constant (dx, dy) step. The (x, y) coordinates are therefore equivalent to a new reference frame (i, j) proportional to the projection’s (x, y) frame.
Transformations between datums and projections is handled by several tools in the python ecosystem, for example GDAL or the more lightweight pyproj, which is the tool used by Salem internally [1].
The concept of Grid added by Salem is useful when transforming data between two structured datasets, or from an unstructured dataset to a structured one.
[1] | Most datasets nowadays are defined in the WGS 84 datum, therefore the concepts of datum and projection are often interchangeable: (lon, lat) coordinates are equivalent to cartesian (x, y) coordinates in the plate carree projection. |
A Grid
is defined by a projection, a reference point in
this projection, a grid spacing and a number of grid points:
In [1]: import numpy as np
In [2]: import salem
In [3]: from salem import wgs84
In [4]: grid = salem.Grid(nxny=(3, 2), dxdy=(1, 1), x0y0=(0.5, 0.5), proj=wgs84)
In [5]: x, y = grid.xy_coordinates
In [6]: x
Out[6]:
array([[0.5, 1.5, 2.5],
[0.5, 1.5, 2.5]])
In [7]: y