Skip to content

CRS & Reprojection Helpers#

CRS-handling helpers in pyramids.base.crs — the single source of truth for osr.SpatialReference construction, WKT/Proj4 → EPSG resolution, and coordinate reprojection. The most-commonly-used ones are also re-exposed as FeatureCollection static methods for ergonomic continuity (e.g. FeatureCollection.reproject_coordinates delegates to pyramids.base.crs.reproject_coordinates).

Functions#

pyramids.base.crs.sr_from_epsg(epsg) #

Build an :class:osr.SpatialReference from an EPSG code.

Parameters:

Name Type Description Default
epsg int

EPSG code; cast to int before being handed to :meth:osr.SpatialReference.ImportFromEPSG.

 required

Returns:

Type Description
SpatialReference

osr.SpatialReference: The constructed SRS.

Raises:

Type Description
ValueError

If GDAL cannot resolve the EPSG code (the non-zero return path from ImportFromEPSG — usually propagates as a GDAL exception when gdal.UseExceptions() is active, which pyramids installs at package import).
Source code in src/pyramids/base/crs.py 
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
def sr_from_epsg(epsg: int) -> osr.SpatialReference:
    """Build an :class:`osr.SpatialReference` from an EPSG code.

    Args:
        epsg: EPSG code; cast to `int` before being handed to
            :meth:`osr.SpatialReference.ImportFromEPSG`.

    Returns:
        osr.SpatialReference: The constructed SRS.

    Raises:
        ValueError: If GDAL cannot resolve the EPSG code (the
            non-zero return path from `ImportFromEPSG` — usually
            propagates as a GDAL exception when
            `gdal.UseExceptions()` is active, which pyramids
            installs at package import).
    """
    sr = osr.SpatialReference()
    err = sr.ImportFromEPSG(int(epsg))
    if err != 0:
        raise ValueError(
            f"Failed to create SRS from EPSG:{epsg} (osr returned error {err})."
        )
    return sr

pyramids.base.crs.sr_from_wkt(wkt) #

Build an :class:osr.SpatialReference from a WKT string.

Thin wrapper around osr.SpatialReference(wkt=wkt) that gives the WKT path a consistent name alongside :func:sr_from_epsg and :func:create_sr_from_proj. Use this when you have a WKT (the most common case in the dataset stack — dataset.crs returns WKT) and want a typed SRS without re-typing the constructor's keyword argument every call site.

Parameters:

Name Type Description Default
wkt str

Well-Known Text representation of the spatial reference.

 required

Returns:

Type Description
SpatialReference

osr.SpatialReference: The constructed SRS.

Examples:

  • Round-trip an EPSG code through WKT: 
    >>> from osgeo import osr
>>> from pyramids.base.crs import sr_from_epsg, sr_from_wkt
>>> wkt = sr_from_epsg(4326).ExportToWkt()
>>> sr = sr_from_wkt(wkt)
>>> sr.IsGeographic()
1
Source code in src/pyramids/base/crs.py 
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
def sr_from_wkt(wkt: str) -> osr.SpatialReference:
    """Build an :class:`osr.SpatialReference` from a WKT string.

    Thin wrapper around `osr.SpatialReference(wkt=wkt)` that gives
    the WKT path a consistent name alongside :func:`sr_from_epsg` and
    :func:`create_sr_from_proj`. Use this when you have a WKT (the
    most common case in the dataset stack — `dataset.crs` returns
    WKT) and want a typed SRS without re-typing the constructor's
    keyword argument every call site.

    Args:
        wkt: Well-Known Text representation of the spatial reference.

    Returns:
        osr.SpatialReference: The constructed SRS.

    Examples:
        - Round-trip an EPSG code through WKT:
            ```python
            >>> from osgeo import osr
            >>> from pyramids.base.crs import sr_from_epsg, sr_from_wkt
            >>> wkt = sr_from_epsg(4326).ExportToWkt()
            >>> sr = sr_from_wkt(wkt)
            >>> sr.IsGeographic()
            1

            ```
    """
    return osr.SpatialReference(wkt=wkt)

pyramids.base.crs.create_sr_from_proj(prj, string_type=None) #

