Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Oblique and Rotated Mercator #5548

Merged
merged 33 commits into from
Oct 24, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
33 commits
Select commit Hold shift + click to select a range
c3e761c
Introduce new coord system classes.
trexfeathers Oct 13, 2023
d60eab9
Add loading code for oblique mercator.
trexfeathers Oct 13, 2023
5ca3873
Fix for azimuth check.
trexfeathers Oct 17, 2023
0176b1d
Add saving code for oblique mercator.
trexfeathers Oct 17, 2023
af074a9
Fix to rotated repr.
trexfeathers Oct 17, 2023
518200a
Scale factor wording fix.
trexfeathers Oct 17, 2023
7a74af9
Tests first pass.
trexfeathers Oct 20, 2023
a81507c
Temp test disable.
trexfeathers Oct 23, 2023
27c486f
Temp RotatedMercator test disable.
trexfeathers Oct 23, 2023
dd24604
Deprecate RotatedMercator.
trexfeathers Oct 23, 2023
b2b6a11
Revert "Temp RotatedMercator test disable."
trexfeathers Oct 23, 2023
687ccd2
First attempted fix for RM test inheritance.
trexfeathers Oct 23, 2023
a01eab2
Revert "Temp test disable."
trexfeathers Oct 23, 2023
2499945
Fix warnings doctests.
trexfeathers Oct 23, 2023
709c92e
Add deprecation test for RotatedMercator.
trexfeathers Oct 23, 2023
e56a0a1
Oblique Mercator loading tests.
trexfeathers Oct 23, 2023
75878f6
Oblique Mercator loading deprecation test.
trexfeathers Oct 23, 2023
c222712
Saving test for Oblique Mercator.
trexfeathers Oct 23, 2023
cdf3ffd
Fix isinstance() check.
trexfeathers Oct 23, 2023
a1005a5
What's New entry.
trexfeathers Oct 23, 2023
77eba55
Temp test disable.
trexfeathers Oct 23, 2023
ff251b7
More temp test disabling.
trexfeathers Oct 23, 2023
8b610e5
WIP testing.
trexfeathers Oct 23, 2023
87a6bce
WIP testing.
trexfeathers Oct 23, 2023
2e700e3
Revert "More temp test disabling."
trexfeathers Oct 23, 2023
f26bbd7
Revert "Temp test disable."
trexfeathers Oct 23, 2023
6a1afae
Use RotatedMercator inheritance for isinstance() check.
trexfeathers Oct 23, 2023
7a63804
Check grid_mapping_name instead of using isinstance().
trexfeathers Oct 24, 2023
55174bf
Better type hinting.
trexfeathers Oct 24, 2023
31aa968
Use return over yield in a fixture.
trexfeathers Oct 24, 2023
4ddd331
Duck typing comment.
trexfeathers Oct 24, 2023
20fca6e
Better grid_mapping_name checking.
trexfeathers Oct 24, 2023
2bcee96
Better structure for test parameterisation.
trexfeathers Oct 24, 2023
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
12 changes: 6 additions & 6 deletions docs/src/further_topics/filtering_warnings.rst
Original file line number Diff line number Diff line change
Expand Up @@ -47,9 +47,9 @@ Warnings:

