From 076d15bcc655ca5a4c5d19bcb1f0496c55d17320 Mon Sep 17 00:00:00 2001
From: Trevor Hillebrand <trhille@lanl.gov>
Date: Thu, 13 Aug 2026 09:00:46 -0700
Subject: [PATCH 1/7] Add ISMIP7 fracture test case for shelf-collapse pathways

Add a 'fracture' test case to the ismip7_forcing test group that
processes the ISMIP7 surface-melt-driven ice shelf collapse forcing
(AIS) via three pathways, each as its own step:

- process_excess_melt (Path A): remaps the excess meltwater field.
  Reconstructs x/y coordinates and corrects the y-axis orientation of
  the source file before conservative remapping.
- process_lake_properties (Path B): remaps supraglacial lake mean depth
  and area fraction (Grau et al. 2025) via bilinear remapping.
- process_shelf_collapse (Path C): remaps the annual ice shelf collapse
  mask via neareststod, rounding to 0/1.

Add per-pathway remap-method config options and documentation.
---
 .../landice/tests/ismip7_forcing/__init__.py  |   2 +
 .../tests/ismip7_forcing/fracture/__init__.py |  51 +++
 .../fracture/process_excess_melt.py           | 333 ++++++++++++++++++
 .../fracture/process_lake_properties.py       | 270 ++++++++++++++
 .../fracture/process_shelf_collapse.py        | 199 +++++++++++
 .../tests/ismip7_forcing/ismip7_forcing.cfg   |  24 ++
 .../ismip7_forcing/ismip7_forcing_test.cfg    |  24 ++
 docs/developers_guide/landice/api.rst         |  12 +
 .../landice/test_groups/ismip7_forcing.rst    |  94 ++++-
 9 files changed, 1006 insertions(+), 3 deletions(-)
 create mode 100644 compass/landice/tests/ismip7_forcing/fracture/__init__.py
 create mode 100644 compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
 create mode 100644 compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
 create mode 100644 compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py

diff --git a/compass/landice/tests/ismip7_forcing/__init__.py b/compass/landice/tests/ismip7_forcing/__init__.py
index 3f621a4c37..cb448350e3 100644
--- a/compass/landice/tests/ismip7_forcing/__init__.py
+++ b/compass/landice/tests/ismip7_forcing/__init__.py
@@ -1,4 +1,5 @@
 from compass.landice.tests.ismip7_forcing.atmosphere import Atmosphere
+from compass.landice.tests.ismip7_forcing.fracture import Fracture
 from compass.landice.tests.ismip7_forcing.ocean_thermal import OceanThermal
 from compass.testgroup import TestGroup
 
@@ -22,3 +23,4 @@ def __init__(self, mpas_core):
 
         self.add_test_case(Atmosphere(test_group=self))
         self.add_test_case(OceanThermal(test_group=self))
