From 521f3a6b4897b56fe697e8b15abd349907b39c09 Mon Sep 17 00:00:00 2001
From: Manuel Schlund <manuel.schlund@dlr.de>
Date: Tue, 24 Mar 2026 09:51:41 +0100
Subject: [PATCH 01/11] Add derivation script for phcint_total

---
 esmvalcore/cmor/_fixes/icon/icon_xpp.py       |  4 ++
 .../tables/cmip6-custom/CMIP6_custom.json     | 12 ++++
 .../preprocessor/_derive/phcint_total.py      | 52 ++++++++++++++
 esmvalcore/preprocessor/_volume.py            | 28 ++------
 .../preprocessor/_derive/test_phcint_total.py | 70 +++++++++++++++++++
 5 files changed, 145 insertions(+), 21 deletions(-)
 create mode 100644 esmvalcore/preprocessor/_derive/phcint_total.py
 create mode 100644 tests/unit/preprocessor/_derive/test_phcint_total.py

diff --git a/esmvalcore/cmor/_fixes/icon/icon_xpp.py b/esmvalcore/cmor/_fixes/icon/icon_xpp.py
index d26e3a15af..5369149096 100644
--- a/esmvalcore/cmor/_fixes/icon/icon_xpp.py
+++ b/esmvalcore/cmor/_fixes/icon/icon_xpp.py
@@ -16,6 +16,10 @@ class AllVars(AllVarsBase):
     DEFAULT_PFULL_VAR_NAME = "pres"
 
 
+class Msftmz(IconFix):
+    """Fixes for ``msftmz``."""
+
+
 class Clwvi(IconFix):
     """Fixes for ``clwvi``."""
 
diff --git a/esmvalcore/cmor/tables/cmip6-custom/CMIP6_custom.json b/esmvalcore/cmor/tables/cmip6-custom/CMIP6_custom.json
index d399838a75..16ade07200 100644
--- a/esmvalcore/cmor/tables/cmip6-custom/CMIP6_custom.json
+++ b/esmvalcore/cmor/tables/cmip6-custom/CMIP6_custom.json
@@ -767,6 +767,18 @@
             "type": "real",
             "units": "kg s-1"
         },
+        "phcint_total": {
+            "cell_measures": "area: areacello",
+            "cell_methods": "area: time: mean where sea",
+            "comment": "This is the total column vertically-integrated heat content derived from potential temperature (thetao).",
+            "dimensions": "longitude latitude time",
+            "long_name": "Total Column Integrated Ocean Heat Content from Potential Temperature",
+            "modeling_realm": "ocean",
+            "out_name": "phcint_total",
+            "standard_name": "integral_wrt_depth_of_sea_water_potential_temperature_expressed_as_heat_content",
+            "type": "real",
+            "units": "J m-2"
+        },
         "ptype": {
             "cell_measures": "area: areacella",
             "cell_methods": "time: mean",
diff --git a/esmvalcore/preprocessor/_derive/phcint_total.py b/esmvalcore/preprocessor/_derive/phcint_total.py
new file mode 100644
index 0000000000..ac6f7c5598
--- /dev/null
+++ b/esmvalcore/preprocessor/_derive/phcint_total.py
@@ -0,0 +1,52 @@
+"""Derivation of variable ``phcint_total``."""
+
+from __future__ import annotations
+
+from typing import TYPE_CHECKING
+
+import iris
+from cf_units import Unit
+from iris import NameConstraint
+
+from esmvalcore.preprocessor._volume import depth_integration
+
+from ._baseclass import DerivedVariableBase
+
+if TYPE_CHECKING:
+    from iris.cube import Cube, CubeList
+
+    from esmvalcore.typing import Facets
+
+RHO_CP = iris.coords.AuxCoord(4.09169e6, units=Unit("kg m-3 J kg-1 K-1"))
+
+
+class DerivedVariable(DerivedVariableBase):
+    """Derivation of variable `ohc`."""
+
+    @staticmethod
+    def required(project: str) -> list[Facets]:  # noqa: ARG004
+        """Declare the variables needed for derivation."""
+        return [{"short_name": "thetao"}]
+
+    @staticmethod
+    def calculate(cubes: CubeList) -> Cube:
+        """Compute total column vertically-integrated heat content.
+
+        Use c_p * rho_0 = 4.09169e+6 J m-3 K-1 (Kuhlbrodt et al., 2015, Clim.
+        Dyn.)
+
+        Arguments
+        ---------
+        cubes:
+           Input cubes.
+
+        Returns
+        -------
+        :
+            Output cube.
+
+        """
+        cube = cubes.extract_cube(NameConstraint(var_name="thetao"))
+        cube.convert_units("K")
+        cube = cube * RHO_CP
+        return depth_integration(cube)
diff --git a/esmvalcore/preprocessor/_volume.py b/esmvalcore/preprocessor/_volume.py
index aefd15c8c1..e7d92f3a32 100644
--- a/esmvalcore/preprocessor/_volume.py
+++ b/esmvalcore/preprocessor/_volume.py
@@ -110,9 +110,7 @@ def extract_volume(
             'Depth extraction bounds can be set to "open", "closed", '
             f'"left_closed", or "right_closed". Got "{interval_bounds}".'
         )
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
 
     z_constraint = iris.Constraint(coord_values=coord_values)
 
@@ -162,18 +160,14 @@ def calculate_volume(cube: Cube) -> np.ndarray | da.Array:
             "Bounds should be 2 in the last dimension to compute the "
             "thickness."
         )
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
 
     # Convert units to get the thickness in meters
     try:
         depth.convert_units("m")
     except ValueError as err:
         msg = f"Cannot compute volume using the Z-axis. {err}"
-        raise ValueError(
-            msg,
-        ) from err
+        raise ValueError(msg) from err
 
     # Calculate Z-direction thickness
     thickness = depth.core_bounds()[..., 1] - depth.core_bounds()[..., 0]
