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 
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
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 
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
97
98
99
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 
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
161
162
163
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 
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
372
373
374
375
376
377
378
379
380
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 not (code and str(code).isdigit()):
        # No usable root EPSG code: either absent (AutoIdentifyEPSG could not tag
        # the root — e.g. a UTM PROJCS whose WKT lacks an AUTHORITY node), or a
        # *non-numeric* authority code — notably OGC:CRS84 (the lon/lat WGS 84 that
        # WMS/WMTS layers report), where GetAuthorityCode returns "CRS84". Try an
        # exact PROJ-database match before giving up. A UTM PROJCS resolves to its
        # numeric code here; CRS84's best match is CRS84 itself (still non-numeric),
        # so it is dropped to None below and raises CRSError — the soft
        # epsg_from_wkt() then supplies its default. The point of this branch is to
        # never crash on int("CRS84"), not to equate CRS84 with EPSG:4326.
        match = _epsg_from_db_match(srs)
        code = match if (match and str(match).isdigit()) else None
    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 
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
454
455
456
457
458
459
460
461
462
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 
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
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()