+        self.add_test_case(Fracture(test_group=self))
diff --git a/compass/landice/tests/ismip7_forcing/fracture/__init__.py b/compass/landice/tests/ismip7_forcing/fracture/__init__.py
new file mode 100644
index 0000000000..b6a7e63c6a
--- /dev/null
+++ b/compass/landice/tests/ismip7_forcing/fracture/__init__.py
@@ -0,0 +1,51 @@
+from compass.landice.tests.ismip7_forcing.configure import (
+    configure as configure_testgroup,
+)
+from compass.landice.tests.ismip7_forcing.fracture.process_excess_melt import (
+    ProcessExcessMelt,
+)
+from compass.landice.tests.ismip7_forcing.fracture.process_lake_properties import (  # noqa: E501
+    ProcessLakeProperties,
+)
+from compass.landice.tests.ismip7_forcing.fracture.process_shelf_collapse import (  # noqa: E501
+    ProcessShelfCollapse,
+)
+from compass.testcase import TestCase
+
+
+class Fracture(TestCase):
+    """
+    A test case for processing ISMIP7 fracture forcing data.
+    Implements the surface-melt-driven ice shelf collapse pathways:
+
+    * Path A (``process_excess_melt``): excess meltwater after firn air
+      content depletion.
+    * Path B (``process_lake_properties``): supraglacial lake mean depth
+      and area fraction from Grau et al. (2025).
+    * Path C (``process_shelf_collapse``): the ice shelf collapse mask, in
+      which a floating grid cell is flagged as collapsed when excess
+      meltwater exceeds 72.5 mm/yr for 10 consecutive years.
+    """
+
+    def __init__(self, test_group):
+        """
+        Create the test case
+
+        Parameters
+        ----------
+        test_group : compass.landice.tests.ismip7_forcing.Ismip7Forcing
+            The test group that this test case belongs to
+        """
+        name = "fracture"
+        subdir = name
+        super().__init__(test_group=test_group, name=name, subdir=subdir)
+
+        self.add_step(ProcessExcessMelt(test_case=self))
+        self.add_step(ProcessLakeProperties(test_case=self))
+        self.add_step(ProcessShelfCollapse(test_case=self))
+
+    def configure(self):
+        """
+        Configures test case
+        """
+        configure_testgroup(config=self.config)
diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py b/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
new file mode 100644
index 0000000000..b7af781031
--- /dev/null
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
@@ -0,0 +1,333 @@
+import glob
+import os
+import shutil
+
+import numpy as np
+import xarray as xr
+from mpas_tools.io import write_netcdf
+from mpas_tools.logging import check_call
+from scipy.ndimage import distance_transform_edt
+
+from compass.landice.tests.ismip7_forcing.create_mapfile import (
+    build_mapping_file,
+)
+from compass.step import Step
+
+
+class ProcessExcessMelt(Step):
+    """
+    A step for processing the ISMIP7 excess meltwater field (Path A).
+    Remaps the annual excess melt (melt + rain after firn air content
+    depletion) from the ISMIP7 polar stereographic grid to the MALI
+    unstructured mesh.
+
+    The excess melt file lacks ``x``/``y`` coordinate variables and its
+    array is flipped along the y axis relative to the other fracture
+    files (it was produced with CDO). This step reconstructs a source
+    grid with ``x``/``y`` coordinates and the correct orientation before
+    remapping.
+    """
+
+    def __init__(self, test_case):
+        """
+        Create the step
+
+        Parameters
+        ----------
+        test_case : compass.landice.tests.ismip7_forcing.fracture.Fracture
+            The test case this step belongs to
+        """
+        super().__init__(test_case=test_case, name="process_excess_melt")
+
+    def setup(self):
+        """
+        Set up this step of the test case
+        """
+        config = self.config
+        section = config["ismip7"]
+        base_path_mali = section.get("base_path_mali")
+        mali_mesh_file = section.get("mali_mesh_file")
+
+        self.add_input_file(filename=mali_mesh_file,
+                            target=os.path.join(base_path_mali,
+                                                mali_mesh_file))
+
+    def run(self):
+        """
+        Run this step of the test case
+        """
+        logger = self.logger
+        config = self.config
+
+        section = config["ismip7"]
+        base_path_ismip7 = section.get("base_path_ismip7")
+        mali_mesh_name = section.get("mali_mesh_name")
+        mali_mesh_file = section.get("mali_mesh_file")
+        model = section.get("model")
+        scenario = section.get("scenario")
+        output_base_path = section.get("output_base_path")
+        ice_sheet = section.get("ice_sheet")
+
+        section = config["ismip7_fracture"]
+        method_remap = section.get("method_remap_excess_melt")
+        version = section.get("version")
+        start_year = section.getint("start_year")
+        end_year = section.getint("end_year")
+
+        # Discover the excess melt file
+        input_path = os.path.join(base_path_ismip7, "fracture", version)
+        file_pattern = "excess_melt_*.nc"
+        all_files = sorted(glob.glob(os.path.join(input_path, file_pattern)))
+
+        if not all_files:
+            raise FileNotFoundError(
+                f"No excess melt file found matching pattern:\n"
+                f"  {os.path.join(input_path, file_pattern)}")
+        if len(all_files) > 1:
+            raise ValueError(
+                f"Expected a single excess melt file but found "
+                f"{len(all_files)}:\n  " + "\n  ".join(all_files))
+
+        input_file = all_files[0]
+        basename = os.path.basename(input_file)
+        logger.info(f"Processing excess melt: {basename}")
+
+        # Build a source file with x/y coordinates and correct orientation
+        gridded_file = f"gridded_{basename}"
+        self._prepare_source_grid(input_file, input_path, gridded_file,
+                                  logger)
+
+        # Build mapping file. Excess melt is a flux, so conservative
+        # remapping is appropriate by default.
+        mapping_file = (f"map_ismip7_{ice_sheet}_fracture_to_"
+                        f"{mali_mesh_name}_{method_remap}.nc")
+
+        if not os.path.exists(mapping_file):
+            logger.info("Building mapping file for the excess melt grid...")
+            build_mapping_file(config, logger,
+                               gridded_file, mapping_file,
+                               mali_mesh_file=mali_mesh_file,
+                               method_remap=method_remap)
+
+        # Extrapolate fill values on the source grid before remapping so
+        # they don't pollute neighboring cells during interpolation
+        extrap_file = f"extrap_{basename}"
+        self._extrapolate_source(gridded_file, extrap_file, "excess_melt",
+                                 logger)
+
+        # Remap the excess melt onto the MALI mesh
+        remapped_file = f"remapped_{basename}"
+        logger.info(f"Remapping: {basename}")
+        args = ["ncremap",
+                "-i", extrap_file,
+                "-o", remapped_file,
+                "-m", mapping_file,
+                "-v", "excess_melt"]
+        check_call(args, logger=logger)
+
+        # Rename to MALI conventions
+        logger.info("Renaming variables to MALI conventions...")
+        output_file = f"{mali_mesh_name}_{basename}"
+        self._rename_to_mali_vars(remapped_file, output_file,
+                                  start_year, end_year)
+
+        # Clean up temporary files
+        for f in [gridded_file, extrap_file, remapped_file]:
+            if os.path.exists(f):
+                os.remove(f)
+
+        # Place output in the appropriate directory
+        output_path = os.path.join(output_base_path, "excess_melt",
+                                   f"{model}_{scenario}")
+        if not os.path.exists(output_path):
+            os.makedirs(output_path)
+
+        dst = os.path.join(output_path, output_file)
+        shutil.copy(output_file, dst)
+
+        logger.info(f"Done. Output: {dst}")
+
+    def _prepare_source_grid(self, input_file, input_path, output_file,
+                             logger):
+        """
+        Build a source file for the excess melt data with ``x``/``y``
+        coordinate variables and the standard ISMIP7 orientation
+        (south-to-north). The excess melt file only has 2D ``lat``/``lon``
+        and its array is flipped along the y axis relative to the other
+        fracture files, so ``x``/``y`` are borrowed from a sibling
+        fracture file and the data are flipped to match.
+
+        Parameters
+        ----------
+        input_file : str
+            Path to the excess melt file
+
+        input_path : str
+            Directory containing the fracture forcing files
+
+        output_file : str
+            Path to write the reconstructed source file
+
+        logger : logging.Logger
+            Logger for status messages
+        """
+        # Find a sibling fracture file that has x/y coordinate variables
+        grid_donor = None
+        for pattern in ["lake_properties_*.nc",
+                        "ice_shelf_collapse_mask_*.nc"]:
+            for candidate in sorted(glob.glob(os.path.join(input_path,
+                                                           pattern))):
+                with xr.open_dataset(candidate) as dtest:
+                    if "x" in dtest.variables and "y" in dtest.variables:
+                        grid_donor = candidate
+                        break
+            if grid_donor is not None:
+                break
+
+        if grid_donor is None:
+            raise FileNotFoundError(
+                "Could not find a sibling fracture file with x/y "
+                "coordinate variables to define the excess melt grid.")
+
+        logger.info(f"Using grid from sibling file: "
+                    f"{os.path.basename(grid_donor)}")
+
+        with xr.open_dataset(grid_donor) as donor:
+            x = donor["x"].values
+            y = donor["y"].values
+            donor_lat = donor["lat"].values
+
+        ds = xr.open_dataset(input_file, decode_times=False)
+
+        # Flip along the y axis to match the sibling grid orientation
+        flipped_lat = ds["lat"].values[::-1, :]
+        if np.nanmax(np.abs(flipped_lat - donor_lat)) > 1.0e-3:
+            raise ValueError(
+                "Excess melt grid does not match the sibling grid after a "
+                "y-axis flip; orientation cannot be determined "
+                "automatically.")
+
+        excess = ds["excess_melt"].values[:, ::-1, :]
+        years = ds["year"].values.astype(int)
+
+        out = xr.Dataset()
+        out["x"] = ("x", x)
+        out["y"] = ("y", y)
+        out["time"] = ("time", years)
+        out["excess_melt"] = (("time", "y", "x"), excess)
+        out["excess_melt"].attrs = dict(ds["excess_melt"].attrs)
+
+        write_netcdf(out, output_file)
+        ds.close()
+
+    def _rename_to_mali_vars(self, remapped_file, output_file,
+                             start_year, end_year):
+        """
+        Rename dimensions/variables of the remapped excess melt to MALI
+        conventions, restrict to the requested year range, and add an
+        ``xtime`` variable.
+
+        Parameters
+        ----------
+        remapped_file : str
+            Excess melt remapped onto the MALI mesh
+
+        output_file : str
+            Output file with MALI variable/dimension names
+
+        start_year : int
+            First year (inclusive) to retain
+
+        end_year : int
+            Last year (inclusive) to retain
+        """
+        # The time coordinate has units="year" (integer years), which is
+        # not CF-compliant, so disable time decoding.
+        ds = xr.open_dataset(remapped_file, decode_times=False)
+
+        rename_dims = {}
+        if "ncol" in ds.dims:
+            rename_dims["ncol"] = "nCells"
+        if "time" in ds.dims:
+            rename_dims["time"] = "Time"
+        if rename_dims:
+            ds = ds.rename(rename_dims)
+
+        if "excess_melt" in ds:
+            ds = ds.rename({"excess_melt": "ismip7ExcessMelt"})
+
+        # Restrict to the requested year range
+        years = ds["time"].values.astype(int)
+        keep = (years >= start_year) & (years <= end_year)
+        ds = ds.isel(Time=keep)
+        years = years[keep]
+
+        # Excess melt is an annual field applied at the start of each year
+        xtime = [f"{int(yr):04d}-01-01_00:00:00".ljust(64) for yr in years]
+        ds["xtime"] = ("Time", xtime)
+        ds["xtime"] = ds.xtime.astype("S")
+
+        ds["ismip7ExcessMelt"].attrs = {
+            "long_name": "excess meltwater after firn air content depletion",
+            "units": "mm w.e. yr-1",
+        }
+
+        vars_to_drop = [v for v in ["lat_vertices", "lon_vertices", "lat",
+                                    "lon", "area", "time"]
+                        if v in ds]
+        if vars_to_drop:
+            ds = ds.drop_vars(vars_to_drop)
+
+        write_netcdf(ds, output_file)
+        ds.close()
+
+    def _extrapolate_source(self, input_file, output_file, varname, logger):
+        """
+        Extrapolate fill/missing values on the source polar stereographic
+        grid using nearest-neighbor via distance_transform_edt. This must
+        be done before remapping so that fill values don't contaminate the
+        interpolation stencil.
+
+        Parameters
+        ----------
+        input_file : str
+            Path to the input NetCDF file on the source grid
+
+        output_file : str
+            Path to write the extrapolated file
+
+        varname : str
+            Name of the variable to extrapolate
+
+        logger : logging.Logger
+            Logger for status messages
+        """
+        logger.info(f"    Extrapolating fill values on source grid: "
+                    f"{os.path.basename(input_file)}")
+
+        ds = xr.open_dataset(input_file, decode_times=False)
+        data = ds[varname]
+
+        values = data.values.copy()
+        non_spatial_shape = values.shape[:-2]
+
+        for idx in np.ndindex(non_spatial_shape):
+            slab = values[idx]
+            valid_mask = np.isfinite(slab)
+            if valid_mask.all() or not valid_mask.any():
+                continue
+            nearest_inds = distance_transform_edt(
+                ~valid_mask, return_distances=False, return_indices=True)
+            invalid = ~valid_mask
+            values[idx][invalid] = slab[
+                nearest_inds[0, invalid],
+                nearest_inds[1, invalid]]
+
+        ds[varname] = (data.dims, values)
+        ds[varname].attrs = data.attrs
+
+        if "_FillValue" in ds[varname].encoding:
+            del ds[varname].encoding["_FillValue"]
+
+        write_netcdf(ds, output_file)
+        ds.close()
diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py b/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
new file mode 100644
index 0000000000..d5e4677ab5
--- /dev/null
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
@@ -0,0 +1,270 @@
+import glob
+import os
+import shutil
+
+import numpy as np
+import xarray as xr
+from mpas_tools.io import write_netcdf
+from mpas_tools.logging import check_call
+from scipy.ndimage import distance_transform_edt
+
+from compass.landice.tests.ismip7_forcing.create_mapfile import (
+    build_mapping_file,
+)
+from compass.step import Step
+
+
+class ProcessLakeProperties(Step):
+    """
+    A step for processing the ISMIP7 supraglacial lake properties (Path B).
+    Remaps the annual mean lake depth and lake area fraction (from the
+    Grau et al. (2025) parameterization) from the ISMIP7 polar
+    stereographic grid to the MALI unstructured mesh.
+    """
+
+    # Source variable name -> MALI output variable name and attributes
+    _variables = {
+        "lake_depth": {
+            "mali_name": "ismip7LakeDepth",
+            "long_name": "mean supraglacial lake depth",
+            "units": "m",
+        },
+        "fraction_lake_area": {
+            "mali_name": "ismip7LakeAreaFraction",
+            "long_name": "supraglacial lake area fraction",
+            "units": "1",
+        },
+    }
+
+    def __init__(self, test_case):
+        """
+        Create the step
+
+        Parameters
+        ----------
+        test_case : compass.landice.tests.ismip7_forcing.fracture.Fracture
+            The test case this step belongs to
+        """
+        super().__init__(test_case=test_case, name="process_lake_properties")
+
+    def setup(self):
+        """
+        Set up this step of the test case
+        """
+        config = self.config
+        section = config["ismip7"]
+        base_path_mali = section.get("base_path_mali")
+        mali_mesh_file = section.get("mali_mesh_file")
+
+        self.add_input_file(filename=mali_mesh_file,
+                            target=os.path.join(base_path_mali,
+                                                mali_mesh_file))
+
+    def run(self):
+        """
+        Run this step of the test case
+        """
+        logger = self.logger
+        config = self.config
+
+        section = config["ismip7"]
+        base_path_ismip7 = section.get("base_path_ismip7")
+        mali_mesh_name = section.get("mali_mesh_name")
+        mali_mesh_file = section.get("mali_mesh_file")
+        model = section.get("model")
+        scenario = section.get("scenario")
+        output_base_path = section.get("output_base_path")
+        ice_sheet = section.get("ice_sheet")
+
+        section = config["ismip7_fracture"]
+        method_remap = section.get("method_remap_lake_properties")
+        version = section.get("version")
+        start_year = section.getint("start_year")
+        end_year = section.getint("end_year")
+
+        # Discover the lake properties file
+        input_path = os.path.join(base_path_ismip7, "fracture", version)
+        file_pattern = "lake_properties_*.nc"
+        all_files = sorted(glob.glob(os.path.join(input_path, file_pattern)))
+
+        if not all_files:
+            raise FileNotFoundError(
+                f"No lake properties file found matching pattern:\n"
+                f"  {os.path.join(input_path, file_pattern)}")
+        if len(all_files) > 1:
+            raise ValueError(
+                f"Expected a single lake properties file but found "
+                f"{len(all_files)}:\n  " + "\n  ".join(all_files))
+
+        input_file = all_files[0]
+        basename = os.path.basename(input_file)
+        logger.info(f"Processing lake properties: {basename}")
+
+        # Build mapping file. Lake properties are continuous fields, so
+        # bilinear remapping is appropriate by default.
+        mapping_file = (f"map_ismip7_{ice_sheet}_fracture_to_"
+                        f"{mali_mesh_name}_{method_remap}.nc")
+
+        if not os.path.exists(mapping_file):
+            logger.info("Building mapping file for the lake properties "
+                        "grid...")
+            build_mapping_file(config, logger,
+                               input_file, mapping_file,
+                               mali_mesh_file=mali_mesh_file,
+                               method_remap=method_remap)
+
+        # Extrapolate fill values on the source grid before remapping so
+        # they don't pollute neighboring cells during interpolation
+        extrap_file = f"extrap_{basename}"
+        self._extrapolate_source(input_file, extrap_file,
+                                 list(self._variables.keys()), logger)
+
+        # Remap both lake property variables onto the MALI mesh
+        remapped_file = f"remapped_{basename}"
+        logger.info(f"Remapping: {basename}")
+        args = ["ncremap",
+                "-i", extrap_file,
+                "-o", remapped_file,
+                "-m", mapping_file,
+                "-v", ",".join(self._variables.keys())]
+        check_call(args, logger=logger)
+
+        # Rename to MALI conventions
+        logger.info("Renaming variables to MALI conventions...")
+        output_file = f"{mali_mesh_name}_{basename}"
+        self._rename_to_mali_vars(remapped_file, output_file,
+                                  start_year, end_year)
+
+        # Clean up temporary files
+        for f in [extrap_file, remapped_file]:
+            if os.path.exists(f):
+                os.remove(f)
+
+        # Place output in the appropriate directory
+        output_path = os.path.join(output_base_path, "lake_properties",
+                                   f"{model}_{scenario}")
+        if not os.path.exists(output_path):
+            os.makedirs(output_path)
+
+        dst = os.path.join(output_path, output_file)
+        shutil.copy(output_file, dst)
+
+        logger.info(f"Done. Output: {dst}")
+
+    def _rename_to_mali_vars(self, remapped_file, output_file,
+                             start_year, end_year):
+        """
+        Rename dimensions/variables of the remapped lake properties to MALI
+        conventions, restrict to the requested year range, and add an
+        ``xtime`` variable.
+
+        Parameters
+        ----------
+        remapped_file : str
+            Lake properties remapped onto the MALI mesh
+
+        output_file : str
+            Output file with MALI variable/dimension names
+
+        start_year : int
+            First year (inclusive) to retain
+
+        end_year : int
+            Last year (inclusive) to retain
+        """
+        # The time coordinate has units="year" (integer years), which is
+        # not CF-compliant, so disable time decoding.
+        ds = xr.open_dataset(remapped_file, decode_times=False)
+
+        rename_dims = {}
+        if "ncol" in ds.dims:
+            rename_dims["ncol"] = "nCells"
+        if "time" in ds.dims:
+            rename_dims["time"] = "Time"
+        if rename_dims:
+            ds = ds.rename(rename_dims)
+
+        rename_vars = {src: info["mali_name"]
+                       for src, info in self._variables.items()
+                       if src in ds}
+        ds = ds.rename(rename_vars)
+
+        # Restrict to the requested year range
+        years = ds["time"].values.astype(int)
+        keep = (years >= start_year) & (years <= end_year)
+        ds = ds.isel(Time=keep)
+        years = years[keep]
+
+        # Lake properties are annual fields applied at the start of each year
+        xtime = [f"{int(yr):04d}-01-01_00:00:00".ljust(64) for yr in years]
+        ds["xtime"] = ("Time", xtime)
+        ds["xtime"] = ds.xtime.astype("S")
+
+        for info in self._variables.values():
+            mali_name = info["mali_name"]
+            if mali_name in ds:
+                ds[mali_name].attrs = {
+                    "long_name": info["long_name"],
+                    "units": info["units"],
+                }
+
+        vars_to_drop = [v for v in ["lat_vertices", "lon_vertices", "lat",
+                                    "lon", "area", "time"]
+                        if v in ds]
+        if vars_to_drop:
+            ds = ds.drop_vars(vars_to_drop)
+
+        write_netcdf(ds, output_file)
+        ds.close()
+
+    def _extrapolate_source(self, input_file, output_file, varnames, logger):
+        """
+        Extrapolate fill/missing values on the source polar stereographic
+        grid using nearest-neighbor via distance_transform_edt. This must
+        be done before remapping so that fill values don't contaminate the
+        interpolation stencil.
+
+        Parameters
+        ----------
+        input_file : str
+            Path to the input NetCDF file on the source grid
+
+        output_file : str
+            Path to write the extrapolated file
+
+        varnames : list of str
+            Names of the variables to extrapolate
+
+        logger : logging.Logger
+            Logger for status messages
+        """
+        logger.info(f"    Extrapolating fill values on source grid: "
+                    f"{os.path.basename(input_file)}")
+
+        ds = xr.open_dataset(input_file, decode_times=False)
+
+        for varname in varnames:
+            data = ds[varname]
+            values = data.values.copy()
+            non_spatial_shape = values.shape[:-2]
+
+            for idx in np.ndindex(non_spatial_shape):
+                slab = values[idx]
+                valid_mask = np.isfinite(slab)
+                if valid_mask.all() or not valid_mask.any():
+                    continue
+                nearest_inds = distance_transform_edt(
+                    ~valid_mask, return_distances=False,
+                    return_indices=True)
+                invalid = ~valid_mask
+                values[idx][invalid] = slab[
+                    nearest_inds[0, invalid],
+                    nearest_inds[1, invalid]]
+
+            ds[varname] = (data.dims, values)
+            ds[varname].attrs = data.attrs
+            if "_FillValue" in ds[varname].encoding:
+                del ds[varname].encoding["_FillValue"]
+
+        write_netcdf(ds, output_file)
+        ds.close()
diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py b/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py
new file mode 100644
index 0000000000..02b2d0312c
--- /dev/null
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py
@@ -0,0 +1,199 @@
+import glob
+import os
+import shutil
+
+import xarray as xr
+from mpas_tools.io import write_netcdf
+from mpas_tools.logging import check_call
+
+from compass.landice.tests.ismip7_forcing.create_mapfile import (
+    build_mapping_file,
+)
+from compass.step import Step
+
+
+class ProcessShelfCollapse(Step):
+    """
+    A step for processing the ISMIP7 ice shelf collapse mask (Path C).
+    Remaps the annual 0/1 collapse mask from the ISMIP7 polar
+    stereographic grid to the MALI unstructured mesh and renames
+    variables to MALI conventions. The mask flags floating grid cells
+    that collapse once excess meltwater (after firn air content
+    depletion) exceeds 72.5 mm/yr for 10 consecutive years.
+    """
+
+    def __init__(self, test_case):
+        """
+        Create the step
+
+        Parameters
+        ----------
+        test_case : compass.landice.tests.ismip7_forcing.fracture.Fracture
+            The test case this step belongs to
+        """
+        super().__init__(test_case=test_case, name="process_shelf_collapse")
+
+    def setup(self):
+        """
+        Set up this step of the test case
+        """
+        config = self.config
+        section = config["ismip7"]
+        base_path_mali = section.get("base_path_mali")
+        mali_mesh_file = section.get("mali_mesh_file")
+
+        self.add_input_file(filename=mali_mesh_file,
+                            target=os.path.join(base_path_mali,
+                                                mali_mesh_file))
+
+    def run(self):
+        """
+        Run this step of the test case
+        """
+        logger = self.logger
+        config = self.config
+
+        section = config["ismip7"]
+        base_path_ismip7 = section.get("base_path_ismip7")
+        mali_mesh_name = section.get("mali_mesh_name")
+        mali_mesh_file = section.get("mali_mesh_file")
+        model = section.get("model")
+        scenario = section.get("scenario")
+        output_base_path = section.get("output_base_path")
+        ice_sheet = section.get("ice_sheet")
+
+        section = config["ismip7_fracture"]
+        method_remap = section.get("method_remap_shelf_collapse")
+        version = section.get("version")
+        start_year = section.getint("start_year")
+        end_year = section.getint("end_year")
+
+        # Discover the ice shelf collapse mask file
+        input_path = os.path.join(base_path_ismip7, "fracture", version)
+        file_pattern = "ice_shelf_collapse_mask_*.nc"
+        all_files = sorted(glob.glob(os.path.join(input_path, file_pattern)))
+
+        if not all_files:
+            raise FileNotFoundError(
+                f"No ice shelf collapse mask file found matching pattern:\n"
+                f"  {os.path.join(input_path, file_pattern)}")
+        if len(all_files) > 1:
+            raise ValueError(
+                f"Expected a single ice shelf collapse mask file but found "
+                f"{len(all_files)}:\n  " + "\n  ".join(all_files))
+
+        input_file = all_files[0]
+        basename = os.path.basename(input_file)
+        logger.info(f"Processing ice shelf collapse mask: {basename}")
+
+        # Build mapping file. neareststod preserves the 0/1 mask values.
+        mapping_file = (f"map_ismip7_{ice_sheet}_fracture_to_"
+                        f"{mali_mesh_name}_{method_remap}.nc")
+
+        if not os.path.exists(mapping_file):
+            logger.info("Building mapping file for the collapse mask grid...")
+            build_mapping_file(config, logger,
+                               input_file, mapping_file,
+                               mali_mesh_file=mali_mesh_file,
+                               method_remap=method_remap)
+
+        # Remap the collapse mask onto the MALI mesh
+        remapped_file = f"remapped_{basename}"
+        if not os.path.exists(remapped_file):
+            logger.info(f"Remapping: {basename}")
+            args = ["ncremap",
+                    "-i", input_file,
+                    "-o", remapped_file,
+                    "-m", mapping_file,
+                    "-v", "mask"]
+            check_call(args, logger=logger)
+
+        # Combine time slice and rename to MALI conventions
+        logger.info("Renaming variables to MALI conventions...")
+        output_file = f"{mali_mesh_name}_{basename}"
+        self._rename_to_mali_vars(remapped_file, output_file,
+                                  start_year, end_year)
+
+        # Clean up temporary remapped file
+        if os.path.exists(remapped_file):
+            os.remove(remapped_file)
+
+        # Place output in the appropriate directory
+        output_path = os.path.join(output_base_path, "shelf_collapse",
+                                   f"{model}_{scenario}")
+        if not os.path.exists(output_path):
+            os.makedirs(output_path)
+
+        dst = os.path.join(output_path, output_file)
+        shutil.copy(output_file, dst)
+
+        logger.info(f"Done. Output: {dst}")
+
+    def _rename_to_mali_vars(self, remapped_file, output_file,
+                             start_year, end_year):
+        """
+        Rename dimensions/variables of the remapped collapse mask to MALI
+        conventions, restrict to the requested year range, round the mask
+        to 0/1, and add an ``xtime`` variable.
+
+        Parameters
+        ----------
+        remapped_file : str
+            Collapse mask remapped onto the MALI mesh
+
+        output_file : str
+            Output file with MALI variable/dimension names
+
+        start_year : int
+            First year (inclusive) to retain
+
+        end_year : int
+            Last year (inclusive) to retain
+        """
+        # The collapse mask time coordinate has units="year" (integer years),
+        # which is not CF-compliant, so disable time decoding.
+        ds = xr.open_dataset(remapped_file, decode_times=False)
+
+        # Rename dimensions to MALI conventions
+        rename_dims = {}
+        if "ncol" in ds.dims:
+            rename_dims["ncol"] = "nCells"
+        if "time" in ds.dims:
+            rename_dims["time"] = "Time"
+        if rename_dims:
+            ds = ds.rename(rename_dims)
+
+        # Rename variable
+        if "mask" in ds:
+            ds = ds.rename({"mask": "calvingMask"})
+
+        # Restrict to the requested year range
+        years = ds["time"].values.astype(int)
+        keep = (years >= start_year) & (years <= end_year)
+        ds = ds.isel(Time=keep)
+        years = years[keep]
+
+        # Round the remapped mask to 0/1 and store as integers. Ice shelves
+        # collapse on January 1st, so the mask is applied at the start of
+        # each year.
+        calving_mask = (ds["calvingMask"] >= 0.5).astype(int)
+        ds["calvingMask"] = calving_mask
+
+        # Add xtime variable, one entry per year at January 1st
+        xtime = [f"{int(yr):04d}-01-01_00:00:00".ljust(64) for yr in years]
+        ds["xtime"] = ("Time", xtime)
+        ds["xtime"] = ds.xtime.astype("S")
+
+        ds["calvingMask"].attrs = {
+            "long_name": "ice shelf collapse mask (1 = collapse)",
+        }
+
+        # Drop auxiliary variables carried over from remapping
+        vars_to_drop = [v for v in ["lat_vertices", "lon_vertices", "lat",
+                                    "lon", "area", "time"]
+                        if v in ds]
+        if vars_to_drop:
+            ds = ds.drop_vars(vars_to_drop)
+
+        write_netcdf(ds, output_file)
+        ds.close()
diff --git a/compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg b/compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg
index b996d2aaba..a6eadbd26e 100644
--- a/compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg
+++ b/compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg
@@ -69,3 +69,27 @@ method_remap = bilinear
 # Base path to observational climatology data
 # (directory containing tf/, so/, thetao/ subdirs)
 base_path_climatology = /path/to/ISMIP7/forcing/AIS/obs/zhou_annual_06_nov
