Skip to content

Commit

Permalink
BUG: Add datum ensemble support to GeographicCRS
Browse files Browse the repository at this point in the history
  • Loading branch information
@snowman2
snowman2 committed Mar 23, 2023
1 parent f82cae5 commit dac3ee8
Show file tree
Hide file tree
Showing 3 changed files with 32 additions and 26 deletions.
13 changes: 8 additions & 5 deletions pyproj/crs/_cf1x8.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -34,6 +34,12 @@


def _horizontal_datum_from_params(cf_params):
datum_name = cf_params.get("horizontal_datum_name")
if datum_name and datum_name not in ("undefined", "unknown"):
try:
return Datum.from_name(datum_name)
except CRSError:
pass
# step 1: build ellipsoid
ellipsoid = None
ellipsoid_name = cf_params.get("reference_ellipsoid_name")
Expand All @@ -46,7 +52,7 @@ def _horizontal_datum_from_params(cf_params):
radius=cf_params.get("earth_radius"),
)
except CRSError:
if ellipsoid_name:
if ellipsoid_name and ellipsoid_name not in ("undefined", "unknown"):
ellipsoid = Ellipsoid.from_name(ellipsoid_name)

# step 2: build prime meridian
Expand All @@ -58,19 +64,16 @@ def _horizontal_datum_from_params(cf_params):
longitude=cf_params["longitude_of_prime_meridian"],
)
except KeyError:
if prime_meridian_name:
if prime_meridian_name and prime_meridian_name not in ("undefined", "unknown"):
prime_meridian = PrimeMeridian.from_name(prime_meridian_name)

# step 3: build datum
datum_name = cf_params.get("horizontal_datum_name")
if ellipsoid or prime_meridian:
return CustomDatum(
name=datum_name or "undefined",
ellipsoid=ellipsoid or "WGS 84",
prime_meridian=prime_meridian or "Greenwich",
)
if datum_name:
return Datum.from_name(datum_name)
return None