@@ -314,9 +308,7 @@ def volume_statistics(
             "This may indicate Z axis depending on other dimension than "
             "space that could provoke invalid aggregation..."
         )
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
 
     (agg, agg_kwargs) = get_iris_aggregator(operator, **operator_kwargs)
     agg_kwargs = update_weights_kwargs(
@@ -391,9 +383,7 @@ def axis_statistics(
         coord = cube.coord(axis=axis)
     except iris.exceptions.CoordinateNotFoundError as err:
         msg = f"Axis {axis} not found in cube {cube.summary(shorten=True)}"
-        raise ValueError(
-            msg,
-        ) from err
+        raise ValueError(msg) from err
 
     # Multidimensional coordinates are currently not supported
     coord_dims = cube.coord_dims(coord)
@@ -401,9 +391,7 @@ def axis_statistics(
         msg = (
             "axis_statistics not implemented for multidimensional coordinates."
         )
-        raise NotImplementedError(
-            msg,
-        )
+        raise NotImplementedError(msg)
 
     # For weighted operations, create a dummy weights coordinate using the
     # bounds of the original coordinate (this handles units properly, e.g., for
@@ -621,9 +609,7 @@ def extract_trajectory(
     """
     if len(latitudes) != len(longitudes):
         msg = "Longitude & Latitude coordinates have different lengths"
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
 
     if len(latitudes) == len(longitudes) == 2:
         minlat, maxlat = np.min(latitudes), np.max(latitudes)
diff --git a/tests/unit/preprocessor/_derive/test_phcint_total.py b/tests/unit/preprocessor/_derive/test_phcint_total.py
new file mode 100644
index 0000000000..5b78861181
--- /dev/null
+++ b/tests/unit/preprocessor/_derive/test_phcint_total.py
@@ -0,0 +1,70 @@
+"""Test derivation of ``phcint_total``."""
+
+from __future__ import annotations
+
+from typing import TYPE_CHECKING
+
+import numpy as np
+import pytest
+from iris.coords import DimCoord
+from iris.cube import CubeList
+
+from esmvalcore.preprocessor._derive import derive, get_required
+
+if TYPE_CHECKING:
+    from iris.cube import Cube
+
+
+@pytest.fixture
+def cubes(realistic_4d_cube: Cube) -> CubeList:
+    depth_coord = DimCoord(
+        [500.0],
+        bounds=[[0.0, 1000.0]],
+        standard_name="depth",
+        units="m",
+        attributes={"positive": "down"},
+    )
+    realistic_4d_cube.remove_coord("air_pressure")
+    realistic_4d_cube.add_dim_coord(depth_coord, 1)
+    realistic_4d_cube.var_name = "thetao"
+    return CubeList([realistic_4d_cube])
+
+
+@pytest.mark.parametrize("project", ["CMIP3", "CMIP5", "CMIP6", "CMIP7"])
+def test_get_required(project: str) -> None:
+    assert get_required("phcint_total", project) == [{"short_name": "thetao"}]
+
+
+def test_derive(cubes: CubeList) -> None:
+    short_name = "phcint_total"
+    long_name = (
+        "Total Column Integrated Ocean Heat Content from Potential Temperature"
+    )
+    units = "J m-2"
+    standard_name = "integral_wrt_depth_of_sea_water_potential_temperature_expressed_as_heat_content"
+    derived_cube = derive(
+        cubes,
+        short_name=short_name,
+        long_name=long_name,
+        units=units,
+        standard_name=standard_name,
+    )
+
+    assert derived_cube.standard_name == standard_name
+    assert derived_cube.long_name == long_name
+    assert derived_cube.var_name == short_name
+    assert derived_cube.units == units
+
+    expected_data = np.ma.masked_invalid(
+        [
+            [
+                [0.0, np.nan, np.nan],
+                [np.nan, 16366760000.0, 20458450000.0],
+            ],
+            [
+                [24550140000.0, 28641830000.0, 32733520000.0],
+                [36825210000.0, 40916900000.0, 45008590000.0],
+            ],
+        ],
+    )
+    np.testing.assert_allclose(derived_cube.data, expected_data)

From 2d5c1d85dff9c5d9512f803f814524272ed2d184 Mon Sep 17 00:00:00 2001
From: Manuel Schlund <manuel.schlund@dlr.de>
Date: Tue, 24 Mar 2026 14:01:22 +0100
Subject: [PATCH 02/11] Remove unnecessary line breaks

---
 esmvalcore/cmor/_fixes/native_datasets.py     |  8 ++---
 esmvalcore/preprocessor/_mask.py              | 20 +++---------
 esmvalcore/preprocessor/_shared.py            | 32 +++++--------------
 .../preprocessor/_supplementary_vars.py       |  4 +--
 esmvalcore/preprocessor/_weighting.py         |  8 ++---
 5 files changed, 18 insertions(+), 54 deletions(-)

diff --git a/esmvalcore/cmor/_fixes/native_datasets.py b/esmvalcore/cmor/_fixes/native_datasets.py
index 2f7f24b763..a4974c80c3 100644
--- a/esmvalcore/cmor/_fixes/native_datasets.py
+++ b/esmvalcore/cmor/_fixes/native_datasets.py
@@ -86,9 +86,7 @@ def fix_var_metadata(self, cube: Cube) -> None:
                     f"Failed to fix invalid units '{invalid_units}' for "
                     f"variable '{self.vardef.short_name}'"
                 )
-                raise ValueError(
-                    msg,
-                ) from exc
+                raise ValueError(msg) from exc
         safe_convert_units(cube, self.vardef.units)
 
         # Fix attributes
@@ -132,9 +130,7 @@ def get_cube(
                 f"Variable '{var_name}' used to extract "
                 f"'{self.vardef.short_name}' is not available in input file"
             )
-            raise ValueError(
-                msg,
-            )
+            raise ValueError(msg)
         return cubes.extract_cube(NameConstraint(var_name=var_name))
 
     def fix_regular_time(
diff --git a/esmvalcore/preprocessor/_mask.py b/esmvalcore/preprocessor/_mask.py
index 790a58b3bd..916a724f5e 100644
--- a/esmvalcore/preprocessor/_mask.py
+++ b/esmvalcore/preprocessor/_mask.py
@@ -153,9 +153,7 @@ def mask_landsea(cube: Cube, mask_out: Literal["land", "sea"]) -> Cube:
             "Use of shapefiles with irregular grids not yet implemented, "
             "land-sea mask not applied."
         )
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
 
     return cube
 
@@ -528,9 +526,7 @@ def _get_shape(cubes):
     shapes = {cube.shape for cube in cubes}
     if len(shapes) > 1:
         msg = f"Expected cubes with identical shapes, got shapes {shapes}"
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
     return next(iter(shapes))
 
 
@@ -616,9 +612,7 @@ def mask_multimodel(products):
         f"iris.cube.Cube or esmvalcore.preprocessor.PreprocessorFile, "
         f"got {product_types}"
     )
-    raise TypeError(
-        msg,
-    )
+    raise TypeError(msg)
 
 
 def mask_fillvalues(
@@ -685,9 +679,7 @@ def mask_fillvalues(
                 valid = ~mask.all(axis=(-2, -1), keepdims=True)
             else:
                 msg = f"Unable to handle {mask.ndim} dimensional data"
-                raise NotImplementedError(
-                    msg,
-                )
+                raise NotImplementedError(msg)
             combined_mask = array_module.where(
                 valid,
                 combined_mask | mask,
@@ -730,9 +722,7 @@ def _get_fillvalues_mask(
             f"Fraction of missing values {threshold_fraction} should be "
             f"between 0 and 1.0"
         )
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
     nr_time_points = len(cube.coord("time").points)
     if time_window > nr_time_points:
         msg = "Time window (in time units) larger than total time span. Stop."
diff --git a/esmvalcore/preprocessor/_shared.py b/esmvalcore/preprocessor/_shared.py
index ecf58028e2..f5709b3653 100644
--- a/esmvalcore/preprocessor/_shared.py
+++ b/esmvalcore/preprocessor/_shared.py
@@ -173,9 +173,7 @@ def update_weights_kwargs(
     kwargs = dict(kwargs)
     if not aggregator_accept_weights(aggregator) and "weights" in kwargs:
         msg = f"Aggregator '{operator}' does not support 'weights' option"
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
     if aggregator_accept_weights(aggregator) and kwargs.get("weights", True):
         kwargs["weights"] = weights
         if cube is not None and callback is not None:
@@ -232,9 +230,7 @@ def get_normalized_cube(
             f"Expected 'subtract' or 'divide' for `normalize`, got "
             f"'{normalize}'"
         )
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
 
     # Keep old metadata except for units
     new_units = normalized_cube.units
@@ -274,9 +270,7 @@ def preserve_float_dtype(func: Callable) -> Callable:
             f"Cannot preserve float dtype during function '{func.__name__}', "
             f"function takes no arguments"
         )
-        raise TypeError(
-            msg,
-        )
+        raise TypeError(msg)
 
     @wraps(func)
     def wrapper(*args: Any, **kwargs: Any) -> DataType:
@@ -298,9 +292,7 @@ def wrapper(*args: Any, **kwargs: Any) -> DataType:
                 f"type {type(result)} do not have the necessary attribute "
                 f"'dtype'"
             )
-            raise TypeError(
-                msg,
-            )
+            raise TypeError(msg)
 
         return result
 
@@ -381,9 +373,7 @@ def get_weights(
                 f"`cell_area` can be given to the cube as supplementary "
                 f"variable)"
             )
-            raise CoordinateNotFoundError(
-                msg,
-            )
+            raise CoordinateNotFoundError(msg)
         try_adding_calculated_cell_area(cube)
         area_weights = cube.cell_measure("cell_area").core_data()
         if cube.has_lazy_data():
@@ -438,18 +428,14 @@ def get_coord_weights(
             f"Cannot calculate weights for coordinate '{coord.name()}' "
             f"without bounds"
         )
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
     if coord.core_bounds().shape[-1] != 2:
         msg = (
             f"Cannot calculate weights for coordinate '{coord.name()}' "
             f"with {coord.core_bounds().shape[-1]} bounds per point, expected "
             f"2 bounds per point"
         )
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
 
     # Calculate weights of same shape as coordinate and make sure to use
     # identical chunks as parent cube for non-scalar lazy data
@@ -573,9 +559,7 @@ def get_all_coords(
                 f"{cube.summary(shorten=True)} must not have unnamed "
                 f"dimensions"
             )
-            raise ValueError(
-                msg,
-            )
+            raise ValueError(msg)
     return coords
 
 
diff --git a/esmvalcore/preprocessor/_supplementary_vars.py b/esmvalcore/preprocessor/_supplementary_vars.py
index f7a1491d09..f1ed2915ff 100644
--- a/esmvalcore/preprocessor/_supplementary_vars.py
+++ b/esmvalcore/preprocessor/_supplementary_vars.py
@@ -82,9 +82,7 @@ def add_cell_measure(
     """
     if measure not in ["area", "volume"]:
         msg = f"measure name must be 'area' or 'volume', got {measure} instead"
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
     coord_dims = tuple(
         range(cube.ndim - len(cell_measure_cube.shape), cube.ndim),
     )
diff --git a/esmvalcore/preprocessor/_weighting.py b/esmvalcore/preprocessor/_weighting.py
index 9089e1d62e..dc42727616 100644
--- a/esmvalcore/preprocessor/_weighting.py
+++ b/esmvalcore/preprocessor/_weighting.py
@@ -73,18 +73,14 @@ def weighting_landsea_fraction(cube, area_type):
     """
     if area_type not in ("land", "sea"):
         msg = f"Expected 'land' or 'sea' for area_type, got '{area_type}'"
-        raise TypeError(
-            msg,
-        )
+        raise TypeError(msg)
     (land_fraction, errors) = _get_land_fraction(cube)
     if land_fraction is None:
         msg = (
             f"Weighting of '{cube.var_name}' with '{area_type}' fraction "
             f"failed because of the following errors: {' '.join(errors)}"
         )
-        raise ValueError(
-            msg,
-        )
+        raise ValueError(msg)
     core_data = cube.core_data()
     if area_type == "land":
         cube.data = core_data * land_fraction

From eee2cd187b8d067b2bfeb81f1a8178581cb74adf Mon Sep 17 00:00:00 2001
From: Manuel Schlund <manuel.schlund@dlr.de>
Date: Tue, 24 Mar 2026 14:01:58 +0100
Subject: [PATCH 03/11] Fix phcint derivation (cell measures were dropped
 before)

---
 esmvalcore/preprocessor/_derive/phcint_total.py | 9 ++++++---
 1 file changed, 6 insertions(+), 3 deletions(-)

diff --git a/esmvalcore/preprocessor/_derive/phcint_total.py b/esmvalcore/preprocessor/_derive/phcint_total.py
index ac6f7c5598..6459abe7c9 100644
--- a/esmvalcore/preprocessor/_derive/phcint_total.py
+++ b/esmvalcore/preprocessor/_derive/phcint_total.py
@@ -4,7 +4,6 @@
 
 from typing import TYPE_CHECKING
 
-import iris
 from cf_units import Unit
 from iris import NameConstraint
 
@@ -17,7 +16,8 @@
 
     from esmvalcore.typing import Facets
 
-RHO_CP = iris.coords.AuxCoord(4.09169e6, units=Unit("kg m-3 J kg-1 K-1"))
+RHO_CP = 4.09169e6
+RHO_CP_UNIT = Unit("kg m-3 J kg-1 K-1")
 
 
 class DerivedVariable(DerivedVariableBase):
@@ -48,5 +48,8 @@ def calculate(cubes: CubeList) -> Cube:
         """
         cube = cubes.extract_cube(NameConstraint(var_name="thetao"))
         cube.convert_units("K")
-        cube = cube * RHO_CP
+        cube.data = (
+            cube.core_data() * RHO_CP
+        )  # https://github.com/SciTools/iris/issues/6990
+        cube.units *= RHO_CP_UNIT
         return depth_integration(cube)

From 548405b777cd8a1074da3b062a6278e47f169d25 Mon Sep 17 00:00:00 2001
From: Manuel Schlund <manuel.schlund@dlr.de>
Date: Tue, 24 Mar 2026 15:01:41 +0100
Subject: [PATCH 04/11] Add ICON-XPP fix for msftmz

---
 esmvalcore/cmor/_fixes/icon/icon_xpp.py       | 50 +++++++++++++--
 .../cmor/_fixes/icon/test_icon_xpp.py         | 62 ++++++++++++++++++-
 2 files changed, 107 insertions(+), 5 deletions(-)

diff --git a/esmvalcore/cmor/_fixes/icon/icon_xpp.py b/esmvalcore/cmor/_fixes/icon/icon_xpp.py
index 5369149096..baeb9bee8d 100644
--- a/esmvalcore/cmor/_fixes/icon/icon_xpp.py
+++ b/esmvalcore/cmor/_fixes/icon/icon_xpp.py
@@ -2,6 +2,7 @@
 
 import logging
 
+from iris.coords import AuxCoord
 from iris.cube import CubeList
 from scipy import constants
 
@@ -16,10 +17,6 @@ class AllVars(AllVarsBase):
     DEFAULT_PFULL_VAR_NAME = "pres"
 
 
-class Msftmz(IconFix):
-    """Fixes for ``msftmz``."""
-
-
 class Clwvi(IconFix):
     """Fixes for ``clwvi``."""
 
@@ -60,6 +57,51 @@ def fix_metadata(self, cubes: CubeList) -> CubeList:
 Hfss = NegateData
 
 
+class Msftmz(IconFix):
+    """Fixes for ``msftmz``."""
+
+    def fix_metadata(self, cubes: CubeList) -> CubeList:
+        """Fix metadata."""
+        preprocessed_cubes = CubeList([])
+        basin_coord = AuxCoord(
+            "placeholder",
+            standard_name="region",
+            long_name="ocean basin",
+            var_name="basin",
+        )
+        var_names = {
+            "atlantic_moc": "atlantic_arctic_ocean",
+            "pacific_moc": "indian_pacific_ocean",
+            "global_moc": "global_ocean",
+        }
+        for var_name, basin in var_names.items():
+            cube = self.get_cube(cubes, var_name=var_name)
+            cube.var_name = "msftmz"
+            cube.long_name = None
+            cube.attributes.locals = {}
+
+            # Remove longitude coordinate (with length 1)
+            cube = cube[..., 0]
+            cube.remove_coord("longitude")
+
+            # Add scalar basin coordinate
+            cube.add_aux_coord(basin_coord.copy(basin), ())
+            preprocessed_cubes.append(cube)
+
+        msftmz_cube = preprocessed_cubes.merge_cube()
+
+        # Swap time and basin coordinates
+        msftmz_cube.transpose([1, 0, 2, 3])
+
+        # By default, merge_cube() sorts the coordinate alphabetically (i.e.,
+        # atlantic_arctic_ocean -> global_ocean -> indian_pacific_ocean). Thus,
+        # we need to restore the desired order (atlantic_arctic_ocean ->
+        # indian_pacific_ocean -> global_ocean).
+        msftmz_cube = msftmz_cube[:, [0, 2, 1], ...]
+
+        return CubeList([msftmz_cube])
+
+
 Rlut = NegateData
 
 
diff --git a/tests/integration/cmor/_fixes/icon/test_icon_xpp.py b/tests/integration/cmor/_fixes/icon/test_icon_xpp.py
index aa50be5ebb..369bd161d7 100644
--- a/tests/integration/cmor/_fixes/icon/test_icon_xpp.py
+++ b/tests/integration/cmor/_fixes/icon/test_icon_xpp.py
@@ -5,8 +5,9 @@
 import pytest
 from cf_units import Unit
 from iris import NameConstraint
-from iris.coords import DimCoord
+from iris.coords import AuxCoord, DimCoord
 from iris.cube import Cube, CubeList
+from iris.util import new_axis
 
 import esmvalcore.cmor._fixes.icon.icon_xpp
 from esmvalcore.cmor._fixes.fix import GenericFix
@@ -18,6 +19,7 @@
     Gpp,
     Hfls,
     Hfss,
+    Msftmz,
     Rlut,
     Rlutcs,
     Rsutcs,
@@ -704,6 +706,64 @@ def test_hfss_fix(cubes_regular_grid):
     np.testing.assert_allclose(fixed_cube.data, [[[0.0, -1.0], [-2.0, -3.0]]])
 
 
+# Test msftmz (for extra fix)
+
+
+def test_get_msftmz_fix():
+    """Test getting of fix."""
+    fix = Fix.get_fixes("ICON", "ICON-XPP", "Omon", "msftmz")
+    assert fix == [Msftmz(None), AllVars(None), GenericFix(None)]
+
+
+def test_msftmz_fix(cubes_regular_grid):
+    """Test fix."""
+    depth_coord = AuxCoord(
+        10.0,
+        standard_name="depth",
+        long_name="depth below sea",
+        units="m",
+        attributes={"positive": "down"},
+    )
+    cube = cubes_regular_grid[0][..., [0]]
+    cube.coord("latitude").var_name = "lat"
+    cube.add_aux_coord(depth_coord, ())
+    cube = new_axis(cube, "depth")
+    cube.transpose([1, 0, 2, 3])
+    cubes = CubeList([cube.copy() * 0.0, cube.copy() * 1.0, cube.copy() * 2.0])
+    cubes[0].var_name = "atlantic_moc"
+    cubes[0].units = "kg s-1"
+    cubes[1].var_name = "pacific_moc"
+    cubes[1].units = "kg s-1"
+    cubes[2].var_name = "global_moc"
+    cubes[2].units = "kg s-1"
+
+    fixed_cubes = fix_metadata(cubes, "Omon", "msftmz")
+
+    assert len(fixed_cubes) == 1
+    cube = fixed_cubes[0]
+    assert cube.var_name == "msftmz"
+    assert (
+        cube.standard_name
+        == "ocean_meridional_overturning_mass_streamfunction"
+    )
+    assert cube.long_name == "Ocean Meridional Overturning Mass Streamfunction"
+
+    assert cube.units == "kg s-1"
+    assert "positive" not in cube.attributes
+    assert "invalid_units" not in cube.attributes
+
+    np.testing.assert_equal(
+        cube.coord("region").points,
+        ["atlantic_arctic_ocean", "indian_pacific_ocean", "global_ocean"],
+    )
+
+    assert cube.shape == (1, 3, 1, 2)
+    np.testing.assert_allclose(
+        cube.data,
+        [[[[0.0, 0.0]], [[0.0, 2.0]], [[0.0, 4.0]]]],
+    )
+
+
 # Test rlut (for extra fix)
 
 

From 3c1044446e4c5de8fd139466ed9e1da35501dcf7 Mon Sep 17 00:00:00 2001
From: Manuel Schlund <manuel.schlund@dlr.de>
Date: Tue, 24 Mar 2026 15:01:53 +0100
Subject: [PATCH 05/11] Expand ICON-XPP extra facets

---
 esmvalcore/config/configurations/defaults/extra_facets_icon.yml | 2 ++
 1 file changed, 2 insertions(+)

diff --git a/esmvalcore/config/configurations/defaults/extra_facets_icon.yml b/esmvalcore/config/configurations/defaults/extra_facets_icon.yml
index ab8418d5bb..4c51fe38c9 100644
--- a/esmvalcore/config/configurations/defaults/extra_facets_icon.yml
+++ b/esmvalcore/config/configurations/defaults/extra_facets_icon.yml
@@ -135,6 +135,7 @@ projects:
           # 2D ocean variables
           hfds: {raw_name: HeatFlux_Total, var_type: oce_dbg}
           mlotst: {raw_name: mld, var_type: oce_dbg}
+          phcint_total: {raw_name: heat_content_total, var_type: oce_def}
           siconc: {raw_name: conc, var_type: oce_ice, raw_units: '1'}
           siconca: {raw_name: fr_seaice, var_type: atm_2d_ml}
           sithick: {raw_name: hi, var_type: oce_ice}
@@ -143,6 +144,7 @@ projects:
           zos: {raw_name: zos, var_type: oce_dbg}
 
           # 3D ocean variables
+          msftmz: {var_type: oce_moc, lat_var: lat}
           so: {var_type: oce_def, raw_units: "0.001"}
           thetao: {raw_name: to, var_type: oce_def, raw_units: degC}
           uo: {raw_name: u, var_type: oce_def}

From 382fadeb5c443b73d9fecf9cdd24cde28beef833 Mon Sep 17 00:00:00 2001
From: Manuel Schlund <manuel.schlund@dlr.de>
Date: Tue, 24 Mar 2026 15:45:26 +0100
Subject: [PATCH 06/11] Support amoc for ICON data

---
 .../config/configurations/defaults/extra_facets_icon.yml       | 3 +++
 esmvalcore/preprocessor/_derive/amoc.py                        | 2 ++
 2 files changed, 5 insertions(+)

diff --git a/esmvalcore/config/configurations/defaults/extra_facets_icon.yml b/esmvalcore/config/configurations/defaults/extra_facets_icon.yml
index 4c51fe38c9..5768293697 100644
--- a/esmvalcore/config/configurations/defaults/extra_facets_icon.yml
+++ b/esmvalcore/config/configurations/defaults/extra_facets_icon.yml
@@ -132,6 +132,9 @@ projects:
           gpp: {raw_name: assimi_gross_assimilation_box, var_type: jsb_2d_ml}
           lai: {raw_name: pheno_lai_box, var_type: jsb_2d_ml, raw_units: '1'}
 
+          # 1D ocean variables
+          amoc: {var_type: oce_moc}
+
           # 2D ocean variables
           hfds: {raw_name: HeatFlux_Total, var_type: oce_dbg}
           mlotst: {raw_name: mld, var_type: oce_dbg}
diff --git a/esmvalcore/preprocessor/_derive/amoc.py b/esmvalcore/preprocessor/_derive/amoc.py
index 3607aa1d62..e6f8fe3387 100644
--- a/esmvalcore/preprocessor/_derive/amoc.py
+++ b/esmvalcore/preprocessor/_derive/amoc.py
@@ -21,6 +21,8 @@ def required(project):
                 {"short_name": "msftmz", "optional": True},
                 {"short_name": "msftyz", "optional": True},
             ]
+        elif project == "ICON":
+            required = [{"short_name": "msftmz", "mip": "Omon"}]
         else:
             msg = f"Project {project} can not be used for Amoc derivation."
             raise ValueError(msg)

From 6496a33b0e8f9832272b5a74f8c4658525828a3c Mon Sep 17 00:00:00 2001
From: Manuel Schlund <manuel.schlund@dlr.de>
Date: Tue, 24 Mar 2026 15:51:13 +0100
Subject: [PATCH 07/11] Add variable phcint_total to custom CMIP5 tables

---
 .../tables/cmip5-custom/CMOR_phcint_total.dat | 21 +++++++++++++++++++
 1 file changed, 21 insertions(+)
 create mode 100644 esmvalcore/cmor/tables/cmip5-custom/CMOR_phcint_total.dat

diff --git a/esmvalcore/cmor/tables/cmip5-custom/CMOR_phcint_total.dat b/esmvalcore/cmor/tables/cmip5-custom/CMOR_phcint_total.dat
new file mode 100644
index 0000000000..1544900a12
--- /dev/null
+++ b/esmvalcore/cmor/tables/cmip5-custom/CMOR_phcint_total.dat
@@ -0,0 +1,21 @@
+SOURCE: CMIP7
+!============
+variable_entry:    phcint_total
+!============
+modeling_realm:    ocean
+!----------------------------------
+! Variable attributes:
+!----------------------------------
+standard_name:     integral_wrt_depth_of_sea_water_potential_temperature_expressed_as_heat_content
+units:             J m-2
+cell_methods:      area: time: mean where sea
+cell_measures:     area: areacello
+long_name:         Total Column Integrated Ocean Heat Content from Potential Temperature
+comment:           This is the total column vertically-integrated heat content derived from potential temperature (thetao).
+!----------------------------------
+! Additional variable information:
+!----------------------------------
+dimensions:        longitude latitude time
+type:              real
+!----------------------------------
+!

From 3b265e7b065461af1bb3b98211b1ff8f2284062b Mon Sep 17 00:00:00 2001
From: Manuel Schlund <manuel.schlund@dlr.de>
Date: Wed, 25 Mar 2026 10:27:17 +0100
Subject: [PATCH 08/11] phcint_total -> phcint

---
 ...{CMOR_phcint_total.dat => CMOR_phcint.dat} |  9 ++++---
 .../tables/cmip6-custom/CMIP6_custom.json     | 10 ++++----
 .../cmor/tables/cmip7/tables/CMIP7_ocean.json |  2 +-
 .../defaults/extra_facets_icon.yml            |  2 +-
 .../_derive/{phcint_total.py => phcint.py}    | 16 +++++++++---
 .../{test_phcint_total.py => test_phcint.py}  | 25 +++++++++++--------
 6 files changed, 38 insertions(+), 26 deletions(-)
 rename esmvalcore/cmor/tables/cmip5-custom/{CMOR_phcint_total.dat => CMOR_phcint.dat} (64%)
 rename esmvalcore/preprocessor/_derive/{phcint_total.py => phcint.py} (71%)
 rename tests/unit/preprocessor/_derive/{test_phcint_total.py => test_phcint.py} (73%)

diff --git a/esmvalcore/cmor/tables/cmip5-custom/CMOR_phcint_total.dat b/esmvalcore/cmor/tables/cmip5-custom/CMOR_phcint.dat
similarity index 64%
rename from esmvalcore/cmor/tables/cmip5-custom/CMOR_phcint_total.dat
rename to esmvalcore/cmor/tables/cmip5-custom/CMOR_phcint.dat
index 1544900a12..b28c63cb6a 100644
--- a/esmvalcore/cmor/tables/cmip5-custom/CMOR_phcint_total.dat
+++ b/esmvalcore/cmor/tables/cmip5-custom/CMOR_phcint.dat
@@ -1,6 +1,7 @@
 SOURCE: CMIP7
+generic_levels:    olevel
 !============
-variable_entry:    phcint_total
+variable_entry:    phcint
 !============
 modeling_realm:    ocean
 !----------------------------------
@@ -10,12 +11,12 @@ standard_name:     integral_wrt_depth_of_sea_water_potential_temperature_express
 units:             J m-2
 cell_methods:      area: time: mean where sea
 cell_measures:     area: areacello
-long_name:         Total Column Integrated Ocean Heat Content from Potential Temperature
-comment:           This is the total column vertically-integrated heat content derived from potential temperature (thetao).
+long_name:         Integrated Ocean Heat Content from Potential Temperature
+comment:           This is the vertically-integrated heat content derived from potential temperature (thetao).
 !----------------------------------
 ! Additional variable information:
 !----------------------------------
-dimensions:        longitude latitude time
+dimensions:        longitude latitude time olevel
 type:              real
 !----------------------------------
 !
diff --git a/esmvalcore/cmor/tables/cmip6-custom/CMIP6_custom.json b/esmvalcore/cmor/tables/cmip6-custom/CMIP6_custom.json
index 16ade07200..fa0b8bdd8c 100644
--- a/esmvalcore/cmor/tables/cmip6-custom/CMIP6_custom.json
+++ b/esmvalcore/cmor/tables/cmip6-custom/CMIP6_custom.json
@@ -767,14 +767,14 @@
             "type": "real",
             "units": "kg s-1"
         },
-        "phcint_total": {
+        "phcint": {
             "cell_measures": "area: areacello",
             "cell_methods": "area: time: mean where sea",
-            "comment": "This is the total column vertically-integrated heat content derived from potential temperature (thetao).",
-            "dimensions": "longitude latitude time",
-            "long_name": "Total Column Integrated Ocean Heat Content from Potential Temperature",
+            "comment": "This is the vertically-integrated heat content derived from potential temperature (thetao).",
+            "dimensions": "longitude latitude time olevel",
+            "long_name": "Integrated Ocean Heat Content from Potential Temperature",
             "modeling_realm": "ocean",
-            "out_name": "phcint_total",
+            "out_name": "phcint",
             "standard_name": "integral_wrt_depth_of_sea_water_potential_temperature_expressed_as_heat_content",
             "type": "real",
             "units": "J m-2"
diff --git a/esmvalcore/cmor/tables/cmip7/tables/CMIP7_ocean.json b/esmvalcore/cmor/tables/cmip7/tables/CMIP7_ocean.json
index bd6259b6d4..bead80ef17 100644
--- a/esmvalcore/cmor/tables/cmip7/tables/CMIP7_ocean.json
+++ b/esmvalcore/cmor/tables/cmip7/tables/CMIP7_ocean.json
@@ -3320,4 +3320,4 @@
             "units": "m"
         }
     }
-}
\ No newline at end of file
+}
diff --git a/esmvalcore/config/configurations/defaults/extra_facets_icon.yml b/esmvalcore/config/configurations/defaults/extra_facets_icon.yml
index 5768293697..bc346b096f 100644
--- a/esmvalcore/config/configurations/defaults/extra_facets_icon.yml
+++ b/esmvalcore/config/configurations/defaults/extra_facets_icon.yml
@@ -138,7 +138,6 @@ projects:
           # 2D ocean variables
           hfds: {raw_name: HeatFlux_Total, var_type: oce_dbg}
           mlotst: {raw_name: mld, var_type: oce_dbg}