+
+# config options for ismip7 fracture (Path C, ice shelf collapse) forcing
+[ismip7_fracture]
+
+# Remapping method for the ice shelf collapse mask (Path C). neareststod
+# preserves the 0/1 mask values. Options: bilinear, neareststod, conserve
+method_remap_shelf_collapse = neareststod
+
+# Remapping method for the excess meltwater field (Path A), a flux.
+# Options: bilinear, neareststod, conserve
+method_remap_excess_melt = conserve
+
+# Remapping method for the supraglacial lake properties (Path B).
+# Options: bilinear, neareststod, conserve
+method_remap_lake_properties = bilinear
+
+# Version subdirectory of the fracture forcing data
+version = v2
+
+# Start year for processing
+start_year = 1850
+
+# End year for processing
+end_year = 2014
diff --git a/compass/landice/tests/ismip7_forcing/ismip7_forcing_test.cfg b/compass/landice/tests/ismip7_forcing/ismip7_forcing_test.cfg
index 270f5d396c..28915fe41d 100644
--- a/compass/landice/tests/ismip7_forcing/ismip7_forcing_test.cfg
+++ b/compass/landice/tests/ismip7_forcing/ismip7_forcing_test.cfg
@@ -69,3 +69,27 @@ method_remap = bilinear
 # Base path to observational climatology data
 # (directory containing tf/, so/, thetao/ subdirs)
 base_path_climatology = /global/cfs/cdirs/m4288/users/trhille/ISMIP7/forcing/AIS/obs/zhou_annual_06_nov
