transform_factory#

pycraf.geospatial.transform_factory(sys_in, sys_out, code_in='epsg', code_out='epsg')[source]#

A factory to produce conversion functions (decorated with ranged_quantity_input for unit handling).

The returned function will automatically have the correct signature (e.g., lon/lat –> x/y, or x/y/z –> lon/lat/height) depending on the chosen input and output projections. Also a proper docstring is produced.

While geospatial comes with several convenience transform functions, the user could also simply built all of the with the transform_factory. For example, the conversion from WGS84 to ETRS89 is simply:

>>> from pycraf import geospatial
>>> import astropy.units as u

>>> wgs84_to_etrs89 = geospatial.transform_factory(
...     geospatial.EPSG.WGS84, geospatial.EPSG.ETRS89
...     )

The transform_factory can be used to work with more uncommon projections. As an example, consider the Gauss-Kruger projection, which is defined for various zones (like UTM) but also with different zone widths. This would make it completely unhandy, to wrap all of the pyproj transforms in pycraf only for adding unit-checking.

A lot (thousands!) of projections are pre-defined and have a so-called EPSG code, e.g. WGS84="EPSG:4326". The website spatialreference.org is very useful to find information on EPSG codes.

Parameters:

sys_in, sys_outstr, int, EPSG, ESRI: Input and output projection for the desired coordinate transformation. If a string, this is directly fed into the Proj class. If an integer, it is treated as an EPSG or ESRI code and the Proj class is instantiated with "epsg:{code}" or "esri:{code}", depending on the code_[in,out] value. For convenience, a couple of common systems are defined in the EPSG and ESRI Enums (which just hold the associated EPSG or ESRI code for each name), e.g., EPSG.WGS84 == 4326.
code_in, code_out{‘epsg’, ‘esri’}: Whether to interpret integer-valued sys_[in,out] arguments as EPSG or ESRI code. (default: ‘epsg’)

Returns:

transform_funcFunction: Transform function with unit decorator (ranged_quantity_input).

Notes

pyproj also allows to preserve the units (preserve_units=True); some EPSG systems use different distance units, such as ft. Per default, pyproj will convert these to meters. The ranged_quantity_input will automatically convert your inputs to the correct unit (i.e., ‘m’ for distances), but the outputs will always be in meters (or degrees). Use units functions to convert to something else such as ‘ft’ if desired.

Examples

To define a “Gauss-Kruger Zone 4 (Germany)” to GPS/WGS84 transform:

>>> from pycraf import geospatial
>>> import astropy.units as u

>>> # see http://spatialreference.org/ref/epsg/31467/
>>> wgs84_to_etrs89 = geospatial.transform_factory(
...     geospatial.EPSG.WGS84, 31467
...     )
>>> wgs84_to_etrs89(6 * u.deg, 50 * u.deg)  
(<Quantity 3285005.658 m>, <Quantity 5544721.322 m>)

If an ESRI system is desired, either use the integer code and set code_[in,out] = 'esri', or use the ESRI enum:

>>> from pycraf import geospatial
>>> import astropy.units as u

>>> # see http://spatialreference.org/ref/esri/54009/
>>> wgs84_to_mollweide = geospatial.transform_factory(
...     geospatial.EPSG.WGS84, geospatial.ESRI.MOLLWEIDE
...     )
>>> wgs84_to_mollweide(6 * u.deg, 50 * u.deg)  
(<Quantity 456379.912 m>, <Quantity 5873471.956 m>)

>>> mollweide_to_wgs84 = geospatial.transform_factory(
...     geospatial.ESRI.MOLLWEIDE, geospatial.EPSG.WGS84
...     )
>>> mollweide_to_wgs84(456379.912 * u.m, 5873471.956 * u.m)  
(<Quantity 6. deg>, <Quantity 50. deg>)