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

Maximum value over a height coordinate #1945

Merged
merged 9 commits into from
Oct 10, 2023
Merged
69 changes: 69 additions & 0 deletions improver/cli/max_in_height.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,69 @@
#!/usr/bin/env python
# -*- coding: utf-8 -*-
# -----------------------------------------------------------------------------
# (C) British Crown copyright. The Met Office.
# All rights reserved.
#
# Redistribution and use in source and binary forms, with or without
# modification, are permitted provided that the following conditions are met:
#
# * Redistributions of source code must retain the above copyright notice, this
# list of conditions and the following disclaimer.
#
# * Redistributions in binary form must reproduce the above copyright notice,
# this list of conditions and the following disclaimer in the documentation
# and/or other materials provided with the distribution.
#
# * Neither the name of the copyright holder nor the names of its
# contributors may be used to endorse or promote products derived from
# this software without specific prior written permission.
#
# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
# ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
# LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
# CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
# SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
# INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
# CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
# ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
# POSSIBILITY OF SUCH DAMAGE.
"""Script to calculate the maximum over the height coordinate"""

from improver import cli


@cli.clizefy
@cli.with_output
def process(
cube: cli.inputcube,
*,
lower_height_bound: float = None,
upper_height_bound: float = None,
):
"""Calculate the maximum value over the height coordinate of a cube. If height bounds are
specified then the maximum value between these height levels is calculated.

Args:
cube (iris.cube.Cube):
A cube with a height coordinate.
lower_height_bound (float):
The lower bound for the height coordinate. This is either a float or None if no lower
bound is desired. Any specified bounds should have the same units as the height
coordinate of cube.
upper_height_bound (float):
The upper bound for the height coordinate. This is either a float or None if no upper
bound is desired. Any specified bounds should have the same units as the height
coordinate of cube.
Returns:
A cube of the maximum value over the height coordinate or maximum value between the provided
height bounds."""

from improver.utilities.cube_manipulation import maximum_in_height