+
+# config options for ismip7 fracture (Path C, ice shelf collapse) forcing
+[ismip7_fracture]
+
+# Remapping method for the ice shelf collapse mask (Path C). neareststod
+# preserves the 0/1 mask values. Options: bilinear, neareststod, conserve
+method_remap_shelf_collapse = neareststod
+
+# Remapping method for the excess meltwater field (Path A), a flux.
+# Options: bilinear, neareststod, conserve
+method_remap_excess_melt = conserve
+
+# Remapping method for the supraglacial lake properties (Path B).
+# Options: bilinear, neareststod, conserve
+method_remap_lake_properties = bilinear
+
+# Version subdirectory of the fracture forcing data
+version = v2
+
+# Start year for processing
+start_year = 2015
+
+# End year for processing
+end_year = 2300
diff --git a/docs/developers_guide/landice/api.rst b/docs/developers_guide/landice/api.rst
index 57514fd526..100a147ad8 100644
--- a/docs/developers_guide/landice/api.rst
+++ b/docs/developers_guide/landice/api.rst
@@ -410,6 +410,18 @@ ismip7_forcing
    ocean_thermal.process_thermal_forcing.ProcessThermalForcing.setup
    ocean_thermal.process_thermal_forcing.ProcessThermalForcing.run
 
+   fracture.Fracture
+   fracture.Fracture.configure
+   fracture.process_excess_melt.ProcessExcessMelt
+   fracture.process_excess_melt.ProcessExcessMelt.setup
+   fracture.process_excess_melt.ProcessExcessMelt.run
+   fracture.process_lake_properties.ProcessLakeProperties
+   fracture.process_lake_properties.ProcessLakeProperties.setup
+   fracture.process_lake_properties.ProcessLakeProperties.run
+   fracture.process_shelf_collapse.ProcessShelfCollapse
+   fracture.process_shelf_collapse.ProcessShelfCollapse.setup
+   fracture.process_shelf_collapse.ProcessShelfCollapse.run
+
 isunnguata_sermia
 ~~~~~~~~~~~~~~~~~
 
diff --git a/docs/users_guide/landice/test_groups/ismip7_forcing.rst b/docs/users_guide/landice/test_groups/ismip7_forcing.rst
index e26d3f88e9..d134093dce 100644
--- a/docs/users_guide/landice/test_groups/ismip7_forcing.rst
+++ b/docs/users_guide/landice/test_groups/ismip7_forcing.rst
@@ -10,7 +10,8 @@ force MALI in its simulations under the ISMIP7 experimental protocol.
 The test group supports both the Antarctic Ice Sheet (AIS) and the Greenland
 Ice Sheet (GrIS), controlled by a single ``ice_sheet`` config option.
 
-The test group includes two test cases: ``atmosphere`` and ``ocean_thermal``.
+The test group includes three test cases: ``atmosphere``, ``ocean_thermal``,
+and ``fracture``.
 
 * The ``atmosphere`` test case has five steps:
   ``process_smb``, ``process_temperature``, ``process_smb_gradient``,
@@ -22,9 +23,15 @@ The test group includes two test cases: ``atmosphere`` and ``ocean_thermal``.
   process the observational ocean thermal forcing climatology (Zhou et al.)
   for AIS, controlled by the ``process_ocean_climatology`` config option.
 
+* The ``fracture`` test case has three steps: ``process_excess_melt``
+  (Path A), ``process_lake_properties`` (Path B), and
+  ``process_shelf_collapse`` (Path C). It processes the ISMIP7
+  surface-melt-driven ice shelf collapse forcing (AIS only).
+
 (For more details on the steps of each test case, see
-:ref:`landice_ismip7_forcing_atmosphere` and
-:ref:`landice_ismip7_forcing_ocean_thermal`.)
+:ref:`landice_ismip7_forcing_atmosphere`,
+:ref:`landice_ismip7_forcing_ocean_thermal`, and
+:ref:`landice_ismip7_forcing_fracture`.)
 
 .. _landice_ismip7_forcing_usage:
 
@@ -43,6 +50,9 @@ To use this test group, users need to:
 
 5. Run the ``ocean_thermal`` test case for each model and scenario combination.
 
+6. Run the ``fracture`` test case (AIS only) for each model and scenario
+   combination to process the Path C ice shelf collapse mask.
+
 Example user config files are provided in the source tree for local testing:
 
 * ``compass/landice/tests/ismip7_forcing/ismip7_forcing_test.cfg``
@@ -88,6 +98,12 @@ For AIS ocean thermal climatology (8km, 30 depth levels, static):
 
    {base_path_climatology}/tf/v3/tf_AIS_obs_ocean_climatology_*.nc
 
+For AIS fracture / ice shelf collapse mask (8km, annual, Path C):
+
+.. code-block:: none
+
+   fracture/v2/ice_shelf_collapse_mask_*.nc
+
 For GrIS atmosphere (1km, polar stereographic EPSG:3413):
 
 .. code-block:: none
@@ -183,6 +199,28 @@ values are:
    # Base path to observational climatology data
    base_path_climatology = /path/to/ISMIP7/forcing/AIS/obs/zhou_annual_06_nov
 
+   # config options for ismip7 fracture (Path C, ice shelf collapse) forcing
+   [ismip7_fracture]
+
+   # Remapping method for the ice shelf collapse mask (Path C).
+   # neareststod preserves the 0/1 mask values
+   method_remap_shelf_collapse = neareststod
+
+   # Remapping method for the excess meltwater field (Path A), a flux
+   method_remap_excess_melt = conserve
+
+   # Remapping method for the supraglacial lake properties (Path B)
+   method_remap_lake_properties = bilinear
+
+   # Version subdirectory of the fracture forcing data
+   version = v2
+
+   # Start year for processing
+   start_year = 1850
+
+   # End year for processing
+   end_year = 2014
+
 All ``NotAvailable`` options must be overridden in a user config file passed
 at setup time (e.g., ``compass setup ... -f my_ismip7.cfg``).
 
@@ -253,3 +291,53 @@ but without a Time dimension, producing a single static file.
 For **GrIS**, thermal forcing is 2D (depth-averaged), with monthly temporal
 resolution and yearly input files. The output variable is
 ``ismip6_2dThermalForcing``.
