diff --git a/compass/landice/tests/ismip7_forcing/atmosphere/process_runoff.py b/compass/landice/tests/ismip7_forcing/atmosphere/process_runoff.py
index e03e85ecd6..801844edf0 100644
--- a/compass/landice/tests/ismip7_forcing/atmosphere/process_runoff.py
+++ b/compass/landice/tests/ismip7_forcing/atmosphere/process_runoff.py
@@ -2,16 +2,15 @@
 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.ice_sheet_params import get_params
+from compass.landice.tests.ismip7_forcing.remap_utils import extrapolate_source
 from compass.step import Step
 
 
@@ -72,6 +71,11 @@ def run(self):
         prefix = params['prefix']
         resolution = params['atm_resolution']
         version = params['atm_version']
+        if params['atm_model'] is not None:
+            forcing_group = scenario
+            model = params['atm_model']
+        else:
+            forcing_group = f"{model}_{scenario}"
         input_path = os.path.join(base_path_ismip7, "mrro", version)
         file_pattern = (f"mrro_{prefix}_{model}_{scenario}_"
                         f"SDBN1-{resolution}_{version}_*.nc")
@@ -85,7 +89,11 @@ def run(self):
         # Filter to requested year range
         input_files = []
         for f in all_files:
-            year = int(os.path.basename(f).split("_")[-1].replace(".nc", ""))
+            # skip non-yearly files such as climatology averages (*_avg.nc)
+            token = os.path.basename(f).split("_")[-1].replace(".nc", "")
+            if not token.isdigit():
+                continue
+            year = int(token)
             if start_year <= year <= end_year:
                 input_files.append(f)
 
@@ -123,8 +131,7 @@ def run(self):
             # so they don't pollute neighboring cells during interpolation
             extrap_file = f"extrap_{basename}"
             if not os.path.exists(extrap_file):
-                self._extrapolate_source(input_file, extrap_file, "mrro",
-                                         logger)
+                extrapolate_source(input_file, extrap_file, "mrro", logger)
 
             logger.info(f"  Remapping: {basename}")
             args = ["ncremap",
@@ -152,8 +159,8 @@ def run(self):
                 os.remove(f)
 
         # Place output in appropriate directory
-        output_path = os.path.join(output_base_path, "atmosphere_forcing",
-                                   f"{model}_{scenario}")
+        output_path = os.path.join(output_base_path, forcing_group,
+                                   "atmosphere")
         if not os.path.exists(output_path):
             os.makedirs(output_path)
 
@@ -176,7 +183,8 @@ def _combine_and_rename(self, remapped_files, output_file):
             Output file path
         """
         ds = xr.open_mfdataset(remapped_files, concat_dim="time",
-                               combine="nested", engine="netcdf4")
+                               combine="nested", engine="netcdf4",
+                               drop_variables="time_bnds")
 
         # Rename dimensions to MALI conventions
         rename_dims = {}
@@ -225,56 +233,3 @@ def _combine_and_rename(self, remapped_files, output_file):
             ds = ds.drop_vars("Time")
 
         write_netcdf(ds, output_file)
-
-    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 (e.g., "mrro")
-
-        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, engine="netcdf4")
-        data = ds[varname]
-
-        # Process each time step
-        # Source files have dims like (time, y, x)
-        values = data.values.copy()
-        non_spatial_shape = values.shape[:-2]  # (time,)
-
-        for idx in np.ndindex(non_spatial_shape):
-            slab = values[idx]  # shape (ny, nx)
-            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
-
-        # Remove _FillValue encoding so output has no masked values
-        if "_FillValue" in ds[varname].encoding:
-            del ds[varname].encoding["_FillValue"]
-
-        write_netcdf(ds, output_file)
diff --git a/compass/landice/tests/ismip7_forcing/atmosphere/process_smb.py b/compass/landice/tests/ismip7_forcing/atmosphere/process_smb.py
index 222fb6e92e..3bb864710e 100644
--- a/compass/landice/tests/ismip7_forcing/atmosphere/process_smb.py
+++ b/compass/landice/tests/ismip7_forcing/atmosphere/process_smb.py
@@ -10,6 +10,7 @@
     build_mapping_file,
 )
 from compass.landice.tests.ismip7_forcing.ice_sheet_params import get_params
+from compass.landice.tests.ismip7_forcing.remap_utils import extrapolate_source
 from compass.step import Step
 
 
@@ -69,6 +70,11 @@ def run(self):
         prefix = params['prefix']
         resolution = params['atm_resolution']
         version = params['atm_version']
+        if params['atm_model'] is not None:
+            forcing_group = scenario
+            model = params['atm_model']
+        else:
+            forcing_group = f"{model}_{scenario}"
         input_path = os.path.join(base_path_ismip7, "acabf", version)
         file_pattern = (f"acabf_{prefix}_{model}_{scenario}_"
                         f"SDBN1-{resolution}_{version}_*.nc")
@@ -82,8 +88,12 @@ def run(self):
         # Filter to requested year range
         input_files = []
         for f in all_files:
-            # Extract year from filename (last part before .nc)
-            year = int(os.path.basename(f).split("_")[-1].replace(".nc", ""))
+            # Extract year from filename (last part before .nc); skip
+            # non-yearly files such as climatology averages (e.g. *_avg.nc)
+            token = os.path.basename(f).split("_")[-1].replace(".nc", "")
+            if not token.isdigit():
+                continue
+            year = int(token)
             if start_year <= year <= end_year:
                 input_files.append(f)
 
@@ -117,15 +127,24 @@ def run(self):
                 logger.info(f"  Remapped file exists, skipping: {basename}")
                 continue
 
+            # Extrapolate fill values on source grid before remapping
+            # so they don't pollute neighboring cells during interpolation
+            extrap_file = f"extrap_{basename}"
+            if not os.path.exists(extrap_file):
+                extrapolate_source(input_file, extrap_file, "acabf", logger)
+
             logger.info(f"  Remapping: {basename}")
             args = ["ncremap",
-                    "-i", input_file,
+                    "-i", extrap_file,
                     "-o", remapped_file,
                     "-m", mapping_file,
                     "-v", "acabf"]
 
             check_call(args, logger=logger)
 
+            # Clean up extrapolated source file
+            os.remove(extrap_file)
+
         # Combine remapped files and rename to MALI conventions
         logger.info("Combining remapped files and renaming variables...")
         output_file = (f"{mali_mesh_name}_SMB_{model}_{scenario}_"
@@ -140,8 +159,8 @@ def run(self):
                 os.remove(f)
 
         # Place output in appropriate directory
-        output_path = os.path.join(output_base_path, "atmosphere_forcing",
-                                   f"{model}_{scenario}")
+        output_path = os.path.join(output_base_path, forcing_group,
+                                   "atmosphere")
         if not os.path.exists(output_path):
             os.makedirs(output_path)
 
@@ -164,7 +183,8 @@ def _combine_and_rename(self, remapped_files, output_file):
             Output file path
         """
         ds = xr.open_mfdataset(remapped_files, concat_dim="time",
-                               combine="nested", engine="netcdf4")
+                               combine="nested", engine="netcdf4",
+                               drop_variables="time_bnds")
 
         # Rename dimensions to MALI conventions
         rename_dims = {}