>>> my_operation()
...
iris/coord_systems.py:454: IrisUserWarning: Setting inverse_flattening does not affect other properties of the GeogCS object. To change other properties set them explicitly or create a new GeogCS instance.
iris/coord_systems.py:456: IrisUserWarning: Setting inverse_flattening does not affect other properties of the GeogCS object. To change other properties set them explicitly or create a new GeogCS instance.
warnings.warn(wmsg, category=iris.exceptions.IrisUserWarning)
iris/coord_systems.py:821: IrisDefaultingWarning: Discarding false_easting and false_northing that are not used by Cartopy.
iris/coord_systems.py:823: IrisDefaultingWarning: Discarding false_easting and false_northing that are not used by Cartopy.
warnings.warn(

Warnings can be suppressed using the Python warnings filter with the ``ignore``
Expand Down Expand Up @@ -110,7 +110,7 @@ You can target specific Warning messages, e.g.
... warnings.filterwarnings("ignore", message="Discarding false_easting")
... my_operation()
...
iris/coord_systems.py:454: IrisUserWarning: Setting inverse_flattening does not affect other properties of the GeogCS object. To change other properties set them explicitly or create a new GeogCS instance.
iris/coord_systems.py:456: IrisUserWarning: Setting inverse_flattening does not affect other properties of the GeogCS object. To change other properties set them explicitly or create a new GeogCS instance.
warnings.warn(wmsg, category=iris.exceptions.IrisUserWarning)

::
Expand All @@ -125,10 +125,10 @@ Or you can target Warnings raised by specific lines of specific modules, e.g.
.. doctest:: filtering_warnings

>>> with warnings.catch_warnings():
... warnings.filterwarnings("ignore", module="iris.coord_systems", lineno=454)
... warnings.filterwarnings("ignore", module="iris.coord_systems", lineno=456)
... my_operation()
...
iris/coord_systems.py:821: IrisDefaultingWarning: Discarding false_easting and false_northing that are not used by Cartopy.
iris/coord_systems.py:823: IrisDefaultingWarning: Discarding false_easting and false_northing that are not used by Cartopy.
warnings.warn(

::
Expand Down Expand Up @@ -188,7 +188,7 @@ module during execution:
... )
... my_operation()
...
iris/coord_systems.py:454: IrisUserWarning: Setting inverse_flattening does not affect other properties of the GeogCS object. To change other properties set them explicitly or create a new GeogCS instance.
iris/coord_systems.py:456: IrisUserWarning: Setting inverse_flattening does not affect other properties of the GeogCS object. To change other properties set them explicitly or create a new GeogCS instance.
warnings.warn(wmsg, category=iris.exceptions.IrisUserWarning)

----
Expand Down
4 changes: 4 additions & 0 deletions docs/src/whatsnew/latest.rst
Original file line number Diff line number Diff line change
Expand Up @@ -34,6 +34,10 @@ This document explains the changes made to Iris for this release
:class:`UserWarning`\s for richer filtering. The full index of
sub-categories can be seen here: :mod:`iris.exceptions` . (:pull:`5498`)

#. `@trexfeathers`_ added the :class:`~iris.coord_systems.ObliqueMercator`
and :class:`~iris.coord_systems.RotatedMercator` coordinate systems,
complete with NetCDF loading and saving. (:pull:`5548`)


🐛 Bugs Fixed
=============
Expand Down
196 changes: 196 additions & 0 deletions lib/iris/coord_systems.py
Original file line number Diff line number Diff line change
Expand Up @@ -10,11 +10,13 @@

from abc import ABCMeta, abstractmethod
from functools import cached_property
import re
import warnings

import cartopy.crs as ccrs
import numpy as np

from iris._deprecation import warn_deprecated
import iris.exceptions


Expand Down Expand Up @@ -1634,3 +1636,197 @@ def as_cartopy_crs(self):

def as_cartopy_projection(self):
return self.as_cartopy_crs()


class ObliqueMercator(CoordSystem):
"""
A cylindrical map projection, with XY coordinates measured in metres.

Designed for regions not well suited to :class:`Mercator` or
:class:`TransverseMercator`, as the positioning of the cylinder is more
customisable.

See Also
--------
:class:`RotatedMercator`

"""

grid_mapping_name = "oblique_mercator"

def __init__(
self,
azimuth_of_central_line,
latitude_of_projection_origin,
longitude_of_projection_origin,
false_easting=None,
false_northing=None,
scale_factor_at_projection_origin=None,
ellipsoid=None,
):
"""
Constructs an ObliqueMercator object.

Parameters
----------
azimuth_of_central_line : float
Azimuth of centerline clockwise from north at the center point of
the centre line.
latitude_of_projection_origin : float
The true longitude of the central meridian in degrees.
longitude_of_projection_origin: float
The true latitude of the planar origin in degrees.
false_easting: float, optional
X offset from the planar origin in metres.
Defaults to 0.0 .
false_northing: float, optional
Y offset from the planar origin in metres.
Defaults to 0.0 .
scale_factor_at_projection_origin: float, optional
Scale factor at the central meridian.
Defaults to 1.0 .
ellipsoid: :class:`GeogCS`, optional
If given, defines the ellipsoid.

Examples
--------
>>> from iris.coord_systems import GeogCS, ObliqueMercator
>>> my_ellipsoid = GeogCS(6371229.0, None, 0.0)
>>> ObliqueMercator(90.0, -22.0, -59.0, -25000.0, -25000.0, 1., my_ellipsoid)
ObliqueMercator(azimuth_of_central_line=90.0, latitude_of_projection_origin=-22.0, longitude_of_projection_origin=-59.0, false_easting=-25000.0, false_northing=-25000.0, scale_factor_at_projection_origin=1.0, ellipsoid=GeogCS(6371229.0))

"""
#: Azimuth of centerline clockwise from north.
self.azimuth_of_central_line = float(azimuth_of_central_line)

#: True latitude of planar origin in degrees.
self.latitude_of_projection_origin = float(
latitude_of_projection_origin
)

#: True longitude of planar origin in degrees.
self.longitude_of_projection_origin = float(
longitude_of_projection_origin
)

#: X offset from planar origin in metres.
self.false_easting = _arg_default(false_easting, 0)

#: Y offset from planar origin in metres.
self.false_northing = _arg_default(false_northing, 0)

#: Scale factor at the central meridian.
self.scale_factor_at_projection_origin = _arg_default(
scale_factor_at_projection_origin, 1.0
)

#: Ellipsoid definition (:class:`GeogCS` or None).
self.ellipsoid = ellipsoid

def __repr__(self):
return (
"{!s}(azimuth_of_central_line={!r}, "
"latitude_of_projection_origin={!r}, "
"longitude_of_projection_origin={!r}, false_easting={!r}, "
"false_northing={!r}, scale_factor_at_projection_origin={!r}, "
"ellipsoid={!r})".format(
self.__class__.__name__,
self.azimuth_of_central_line,
self.latitude_of_projection_origin,
self.longitude_of_projection_origin,
self.false_easting,
self.false_northing,
self.scale_factor_at_projection_origin,
self.ellipsoid,
)
)

def as_cartopy_crs(self):
globe = self._ellipsoid_to_globe(self.ellipsoid, None)

return ccrs.ObliqueMercator(
central_longitude=self.longitude_of_projection_origin,
central_latitude=self.latitude_of_projection_origin,
false_easting=self.false_easting,
false_northing=self.false_northing,
scale_factor=self.scale_factor_at_projection_origin,
azimuth=self.azimuth_of_central_line,
globe=globe,
)

def as_cartopy_projection(self):
return self.as_cartopy_crs()


class RotatedMercator(ObliqueMercator):
"""
:class:`ObliqueMercator` with ``azimuth_of_central_line=90``.

As noted in CF versions 1.10 and earlier:

The Rotated Mercator projection is an Oblique Mercator projection
with azimuth = +90.

.. deprecated:: 3.8.0
This coordinate system was introduced as already scheduled for removal
in a future release, since CF version 1.11 onwards now requires use of
:class:`ObliqueMercator` with ``azimuth_of_central_line=90.`` .
Any :class:`RotatedMercator` instances will always be saved to NetCDF
as the ``oblique_mercator`` grid mapping.

"""

def __init__(
self,
latitude_of_projection_origin,
longitude_of_projection_origin,
false_easting=None,
false_northing=None,
scale_factor_at_projection_origin=None,
ellipsoid=None,
):
"""
Constructs a RotatedMercator object.

Parameters
----------
latitude_of_projection_origin : float
The true longitude of the central meridian in degrees.
longitude_of_projection_origin: float
The true latitude of the planar origin in degrees.
false_easting: float, optional
X offset from the planar origin in metres.
Defaults to 0.0 .
false_northing: float, optional
Y offset from the planar origin in metres.
Defaults to 0.0 .
scale_factor_at_projection_origin: float, optional
Scale factor at the central meridian.
Defaults to 1.0 .
ellipsoid: :class:`GeogCS`, optional
If given, defines the ellipsoid.

"""
message = (
"iris.coord_systems.RotatedMercator is deprecated, and will be "
"removed in a future release. Instead please use "
"iris.coord_systems.ObliqueMercator with "
"azimuth_of_central_line=90 ."
)
warn_deprecated(message)

super().__init__(
90.0,
latitude_of_projection_origin,
longitude_of_projection_origin,
false_easting,
false_northing,
scale_factor_at_projection_origin,
ellipsoid,
)

def __repr__(self):
# Remove the azimuth argument from the parent repr.
result = super().__repr__()
result = re.sub(r"azimuth_of_central_line=\d*\.?\d*, ", "", result)
return result
8 changes: 8 additions & 0 deletions lib/iris/fileformats/_nc_load_rules/actions.py
Original file line number Diff line number Diff line change
Expand Up @@ -156,6 +156,14 @@ def action_default(engine):
None,
hh.build_geostationary_coordinate_system,
),
hh.CF_GRID_MAPPING_OBLIQUE: (
None,
hh.build_oblique_mercator_coordinate_system,
),
hh.CF_GRID_MAPPING_ROTATED_MERCATOR: (
None,
hh.build_oblique_mercator_coordinate_system,
),
}