+
+.. _landice_ismip7_forcing_fracture:
+
+fracture
+--------
+
+The ``landice/ismip7_forcing/fracture`` test case processes the ISMIP7
+surface-melt-driven ice shelf collapse forcing (AIS only). It implements the
+three ISMIP7 pathways as separate steps, each remapping annual fields from
+the native 8km polar stereographic grid onto the MALI unstructured mesh.
+
+All three source files are discovered from the ``fracture/{version}/``
+subdirectory of ``base_path_ismip7``.
+
+* **process_excess_melt** (Path A): Remaps the excess meltwater field
+  (melt + rain after firn air content depletion), matching
+  ``excess_melt_*.nc``. The output variable is ``ismip7ExcessMelt``
+  (mm w.e. yr-1) and is written to
+  ``{output_base_path}/excess_melt/{model}_{scenario}/``. Conservative
+  remapping is used by default since this is a flux. This source file has no
+  ``x``/``y`` coordinate variables and its array is flipped along the y axis
+  relative to the other fracture files, so the step reconstructs the source
+  grid (borrowing ``x``/``y`` from a sibling fracture file and flipping the
+  data to match) before remapping.
+
+* **process_lake_properties** (Path B): Remaps the supraglacial lake mean
+  depth and area fraction from the Grau et al. (2025) parameterization,
+  matching ``lake_properties_*.nc``. The output variables are
+  ``ismip7LakeDepth`` (m) and ``ismip7LakeAreaFraction`` (unitless), written
+  to ``{output_base_path}/lake_properties/{model}_{scenario}/``. Bilinear
+  remapping is used by default.
+
+* **process_shelf_collapse** (Path C): Remaps the annual ice shelf collapse
+  mask, matching ``ice_shelf_collapse_mask_*.nc``.
+  An ice shelf grid cell is flagged as collapsed (mask value 1) when excess
+  meltwater, computed after firn air content depletion, exceeds 72.5 mm/yr for
+  10 consecutive years; otherwise the mask value is 0. The mask is applied on
+  floating areas only, similar to ISMIP6, and ice shelves collapse on
+  January 1st of each year. The remapping uses ``neareststod`` by default so
+  that the 0/1 mask values are preserved, and the remapped mask is rounded to
+  0/1. The output variable is ``calvingMask`` with an accompanying ``xtime``
+  variable, and the result is written to
+  ``{output_base_path}/shelf_collapse/{model}_{scenario}/``.
+
+All three pathways produce continuous fields (Paths A and B) or a discrete
+mask (Path C) with an accompanying ``xtime`` variable. The output variable
+names for Paths A and B (``ismip7ExcessMelt``, ``ismip7LakeDepth``,
+``ismip7LakeAreaFraction``) are descriptive placeholders and may need to be
+aligned with the MALI Registry once the corresponding model input fields are
+defined.

From 73cab6a78abe08cbd6dd9d9366c2dc2ef998586c Mon Sep 17 00:00:00 2001
From: Trevor Hillebrand <trhille@lanl.gov>
Date: Thu, 13 Aug 2026 09:07:01 -0700
Subject: [PATCH 2/7] Document fracture test case in developer's guide

---
 .../landice/test_groups/ismip7_forcing.rst    | 45 ++++++++++++++++++-
 1 file changed, 44 insertions(+), 1 deletion(-)

diff --git a/docs/developers_guide/landice/test_groups/ismip7_forcing.rst b/docs/developers_guide/landice/test_groups/ismip7_forcing.rst
index 3191ce9a6b..f59a3c7b20 100644
--- a/docs/developers_guide/landice/test_groups/ismip7_forcing.rst
+++ b/docs/developers_guide/landice/test_groups/ismip7_forcing.rst
@@ -9,7 +9,7 @@ The ``ismip7_forcing`` test group
 the Ice Sheet Model Intercomparison for CMIP7 (ISMIP7) protocol from its
 native polar stereographic grid to the MALI unstructured mesh. The test group
 supports both AIS and GrIS via the ``ice_sheet`` config option. It includes
-two test cases: ``atmosphere`` and ``ocean_thermal``.
+three test cases: ``atmosphere``, ``ocean_thermal``, and ``fracture``.
 
 .. _dev_landice_ismip7_forcing_framework:
 
@@ -121,3 +121,46 @@ For GrIS, the step:
 
 * Remaps 2D monthly thermal forcing
 * Produces ``ismip6_2dThermalForcing`` (dims: Time × nCells)
+
+.. _dev_landice_ismip7_forcing_fracture:
+
+fracture
+~~~~~~~~
+
+The :py:class:`compass.landice.tests.ismip7_forcing.fracture.Fracture`
+test case processes the ISMIP7 surface-melt-driven ice shelf collapse
+forcing (AIS only). It implements the three ISMIP7 pathways as independent
+steps, each discovering its source file from the ``fracture/{version}/``
+subdirectory of ``base_path_ismip7``, building or reusing a mapping file,
+remapping with ``ncremap``, and renaming the result to MALI conventions with
+an accompanying ``xtime`` variable. Per-pathway remapping methods are set in
+the ``[ismip7_fracture]`` config section.
+
+Steps:
+
+* :py:class:`~compass.landice.tests.ismip7_forcing.fracture.process_excess_melt.ProcessExcessMelt`
+  (Path A) — ``excess_melt`` → ``ismip7ExcessMelt``. The excess melt file
+  lacks ``x``/``y`` coordinate variables and its array is flipped along the
+  y axis relative to the other fracture files (it was produced with CDO).
+  The ``_prepare_source_grid()`` method borrows ``x``/``y`` from a sibling
+  fracture file, flips the data to match (raising if the flipped ``lat`` does
+  not match the sibling grid), and writes a reconstructed source file. The
+  field is then extrapolated (nearest neighbor, filling NaNs) and remapped
+  conservatively by default (it is a flux).
+* :py:class:`~compass.landice.tests.ismip7_forcing.fracture.process_lake_properties.ProcessLakeProperties`
+  (Path B) — ``lake_depth`` → ``ismip7LakeDepth`` and
+  ``fraction_lake_area`` → ``ismip7LakeAreaFraction``. Both variables are
+  extrapolated and remapped in a single ``ncremap`` call (bilinear by
+  default).
+* :py:class:`~compass.landice.tests.ismip7_forcing.fracture.process_shelf_collapse.ProcessShelfCollapse`
+  (Path C) — ``mask`` → ``calvingMask``. Remapped with ``neareststod`` by
+  default and rounded to 0/1 so the discrete collapse mask is preserved.
+
+The annual source fields use an integer ``year``/``time`` coordinate with
+``units="year"`` (not CF-compliant), so each step opens the data with
+``decode_times=False`` and constructs ``xtime`` at January 1st of each year.
+
+The output variable names for Paths A and B (``ismip7ExcessMelt``,
+``ismip7LakeDepth``, ``ismip7LakeAreaFraction``) are descriptive placeholders
+and may need to be aligned with the MALI Registry once the corresponding
+model input fields are defined.

From d8b1c111df9d392f1db46267201e1fe28051a9e4 Mon Sep 17 00:00:00 2001
From: Trevor Hillebrand <trhille@lanl.gov>
Date: Thu, 13 Aug 2026 09:26:53 -0700
Subject: [PATCH 3/7] Allow skipping fracture pathways via None remap method

Each fracture step returns early without processing its file when its
[ismip7_fracture] remapping-method option is set to None. This lets a
user process only the pathways whose source files are available,
without adding separate enable/disable flags.
---
 .../tests/ismip7_forcing/fracture/process_excess_melt.py  | 6 ++++++
 .../ismip7_forcing/fracture/process_lake_properties.py    | 6 ++++++
 .../ismip7_forcing/fracture/process_shelf_collapse.py     | 6 ++++++
 compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg   | 7 ++++---
 .../landice/tests/ismip7_forcing/ismip7_forcing_test.cfg  | 7 ++++---
 .../landice/test_groups/ismip7_forcing.rst                | 4 +++-
 docs/users_guide/landice/test_groups/ismip7_forcing.rst   | 8 ++++++++
 7 files changed, 37 insertions(+), 7 deletions(-)

diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py b/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
index b7af781031..bb23428574 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
@@ -74,6 +74,12 @@ def run(self):
         start_year = section.getint("start_year")
         end_year = section.getint("end_year")
 
+        # Skip this pathway if no remapping method is requested
+        if method_remap.lower() == "none":
+            logger.info("method_remap_excess_melt is None; skipping excess "
+                        "melt (Path A) processing.")
+            return
+
         # Discover the excess melt file
         input_path = os.path.join(base_path_ismip7, "fracture", version)
         file_pattern = "excess_melt_*.nc"
diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py b/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
index d5e4677ab5..2a60443c5c 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
@@ -82,6 +82,12 @@ def run(self):
         start_year = section.getint("start_year")
         end_year = section.getint("end_year")
 
+        # Skip this pathway if no remapping method is requested
+        if method_remap.lower() == "none":
+            logger.info("method_remap_lake_properties is None; skipping lake "
+                        "properties (Path B) processing.")
+            return
+
         # Discover the lake properties file
         input_path = os.path.join(base_path_ismip7, "fracture", version)
         file_pattern = "lake_properties_*.nc"
diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py b/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py
index 02b2d0312c..5924f2a972 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py
@@ -68,6 +68,12 @@ def run(self):
         start_year = section.getint("start_year")
         end_year = section.getint("end_year")
 
+        # Skip this pathway if no remapping method is requested
+        if method_remap.lower() == "none":
+            logger.info("method_remap_shelf_collapse is None; skipping ice "
+                        "shelf collapse mask (Path C) processing.")
+            return
+
         # Discover the ice shelf collapse mask file
         input_path = os.path.join(base_path_ismip7, "fracture", version)
         file_pattern = "ice_shelf_collapse_mask_*.nc"
diff --git a/compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg b/compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg
index a6eadbd26e..f02b79067a 100644
--- a/compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg
+++ b/compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg
@@ -73,16 +73,17 @@ base_path_climatology = /path/to/ISMIP7/forcing/AIS/obs/zhou_annual_06_nov
 # config options for ismip7 fracture (Path C, ice shelf collapse) forcing
 [ismip7_fracture]
 
+# Remapping method for each pathway. Set a method to None to skip processing
+# that pathway's file. Options: bilinear, neareststod, conserve, None
+
 # Remapping method for the ice shelf collapse mask (Path C). neareststod
-# preserves the 0/1 mask values. Options: bilinear, neareststod, conserve
+# preserves the 0/1 mask values.
 method_remap_shelf_collapse = neareststod
 
 # Remapping method for the excess meltwater field (Path A), a flux.
-# Options: bilinear, neareststod, conserve
 method_remap_excess_melt = conserve
 
 # Remapping method for the supraglacial lake properties (Path B).
-# Options: bilinear, neareststod, conserve
 method_remap_lake_properties = bilinear
 
 # Version subdirectory of the fracture forcing data
diff --git a/compass/landice/tests/ismip7_forcing/ismip7_forcing_test.cfg b/compass/landice/tests/ismip7_forcing/ismip7_forcing_test.cfg
index 28915fe41d..ab282b2e06 100644
--- a/compass/landice/tests/ismip7_forcing/ismip7_forcing_test.cfg
+++ b/compass/landice/tests/ismip7_forcing/ismip7_forcing_test.cfg
@@ -73,16 +73,17 @@ base_path_climatology = /global/cfs/cdirs/m4288/users/trhille/ISMIP7/forcing/AIS
 # config options for ismip7 fracture (Path C, ice shelf collapse) forcing
 [ismip7_fracture]
 