Create an :class:osr.SpatialReference from a projection string.

Parameters:

Name Type Description Default
prj str

The projection string (WKT, ESRI WKT, or Proj4).

 required
string_type str | None

One of "WKT", "ESRI wkt", "PROj4", or None for auto-detect (default). Auto-detect uses WKT import and falls back to ESRI WKT or Proj4 based on the prefix.

 None

Returns:

Type Description
SpatialReference

osr.SpatialReference: The constructed spatial reference.

Examples:

  • Parse a standard EPSG:4326 WKT string and inspect the result: 
    >>> from osgeo import osr
>>> ref = osr.SpatialReference()
>>> _ = ref.ImportFromEPSG(4326)
>>> wkt = ref.ExportToWkt()
>>> srs = create_sr_from_proj(wkt)
>>> srs.IsGeographic()
1
>>> srs.GetName()
'WGS 84'
  • Parse a Proj4 string by passing string_type="PROJ4": 
    >>> srs = create_sr_from_proj(
...     "+proj=longlat +datum=WGS84 +no_defs", string_type="PROJ4"
... )
>>> srs.IsGeographic()
1
>>> srs.IsProjected()
0
  • Parse an EPSG:3857 WKT and confirm the axis order is projected: 
    >>> from osgeo import osr
>>> ref = osr.SpatialReference()
>>> _ = ref.ImportFromEPSG(3857)
>>> srs = create_sr_from_proj(ref.ExportToWkt())
>>> srs.IsProjected()
1
>>> srs.GetName()
'WGS 84 / Pseudo-Mercator'
Source code in src/pyramids/base/crs.py 
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
def create_sr_from_proj(
    prj: str, string_type: str | None = None
) -> osr.SpatialReference:
    """Create an :class:`osr.SpatialReference` from a projection string.

    Args:
        prj (str):
            The projection string (WKT, ESRI WKT, or Proj4).
        string_type (str | None):
            One of `"WKT"`, `"ESRI wkt"`, `"PROj4"`, or `None`
            for auto-detect (default). Auto-detect uses WKT import and
            falls back to ESRI WKT or Proj4 based on the prefix.

    Returns:
        osr.SpatialReference: The constructed spatial reference.

    Examples:
        - Parse a standard EPSG:4326 WKT string and inspect the result:
            ```python
            >>> from osgeo import osr
            >>> ref = osr.SpatialReference()
            >>> _ = ref.ImportFromEPSG(4326)
            >>> wkt = ref.ExportToWkt()
            >>> srs = create_sr_from_proj(wkt)
            >>> srs.IsGeographic()
            1
            >>> srs.GetName()
            'WGS 84'

            ```
        - Parse a Proj4 string by passing `string_type="PROJ4"`:
            ```python
            >>> srs = create_sr_from_proj(
            ...     "+proj=longlat +datum=WGS84 +no_defs", string_type="PROJ4"
            ... )
            >>> srs.IsGeographic()
            1
            >>> srs.IsProjected()
            0

            ```
        - Parse an EPSG:3857 WKT and confirm the axis order is projected:
            ```python
            >>> from osgeo import osr
            >>> ref = osr.SpatialReference()
            >>> _ = ref.ImportFromEPSG(3857)
            >>> srs = create_sr_from_proj(ref.ExportToWkt())
            >>> srs.IsProjected()
            1
            >>> srs.GetName()
            'WGS 84 / Pseudo-Mercator'

            ```
    """
    srs = osr.SpatialReference()
    if string_type is None:
        srs.ImportFromWkt(prj)
    elif prj.startswith("PROJCS") or prj.startswith("GEOGCS"):
        srs.ImportFromESRI([prj])
    else:
        srs.ImportFromProj4(prj)
    return srs

pyramids.base.crs.get_epsg_from_prj(prj) #

Return the EPSG code identified by a projection string.