diff --git a/compass/landice/tests/ismip7_forcing/atmosphere/process_smb_gradient.py b/compass/landice/tests/ismip7_forcing/atmosphere/process_smb_gradient.py
index 2f3e4c2cb4..d44192d14c 100644
--- a/compass/landice/tests/ismip7_forcing/atmosphere/process_smb_gradient.py
+++ b/compass/landice/tests/ismip7_forcing/atmosphere/process_smb_gradient.py
@@ -10,6 +10,7 @@
     build_mapping_file,
 )
 from compass.landice.tests.ismip7_forcing.ice_sheet_params import get_params
+from compass.landice.tests.ismip7_forcing.remap_utils import extrapolate_source
 from compass.step import Step
 
 
@@ -70,6 +71,11 @@ def run(self):
         prefix = params['prefix']
         resolution = params['atm_resolution']
         version = params['atm_version']
+        if params['atm_model'] is not None:
+            forcing_group = scenario
+            model = params['atm_model']
+        else:
+            forcing_group = f"{model}_{scenario}"
         input_path = os.path.join(base_path_ismip7, "dacabfdz", version)
         file_pattern = (f"dacabfdz_{prefix}_{model}_{scenario}_"
                         f"SDBN1-{resolution}_{version}_*.nc")
@@ -83,7 +89,11 @@ def run(self):
         # Filter to requested year range
         input_files = []
         for f in all_files:
-            year = int(os.path.basename(f).split("_")[-1].replace(".nc", ""))
+            # skip non-yearly files such as climatology averages (*_avg.nc)
+            token = os.path.basename(f).split("_")[-1].replace(".nc", "")
+            if not token.isdigit():
+                continue
+            year = int(token)
             if start_year <= year <= end_year:
                 input_files.append(f)
 
@@ -118,15 +128,25 @@ def run(self):
                 logger.info(f"  Remapped file exists, skipping: {basename}")
                 continue
 
+            # Extrapolate fill values on source grid before remapping
+            # so they don't pollute neighboring cells during interpolation
+            extrap_file = f"extrap_{basename}"
+            if not os.path.exists(extrap_file):
+                extrapolate_source(input_file, extrap_file, "dacabfdz",
+                                   logger)
+
             logger.info(f"  Remapping: {basename}")
             args = ["ncremap",
-                    "-i", input_file,
+                    "-i", extrap_file,
                     "-o", remapped_file,
                     "-m", mapping_file,
                     "-v", "dacabfdz"]
 
             check_call(args, logger=logger)
 
+            # Clean up extrapolated source file
+            os.remove(extrap_file)
+
         # Combine remapped files and rename to MALI conventions
         logger.info("Combining remapped files and renaming variables...")
         output_file = (f"{mali_mesh_name}_SMB_gradient_{model}_{scenario}_"
@@ -141,8 +161,8 @@ def run(self):
                 os.remove(f)
 
         # Place output in appropriate directory
-        output_path = os.path.join(output_base_path, "atmosphere_forcing",
-                                   f"{model}_{scenario}")
+        output_path = os.path.join(output_base_path, forcing_group,
+                                   "atmosphere")
         if not os.path.exists(output_path):
             os.makedirs(output_path)
 
@@ -165,7 +185,8 @@ def _combine_and_rename(self, remapped_files, output_file):
             Output file path
         """
         ds = xr.open_mfdataset(remapped_files, concat_dim="time",
-                               combine="nested", engine="netcdf4")
+                               combine="nested", engine="netcdf4",
+                               drop_variables="time_bnds")
 
         # Rename dimensions to MALI conventions
         rename_dims = {}
diff --git a/compass/landice/tests/ismip7_forcing/atmosphere/process_temperature.py b/compass/landice/tests/ismip7_forcing/atmosphere/process_temperature.py
index caa8f6be9f..24c9529342 100644
--- a/compass/landice/tests/ismip7_forcing/atmosphere/process_temperature.py
+++ b/compass/landice/tests/ismip7_forcing/atmosphere/process_temperature.py
@@ -10,6 +10,7 @@
     build_mapping_file,
 )
 from compass.landice.tests.ismip7_forcing.ice_sheet_params import get_params
+from compass.landice.tests.ismip7_forcing.remap_utils import extrapolate_source
 from compass.step import Step
 
 
@@ -69,6 +70,11 @@ def run(self):
         prefix = params['prefix']
         resolution = params['atm_resolution']
         version = params['atm_version']
+        if params['atm_model'] is not None:
+            forcing_group = scenario
+            model = params['atm_model']
+        else:
+            forcing_group = f"{model}_{scenario}"
         input_path = os.path.join(base_path_ismip7, "ts", version)
         file_pattern = (f"ts_{prefix}_{model}_{scenario}_"
                         f"SDBN1-{resolution}_{version}_*.nc")
@@ -82,7 +88,11 @@ def run(self):
         # Filter to requested year range
         input_files = []
         for f in all_files:
-            year = int(os.path.basename(f).split("_")[-1].replace(".nc", ""))
+            # skip non-yearly files such as climatology averages (*_avg.nc)
+            token = os.path.basename(f).split("_")[-1].replace(".nc", "")
+            if not token.isdigit():
+                continue
+            year = int(token)
             if start_year <= year <= end_year:
                 input_files.append(f)
 
@@ -117,15 +127,24 @@ def run(self):
                 logger.info(f"  Remapped file exists, skipping: {basename}")
                 continue
 
+            # Extrapolate fill values on source grid before remapping
+            # so they don't pollute neighboring cells during interpolation
+            extrap_file = f"extrap_{basename}"
+            if not os.path.exists(extrap_file):
+                extrapolate_source(input_file, extrap_file, "ts", logger)
+
             logger.info(f"  Remapping: {basename}")
             args = ["ncremap",
-                    "-i", input_file,
+                    "-i", extrap_file,
                     "-o", remapped_file,
                     "-m", mapping_file,
                     "-v", "ts"]
 
             check_call(args, logger=logger)
 
+            # Clean up extrapolated source file
+            os.remove(extrap_file)
+
         # Combine remapped files and rename to MALI conventions
         logger.info("Combining remapped files and renaming variables...")
         output_file = (f"{mali_mesh_name}_temperature_{model}_{scenario}_"
@@ -140,8 +159,8 @@ def run(self):
                 os.remove(f)
 
         # Place output in appropriate directory
-        output_path = os.path.join(output_base_path, "atmosphere_forcing",
-                                   f"{model}_{scenario}")
+        output_path = os.path.join(output_base_path, forcing_group,
+                                   "atmosphere")
         if not os.path.exists(output_path):
             os.makedirs(output_path)
 
@@ -164,7 +183,8 @@ def _combine_and_rename(self, remapped_files, output_file):
             Output file path
         """
         ds = xr.open_mfdataset(remapped_files, concat_dim="time",
-                               combine="nested", engine="netcdf4")
+                               combine="nested", engine="netcdf4",
+                               drop_variables="time_bnds")
 
         # Rename dimensions to MALI conventions
         rename_dims = {}