+# Remapping method for each pathway. Set a method to None to skip processing
+# that pathway's file. Options: bilinear, neareststod, conserve, None
+
 # Remapping method for the ice shelf collapse mask (Path C). neareststod
-# preserves the 0/1 mask values. Options: bilinear, neareststod, conserve
+# preserves the 0/1 mask values.
 method_remap_shelf_collapse = neareststod
 
 # Remapping method for the excess meltwater field (Path A), a flux.
-# Options: bilinear, neareststod, conserve
 method_remap_excess_melt = conserve
 
 # Remapping method for the supraglacial lake properties (Path B).
-# Options: bilinear, neareststod, conserve
 method_remap_lake_properties = bilinear
 
 # Version subdirectory of the fracture forcing data
diff --git a/docs/developers_guide/landice/test_groups/ismip7_forcing.rst b/docs/developers_guide/landice/test_groups/ismip7_forcing.rst
index f59a3c7b20..34648255fa 100644
--- a/docs/developers_guide/landice/test_groups/ismip7_forcing.rst
+++ b/docs/developers_guide/landice/test_groups/ismip7_forcing.rst
@@ -134,7 +134,9 @@ steps, each discovering its source file from the ``fracture/{version}/``
 subdirectory of ``base_path_ismip7``, building or reusing a mapping file,
 remapping with ``ncremap``, and renaming the result to MALI conventions with
 an accompanying ``xtime`` variable. Per-pathway remapping methods are set in
-the ``[ismip7_fracture]`` config section.
+the ``[ismip7_fracture]`` config section. Setting a pathway's remapping-method
+option to ``None`` causes that step to return early without processing its
+file, which is useful when only some pathway source files are available.
 
 Steps:
 
diff --git a/docs/users_guide/landice/test_groups/ismip7_forcing.rst b/docs/users_guide/landice/test_groups/ismip7_forcing.rst
index d134093dce..efb738d86a 100644
--- a/docs/users_guide/landice/test_groups/ismip7_forcing.rst
+++ b/docs/users_guide/landice/test_groups/ismip7_forcing.rst
@@ -202,6 +202,9 @@ values are:
    # config options for ismip7 fracture (Path C, ice shelf collapse) forcing
    [ismip7_fracture]
 
+   # Remapping method for each pathway. Set a method to None to skip
+   # processing that pathway's file.
+
    # Remapping method for the ice shelf collapse mask (Path C).
    # neareststod preserves the 0/1 mask values
    method_remap_shelf_collapse = neareststod
@@ -305,6 +308,11 @@ the native 8km polar stereographic grid onto the MALI unstructured mesh.
 All three source files are discovered from the ``fracture/{version}/``
 subdirectory of ``base_path_ismip7``.
 
+Each pathway is run independently and can be skipped by setting its
+remapping-method config option to ``None`` in the ``[ismip7_fracture]``
+section (for example, ``method_remap_excess_melt = None`` skips Path A). This
+is useful when only some of the pathway source files are available.
+
 * **process_excess_melt** (Path A): Remaps the excess meltwater field
   (melt + rain after firn air content depletion), matching
   ``excess_melt_*.nc``. The output variable is ``ismip7ExcessMelt``

From 6f569c1d223da382bd8915c73b6c7c36af289118 Mon Sep 17 00:00:00 2001
From: Trevor Hillebrand <trhille@lanl.gov>
Date: Thu, 13 Aug 2026 09:27:15 -0700
Subject: [PATCH 4/7] Fix fracture usage description to cover all three
 pathways

---
 docs/users_guide/landice/test_groups/ismip7_forcing.rst | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/docs/users_guide/landice/test_groups/ismip7_forcing.rst b/docs/users_guide/landice/test_groups/ismip7_forcing.rst
index efb738d86a..fbd5289917 100644
--- a/docs/users_guide/landice/test_groups/ismip7_forcing.rst
+++ b/docs/users_guide/landice/test_groups/ismip7_forcing.rst
@@ -51,7 +51,8 @@ To use this test group, users need to:
 5. Run the ``ocean_thermal`` test case for each model and scenario combination.
 
 6. Run the ``fracture`` test case (AIS only) for each model and scenario
-   combination to process the Path C ice shelf collapse mask.
+   combination to process the surface-melt-driven ice shelf collapse
+   pathways (excess melt, lake properties, and the ice shelf collapse mask).
 
 Example user config files are provided in the source tree for local testing:
 

From 6780a2c511a31d807747193d4f3cb63772edb382 Mon Sep 17 00:00:00 2001
From: Trevor Hillebrand <trhille@lanl.gov>
Date: Thu, 13 Aug 2026 10:18:12 -0700
Subject: [PATCH 5/7] Fix time coordinate KeyError in fracture rename step

Capture the integer years from the time coordinate before renaming the
time dimension to Time, since the rename also renames the coordinate
variable. Drop the leftover Time coordinate from the output.
---
 .../tests/ismip7_forcing/fracture/process_excess_melt.py    | 6 ++++--
 .../ismip7_forcing/fracture/process_lake_properties.py      | 6 ++++--
 .../tests/ismip7_forcing/fracture/process_shelf_collapse.py | 6 ++++--
 3 files changed, 12 insertions(+), 6 deletions(-)

diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py b/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
index bb23428574..95334f082e 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
@@ -251,6 +251,9 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
         # not CF-compliant, so disable time decoding.
         ds = xr.open_dataset(remapped_file, decode_times=False)
 
+        # Capture integer years before the time coordinate is renamed
+        years = ds["time"].values.astype(int)
+
         rename_dims = {}
         if "ncol" in ds.dims:
             rename_dims["ncol"] = "nCells"
@@ -263,7 +266,6 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
             ds = ds.rename({"excess_melt": "ismip7ExcessMelt"})
 
         # Restrict to the requested year range
-        years = ds["time"].values.astype(int)
         keep = (years >= start_year) & (years <= end_year)
         ds = ds.isel(Time=keep)
         years = years[keep]
@@ -279,7 +281,7 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
         }
 
         vars_to_drop = [v for v in ["lat_vertices", "lon_vertices", "lat",
-                                    "lon", "area", "time"]
+                                    "lon", "area", "Time"]
                         if v in ds]
         if vars_to_drop:
             ds = ds.drop_vars(vars_to_drop)
diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py b/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
index 2a60443c5c..7fdff49f16 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
@@ -182,6 +182,9 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
         # not CF-compliant, so disable time decoding.
         ds = xr.open_dataset(remapped_file, decode_times=False)
 
+        # Capture integer years before the time coordinate is renamed
+        years = ds["time"].values.astype(int)
+
         rename_dims = {}
         if "ncol" in ds.dims:
             rename_dims["ncol"] = "nCells"
@@ -196,7 +199,6 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
         ds = ds.rename(rename_vars)
 
         # Restrict to the requested year range
-        years = ds["time"].values.astype(int)
         keep = (years >= start_year) & (years <= end_year)
         ds = ds.isel(Time=keep)
         years = years[keep]
@@ -215,7 +217,7 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
                 }
 
         vars_to_drop = [v for v in ["lat_vertices", "lon_vertices", "lat",
-                                    "lon", "area", "time"]
+                                    "lon", "area", "Time"]
                         if v in ds]
         if vars_to_drop:
             ds = ds.drop_vars(vars_to_drop)
diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py b/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py
index 5924f2a972..5ffc5a65c6 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py
@@ -160,6 +160,9 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
         # which is not CF-compliant, so disable time decoding.
         ds = xr.open_dataset(remapped_file, decode_times=False)
 
+        # Capture integer years before the time coordinate is renamed
+        years = ds["time"].values.astype(int)
+
         # Rename dimensions to MALI conventions
         rename_dims = {}
         if "ncol" in ds.dims:
@@ -174,7 +177,6 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
             ds = ds.rename({"mask": "calvingMask"})
 
         # Restrict to the requested year range
-        years = ds["time"].values.astype(int)
         keep = (years >= start_year) & (years <= end_year)
         ds = ds.isel(Time=keep)
         years = years[keep]
@@ -196,7 +198,7 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
 
         # Drop auxiliary variables carried over from remapping
         vars_to_drop = [v for v in ["lat_vertices", "lon_vertices", "lat",
-                                    "lon", "area", "time"]
+                                    "lon", "area", "Time"]
                         if v in ds]
         if vars_to_drop:
             ds = ds.drop_vars(vars_to_drop)

From 55cfdbb4e7c8673ef589fa0d50b83df667444b58 Mon Sep 17 00:00:00 2001
From: Trevor Hillebrand <trhille@lanl.gov>
Date: Thu, 13 Aug 2026 10:43:19 -0700
Subject: [PATCH 6/7] Convert ismip7ExcessMelt to SI units (kg m-2 s-1)

Convert the excess melt field from mm w.e. yr-1 to kg m-2 s-1 using
1 mm w.e. = 1 kg m-2 and a 365-day year, consistent with other MALI
mass fluxes.
---
 .../tests/ismip7_forcing/fracture/process_excess_melt.py   | 7 ++++++-
 docs/users_guide/landice/test_groups/ismip7_forcing.rst    | 2 +-
 2 files changed, 7 insertions(+), 2 deletions(-)

diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py b/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
index 95334f082e..f338695b00 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
@@ -275,9 +275,14 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
         ds["xtime"] = ("Time", xtime)
         ds["xtime"] = ds.xtime.astype("S")
 