-          phcint_total: {raw_name: heat_content_total, var_type: oce_def}
           siconc: {raw_name: conc, var_type: oce_ice, raw_units: '1'}
           siconca: {raw_name: fr_seaice, var_type: atm_2d_ml}
           sithick: {raw_name: hi, var_type: oce_ice}
@@ -148,6 +147,7 @@ projects:
 
           # 3D ocean variables
           msftmz: {var_type: oce_moc, lat_var: lat}
+          phcint: {var_type: oce_def}
           so: {var_type: oce_def, raw_units: "0.001"}
           thetao: {raw_name: to, var_type: oce_def, raw_units: degC}
           uo: {raw_name: u, var_type: oce_def}
diff --git a/esmvalcore/preprocessor/_derive/phcint_total.py b/esmvalcore/preprocessor/_derive/phcint.py
similarity index 71%
rename from esmvalcore/preprocessor/_derive/phcint_total.py
rename to esmvalcore/preprocessor/_derive/phcint.py
index 6459abe7c9..c7e9de623a 100644
--- a/esmvalcore/preprocessor/_derive/phcint_total.py
+++ b/esmvalcore/preprocessor/_derive/phcint.py
@@ -1,4 +1,4 @@
-"""Derivation of variable ``phcint_total``."""
+"""Derivation of variable ``phcint``."""
 
 from __future__ import annotations
 