diff --git a/compass/landice/tests/ismip7_forcing/atmosphere/process_temperature_gradient.py b/compass/landice/tests/ismip7_forcing/atmosphere/process_temperature_gradient.py
index 5dec7b78bc..45d39b9a1d 100644
--- a/compass/landice/tests/ismip7_forcing/atmosphere/process_temperature_gradient.py
+++ b/compass/landice/tests/ismip7_forcing/atmosphere/process_temperature_gradient.py
@@ -10,6 +10,7 @@
     build_mapping_file,
 )
 from compass.landice.tests.ismip7_forcing.ice_sheet_params import get_params
+from compass.landice.tests.ismip7_forcing.remap_utils import extrapolate_source
 from compass.step import Step
 
 
@@ -71,6 +72,11 @@ def run(self):
         prefix = params['prefix']
         resolution = params['atm_resolution']
         version = params['atm_version']
+        if params['atm_model'] is not None:
+            forcing_group = scenario
+            model = params['atm_model']
+        else:
+            forcing_group = f"{model}_{scenario}"
         input_path = os.path.join(base_path_ismip7, "dtsdz", version)
         file_pattern = (f"dtsdz_{prefix}_{model}_{scenario}_"
                         f"SDBN1-{resolution}_{version}_*.nc")
@@ -84,7 +90,11 @@ def run(self):
         # Filter to requested year range
         input_files = []
         for f in all_files:
-            year = int(os.path.basename(f).split("_")[-1].replace(".nc", ""))
+            # skip non-yearly files such as climatology averages (*_avg.nc)
+            token = os.path.basename(f).split("_")[-1].replace(".nc", "")
+            if not token.isdigit():
+                continue
+            year = int(token)
             if start_year <= year <= end_year:
                 input_files.append(f)
 
@@ -119,15 +129,24 @@ def run(self):
                 logger.info(f"  Remapped file exists, skipping: {basename}")
                 continue
 
+            # Extrapolate fill values on source grid before remapping
+            # so they don't pollute neighboring cells during interpolation
+            extrap_file = f"extrap_{basename}"
+            if not os.path.exists(extrap_file):
+                extrapolate_source(input_file, extrap_file, "dtsdz", logger)
+
             logger.info(f"  Remapping: {basename}")
             args = ["ncremap",
-                    "-i", input_file,
+                    "-i", extrap_file,
                     "-o", remapped_file,
                     "-m", mapping_file,
                     "-v", "dtsdz"]
 
             check_call(args, logger=logger)
 
+            # Clean up extrapolated source file
+            os.remove(extrap_file)
+
         # Combine remapped files and rename to MALI conventions
         logger.info("Combining remapped files and renaming variables...")
         output_file = (f"{mali_mesh_name}_temperature_gradient_{model}_"
@@ -142,8 +161,8 @@ def run(self):
                 os.remove(f)
 
         # Place output in appropriate directory
-        output_path = os.path.join(output_base_path, "atmosphere_forcing",
-                                   f"{model}_{scenario}")
+        output_path = os.path.join(output_base_path, forcing_group,
+                                   "atmosphere")
         if not os.path.exists(output_path):
             os.makedirs(output_path)
 
@@ -166,7 +185,8 @@ def _combine_and_rename(self, remapped_files, output_file):
             Output file path
         """
         ds = xr.open_mfdataset(remapped_files, concat_dim="time",
-                               combine="nested", engine="netcdf4")
+                               combine="nested", engine="netcdf4",
+                               drop_variables="time_bnds")
 
         # Rename dimensions to MALI conventions
         rename_dims = {}
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 fbf3f75b83..8a9275c962 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_excess_melt.py
@@ -12,9 +12,9 @@
 )
 from compass.landice.tests.ismip7_forcing.fracture.remap_utils import (
     add_xtime_and_write,
-    extrapolate_source,
     open_rename_and_trim,
 )
+from compass.landice.tests.ismip7_forcing.remap_utils import extrapolate_source
 from compass.step import Step
 
 
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 10df254c69..6867f20735 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/process_lake_properties.py
@@ -9,9 +9,9 @@
 )
 from compass.landice.tests.ismip7_forcing.fracture.remap_utils import (
     add_xtime_and_write,
-    extrapolate_source,
     open_rename_and_trim,
 )
+from compass.landice.tests.ismip7_forcing.remap_utils import extrapolate_source
 from compass.step import Step
 
 
diff --git a/compass/landice/tests/ismip7_forcing/fracture/remap_utils.py b/compass/landice/tests/ismip7_forcing/fracture/remap_utils.py
index 8da5de4ae4..14518c7c79 100644
--- a/compass/landice/tests/ismip7_forcing/fracture/remap_utils.py
+++ b/compass/landice/tests/ismip7_forcing/fracture/remap_utils.py
@@ -1,67 +1,8 @@
 """
 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):