+        # Convert from mm w.e. yr-1 to SI units of kg m-2 s-1
+        # (1 mm w.e. = 1 kg m-2; 365-day year, as used elsewhere in MALI)
+        seconds_per_year = 365.0 * 24.0 * 3600.0
+        ds["ismip7ExcessMelt"] = ds["ismip7ExcessMelt"] / seconds_per_year
+
         ds["ismip7ExcessMelt"].attrs = {
             "long_name": "excess meltwater after firn air content depletion",
-            "units": "mm w.e. yr-1",
+            "units": "kg m-2 s-1",
         }
 
         vars_to_drop = [v for v in ["lat_vertices", "lon_vertices", "lat",
diff --git a/docs/users_guide/landice/test_groups/ismip7_forcing.rst b/docs/users_guide/landice/test_groups/ismip7_forcing.rst
index fbd5289917..fdec1eee14 100644
--- a/docs/users_guide/landice/test_groups/ismip7_forcing.rst
+++ b/docs/users_guide/landice/test_groups/ismip7_forcing.rst
@@ -317,7 +317,7 @@ is useful when only some of the pathway source files are available.
 * **process_excess_melt** (Path A): Remaps the excess meltwater field
   (melt + rain after firn air content depletion), matching
   ``excess_melt_*.nc``. The output variable is ``ismip7ExcessMelt``
-  (mm w.e. yr-1) and is written to
+  (converted from mm w.e. yr-1 to SI units of kg m-2 s-1) and is written to
   ``{output_base_path}/excess_melt/{model}_{scenario}/``. Conservative
   remapping is used by default since this is a flux. This source file has no
   ``x``/``y`` coordinate variables and its array is flipped along the y axis

From 87d40e361000841aca4fd79a4dd1fd0428de68be Mon Sep 17 00:00:00 2001
From: Trevor Hillebrand <trhille@lanl.gov>
Date: Thu, 13 Aug 2026 11:36:34 -0700
Subject: [PATCH 7/7] Consolidate shared fracture remap helpers into
 remap_utils

Move the duplicated extrapolation and rename/xtime boilerplate from the
three fracture steps into a shared remap_utils module (extrapolate_source,
open_rename_and_trim, add_xtime_and_write). The step-specific logic
(unit conversion, mask rounding, per-variable attrs) stays in each step.
---
 .../fracture/process_excess_melt.py           |  99 ++----------
 .../fracture/process_lake_properties.py       | 105 ++----------
 .../fracture/process_shelf_collapse.py        |  51 ++----
 .../ismip7_forcing/fracture/remap_utils.py    | 149 ++++++++++++++++++
 docs/developers_guide/landice/api.rst         |   3 +
 .../landice/test_groups/ismip7_forcing.rst    |   8 +
 6 files changed, 190 insertions(+), 225 deletions(-)
 create mode 100644 compass/landice/tests/ismip7_forcing/fracture/remap_utils.py

diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py b/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
index f338695b00..fbf3f75b83 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
@@ -6,11 +6,15 @@
 import xarray as xr
 from mpas_tools.io import write_netcdf
 from mpas_tools.logging import check_call
-from scipy.ndimage import distance_transform_edt
 
 from compass.landice.tests.ismip7_forcing.create_mapfile import (
     build_mapping_file,
 )
+from compass.landice.tests.ismip7_forcing.fracture.remap_utils import (
+    add_xtime_and_write,
+    extrapolate_source,
+    open_rename_and_trim,
+)
 from compass.step import Step
 
 
@@ -118,8 +122,7 @@ def run(self):
         # Extrapolate fill values on the source grid before remapping so
         # they don't pollute neighboring cells during interpolation
         extrap_file = f"extrap_{basename}"
-        self._extrapolate_source(gridded_file, extrap_file, "excess_melt",
-                                 logger)
+        extrapolate_source(gridded_file, extrap_file, "excess_melt", logger)
 
         # Remap the excess melt onto the MALI mesh
         remapped_file = f"remapped_{basename}"
@@ -247,33 +250,9 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
         end_year : int
             Last year (inclusive) to retain
         """
-        # The time coordinate has units="year" (integer years), which is
-        # not CF-compliant, so disable time decoding.
-        ds = xr.open_dataset(remapped_file, decode_times=False)
-
-        # Capture integer years before the time coordinate is renamed
-        years = ds["time"].values.astype(int)
-
-        rename_dims = {}
-        if "ncol" in ds.dims:
-            rename_dims["ncol"] = "nCells"
-        if "time" in ds.dims:
-            rename_dims["time"] = "Time"
-        if rename_dims:
-            ds = ds.rename(rename_dims)
-
-        if "excess_melt" in ds:
-            ds = ds.rename({"excess_melt": "ismip7ExcessMelt"})
-
-        # Restrict to the requested year range
-        keep = (years >= start_year) & (years <= end_year)
-        ds = ds.isel(Time=keep)
-        years = years[keep]
-
-        # Excess melt is an annual field applied at the start of each year
-        xtime = [f"{int(yr):04d}-01-01_00:00:00".ljust(64) for yr in years]
-        ds["xtime"] = ("Time", xtime)
-        ds["xtime"] = ds.xtime.astype("S")
+        ds, years = open_rename_and_trim(
+            remapped_file, {"excess_melt": "ismip7ExcessMelt"},
+            start_year, end_year)
 
         # Convert from mm w.e. yr-1 to SI units of kg m-2 s-1
         # (1 mm w.e. = 1 kg m-2; 365-day year, as used elsewhere in MALI)
@@ -285,62 +264,4 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
             "units": "kg m-2 s-1",
         }
 
-        vars_to_drop = [v for v in ["lat_vertices", "lon_vertices", "lat",
-                                    "lon", "area", "Time"]
-                        if v in ds]
-        if vars_to_drop:
-            ds = ds.drop_vars(vars_to_drop)
-
-        write_netcdf(ds, output_file)
-        ds.close()
-
-    def _extrapolate_source(self, input_file, output_file, varname, logger):
-        """
-        Extrapolate fill/missing values on the source polar stereographic
-        grid using nearest-neighbor via distance_transform_edt. This must
-        be done before remapping so that fill values don't contaminate the
-        interpolation stencil.
-
-        Parameters
-        ----------
-        input_file : str
-            Path to the input NetCDF file on the source grid
-
-        output_file : str
-            Path to write the extrapolated file
-
-        varname : str
-            Name of the variable to extrapolate
-
-        logger : logging.Logger
-            Logger for status messages
-        """
-        logger.info(f"    Extrapolating fill values on source grid: "
-                    f"{os.path.basename(input_file)}")
-
-        ds = xr.open_dataset(input_file, decode_times=False)
-        data = ds[varname]
-
-        values = data.values.copy()
-        non_spatial_shape = values.shape[:-2]
-
-        for idx in np.ndindex(non_spatial_shape):
-            slab = values[idx]
-            valid_mask = np.isfinite(slab)
-            if valid_mask.all() or not valid_mask.any():
-                continue
-            nearest_inds = distance_transform_edt(
-                ~valid_mask, return_distances=False, return_indices=True)
-            invalid = ~valid_mask
-            values[idx][invalid] = slab[
-                nearest_inds[0, invalid],
-                nearest_inds[1, invalid]]
-
-        ds[varname] = (data.dims, values)
-        ds[varname].attrs = data.attrs
-
-        if "_FillValue" in ds[varname].encoding:
-            del ds[varname].encoding["_FillValue"]
-
-        write_netcdf(ds, output_file)
-        ds.close()
+        add_xtime_and_write(ds, years, output_file)
diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py b/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
index 7fdff49f16..10df254c69 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
@@ -2,15 +2,16 @@
 import os
 import shutil
 
-import numpy as np
-import xarray as xr
-from mpas_tools.io import write_netcdf
 from mpas_tools.logging import check_call
-from scipy.ndimage import distance_transform_edt
 
 from compass.landice.tests.ismip7_forcing.create_mapfile import (
     build_mapping_file,
 )
+from compass.landice.tests.ismip7_forcing.fracture.remap_utils import (
+    add_xtime_and_write,
+    extrapolate_source,
+    open_rename_and_trim,
+)
 from compass.step import Step
 
 
@@ -122,8 +123,8 @@ def run(self):
         # Extrapolate fill values on the source grid before remapping so
         # they don't pollute neighboring cells during interpolation
         extrap_file = f"extrap_{basename}"
-        self._extrapolate_source(input_file, extrap_file,
-                                 list(self._variables.keys()), logger)
+        extrapolate_source(input_file, extrap_file,
+                           list(self._variables.keys()), logger)
 
         # Remap both lake property variables onto the MALI mesh
         remapped_file = f"remapped_{basename}"
@@ -178,35 +179,10 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
         end_year : int
             Last year (inclusive) to retain
         """
-        # The time coordinate has units="year" (integer years), which is
-        # not CF-compliant, so disable time decoding.
-        ds = xr.open_dataset(remapped_file, decode_times=False)
-
-        # Capture integer years before the time coordinate is renamed
-        years = ds["time"].values.astype(int)
-
-        rename_dims = {}
-        if "ncol" in ds.dims:
-            rename_dims["ncol"] = "nCells"
-        if "time" in ds.dims:
-            rename_dims["time"] = "Time"
-        if rename_dims:
-            ds = ds.rename(rename_dims)
-
         rename_vars = {src: info["mali_name"]
-                       for src, info in self._variables.items()
-                       if src in ds}
-        ds = ds.rename(rename_vars)
-
-        # Restrict to the requested year range
-        keep = (years >= start_year) & (years <= end_year)
-        ds = ds.isel(Time=keep)
-        years = years[keep]
-
-        # Lake properties are annual fields applied at the start of each year
-        xtime = [f"{int(yr):04d}-01-01_00:00:00".ljust(64) for yr in years]
-        ds["xtime"] = ("Time", xtime)
-        ds["xtime"] = ds.xtime.astype("S")
+                       for src, info in self._variables.items()}
+        ds, years = open_rename_and_trim(remapped_file, rename_vars,
+                                         start_year, end_year)
 
         for info in self._variables.values():
             mali_name = info["mali_name"]
@@ -216,63 +192,4 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
                     "units": info["units"],
                 }
 
-        vars_to_drop = [v for v in ["lat_vertices", "lon_vertices", "lat",
-                                    "lon", "area", "Time"]
-                        if v in ds]
-        if vars_to_drop:
-            ds = ds.drop_vars(vars_to_drop)
-
-        write_netcdf(ds, output_file)
-        ds.close()
-
-    def _extrapolate_source(self, input_file, output_file, varnames, logger):
-        """
-        Extrapolate fill/missing values on the source polar stereographic
-        grid using nearest-neighbor via distance_transform_edt. This must
-        be done before remapping so that fill values don't contaminate the
-        interpolation stencil.
-
-        Parameters
-        ----------
-        input_file : str
-            Path to the input NetCDF file on the source grid
-
-        output_file : str
-            Path to write the extrapolated file
-
-        varnames : list of str
-            Names of the variables to extrapolate
-
-        logger : logging.Logger
-            Logger for status messages
-        """
-        logger.info(f"    Extrapolating fill values on source grid: "
-                    f"{os.path.basename(input_file)}")
-
-        ds = xr.open_dataset(input_file, decode_times=False)
-
-        for varname in varnames:
-            data = ds[varname]
-            values = data.values.copy()
-            non_spatial_shape = values.shape[:-2]
-
-            for idx in np.ndindex(non_spatial_shape):
-                slab = values[idx]
-                valid_mask = np.isfinite(slab)
-                if valid_mask.all() or not valid_mask.any():
-                    continue
-                nearest_inds = distance_transform_edt(
-                    ~valid_mask, return_distances=False,
-                    return_indices=True)
-                invalid = ~valid_mask
-                values[idx][invalid] = slab[
-                    nearest_inds[0, invalid],
-                    nearest_inds[1, invalid]]
-
-            ds[varname] = (data.dims, values)
-            ds[varname].attrs = data.attrs
-            if "_FillValue" in ds[varname].encoding:
-                del ds[varname].encoding["_FillValue"]
-
-        write_netcdf(ds, output_file)
-        ds.close()
+        add_xtime_and_write(ds, years, output_file)
diff --git a/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py b/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py
index 5ffc5a65c6..e9a8e1cfb5 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_shelf_collapse.py
@@ -2,13 +2,15 @@
 import os
 import shutil
 