Expand Down
56 changes: 56 additions & 0 deletions lib/iris/fileformats/_nc_load_rules/helpers.py
Original file line number Diff line number Diff line change
Expand Up @@ -23,6 +23,7 @@
import pyproj

import iris
from iris._deprecation import warn_deprecated
import iris.aux_factory
from iris.common.mixin import _get_valid_standard_name
import iris.coord_systems
Expand Down Expand Up @@ -124,6 +125,8 @@
CF_GRID_MAPPING_TRANSVERSE = "transverse_mercator"
CF_GRID_MAPPING_VERTICAL = "vertical_perspective"
CF_GRID_MAPPING_GEOSTATIONARY = "geostationary"
CF_GRID_MAPPING_OBLIQUE = "oblique_mercator"
CF_GRID_MAPPING_ROTATED_MERCATOR = "rotated_mercator"

#
# CF Attribute Names.
Expand Down Expand Up @@ -154,6 +157,7 @@
CF_ATTR_GRID_STANDARD_PARALLEL = "standard_parallel"
CF_ATTR_GRID_PERSPECTIVE_HEIGHT = "perspective_point_height"
CF_ATTR_GRID_SWEEP_ANGLE_AXIS = "sweep_angle_axis"
CF_ATTR_GRID_AZIMUTH_CENT_LINE = "azimuth_of_central_line"
CF_ATTR_POSITIVE = "positive"
CF_ATTR_STD_NAME = "standard_name"
CF_ATTR_LONG_NAME = "long_name"
Expand Down Expand Up @@ -893,6 +897,58 @@ def build_geostationary_coordinate_system(engine, cf_grid_var):
return cs