diff --git a/compass/landice/tests/ismip7_forcing/ice_sheet_params.py b/compass/landice/tests/ismip7_forcing/ice_sheet_params.py
index 8324236c56..d991eae45c 100644
--- a/compass/landice/tests/ismip7_forcing/ice_sheet_params.py
+++ b/compass/landice/tests/ismip7_forcing/ice_sheet_params.py
@@ -10,8 +10,11 @@
         'atm_resolution': '2000m',
         'atm_version': 'v2',
         'ocean_version': 'v3',
+        'ocean_grid': 'ocean',
         'ocean_3d': True,
         'ocean_temporal': 'decade',
+        'atm_model': None,
+        'ocean_model': None,
     },
     'gis': {
         'projection': 'gis-bamber',
@@ -19,8 +22,34 @@
         'atm_resolution': '1000m',
         'atm_version': 'v2',
         'ocean_version': 'v2',
+        'ocean_grid': 'ocean',
         'ocean_3d': False,
         'ocean_temporal': 'yearly',
+        'atm_model': None,
+        'ocean_model': None,
+    },
+}
+
+# Overrides applied for the OCX (reanalysis) scenario. OCX has no distinct
+# ESM model: it uses fixed reanalysis products (RACMO for the atmosphere and
+# EN4 for the ocean) at data version v1. When scenario is 'OCX' the [ismip7]
+# model option is ignored and these sources are used instead.
+_OCX_OVERRIDES = {
+    'gis': {
+        'atm_version': 'v1',
+        'ocean_version': 'v1',
+        'ocean_grid': 'ocean-1000m',
+        'atm_model': 'RACMO2.3p2-ERA',
+        'ocean_model': 'EN4',
+    },
+    'ais': {
+        'atm_version': 'v1',
+        'ocean_version': 'v1',
+        'atm_model': 'RACMO2.3p2-ERA',
+        'ocean_model': None,
+        # AIS OCX ocean files have no model token and live in per-choice
+        # subdirectories (main/cold/warm/vary); see process_thermal_forcing.
+        'ocean_choice_layout': True,
     },
 }
 
@@ -44,4 +73,14 @@ def get_params(config):
         raise ValueError(
             f"Unknown ice_sheet '{ice_sheet}'. "
             f"Must be one of: {list(_PARAMS.keys())}")
-    return _PARAMS[ice_sheet]
+    params = dict(_PARAMS[ice_sheet])
+
+    scenario = config.get("ismip7", "scenario")
+    if scenario == "OCX":
+        if ice_sheet not in _OCX_OVERRIDES:
+            raise ValueError(
+                f"The OCX scenario is not yet supported for ice_sheet "
+                f"'{ice_sheet}'.")
+        params.update(_OCX_OVERRIDES[ice_sheet])
+
+    return params
diff --git a/compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg b/compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg
index f02b79067a..be7d195aa4 100644
--- a/compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg
+++ b/compass/landice/tests/ismip7_forcing/ismip7_forcing.cfg
@@ -54,6 +54,10 @@ end_year = 2014
 # Remapping method. Options: bilinear, neareststod, conserve
 method_remap = bilinear
 
+# Ocean forcing choice(s), only used for the AIS OCX scenario. Comma-separated
+# list of any of: main, cold, warm, vary. Use 'all' to process every choice.
+ocean_choice = main
+
 # Start year for processing
 start_year = 1850
 