-import xarray as xr
-from mpas_tools.io import write_netcdf
 from mpas_tools.logging import check_call
 
 from compass.landice.tests.ismip7_forcing.create_mapfile import (
     build_mapping_file,
 )
+from compass.landice.tests.ismip7_forcing.fracture.remap_utils import (
+    add_xtime_and_write,
+    open_rename_and_trim,
+)
 from compass.step import Step
 
 
@@ -156,52 +158,17 @@ def _rename_to_mali_vars(self, remapped_file, output_file,
         end_year : int
             Last year (inclusive) to retain
         """
-        # The collapse mask time coordinate has units="year" (integer years),
-        # which is not CF-compliant, so disable time decoding.
-        ds = xr.open_dataset(remapped_file, decode_times=False)
-
-        # Capture integer years before the time coordinate is renamed
-        years = ds["time"].values.astype(int)
-
-        # Rename dimensions to MALI conventions
-        rename_dims = {}
-        if "ncol" in ds.dims:
-            rename_dims["ncol"] = "nCells"
-        if "time" in ds.dims:
-            rename_dims["time"] = "Time"
-        if rename_dims:
-            ds = ds.rename(rename_dims)
-
-        # Rename variable
-        if "mask" in ds:
-            ds = ds.rename({"mask": "calvingMask"})
-
-        # Restrict to the requested year range
-        keep = (years >= start_year) & (years <= end_year)
-        ds = ds.isel(Time=keep)
-        years = years[keep]
+        ds, years = open_rename_and_trim(remapped_file,
+                                         {"mask": "calvingMask"},
+                                         start_year, end_year)
 
         # Round the remapped mask to 0/1 and store as integers. Ice shelves
         # collapse on January 1st, so the mask is applied at the start of
         # each year.
-        calving_mask = (ds["calvingMask"] >= 0.5).astype(int)
-        ds["calvingMask"] = calving_mask
-
-        # Add xtime variable, one entry per year at January 1st
-        xtime = [f"{int(yr):04d}-01-01_00:00:00".ljust(64) for yr in years]
-        ds["xtime"] = ("Time", xtime)
-        ds["xtime"] = ds.xtime.astype("S")
+        ds["calvingMask"] = (ds["calvingMask"] >= 0.5).astype(int)
 
         ds["calvingMask"].attrs = {
             "long_name": "ice shelf collapse mask (1 = collapse)",
         }
 
-        # Drop auxiliary variables carried over from remapping
-        vars_to_drop = [v for v in ["lat_vertices", "lon_vertices", "lat",
-                                    "lon", "area", "Time"]
-                        if v in ds]
-        if vars_to_drop:
-            ds = ds.drop_vars(vars_to_drop)
-
-        write_netcdf(ds, output_file)
-        ds.close()
+        add_xtime_and_write(ds, years, output_file)
diff --git a/compass/landice/tests/ismip7_forcing/fracture/remap_utils.py b/compass/landice/tests/ismip7_forcing/fracture/remap_utils.py
new file mode 100644
index 0000000000..8da5de4ae4
--- /dev/null
+++ b/compass/landice/tests/ismip7_forcing/fracture/remap_utils.py
@@ -0,0 +1,149 @@
+"""
+Shared helpers for remapping ISMIP7 fracture forcing data to the MALI mesh.
+"""
+import os
+
+import numpy as np
+import xarray as xr
+from mpas_tools.io import write_netcdf
+from scipy.ndimage import distance_transform_edt
+
+
+def extrapolate_source(input_file, output_file, varnames, logger):
+    """
+    Extrapolate fill/missing values on the source polar stereographic grid
+    using nearest-neighbor via ``distance_transform_edt``. This must be done
+    before remapping so that fill values don't contaminate the interpolation
+    stencil.
+
+    Parameters
+    ----------
+    input_file : str
+        Path to the input NetCDF file on the source grid
+
+    output_file : str
+        Path to write the extrapolated file
+
+    varnames : str or list of str
+        Name(s) of the variable(s) to extrapolate
+
+    logger : logging.Logger
+        Logger for status messages
+    """
+    if isinstance(varnames, str):
+        varnames = [varnames]
+
+    logger.info(f"    Extrapolating fill values on source grid: "
+                f"{os.path.basename(input_file)}")
+
+    ds = xr.open_dataset(input_file, decode_times=False)
+
+    for varname in varnames:
+        data = ds[varname]
+        values = data.values.copy()
+        non_spatial_shape = values.shape[:-2]
+
+        for idx in np.ndindex(non_spatial_shape):
+            slab = values[idx]
+            valid_mask = np.isfinite(slab)
+            if valid_mask.all() or not valid_mask.any():
+                continue
+            nearest_inds = distance_transform_edt(
+                ~valid_mask, return_distances=False, return_indices=True)
+            invalid = ~valid_mask
+            values[idx][invalid] = slab[
+                nearest_inds[0, invalid],
+                nearest_inds[1, invalid]]
+
+        ds[varname] = (data.dims, values)
+        ds[varname].attrs = data.attrs
+        if "_FillValue" in ds[varname].encoding:
+            del ds[varname].encoding["_FillValue"]
+
+    write_netcdf(ds, output_file)
+    ds.close()
+
+
+def open_rename_and_trim(remapped_file, rename_vars, start_year, end_year):
+    """
+    Open a remapped file, rename dimensions/variables to MALI conventions,
+    and restrict to the requested year range.
+
+    Parameters
+    ----------
+    remapped_file : str
+        Data remapped onto the MALI mesh
+
+    rename_vars : dict
+        Mapping of source variable names to MALI variable names
+
+    start_year : int
+        First year (inclusive) to retain
+
+    end_year : int
+        Last year (inclusive) to retain
+
+    Returns
+    -------
+    ds : xarray.Dataset
+        The renamed and trimmed dataset
+
+    years : numpy.ndarray
+        The integer years retained
+    """
+    # The time coordinate has units="year" (integer years), which is not
+    # CF-compliant, so disable time decoding.
+    ds = xr.open_dataset(remapped_file, decode_times=False)
+
+    # Capture integer years before the time coordinate is renamed
+    years = ds["time"].values.astype(int)
+
+    rename_dims = {}
+    if "ncol" in ds.dims:
+        rename_dims["ncol"] = "nCells"
+    if "time" in ds.dims:
+        rename_dims["time"] = "Time"
+    if rename_dims:
+        ds = ds.rename(rename_dims)
+
+    rename_vars = {src: dst for src, dst in rename_vars.items() if src in ds}
+    if rename_vars:
+        ds = ds.rename(rename_vars)
+
+    # Restrict to the requested year range
+    keep = (years >= start_year) & (years <= end_year)
+    ds = ds.isel(Time=keep)
+    years = years[keep]
+
+    return ds, years
+
+
+def add_xtime_and_write(ds, years, output_file):
+    """
+    Add an ``xtime`` variable (January 1st of each year), drop auxiliary
+    remapping variables, and write the dataset.
+
+    Parameters
+    ----------
+    ds : xarray.Dataset
+        The dataset to finalize (annual fields applied at the start of the
+        year)
+
+    years : numpy.ndarray
+        The integer years, one per Time index
+
+    output_file : str
+        Output file path
+    """
+    xtime = [f"{int(yr):04d}-01-01_00:00:00".ljust(64) for yr in years]
+    ds["xtime"] = ("Time", xtime)
+    ds["xtime"] = ds.xtime.astype("S")
+
+    vars_to_drop = [v for v in ["lat_vertices", "lon_vertices", "lat",
+                                "lon", "area", "Time"]
+                    if v in ds]
+    if vars_to_drop:
+        ds = ds.drop_vars(vars_to_drop)
+
+    write_netcdf(ds, output_file)
+    ds.close()
diff --git a/docs/developers_guide/landice/api.rst b/docs/developers_guide/landice/api.rst
index 100a147ad8..d2cd995145 100644
--- a/docs/developers_guide/landice/api.rst
+++ b/docs/developers_guide/landice/api.rst
@@ -412,6 +412,9 @@ ismip7_forcing
 
    fracture.Fracture
    fracture.Fracture.configure
+   fracture.remap_utils.extrapolate_source
+   fracture.remap_utils.open_rename_and_trim
+   fracture.remap_utils.add_xtime_and_write
    fracture.process_excess_melt.ProcessExcessMelt
    fracture.process_excess_melt.ProcessExcessMelt.setup
    fracture.process_excess_melt.ProcessExcessMelt.run
diff --git a/docs/developers_guide/landice/test_groups/ismip7_forcing.rst b/docs/developers_guide/landice/test_groups/ismip7_forcing.rst
index 34648255fa..b535bf081a 100644
--- a/docs/developers_guide/landice/test_groups/ismip7_forcing.rst
+++ b/docs/developers_guide/landice/test_groups/ismip7_forcing.rst
@@ -162,6 +162,14 @@ The annual source fields use an integer ``year``/``time`` coordinate with
 ``units="year"`` (not CF-compliant), so each step opens the data with
 ``decode_times=False`` and constructs ``xtime`` at January 1st of each year.
 
+Shared remapping helpers used by the fracture steps live in
+:py:mod:`compass.landice.tests.ismip7_forcing.fracture.remap_utils`:
+``extrapolate_source`` (nearest-neighbor fill of NaNs on the source grid),
+``open_rename_and_trim`` (open a remapped file, rename dimensions/variables
+to MALI conventions, and restrict to the requested year range), and
+``add_xtime_and_write`` (add the ``xtime`` variable, drop auxiliary remapping
+variables, and write the output).
+
 The output variable names for Paths A and B (``ismip7ExcessMelt``,
 ``ismip7LakeDepth``, ``ismip7LakeAreaFraction``) are descriptive placeholders
 and may need to be aligned with the MALI Registry once the corresponding