################################################################################
def build_oblique_mercator_coordinate_system(engine, cf_grid_var):
"""
Create an oblique mercator coordinate system from the CF-netCDF
grid mapping variable.

"""
ellipsoid = _get_ellipsoid(cf_grid_var)

azimuth_of_central_line = getattr(
cf_grid_var, CF_ATTR_GRID_AZIMUTH_CENT_LINE, None
)
latitude_of_projection_origin = getattr(
cf_grid_var, CF_ATTR_GRID_LAT_OF_PROJ_ORIGIN, None
)
longitude_of_projection_origin = getattr(
cf_grid_var, CF_ATTR_GRID_LON_OF_PROJ_ORIGIN, None
)
scale_factor_at_projection_origin = getattr(
cf_grid_var, CF_ATTR_GRID_SCALE_FACTOR_AT_PROJ_ORIGIN, None
)
false_easting = getattr(cf_grid_var, CF_ATTR_GRID_FALSE_EASTING, None)
false_northing = getattr(cf_grid_var, CF_ATTR_GRID_FALSE_NORTHING, None)
kwargs = dict(
azimuth_of_central_line=azimuth_of_central_line,
latitude_of_projection_origin=latitude_of_projection_origin,
longitude_of_projection_origin=longitude_of_projection_origin,
scale_factor_at_projection_origin=scale_factor_at_projection_origin,
false_easting=false_easting,
false_northing=false_northing,
ellipsoid=ellipsoid,
)

# Handle the alternative form noted in CF: rotated mercator.
grid_mapping_name = getattr(cf_grid_var, CF_ATTR_GRID_MAPPING_NAME)
candidate_systems = dict(
oblique_mercator=iris.coord_systems.ObliqueMercator,
rotated_mercator=iris.coord_systems.RotatedMercator,
)
if grid_mapping_name == "rotated_mercator":
message = (
"Iris will stop loading the rotated_mercator grid mapping name in "
"a future release, in accordance with CF version 1.11 . Instead "
"please use oblique_mercator with azimuth_of_central_line = 90 ."
)
warn_deprecated(message)
del kwargs[CF_ATTR_GRID_AZIMUTH_CENT_LINE]

cs = candidate_systems[grid_mapping_name](**kwargs)
return cs


################################################################################
def get_attr_units(cf_var, attributes):
attr_units = getattr(cf_var, CF_ATTR_UNITS, UNKNOWN_UNIT_STRING)
Expand Down
Loading
Loading