diff --git a/compass/landice/tests/ismip7_forcing/ismip7_forcing_ocx_ais.cfg b/compass/landice/tests/ismip7_forcing/ismip7_forcing_ocx_ais.cfg
new file mode 100644
index 0000000000..975cb89b85
--- /dev/null
+++ b/compass/landice/tests/ismip7_forcing/ismip7_forcing_ocx_ais.cfg
@@ -0,0 +1,84 @@
+# Example config for processing the AIS ISMIP7 OCX (reanalysis) forcing.
+# OCX has no distinct ESM model: it uses RACMO2.3p2-ERA for the atmosphere and
+# reanalysis ocean products, selected automatically when scenario = OCX. The
+# [ismip7] model option is ignored for OCX (set it to None). Unlike GrIS, the
+# AIS OCX ocean provides several forcing choices (main/cold/warm/vary) selected
+# with the ocean_choice option below. The same config file drives both test
+# cases; set them up individually, e.g.:
+#   compass setup -t landice/ismip7_forcing/atmosphere -w WORKDIR -f ismip7_forcing_ocx_ais.cfg
+#   compass setup -t landice/ismip7_forcing/ocean_thermal -w WORKDIR -f ismip7_forcing_ocx_ais.cfg
+
+# config options for ismip7 forcing data
+[ismip7]
+
+# Ice sheet: ais (Antarctic) or gis (Greenland)
+ice_sheet = ais
+
+# Base path to the input ISMIP7 forcing files. User has to supply.
+base_path_ismip7 = /global/cfs/cdirs/m4288/users/trhille/ISMIP7/forcing/AIS/OCX
+
+# Base path to the MALI mesh. User has to supply.
+base_path_mali = /global/cfs/cdirs/fanssie/MALI_input_files/AIS_4to20km_r01
+
+# Base path to which output forcing files are saved.
+output_base_path = /global/cfs/cdirs/m4288/users/trhille/ISMIP7/test_processing/AIS
+
+# OCX has no distinct model (RACMO for atmosphere, reanalysis for ocean are
+# selected automatically), so this option is ignored.
+model = None
+
+# Scenario for forcing data.
+scenario = OCX
+
+# Name of the MALI mesh. Used to name mapping and output files.
+mali_mesh_name = AIS_4to20km_r01_20220907
+
+# MALI mesh file. User has to supply.
+mali_mesh_file = AIS_4to20km_r01_20220907.nc
+
+# Number of MPI tasks for ESMF_RegridWeightGen
+esmf_ntasks = 512
+
+# Whether to process time-varying ocean thermal forcing (ESM scenario data)
+process_ocean_thermal = true
+
+# Whether to process observational ocean thermal forcing climatology
+process_ocean_climatology = false
+
+# config options for ismip7 atmosphere forcing
+[ismip7_atmosphere]
+
+# Remapping method. Options: bilinear, neareststod, conserve
+method_remap = conserve
+
+# Start year for processing
+start_year = 1979
+
+# End year for processing
+end_year = 2025
+
+# config options for ismip7 ocean thermal forcing
+[ismip7_ocean_thermal]
+
+# Remapping method. Options: bilinear, neareststod, conserve
+method_remap = bilinear
+
+# Ocean forcing choice(s) for the AIS OCX scenario. Comma-separated list of any
+# of: main, cold, warm, vary. Use 'all' to process every choice. Each choice is
+# written to its own OCX_<choice> output directory.
+ocean_choice = main
+
+# Start year for processing
+start_year = 1950
+
+# End year for processing
+end_year = 2025
+
+# config options for ismip7 ocean thermal forcing climatology
+[ismip7_ocean_climatology]
+
+# Remapping method. Options: bilinear, neareststod, conserve
+method_remap = bilinear
+
+# Base path to observational climatology data
+base_path_climatology = None
diff --git a/compass/landice/tests/ismip7_forcing/ismip7_forcing_ocx_gis.cfg b/compass/landice/tests/ismip7_forcing/ismip7_forcing_ocx_gis.cfg
new file mode 100644
index 0000000000..d3658027db
--- /dev/null
+++ b/compass/landice/tests/ismip7_forcing/ismip7_forcing_ocx_gis.cfg
@@ -0,0 +1,77 @@
+# Example config for processing the GrIS ISMIP7 OCX (reanalysis) forcing.
+# OCX has no distinct ESM model: it uses RACMO2.3p2-ERA for the atmosphere and
+# EN4 for the ocean, selected automatically when scenario = OCX. The [ismip7]
+# model option is ignored for OCX (set it to None). The same config file drives
+# both test cases; set them up individually, e.g.:
+#   compass setup -t landice/ismip7_forcing/atmosphere -w WORKDIR -f ismip7_forcing_ocx_gis.cfg
+#   compass setup -t landice/ismip7_forcing/ocean_thermal -w WORKDIR -f ismip7_forcing_ocx_gis.cfg
+
+# config options for ismip7 forcing data
+[ismip7]
+
+# Ice sheet: ais (Antarctic) or gis (Greenland)
+ice_sheet = gis
+
+# Base path to the input ISMIP7 forcing files. User has to supply.
+base_path_ismip7 = /global/cfs/cdirs/m4288/users/trhille/ISMIP7/forcing/GIS/OCX
+
+# Base path to the MALI mesh. User has to supply.
+base_path_mali = /global/cfs/cdirs/fanssie/MALI_input_files/GIS_1to10km_r02/
+
+# Base path to which output forcing files are saved.
+output_base_path = /global/cfs/cdirs/m4288/users/trhille/ISMIP7/test_processing/GIS
+
+# OCX has no distinct model (RACMO for atmosphere, EN4 for ocean are selected
+# automatically), so this option is ignored.
+model = None
+
+# Scenario for forcing data.
+scenario = OCX
+
+# Name of the MALI mesh. Used to name mapping and output files.
+mali_mesh_name = GIS_1to10km_r02_20230202
+
+# MALI mesh file. User has to supply.
+mali_mesh_file = GIS_1to10km_r02_20230202.nc
+
+# Number of MPI tasks for ESMF_RegridWeightGen
+esmf_ntasks = 512
+
+# Whether to process time-varying ocean thermal forcing (ESM scenario data)
+process_ocean_thermal = true
+
+# Whether to process observational ocean thermal forcing climatology
+process_ocean_climatology = false
+
+# config options for ismip7 atmosphere forcing
+[ismip7_atmosphere]
+
+# Remapping method. Options: bilinear, neareststod, conserve
+method_remap = conserve
+
+# Start year for processing
+start_year = 1990
+
+# End year for processing
+end_year = 2025
+
+# config options for ismip7 ocean thermal forcing
+[ismip7_ocean_thermal]
+
+# Remapping method. Options: bilinear, neareststod, conserve
+method_remap = bilinear
+
+# Start year for processing
+start_year = 1990
+
+# End year for processing
+end_year = 2025
+
+# config options for ismip7 ocean thermal forcing climatology
+[ismip7_ocean_climatology]
+
+# Remapping method. Options: bilinear, neareststod, conserve
+method_remap = bilinear
+
+# Base path to observational climatology data
+base_path_climatology = None
diff --git a/compass/landice/tests/ismip7_forcing/ocean_thermal/process_thermal_forcing.py b/compass/landice/tests/ismip7_forcing/ocean_thermal/process_thermal_forcing.py
index a1469c8f2d..cc71341cbb 100644
--- a/compass/landice/tests/ismip7_forcing/ocean_thermal/process_thermal_forcing.py
+++ b/compass/landice/tests/ismip7_forcing/ocean_thermal/process_thermal_forcing.py
@@ -2,16 +2,15 @@
 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.ice_sheet_params import get_params
+from compass.landice.tests.ismip7_forcing.remap_utils import extrapolate_source
 from compass.step import Step
 
 