Resolves the EPSG of the root CRS object in three steps: :meth:osr.SpatialReference.AutoIdentifyEPSG (tags recognisable CRSes), then the root AUTHORITY code, then a confident, unambiguous :meth:osr.SpatialReference.FindMatches PROJ-database lookup for well-known CRSes whose WKT lacks a root authority (e.g. a UTM PROJCS from GDAL's AAIGrid driver). The code of a child unit/datum node is never returned as if it were a CRS — that bug (issue #403) made GRIB rasters resolve to the degree-unit EPSG:9122 and UTM ASCII grids to the WGS_1984 datum EPSG:6326.

An empty input string is no longer silently mapped to 4326; that legacy default masked real configuration errors. Callers that genuinely want a fallback should handle the CRSError themselves, or use :func:epsg_from_wkt which accepts an explicit default.

Parameters:

Name Type Description Default
prj str

Projection string.

 required

Returns:

Name Type Description
int int

The resolved EPSG code.

Raises:

Type Description
CRSError

If prj is an empty string, or if its root CRS carries no EPSG authority and matches no PROJ-database entry (e.g. a custom spherical-earth GRIB GEOGCS). The unit/datum codes of child nodes are never returned as a CRS.

Examples:

  • Resolve EPSG:4326 from its standard WKT representation: 
    >>> from osgeo import osr
>>> ref = osr.SpatialReference()
>>> _ = ref.ImportFromEPSG(4326)
>>> get_epsg_from_prj(ref.ExportToWkt())
4326
  • Resolve EPSG:3857 (Web Mercator) from its WKT representation: 
    >>> from osgeo import osr
>>> ref = osr.SpatialReference()
>>> _ = ref.ImportFromEPSG(3857)
>>> get_epsg_from_prj(ref.ExportToWkt())
3857
  • A well-known CRS whose WKT carries no root authority still resolves, via a confident PROJ-database match (here a "WGS 84 / UTM zone 18N" PROJCS with its root authority stripped): 
    >>> from osgeo import osr
>>> ref = osr.SpatialReference()
>>> _ = ref.ImportFromEPSG(32618)
>>> wkt = ref.ExportToWkt()
>>> wkt = wkt[: wkt.rfind(",AUTHORITY")] + "]"
>>> osr.SpatialReference(wkt=wkt).GetAuthorityCode(None) is None
True
>>> get_epsg_from_prj(wkt)
32618
  • A genuinely custom CRS that matches no database entry raises CRSError rather than returning a child unit/datum code (here GDAL's spherical-earth GRIB GEOGCS, whose only authority node is the degree unit EPSG:9122): 
    >>> grib_wkt = (
...     'GEOGCS["Coordinate System imported from GRIB file",'
...     'DATUM["unnamed",SPHEROID["Sphere",6371229,0]],'
...     'PRIMEM["Greenwich",0],'
...     'UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],'
...     'AXIS["Latitude",NORTH],AXIS["Longitude",EAST]]'
... )
>>> get_epsg_from_prj(grib_wkt)
Traceback (most recent call last):
    ...
pyramids.base._errors.CRSError: get_epsg_from_prj could not resolve an EPSG code ...
  • An empty projection string raises CRSError (a ValueError subclass): 
    >>> get_epsg_from_prj("")
Traceback (most recent call last):
    ...
pyramids.base._errors.CRSError: get_epsg_from_prj received an empty projection string. ...
Source code in src/pyramids/base/crs.py 
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
def get_epsg_from_prj(prj: str) -> int:
    """Return the EPSG code identified by a projection string.

    Resolves the EPSG of the *root* CRS object in three steps:
    :meth:`osr.SpatialReference.AutoIdentifyEPSG` (tags recognisable
    CRSes), then the root ``AUTHORITY`` code, then a confident,
    unambiguous :meth:`osr.SpatialReference.FindMatches` PROJ-database
    lookup for well-known CRSes whose WKT lacks a root authority (e.g. a
    UTM ``PROJCS`` from GDAL's AAIGrid driver). The code of a child
    unit/datum node is never returned as if it were a CRS — that bug
    (issue #403) made GRIB rasters resolve to the degree-unit EPSG:9122
    and UTM ASCII grids to the WGS_1984 datum EPSG:6326.

    An empty input string is no longer silently mapped to `4326`; that
    legacy default masked real configuration errors. Callers that
    genuinely want a fallback should handle the `CRSError` themselves,
    or use :func:`epsg_from_wkt` which accepts an explicit `default`.

    Args:
        prj (str): Projection string.

    Returns:
        int: The resolved EPSG code.

    Raises:
        CRSError: If `prj` is an empty string, or if its root CRS carries
            no EPSG authority *and* matches no PROJ-database entry (e.g. a
            custom spherical-earth GRIB GEOGCS). The unit/datum codes of
            child nodes are never returned as a CRS.

    Examples:
        - Resolve EPSG:4326 from its standard WKT representation:
            ```python
            >>> from osgeo import osr
            >>> ref = osr.SpatialReference()
            >>> _ = ref.ImportFromEPSG(4326)
            >>> get_epsg_from_prj(ref.ExportToWkt())
            4326

            ```
        - Resolve EPSG:3857 (Web Mercator) from its WKT representation:
            ```python
            >>> from osgeo import osr
            >>> ref = osr.SpatialReference()
            >>> _ = ref.ImportFromEPSG(3857)
            >>> get_epsg_from_prj(ref.ExportToWkt())
            3857

            ```
        - A well-known CRS whose WKT carries no root authority still
          resolves, via a confident PROJ-database match (here a
          "WGS 84 / UTM zone 18N" PROJCS with its root authority stripped):
            ```python
            >>> from osgeo import osr
            >>> ref = osr.SpatialReference()
            >>> _ = ref.ImportFromEPSG(32618)
            >>> wkt = ref.ExportToWkt()
            >>> wkt = wkt[: wkt.rfind(",AUTHORITY")] + "]"
            >>> osr.SpatialReference(wkt=wkt).GetAuthorityCode(None) is None
            True
            >>> get_epsg_from_prj(wkt)
            32618

            ```
        - A genuinely custom CRS that matches no database entry raises
          `CRSError` rather than returning a child unit/datum code (here
          GDAL's spherical-earth GRIB GEOGCS, whose only authority node is
          the degree unit EPSG:9122):
            ```python
            >>> grib_wkt = (
            ...     'GEOGCS["Coordinate System imported from GRIB file",'
            ...     'DATUM["unnamed",SPHEROID["Sphere",6371229,0]],'
            ...     'PRIMEM["Greenwich",0],'
            ...     'UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],'
            ...     'AXIS["Latitude",NORTH],AXIS["Longitude",EAST]]'
            ... )
            >>> get_epsg_from_prj(grib_wkt)
            Traceback (most recent call last):
                ...
            pyramids.base._errors.CRSError: get_epsg_from_prj could not resolve an EPSG code ...

            ```
        - An empty projection string raises `CRSError` (a `ValueError` subclass):
            ```python
            >>> get_epsg_from_prj("")
            Traceback (most recent call last):
                ...
            pyramids.base._errors.CRSError: get_epsg_from_prj received an empty projection string. ...

            ```
    """
    if prj == "":
        raise CRSError(
            "get_epsg_from_prj received an empty projection string. "
            "An empty projection is ambiguous and is no longer "
            "silently defaulted to EPSG:4326. If you want "
            "a fallback EPSG, catch CRSError (also a ValueError) "
            "and supply it at the call site, or call "
            "epsg_from_wkt(prj, default=...)."
        )
    srs = create_sr_from_proj(prj)
    try:
        # AutoIdentifyEPSG attaches a root EPSG authority when it can
        # recognise the CRS; we ignore its return code and read the root
        # authority below. It raises "Unsupported SRS" for custom CRSes
        # it cannot identify (e.g. GDAL's spherical-earth GRIB GEOGCS).
        srs.AutoIdentifyEPSG()
    except RuntimeError:
        pass

    # Resolve the EPSG of the *root* CRS object only. Do NOT fall back to
    # GetAttrValue("AUTHORITY", 1): that walks the WKT tree depth-first and
    # returns the first AUTHORITY node, which for a CRS whose root carries
    # no authority is a child unit/datum code (e.g. the degree-unit
    # EPSG:9122 inside a GRIB GEOGCS, or the WGS_1984 datum EPSG:6326 inside
    # a UTM PROJCS) — a non-CRS code that breaks every downstream
    # sr_from_epsg() call. See issue #403.
    code = srs.GetAuthorityCode(None)
    if code is None:
        # AutoIdentifyEPSG could not tag the root, but the CRS may still be a
        # well-known database entry whose WKT simply lacks an AUTHORITY node
        # (e.g. a UTM PROJCS). Try an exact PROJ-database match before giving
        # up, so identifiable CRSes resolve to their true CRS code.
        code = _epsg_from_db_match(srs)
    if code is None:
        raise CRSError(
            "get_epsg_from_prj could not resolve an EPSG code from the "
            "projection: its root CRS carries no EPSG authority and matches "
            "no PROJ-database entry. This is expected for genuinely custom "
            "CRSes such as GDAL's spherical-earth GRIB GEOGCS. Catch CRSError "
            "(also a ValueError) and supply a fallback, or call "
            "epsg_from_wkt(prj, default=...)."
        )
    return int(code)

pyramids.base.crs.epsg_from_wkt(wkt, default=4326) #

Resolve an EPSG code from a WKT / Proj string with a fallback.

Wraps :func:get_epsg_from_prj to absorb the get_epsg_from_prj(wkt) if wkt else default idiom that was previously open-coded in four places across the dataset stack. Returns default when wkt is empty (or None), and also when get_epsg_from_prj cannot resolve an EPSG from a non-empty wkt (it raises :class:CRSError for a custom CRS whose root carries no EPSG authority — e.g. a spherical-earth GRIB GEOGCS); otherwise delegates to :func:get_epsg_from_prj.

Use this in places where an empty projection should be treated as a soft "unknown CRS, assume WGS84" rather than a hard error — for example the Dataset.epsg property on a freshly-built in-memory raster that has no projection metadata yet. Use :func:get_epsg_from_prj directly when you want the strict behaviour where an empty projection raises.

Parameters:

Name Type Description Default
wkt str

Projection string (WKT, ESRI WKT, or Proj4). An empty string or None returns default.

 required
default int

EPSG code to return when wkt is empty / None, or when its CRS cannot be resolved to an EPSG. Defaults to 4326 (the historical pyramids default).

 4326

Returns:

Name Type Description
int int

EPSG code resolved from wkt, or default when wkt is
int

empty or its CRS carries no resolvable EPSG.

Examples:

  • Empty input falls back to the supplied default: 
    >>> from pyramids.base.crs import epsg_from_wkt
>>> epsg_from_wkt("")
4326
>>> epsg_from_wkt("", default=3857)
3857
  • Non-empty WKT delegates to :func:get_epsg_from_prj: 
    >>> from osgeo import osr
>>> from pyramids.base.crs import epsg_from_wkt
>>> ref = osr.SpatialReference()
>>> _ = ref.ImportFromEPSG(3857)
>>> epsg_from_wkt(ref.ExportToWkt())
3857
  • An unresolvable custom CRS falls back to default instead of raising (here GDAL's spherical-earth GRIB GEOGCS): 
    >>> from pyramids.base.crs import epsg_from_wkt
>>> grib_wkt = (
...     'GEOGCS["Coordinate System imported from GRIB file",'
...     'DATUM["unnamed",SPHEROID["Sphere",6371229,0]],'
...     'PRIMEM["Greenwich",0],'
...     'UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],'
...     'AXIS["Latitude",NORTH],AXIS["Longitude",EAST]]'
... )
>>> epsg_from_wkt(grib_wkt)
4326
>>> epsg_from_wkt(grib_wkt, default=3857)
3857
Source code in src/pyramids/base/crs.py 
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
def epsg_from_wkt(wkt: str, default: int = 4326) -> int:
    """Resolve an EPSG code from a WKT / Proj string with a fallback.

    Wraps :func:`get_epsg_from_prj` to absorb the
    `get_epsg_from_prj(wkt) if wkt else default` idiom that was
    previously open-coded in four places across the dataset stack.
    Returns `default` when `wkt` is empty (or `None`), and also when
    `get_epsg_from_prj` cannot resolve an EPSG from a non-empty `wkt`
    (it raises :class:`CRSError` for a custom CRS whose root carries no
    EPSG authority — e.g. a spherical-earth GRIB GEOGCS); otherwise
    delegates to :func:`get_epsg_from_prj`.

    Use this in places where an empty projection should be treated as
    a soft "unknown CRS, assume WGS84" rather than a hard error — for
    example the `Dataset.epsg` property on a freshly-built
    in-memory raster that has no projection metadata yet. Use
    :func:`get_epsg_from_prj` directly when you want the strict
    behaviour where an empty projection raises.

    Args:
        wkt: Projection string (WKT, ESRI WKT, or Proj4). An empty
            string or `None` returns `default`.
        default: EPSG code to return when `wkt` is empty / `None`, or
            when its CRS cannot be resolved to an EPSG. Defaults to
            `4326` (the historical pyramids default).

    Returns:
        int: EPSG code resolved from `wkt`, or `default` when `wkt` is
        empty or its CRS carries no resolvable EPSG.

    Examples:
        - Empty input falls back to the supplied default:
            ```python
            >>> from pyramids.base.crs import epsg_from_wkt
            >>> epsg_from_wkt("")
            4326
            >>> epsg_from_wkt("", default=3857)
            3857

            ```
        - Non-empty WKT delegates to :func:`get_epsg_from_prj`:
            ```python
            >>> from osgeo import osr
            >>> from pyramids.base.crs import epsg_from_wkt
            >>> ref = osr.SpatialReference()
            >>> _ = ref.ImportFromEPSG(3857)
            >>> epsg_from_wkt(ref.ExportToWkt())
            3857

            ```
        - An unresolvable custom CRS falls back to `default` instead of
          raising (here GDAL's spherical-earth GRIB GEOGCS):
            ```python
            >>> from pyramids.base.crs import epsg_from_wkt
            >>> grib_wkt = (
            ...     'GEOGCS["Coordinate System imported from GRIB file",'
            ...     'DATUM["unnamed",SPHEROID["Sphere",6371229,0]],'
            ...     'PRIMEM["Greenwich",0],'
            ...     'UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],'
            ...     'AXIS["Latitude",NORTH],AXIS["Longitude",EAST]]'
            ... )
            >>> epsg_from_wkt(grib_wkt)
            4326
            >>> epsg_from_wkt(grib_wkt, default=3857)
            3857

            ```
    """
    if not wkt:
        result = default
    else:
        try:
            result = get_epsg_from_prj(wkt)
        except CRSError:
            # Non-empty but unresolvable CRS (e.g. a custom spherical-earth
            # GRIB GEOGCS that carries no root EPSG authority). Treat it as
            # the same soft "unknown CRS" case as empty input rather than
            # propagating the hard error to property reads like Dataset.epsg.
            result = default
    return result

pyramids.base.crs.reproject_coordinates(x, y, *, from_crs=4326, to_crs=3857, precision=6) #

Reproject parallel x / y coordinate lists between CRSes.

Argument and return order is (x, y) throughout; accepts any CRS form :meth:pyproj.Transformer.from_crs understands (EPSG int, EPSG string, WKT, Proj4, :class:pyproj.CRS).

Parameters:

Name Type Description Default
x list[float]

X-coordinates in the source CRS (longitudes when from_crs is geographic).

 required
y list[float]

Y-coordinates in the source CRS (latitudes when from_crs is geographic).

 required
from_crs Any

Source CRS. Accepts anything :meth:pyproj.Transformer.from_crs accepts: EPSG integer (4326), authority string ("EPSG:4326"), WKT, Proj4, or a :class:pyproj.CRS instance. Default 4326.

 4326
to_crs Any

Target CRS, same forms as from_crs. Default 3857.

 3857
precision int | None

Decimal places to round each returned coordinate to. Pass None to disable rounding. Default 6.

 6

Returns:

Type Description
tuple[list[float], list[float]]

tuple[list[float], list[float]]: (x, y) in the target CRS.

Raises:

Type Description
ValueError

If len(x)!= len(y).
CRSError

If :meth:pyproj.Transformer.from_crs raises one of pyproj.exceptions.CRSError (malformed WKT / proj string), TypeError (input is not CRS-like — e.g. a bare object()), or ValueError (out-of-range EPSG integer). The wrapper converts each into pyramids' :class:pyramids.base._errors.CRSError so callers do not need to import pyproj to catch bad-CRS failures, and the message names both CRSes plus the underlying explanation. Other exception types (AttributeError, ImportError, …) propagate unchanged — they signal a real bug, not a bad user input.

Examples:

  • Reproject a WGS84 point into Web Mercator: 
    >>> from pyramids.base.crs import reproject_coordinates
>>> x, y = reproject_coordinates(
...     [31.0], [30.0], from_crs=4326, to_crs=3857
... )
>>> round(x[0])
3450904
>>> round(y[0])
3503550
Source code in src/pyramids/base/crs.py 
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
def reproject_coordinates(
    x: list[float],
    y: list[float],
    *,
    from_crs: Any = 4326,
    to_crs: Any = 3857,
    precision: int | None = 6,
) -> tuple[list[float], list[float]]:
    """Reproject parallel x / y coordinate lists between CRSes.

    Argument and return order is `(x, y)` throughout; accepts any
    CRS form :meth:`pyproj.Transformer.from_crs` understands (EPSG
    int, EPSG string, WKT, Proj4, :class:`pyproj.CRS`).

    Args:
        x (list[float]):
            X-coordinates in the source CRS (longitudes when
            `from_crs` is geographic).
        y (list[float]):
            Y-coordinates in the source CRS (latitudes when
            `from_crs` is geographic).
        from_crs:
            Source CRS. Accepts anything
            :meth:`pyproj.Transformer.from_crs` accepts: EPSG integer
            (`4326`), authority string (`"EPSG:4326"`), WKT, Proj4,
            or a :class:`pyproj.CRS` instance. Default `4326`.
        to_crs:
            Target CRS, same forms as `from_crs`. Default `3857`.
        precision (int | None):
            Decimal places to round each returned coordinate to. Pass
            `None` to disable rounding. Default `6`.

    Returns:
        tuple[list[float], list[float]]: `(x, y)` in the target CRS.

    Raises:
        ValueError: If `len(x)!= len(y)`.
        CRSError: If :meth:`pyproj.Transformer.from_crs` raises one
            of `pyproj.exceptions.CRSError` (malformed WKT / proj
            string), `TypeError` (input is not CRS-like — e.g. a
            bare `object()`), or `ValueError` (out-of-range EPSG
            integer). The wrapper converts each into pyramids'
            :class:`pyramids.base._errors.CRSError` so callers do not
            need to import pyproj to catch bad-CRS failures, and the
            message names both CRSes plus the underlying explanation.
            Other exception types (`AttributeError`, `ImportError`,
            …) propagate unchanged — they signal a real bug, not a bad
            user input.

    Examples:
        - Reproject a WGS84 point into Web Mercator:
            ```python
            >>> from pyramids.base.crs import reproject_coordinates
            >>> x, y = reproject_coordinates(
            ...     [31.0], [30.0], from_crs=4326, to_crs=3857
            ... )
            >>> round(x[0])
            3450904
            >>> round(y[0])
            3503550

            ```
    """
    if len(x) != len(y):
        raise ValueError(
            f"x and y must have equal length; got len(x)={len(x)} "
            f"vs. len(y)={len(y)}."
        )
    try:
        transformer = Transformer.from_crs(from_crs, to_crs, always_xy=True)
    except (pyproj.exceptions.CRSError, TypeError, ValueError) as exc:
        raise CRSError(
            f"reproject_coordinates failed to parse CRS "
            f"(from_crs={from_crs!r}, to_crs={to_crs!r}): {exc}"
        ) from exc
    xs = np.full(len(x), np.nan)
    ys = np.full(len(x), np.nan)
    for i in range(len(x)):
        nx, ny = transformer.transform(x[i], y[i])
        if precision is not None:
            nx = round(nx, precision)
            ny = round(ny, precision)
        xs[i] = nx
        ys[i] = ny
    return xs.tolist(), ys.tolist()