Expand Down
12 changes: 7 additions & 5 deletions pyproj/crs/crs.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -718,8 +718,6 @@ def to_cf(
return cf_dict

# handle projected CRS
if self.is_projected and self.datum:
cf_dict["horizontal_datum_name"] = self.datum.name
coordinate_operation = None
if not self.is_bound and self.is_projected:
coordinate_operation = self.coordinate_operation
Expand Down Expand Up @@ -1758,31 +1756,35 @@ class GeographicCRS(CustomConstructorCRS):
def __init__(
self,
name: str = "undefined",
datum: Any = "urn:ogc:def:datum:EPSG::6326",
datum: Any = "urn:ogc:def:ensemble:EPSG::6326",
ellipsoidal_cs: Optional[Any] = None,
) -> None:
"""
Parameters
----------
name: str, default="undefined"
Name of the CRS.
datum: Any, default="urn:ogc:def:datum:EPSG::6326"
datum: Any, default="urn:ogc:def:ensemble:EPSG::6326"
Anything accepted by :meth:`pyproj.crs.Datum.from_user_input` or
a :class:`pyproj.crs.datum.CustomDatum`.
ellipsoidal_cs: Any, optional
Input to create an Ellipsoidal Coordinate System.
Anything accepted by :meth:`pyproj.crs.CoordinateSystem.from_user_input`
or an Ellipsoidal Coordinate System created from :ref:`coordinate_system`.
"""
datum = Datum.from_user_input(datum).to_json_dict()
geographic_crs_json = {
"$schema": "https://proj.org/schemas/v0.2/projjson.schema.json",
"type": "GeographicCRS",
"name": name,
"datum": Datum.from_user_input(datum).to_json_dict(),
"coordinate_system": CoordinateSystem.from_user_input(
ellipsoidal_cs or Ellipsoidal2DCS()
).to_json_dict(),
}
if datum["type"] == "DatumEnsemble":
geographic_crs_json["datum_ensemble"] = datum
else:
geographic_crs_json["datum"] = datum
super().__init__(geographic_crs_json)


Expand Down
33 changes: 17 additions & 16 deletions test/crs/test_crs_cf.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -362,7 +362,7 @@ def test_cf_rotated_latlon():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "rotated_latitude_longitude",
"grid_north_pole_latitude": 32.5,
"grid_north_pole_longitude": 170.0,
Expand Down Expand Up @@ -534,7 +534,7 @@ def test_cf_lambert_conformal_conic_1sp():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "lambert_conformal_conic",
"longitude_of_central_meridian": 265.0,
"false_easting": 0.0,
Expand Down Expand Up @@ -598,7 +598,7 @@ def test_cf_lambert_conformal_conic_2sp(standard_parallel):
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "lambert_conformal_conic",
"standard_parallel": (25.0, 30.0),
"latitude_of_projection_origin": 25.0,
Expand Down Expand Up @@ -753,7 +753,7 @@ def test_geos_crs_sweep():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "geostationary",
"sweep_angle_axis": "x",
"perspective_point_height": 1.0,
Expand Down Expand Up @@ -801,7 +801,7 @@ def test_geos_crs_fixed_angle_axis():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "geostationary",
"sweep_angle_axis": "x",
"perspective_point_height": 1.0,
Expand Down Expand Up @@ -901,7 +901,7 @@ def test_mercator_b():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "mercator",
"standard_parallel": 21.354,
"longitude_of_projection_origin": 10.0,
Expand Down Expand Up @@ -1205,6 +1205,7 @@ def test_azimuthal_equidistant():
assert cf_dict.pop("crs_wkt").startswith("PROJCRS[")
assert cf_dict == expected_cf
# test roundtrip
expected_cf["horizontal_datum_name"] = "World Geodetic System 1984 ensemble"
_test_roundtrip(expected_cf, "PROJCRS[")
# test coordinate system
assert crs.cs_to_cf() == [
Expand Down Expand Up @@ -1232,7 +1233,7 @@ def test_lambert_azimuthal_equal_area():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "lambert_azimuthal_equal_area",
"latitude_of_projection_origin": 1.0,
"longitude_of_projection_origin": 2.0,
Expand Down Expand Up @@ -1272,7 +1273,7 @@ def test_lambert_cylindrical_equal_area():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "lambert_cylindrical_equal_area",
"standard_parallel": 1.0,
"longitude_of_central_meridian": 2.0,
Expand Down Expand Up @@ -1312,7 +1313,7 @@ def test_mercator_a():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "mercator",
"standard_parallel": 0.0,
"longitude_of_projection_origin": 2.0,
Expand Down Expand Up @@ -1353,7 +1354,7 @@ def test_orthographic():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "orthographic",
"latitude_of_projection_origin": 1.0,
"longitude_of_projection_origin": 2.0,
Expand Down Expand Up @@ -1393,7 +1394,7 @@ def test_polar_stereographic_a():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "polar_stereographic",
"latitude_of_projection_origin": 90.0,
"straight_vertical_longitude_from_pole": 1.0,
Expand Down Expand Up @@ -1434,7 +1435,7 @@ def test_polar_stereographic_b():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "polar_stereographic",
"standard_parallel": 0.0,
"straight_vertical_longitude_from_pole": 1.0,
Expand Down Expand Up @@ -1474,7 +1475,7 @@ def test_stereographic():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "stereographic",
"latitude_of_projection_origin": 0.0,
"longitude_of_projection_origin": 1.0,
Expand Down Expand Up @@ -1515,7 +1516,7 @@ def test_sinusoidal():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "sinusoidal",
"longitude_of_projection_origin": 0.0,
"false_easting": 1.0,
Expand Down Expand Up @@ -1554,7 +1555,7 @@ def test_vertical_perspective():
"reference_ellipsoid_name": "WGS 84",
"longitude_of_prime_meridian": 0.0,
"prime_meridian_name": "Greenwich",
"horizontal_datum_name": "World Geodetic System 1984",
"horizontal_datum_name": "World Geodetic System 1984 ensemble",
"grid_mapping_name": "vertical_perspective",
"perspective_point_height": 50.0,
"latitude_of_projection_origin": 0.0,
Expand Down Expand Up @@ -1891,7 +1892,7 @@ def test_3d_ellipsoidal_cs_depth():
"name": "WGS 84 (geographic 3D)",
"datum": {
"type": "GeodeticReferenceFrame",
"name": "World Geodetic System 1984",
"name": "World Geodetic System 1984 ensemble",
"ellipsoid": {
"name": "WGS 84",
"semi_major_axis": 6378137,
Expand Down

0 comments on commit dac3ee8

Please sign in to comment.