@@ -7,7 +7,7 @@
 from cf_units import Unit
 from iris import NameConstraint
 
-from esmvalcore.preprocessor._volume import depth_integration
+from esmvalcore.preprocessor._volume import _add_axis_stats_weights_coord
 
 from ._baseclass import DerivedVariableBase
 
@@ -30,7 +30,7 @@ def required(project: str) -> list[Facets]:  # noqa: ARG004
 
     @staticmethod
     def calculate(cubes: CubeList) -> Cube:
-        """Compute total column vertically-integrated heat content.
+        """Compute vertically-integrated heat content.
 
         Use c_p * rho_0 = 4.09169e+6 J m-3 K-1 (Kuhlbrodt et al., 2015, Clim.
         Dyn.)
@@ -48,8 +48,16 @@ def calculate(cubes: CubeList) -> Cube:
         """
         cube = cubes.extract_cube(NameConstraint(var_name="thetao"))
         cube.convert_units("K")
+
+        # J m-3 (multiply by c_p * rho_0)
         cube.data = (
             cube.core_data() * RHO_CP
         )  # https://github.com/SciTools/iris/issues/6990
         cube.units *= RHO_CP_UNIT
-        return depth_integration(cube)
+
+        # J m-2 (multiply by layer depth)
+        _add_axis_stats_weights_coord(cube, cube.coord(axis="z"))
+        cube = cube * cube.coord("_axis_statistics_weights_")
+        cube.remove_coord("_axis_statistics_weights_")
+
+        return cube
diff --git a/tests/unit/preprocessor/_derive/test_phcint_total.py b/tests/unit/preprocessor/_derive/test_phcint.py
similarity index 73%
rename from tests/unit/preprocessor/_derive/test_phcint_total.py
rename to tests/unit/preprocessor/_derive/test_phcint.py
index 5b78861181..75e118144a 100644
--- a/tests/unit/preprocessor/_derive/test_phcint_total.py
+++ b/tests/unit/preprocessor/_derive/test_phcint.py
@@ -1,4 +1,4 @@
-"""Test derivation of ``phcint_total``."""
+"""Test derivation of ``phcint``."""
 
 from __future__ import annotations
 
@@ -32,16 +32,15 @@ def cubes(realistic_4d_cube: Cube) -> CubeList:
 
 @pytest.mark.parametrize("project", ["CMIP3", "CMIP5", "CMIP6", "CMIP7"])
 def test_get_required(project: str) -> None:
-    assert get_required("phcint_total", project) == [{"short_name": "thetao"}]
+    assert get_required("phcint", project) == [{"short_name": "thetao"}]
 
 
 def test_derive(cubes: CubeList) -> None:
-    short_name = "phcint_total"
-    long_name = (
-        "Total Column Integrated Ocean Heat Content from Potential Temperature"
-    )
+    short_name = "phcint"
+    long_name = "Integrated Ocean Heat Content from Potential Temperature"
     units = "J m-2"
     standard_name = "integral_wrt_depth_of_sea_water_potential_temperature_expressed_as_heat_content"
+
     derived_cube = derive(
         cubes,
         short_name=short_name,
@@ -54,16 +53,20 @@ def test_derive(cubes: CubeList) -> None:
     assert derived_cube.long_name == long_name
     assert derived_cube.var_name == short_name
     assert derived_cube.units == units
-
+    assert derived_cube.shape == (2, 1, 2, 3)
     expected_data = np.ma.masked_invalid(
         [
             [
-                [0.0, np.nan, np.nan],
-                [np.nan, 16366760000.0, 20458450000.0],
+                [
+                    [0.0, np.nan, np.nan],
+                    [np.nan, 16366760000.0, 20458450000.0],
+                ],
             ],
             [
-                [24550140000.0, 28641830000.0, 32733520000.0],
-                [36825210000.0, 40916900000.0, 45008590000.0],
+                [
+                    [24550140000.0, 28641830000.0, 32733520000.0],
+                    [36825210000.0, 40916900000.0, 45008590000.0],
+                ],
             ],
         ],
     )

From 29924804cc358152909afab9a8198943464bdd22 Mon Sep 17 00:00:00 2001
From: Manuel Schlund <manuel.schlund@dlr.de>
Date: Wed, 25 Mar 2026 10:29:44 +0100
Subject: [PATCH 09/11] Undo changes in CMOR tables

---
 esmvalcore/cmor/tables/cmip7/tables/CMIP7_ocean.json | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/esmvalcore/cmor/tables/cmip7/tables/CMIP7_ocean.json b/esmvalcore/cmor/tables/cmip7/tables/CMIP7_ocean.json
index bead80ef17..bd6259b6d4 100644
--- a/esmvalcore/cmor/tables/cmip7/tables/CMIP7_ocean.json
+++ b/esmvalcore/cmor/tables/cmip7/tables/CMIP7_ocean.json
@@ -3320,4 +3320,4 @@
             "units": "m"
         }
     }
-}
+}
\ No newline at end of file

From ed2a3aa58127cacf5be9dcee61455fb34466d09e Mon Sep 17 00:00:00 2001
From: Manuel Schlund <manuel.schlund@dlr.de>
Date: Thu, 26 Mar 2026 09:10:34 +0100
Subject: [PATCH 10/11] Avoid dropping cell measures and ancillary variables in
 the derivation of phcint

---
 esmvalcore/preprocessor/_derive/phcint.py | 21 ++++++++++++---------
 1 file changed, 12 insertions(+), 9 deletions(-)

diff --git a/esmvalcore/preprocessor/_derive/phcint.py b/esmvalcore/preprocessor/_derive/phcint.py
index c7e9de623a..f19e7f3416 100644
--- a/esmvalcore/preprocessor/_derive/phcint.py
+++ b/esmvalcore/preprocessor/_derive/phcint.py
@@ -7,7 +7,7 @@
 from cf_units import Unit
 from iris import NameConstraint
 
-from esmvalcore.preprocessor._volume import _add_axis_stats_weights_coord
+from esmvalcore.preprocessor._shared import get_coord_weights
 
 from ._baseclass import DerivedVariableBase
 
@@ -49,15 +49,18 @@ def calculate(cubes: CubeList) -> Cube:
         cube = cubes.extract_cube(NameConstraint(var_name="thetao"))
         cube.convert_units("K")
 
-        # J m-3 (multiply by c_p * rho_0)
-        cube.data = (
-            cube.core_data() * RHO_CP
-        )  # https://github.com/SciTools/iris/issues/6990
+        # In the following, we modify the cube's data and units instead of the
+        # cube directly to avoid dropping cell measures and ancillary variables
+        # (https://scitools-iris.readthedocs.io/en/stable/further_topics/lenient_maths.html#finer-detail)
+
+        # Multiply by c_p * rho_0 -> J m-3
+        cube.data = cube.core_data() * RHO_CP
         cube.units *= RHO_CP_UNIT
 
-        # J m-2 (multiply by layer depth)
-        _add_axis_stats_weights_coord(cube, cube.coord(axis="z"))
-        cube = cube * cube.coord("_axis_statistics_weights_")
-        cube.remove_coord("_axis_statistics_weights_")
+        # Multiply by layer depth -> J m-2
+        z_coord = cube.coord(axis="z")
+        layer_depth = get_coord_weights(cube, z_coord, broadcast=True)
+        cube.data = cube.core_data() * layer_depth
+        cube.units *= z_coord.units
 
         return cube

From f1d1f8d6f0aa7cb2067b4687adee6027f8d44e43 Mon Sep 17 00:00:00 2001
From: Manuel Schlund <manuel.schlund@dlr.de>
Date: Thu, 26 Mar 2026 14:32:14 +0100
Subject: [PATCH 11/11] 100% diff test coverage

---
 tests/unit/preprocessor/_derive/test_amoc.py | 15 ++++++++++-----
 1 file changed, 10 insertions(+), 5 deletions(-)

diff --git a/tests/unit/preprocessor/_derive/test_amoc.py b/tests/unit/preprocessor/_derive/test_amoc.py
index 3b2150d830..a6fb43b9a4 100644
--- a/tests/unit/preprocessor/_derive/test_amoc.py
+++ b/tests/unit/preprocessor/_derive/test_amoc.py
@@ -48,11 +48,16 @@ def cubes():
 def test_amoc_preamble(cubes):
     derived_var = amoc.DerivedVariable()
 
-    cmip5_required = derived_var.required("CMIP5")
-    assert cmip5_required[0]["short_name"] == "msftmyz"
-    cmip6_required = derived_var.required("CMIP6")
-    assert cmip6_required[0]["short_name"] == "msftmz"
-    assert cmip6_required[1]["short_name"] == "msftyz"
+    assert derived_var.required("CMIP5") == [
+        {"short_name": "msftmyz", "mip": "Omon"},
+    ]
+    assert derived_var.required("CMIP6") == [
+        {"short_name": "msftmz", "optional": True},
+        {"short_name": "msftyz", "optional": True},
+    ]
+    assert derived_var.required("ICON") == [
+        {"short_name": "msftmz", "mip": "Omon"},
+    ]
 
     # if project s neither CMIP5 nor CMIP6
     with pytest.raises(ValueError, match="Project CMIPX can not be used"):