@@ -68,9 +67,9 @@ def run(self):
     def _run_scenario(self):
         """
         Process time-varying ocean thermal forcing from an ESM
-        (e.g., CESM2-WACCM historical or ssp585).
+        (e.g., CESM2-WACCM historical or ssp585) or, for AIS OCX, from one
+        or more reanalysis ocean "choices" (main/cold/warm/vary).
         """
-        logger = self.logger
         config = self.config
         params = get_params(config)
 
@@ -88,14 +87,123 @@ def _run_scenario(self):
         start_year = section.getint("start_year")
         end_year = section.getint("end_year")
 
-        # Discover input files
         prefix = params['prefix']
         ocean_version = params['ocean_version']
+        ocean_grid = params['ocean_grid']
         ocean_3d = params['ocean_3d']
-        input_path = os.path.join(base_path_ismip7, "ocean", "tf",
-                                  ocean_version)
-        file_pattern = (f"tf_{prefix}_{model}_{scenario}_"
-                        f"ocean_{ocean_version}_*.nc")
+
+        # Assemble the list of forcing sources to process. AIS OCX has several
+        # ocean "choices" (main/cold/warm/vary) in per-choice subdirectories
+        # with no model token in the filename; all other cases (ESM scenarios
+        # and GrIS OCX) have a single source.
+        jobs = []
+        if params.get('ocean_choice_layout', False):
+            for choice in self._get_ocean_choices(config):
+                input_path = os.path.join(base_path_ismip7, "ocean", choice,
+                                          ocean_version)
+                file_pattern = (f"tf_{prefix}_{scenario}_ocean_{choice}_"
+                                f"{ocean_version}_*.nc")
+                jobs.append({
+                    'input_path': input_path,
+                    'file_pattern': file_pattern,
+                    'forcing_group': f"{scenario}_{choice}",
+                    'label': f"{scenario}_{choice}",
+                })
+        else:
+            if params['ocean_model'] is not None:
+                forcing_group = scenario
+                model = params['ocean_model']
+            else:
+                forcing_group = f"{model}_{scenario}"
+            input_path = os.path.join(base_path_ismip7, "ocean", "tf",
+                                      ocean_version)
+            file_pattern = (f"tf_{prefix}_{model}_{scenario}_"
+                            f"{ocean_grid}_{ocean_version}_*.nc")
+            jobs.append({
+                'input_path': input_path,
+                'file_pattern': file_pattern,
+                'forcing_group': forcing_group,
+                'label': f"{model}_{scenario}",
+            })
+
+        # Mapping file is shared across sources (identical ocean source grid).
+        mapping_file = (f"map_ismip7_{ice_sheet}_ocean_to_"
+                        f"{mali_mesh_name}_{method_remap}.nc")
+
+        for job in jobs:
+            self._process_ocean_forcing(
+                job, mapping_file, ocean_3d, method_remap,
+                start_year, end_year, mali_mesh_name, mali_mesh_file,
+                output_base_path)
+
+    def _get_ocean_choices(self, config):
+        """
+        Parse the comma-separated ``ocean_choice`` option (AIS OCX only) into a
+        list of ocean choices, expanding ``all`` to every available choice.
+
+        Parameters
+        ----------
+        config : compass.config.CompassConfigParser
+            Configuration options for the test case
+
+        Returns
+        -------
+        choices : list of str
+            Ordered, de-duplicated list of ocean choices to process
+        """
+        valid = ["main", "cold", "warm", "vary"]
+        raw = config.get("ismip7_ocean_thermal", "ocean_choice")
+        tokens = [t.strip() for t in raw.split(",") if t.strip()]
+        if any(t.lower() == "all" for t in tokens):
+            return valid
+
+        choices = []
+        for token in tokens:
+            if token not in valid:
+                raise ValueError(
+                    f"Invalid ocean_choice '{token}'. Must be a "
+                    f"comma-separated list of {valid} or 'all'.")
+            if token not in choices:
+                choices.append(token)
+        if not choices:
+            raise ValueError(
+                "No ocean_choice specified for AIS OCX ocean forcing.")
+        return choices
+
+    def _process_ocean_forcing(self, job, mapping_file, ocean_3d,
+                               method_remap, start_year, end_year,
+                               mali_mesh_name, mali_mesh_file,
+                               output_base_path):
+        """
+        Discover, remap, combine, and save the thermal forcing for a single
+        forcing source described by ``job``.
+
+        Parameters
+        ----------
+        job : dict
+            Source description with keys ``input_path``, ``file_pattern``,
+            ``forcing_group``, and ``label``
+        mapping_file : str
+            Path of the shared ocean-to-MALI mapping file (built on demand)
+        ocean_3d : bool
+            Whether the thermal forcing is 3D (AIS) or 2D (GrIS)
+        method_remap : str
+            Remapping method passed to the mapping-file builder
+        start_year, end_year : int
+            Inclusive year range to process
+        mali_mesh_name, mali_mesh_file : str
+            MALI mesh name and file
+        output_base_path : str
+            Base path under which output is written
+        """
+        logger = self.logger
+        config = self.config
+
+        input_path = job['input_path']
+        file_pattern = job['file_pattern']
+        forcing_group = job['forcing_group']
+        label = job['label']
+
         all_files = sorted(glob.glob(os.path.join(input_path, file_pattern)))
 
         if not all_files:
@@ -104,8 +212,8 @@ def _run_scenario(self):
                 f"  {os.path.join(input_path, file_pattern)}")
 
         # Filter to files that overlap with the requested year range.
-        # AIS files are named with decade ranges (e.g., 1850-1859).
-        # GrIS files are named with single years (e.g., 2015).
+        # AIS files are named with decade or multi-decade ranges (e.g.,
+        # 1850-1859, 1950-2025). GrIS files are named with single years.
         input_files = []
         for f in all_files:
             # Extract year range from filename (last part before .nc)
@@ -125,9 +233,6 @@ def _run_scenario(self):
                     f"overlapping years {start_year}-{end_year}")
 
         # Build mapping file using the first input file as grid template.
-        mapping_file = (f"map_ismip7_{ice_sheet}_ocean_to_"
-                        f"{mali_mesh_name}_{method_remap}.nc")
-
         if not os.path.exists(mapping_file):
             logger.info("Building mapping file for ocean grid...")
             build_mapping_file(config, logger,
@@ -135,7 +240,7 @@ def _run_scenario(self):
                                mali_mesh_file=mali_mesh_file,
                                method_remap=method_remap)
 