return maximum_in_height(
cube,
lower_height_bound=lower_height_bound,
upper_height_bound=upper_height_bound,
)
58 changes: 58 additions & 0 deletions improver/utilities/cube_manipulation.py
Original file line number Diff line number Diff line change
Expand Up @@ -705,3 +705,61 @@ def add_coordinate_to_cube(
output_cube.transpose(final_dim_order)

return output_cube


def maximum_in_height(
cube: Cube, lower_height_bound: float = None, upper_height_bound: float = None
) -> Cube:
"""Calculate the maximum value over the height coordinate. If bounds are specified
then the maximum value between the lower_height_bound and upper_height_bound is calculated.

If either the upper or lower bound is None then no bound is applied. For example if no
lower bound is provided but an upper bound of 300m is provided then the maximum is
calculated for all height levels less than 300m.

Args:
cube:
A cube with a height coordinate.
lower_height_bound:
The lower bound for the height coordinate. This is either a float or None if no
lower bound is desired. Any specified bounds should have the same units as the
height coordinate of cube.
upper_height_bound:
The upper bound for the height coordinate. This is either a float or None if no
upper bound is desired. Any specified bounds should have the same units as the
height coordinate of cube.
Returns:
A cube of the maximum value over the height coordinate or maximum value between the desired
height values. This cube inherits Iris' meta-data updates to the height coordinate and to
the cell methods.

Raises:
ValueError:
If the cube has no height levels between the lower_height_bound and upper_height_bound
"""
height_levels = cube.coord("height").points

# replace None in bounds with a numerical value either below or above the range of height
# levels in the cube so it can be used as a constraint.
if lower_height_bound is None:
lower_height_bound = min(height_levels)
if upper_height_bound is None:
upper_height_bound = max(height_levels)

height_constraint = iris.Constraint(
height=lambda height: lower_height_bound <= height <= upper_height_bound
)
cube_subsetted = cube.extract(height_constraint)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think it is possible for cube_subsetted to be None at this point if certain height bounds are input, which will cause an Attribute_Error on line 750.

I think handling of this issue should be added, with an error message explaining that no heights on the input cube are within the given bounds, so the max over height levels cannot be calculated. Tests should also be added for this error.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You are right this could cause a problem. I've added in an error message and a test for this.


if cube_subsetted is None:
raise ValueError(
f"""The provided cube doesn't have any height levels between the provided bounds.
The provided bounds were {lower_height_bound},{upper_height_bound}."""
)

if len(cube_subsetted.coord("height").points) > 1:
max_cube = cube_subsetted.collapsed("height", iris.analysis.MAX)
else:
max_cube = cube_subsetted

return max_cube
3 changes: 3 additions & 0 deletions improver_tests/acceptance/SHA256SUMS
Original file line number Diff line number Diff line change
Expand Up @@ -424,6 +424,9 @@ a709db26d352457bf4f1bddf5dbade499fe89d455669e567af8965eccbebe9c4 ./manipulate-r
5bb4bb6ac2ce9ab29277275b26964a7f3633ee4c7cfa4aa4f48769e091379188 ./manipulate-reliability-table/basic/reliability_table_precip.nc
7f1093fc474320887110f28cbce1881ca68f3ed30e0fa6b2a265633e31fdd269 ./manipulate-reliability-table/point_by_point/kgo_point_by_point.nc
a162ff6c31dd3f0d84e1633c0cca28083e731876a69d40b7f14c6ff4a3431f25 ./manipulate-reliability-table/point_by_point/reliability_table_point_by_point.nc
0f35c52998ea93be4d8ed08f0ed1164b68bd6a197557afb95ae466fbec81c22e ./max-in-height/input.nc
2856432e02159afc04dadf37a0f8033cba25a9b29dfd73baea10369ccadfb645 ./max-in-height/kgo_with_bounds.nc
823232b882186fd7b7bd6172f5cedcaf5fb980b86a271553d2e4313e5e9b2b4b ./max-in-height/kgo_without_bounds.nc
a2de3ea5608d30d4ac2759e9ff4c4a6573e2c0e8e655595682a2ec2691c2de74 ./max-in-time-window/input_PT0029H00M.nc
ac81531fa507a2a3a12d4adb542ed014594eff4a38137f947d3a68a2063fab49 ./max-in-time-window/input_PT0032H00M.nc
fb02306f960fa36cf9a86921ec6599541e0a0823dfa546856000f810cdd96d73 ./max-in-time-window/kgo.nc
Expand Down
73 changes: 73 additions & 0 deletions improver_tests/acceptance/test_max_in_height.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,73 @@
# -*- coding: utf-8 -*-
# -----------------------------------------------------------------------------
# (C) British Crown copyright. The Met Office.
# All rights reserved.
#
# Redistribution and use in source and binary forms, with or without
# modification, are permitted provided that the following conditions are met:
#
# * Redistributions of source code must retain the above copyright notice, this
# list of conditions and the following disclaimer.
#
# * Redistributions in binary form must reproduce the above copyright notice,
# this list of conditions and the following disclaimer in the documentation
# and/or other materials provided with the distribution.
#
# * Neither the name of the copyright holder nor the names of its
# contributors may be used to endorse or promote products derived from
# this software without specific prior written permission.
#
# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
# ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
# LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
# CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
# SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
# INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
# CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
# ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
# POSSIBILITY OF SUCH DAMAGE.
"""Tests the max_in_height CLI"""

import pytest

from . import acceptance as acc

pytestmark = [pytest.mark.acc, acc.skip_if_kgo_missing]
CLI = acc.cli_name_with_dashes(__file__)
run_cli = acc.run_cli(CLI)


def test_with_bounds(tmp_path):
"""Test max_in_height computation with specified bounds"""

kgo_dir = acc.kgo_root() / "max-in-height"
input_path = kgo_dir / "input.nc"
output_path = tmp_path / "output.nc"
args = [
input_path,
"--upper-height-bound",
"3000",
"--lower-height-bound",
"500",
"--output",
f"{output_path}",
]

kgo_path = kgo_dir / "kgo_with_bounds.nc"
run_cli(args)
acc.compare(output_path, kgo_path)


def test_without_bounds(tmp_path):
"""Test max_in_height computation without bounds."""

kgo_dir = acc.kgo_root() / "max-in-height"
input_path = kgo_dir / "input.nc"
output_path = tmp_path / "output.nc"
args = [input_path, "--output", f"{output_path}"]

kgo_path = kgo_dir / "kgo_without_bounds.nc"
run_cli(args)
acc.compare(output_path, kgo_path)
Original file line number Diff line number Diff line change
@@ -0,0 +1,91 @@
# -*- coding: utf-8 -*-
# -----------------------------------------------------------------------------
# (C) British Crown copyright. The Met Office.
# All rights reserved.
#
# Redistribution and use in source and binary forms, with or without
# modification, are permitted provided that the following conditions are met:
#
# * Redistributions of source code must retain the above copyright notice, this
# list of conditions and the following disclaimer.
#
# * Redistributions in binary form must reproduce the above copyright notice,
# this list of conditions and the following disclaimer in the documentation
# and/or other materials provided with the distribution.
#
# * Neither the name of the copyright holder nor the names of its
# contributors may be used to endorse or promote products derived from
# this software without specific prior written permission.
#
# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
# AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
# ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
# LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
# CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
# SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
# INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
# CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
# ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
# POSSIBILITY OF SUCH DAMAGE.
"""
Unit tests for the function "cube_manipulation.maximum_in_height".
"""

import numpy as np
import pytest
from iris.cube import Cube

from improver.synthetic_data.set_up_test_cubes import set_up_variable_cube
from improver.utilities.cube_manipulation import maximum_in_height


@pytest.fixture()
def wet_bulb_temperature() -> Cube:
"Generate a cube of wet bulb temperature on height levels"
data = np.array(
[
[[100, 200, 100], [100, 200, 100]],
[[300, 400, 100], [300, 400, 100]],
[[200, 300, 300], [200, 300, 300]],
]
)
cube = set_up_variable_cube(
data=data, name="wet_bulb_temperature", height_levels=[100, 200, 300]
)
return cube


@pytest.mark.parametrize(
"lower_bound,upper_bound,expected",
(
(None, None, [300, 400, 300]),
(None, 200, [300, 400, 100]),
(250, None, [200, 300, 300]),
(50, 1000, [300, 400, 300]),
),
)
def test_maximum_in_height(lower_bound, upper_bound, expected, wet_bulb_temperature):
"""Test that the maximum over the height coordinate is correctly calculated for
different combinations of upper and lower bounds."""

result = maximum_in_height(
wet_bulb_temperature,
lower_height_bound=lower_bound,
upper_height_bound=upper_bound,
)

assert np.allclose(result.data, [expected] * 2)
assert "wet_bulb_temperature" == result.name()


def test_height_bounds_error(wet_bulb_temperature):
"""Test an error is raised if the input cube doesn't have any height levels
between the height bounds."""

with pytest.raises(
ValueError, match="any height levels between the provided bounds"
):
maximum_in_height(
wet_bulb_temperature, lower_height_bound=50, upper_height_bound=75
)