-        # Remap each decade file
+        # Remap each file
         remapped_files = []
         for input_file in input_files:
             basename = os.path.basename(input_file)
@@ -150,8 +255,7 @@ def _run_scenario(self):
             # so they don't pollute neighboring cells during interpolation
             extrap_file = f"extrap_{basename}"
             if not os.path.exists(extrap_file):
-                self._extrapolate_source(input_file, extrap_file, "tf",
-                                         logger)
+                extrapolate_source(input_file, extrap_file, "tf", logger)
 
             logger.info(f"  Remapping: {basename}")
             args = ["ncremap",
@@ -167,7 +271,7 @@ def _run_scenario(self):
 
         # Combine remapped files and rename to MALI conventions
         logger.info("Combining remapped files and renaming variables...")
-        output_file = (f"{mali_mesh_name}_thermal_forcing_{model}_{scenario}_"
+        output_file = (f"{mali_mesh_name}_thermal_forcing_{label}_"
                        f"{start_year}-{end_year}.nc")
 
         if ocean_3d:
@@ -184,8 +288,8 @@ def _run_scenario(self):
                 os.remove(f)
 
         # Place output in appropriate directory
-        output_path = os.path.join(output_base_path, "ocean_thermal_forcing",
-                                   f"{model}_{scenario}")
+        output_path = os.path.join(output_base_path, forcing_group,
+                                   "ocean_thermal_forcing")
         if not os.path.exists(output_path):
             os.makedirs(output_path)
 
@@ -246,8 +350,7 @@ def _run_climatology(self):
         if not os.path.exists(remapped_file):
             extrap_file = f"extrap_{basename}"
             if not os.path.exists(extrap_file):
-                self._extrapolate_source(input_file, extrap_file, "tf",
-                                         logger)
+                extrapolate_source(input_file, extrap_file, "tf", logger)
 
             logger.info(f"  Remapping: {basename}")
             args = ["ncremap",
@@ -305,7 +408,8 @@ def _combine_and_rename_3d(self, remapped_files, output_file,
             Last year to include in output
         """
         ds = xr.open_mfdataset(remapped_files, concat_dim="time",
-                               combine="nested", engine="netcdf4")
+                               combine="nested", engine="netcdf4",
+                               drop_variables="time_bnds")
 
         # Subset to requested year range
         years = ds.time.dt.year
@@ -419,7 +523,8 @@ def _combine_and_rename_2d(self, remapped_files, output_file,
             Last year to include in output
         """
         ds = xr.open_mfdataset(remapped_files, concat_dim="time",
-                               combine="nested", engine="netcdf4")
+                               combine="nested", engine="netcdf4",
+                               drop_variables="time_bnds")
 
         # Subset to requested year range
         years = ds.time.dt.year
@@ -563,59 +668,3 @@ def _rename_climatology_3d(self, remapped_file, output_file):
             ds = ds.drop_vars("nISMIP6OceanLayers")
 
         write_netcdf(ds, output_file)
-
-    def _extrapolate_source(self, input_file, output_file, varname, logger):
-        """
-        Extrapolate fill/missing values on the source polar stereographic
-        grid using nearest-neighbor interpolation from valid cells. 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 (e.g., "tf")
-
-        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, engine="netcdf4")
-        data = ds[varname]
-
-        # Process each time step (and z level if 3D)
-        # Source files have dims like (time, z, y, x) or (time, y, x)
-        values = data.values.copy()
-        non_spatial_shape = values.shape[:-2]  # (time,) or (time, z)
-
-        # Use distance_transform_edt with return_indices to find the
-        # nearest valid cell index for each invalid cell. This is O(n)
-        # on the grid and much faster than KD-tree approaches.
-        for idx in np.ndindex(non_spatial_shape):
-            slab = values[idx]  # shape (ny, nx)
-            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
-
-        # Remove _FillValue encoding so output has no masked values
-        if "_FillValue" in ds[varname].encoding:
-            del ds[varname].encoding["_FillValue"]
-
-        write_netcdf(ds, output_file)
diff --git a/compass/landice/tests/ismip7_forcing/remap_utils.py b/compass/landice/tests/ismip7_forcing/remap_utils.py
new file mode 100644
index 0000000000..890a979ab1
--- /dev/null
+++ b/compass/landice/tests/ismip7_forcing/remap_utils.py
@@ -0,0 +1,79 @@
+"""
+Shared helpers for remapping ISMIP7 forcing data to the MALI mesh.
+"""
+import os
+
+import netCDF4
+import numpy as np
+import xarray as xr
+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"]
+
+    # Preserve a fill value for slabs that remain fully invalid after
+    # extrapolation so ncremap ignores them during interpolation.
+    encoding = {}
+    for varname in varnames:
+        dtype = ds[varname].dtype
+        if np.issubdtype(dtype, np.floating) and \
+                bool(np.any(np.isnan(ds[varname].values))):
+            fill = netCDF4.default_fillvals[dtype.str[1:]]
+            encoding[varname] = {"_FillValue": dtype.type(fill)}
+
+    # Write CDF-5 (NETCDF3_64BIT_DATA): a classic-model format with 64-bit
+    # sizes that supports very large variables (e.g. AIS 3D ocean thermal
+    # forcing) without the HDF5 chunk-size limits that make ncremap unable to
+    # open large NETCDF4 files. This is a temporary file consumed by ncremap.
+    ds.to_netcdf(output_file, format="NETCDF3_64BIT_DATA", engine="netcdf4",
+                 encoding=encoding)
+    ds.close()
diff --git a/docs/developers_guide/landice/api.rst b/docs/developers_guide/landice/api.rst
index d2cd995145..75bb6bab8c 100644
--- a/docs/developers_guide/landice/api.rst
+++ b/docs/developers_guide/landice/api.rst
@@ -385,6 +385,7 @@ ismip7_forcing
    configure.configure
    ice_sheet_params.get_params
    create_mapfile.build_mapping_file
+   remap_utils.extrapolate_source
 
    atmosphere.Atmosphere
    atmosphere.Atmosphere.configure
@@ -412,7 +413,6 @@ 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
diff --git a/docs/developers_guide/landice/test_groups/ismip7_forcing.rst b/docs/developers_guide/landice/test_groups/ismip7_forcing.rst
index b535bf081a..0a37029290 100644
--- a/docs/developers_guide/landice/test_groups/ismip7_forcing.rst
+++ b/docs/developers_guide/landice/test_groups/ismip7_forcing.rst
@@ -24,11 +24,19 @@ ice_sheet_params
 
 The module :py:mod:`compass.landice.tests.ismip7_forcing.ice_sheet_params`
 defines a dictionary of ice-sheet-specific parameters (projection, file
-naming prefix, grid resolution, data version, ocean dimensionality) and
-provides the function
+naming prefix, grid resolution, data version, ocean dimensionality, and the
+atmosphere/ocean source names) and provides the function
 :py:func:`compass.landice.tests.ismip7_forcing.ice_sheet_params.get_params`
 to retrieve them based on the ``ice_sheet`` config option.
 
+When ``scenario = OCX``, ``get_params`` applies a set of OCX overrides on top
+of the ice-sheet defaults: data version ``v1``, the ocean file-name grid token
+(e.g. ``ocean-1000m``), and the fixed reanalysis sources (``atm_model`` =
+``RACMO2.3p2-ERA`` and ``ocean_model`` = ``EN4``). The processing steps use
+``atm_model`` / ``ocean_model`` in place of the ``[ismip7] model`` option when
+they are set, so the OCX ``model`` option is ignored. This keeps OCX handling
+centralized and lets a single config file drive both test cases.
+
 configure
 ~~~~~~~~~
 
@@ -39,6 +47,8 @@ the user (i.e., are not ``NotAvailable``).
 Repository-local example user configs are available at
 ``compass/landice/tests/ismip7_forcing/ismip7_forcing_test.cfg`` (AIS) and
 ``compass/landice/tests/ismip7_forcing/ismip7_forcing_test_gis.cfg`` (GrIS).
+The GrIS OCX scenario has a dedicated example
+``ismip7_forcing_ocx_gis.cfg`` in the same directory.
 These are intended for development/testing and include environment-specific
 paths.
 
diff --git a/docs/users_guide/landice/test_groups/ismip7_forcing.rst b/docs/users_guide/landice/test_groups/ismip7_forcing.rst
index fdec1eee14..6afb6b7d92 100644
--- a/docs/users_guide/landice/test_groups/ismip7_forcing.rst
+++ b/docs/users_guide/landice/test_groups/ismip7_forcing.rst
@@ -69,6 +69,39 @@ The AIS example enables both ocean processing modes and uses a 2015-2300
 processing window for both atmosphere and ocean thermal scenario forcing.
 The GrIS example enables scenario ocean processing only and uses 1980-2015.
 
+For the GrIS OCX (reanalysis) scenario, use the dedicated example config
+``compass/landice/tests/ismip7_forcing/ismip7_forcing_ocx_gis.cfg``. OCX has
+no distinct ESM model: it uses ``RACMO2.3p2-ERA`` for the atmosphere and
+``EN4`` for the ocean, both selected automatically when ``scenario = OCX``.
+The ``[ismip7] model`` option is ignored for OCX (set it to ``None``), so a
+single config file processes both the ``atmosphere`` and ``ocean_thermal``
+test cases, just like the ESM scenarios.
+
+.. _landice_ismip7_forcing_output:
+
+Output Layout
+-------------
+
+Processed forcing is written under ``output_base_path`` in a layout that the
+:ref:`landice_ismip7_run` test group can ingest directly:
+
+.. code-block:: none
+
+   {output_base_path}/{group}/atmosphere/{mesh}_SMB_{source}_{scenario}_{years}.nc
+   {output_base_path}/{group}/atmosphere/{mesh}_temperature_{source}_{scenario}_{years}.nc
+   {output_base_path}/{group}/atmosphere/{mesh}_runoff_...  (and the two gradients)
+   {output_base_path}/{group}/ocean_thermal_forcing/{mesh}_thermal_forcing_{source}_{scenario}_{years}.nc
+
+The ``group`` directory is ``{model}_{scenario}`` for the ESM scenarios and
+``{scenario}`` (i.e. ``OCX``) for OCX, whose atmosphere and ocean use
+different sources (``RACMO2.3p2-ERA`` and ``EN4``) but must share one
+directory. To feed :ref:`landice_ismip7_run`, point its ``forcing_basepath``
+at ``output_base_path`` for ESM scenarios, or its ``ocx_forcing_path`` at
+``{output_base_path}/OCX`` for OCX.
+
+Process only one year range into a given ``group`` directory: the run setup
+expects exactly one file per forcing field there.
+
 .. _landice_ismip7_forcing_input_data:
 
 Input Data
@@ -121,6 +154,20 @@ For GrIS ocean thermal (same 1km grid, 2D, yearly files):
 
    ocean/tf/v2/tf_GrIS_{model}_{scenario}_ocean_v2_{year}.nc
 
+The OCX (reanalysis) scenario follows the same directory layout as the ESM
+scenarios, but uses fixed reanalysis sources, data version ``v1``, and a
+named grid resolution in the ocean file names. For GrIS OCX the atmosphere
+source is ``RACMO2.3p2-ERA`` and the ocean source is ``EN4``:
+
+.. code-block:: none
+
+   acabf/v1/acabf_GrIS_RACMO2.3p2-ERA_OCX_SDBN1-1000m_v1_{year}.nc
+   ocean/tf/v1/tf_GrIS_EN4_OCX_ocean-1000m_v1_{year}.nc
+
+Set ``base_path_ismip7`` to the ``OCX`` directory and ``scenario = OCX``.
+The sources, version, and ocean grid token are selected automatically for
+OCX, and the ``[ismip7] model`` option is ignored (set it to ``None``).
+
 .. _landice_ismip7_forcing_config:
 
 config options