diff --git a/conda-env/ci.yml b/conda-env/ci.yml
index b84290d6f..87b384a23 100644
--- a/conda-env/ci.yml
+++ b/conda-env/ci.yml
@@ -25,7 +25,7 @@ dependencies:
   - pywavelets
   - scipy
   - xarray >=2024.3.0
-  - xcdat >=0.7.3,<1.0.0
+  - xcdat >=0.9.1,<1.0.0
   - xesmf >=0.8.7
   - xskillscore >=0.0.20
   # Testing
diff --git a/conda-env/dev.yml b/conda-env/dev.yml
index 69d304a0e..c8e756bb5 100644
--- a/conda-env/dev.yml
+++ b/conda-env/dev.yml
@@ -22,7 +22,7 @@ dependencies:
   - numpy >=2.0.0,<3.0.0
   - pywavelets
   - scipy
-  - xcdat >=0.7.3,<1.0.0
+  - xcdat >=0.9.1,<1.0.0
   - xesmf >=0.8.7
   - xskillscore >=0.0.20
   # Testing
diff --git a/e3sm_diags/driver/arm_diags_driver.py b/e3sm_diags/driver/arm_diags_driver.py
index 01aedf3dc..0f3553a68 100644
--- a/e3sm_diags/driver/arm_diags_driver.py
+++ b/e3sm_diags/driver/arm_diags_driver.py
@@ -279,11 +279,6 @@ def _run_diag_annual_cycle(parameter: ARMDiagsParameter) -> ARMDiagsParameter:
                     vars_funcs = _get_vars_funcs_for_derived_var(ds_ref, var)
                     target_var = list(vars_funcs.keys())[0][0]
 
-                    # NOTE: The bounds dimension can be "nv", which is not
-                    # currently recognized as a valid bounds dimension
-                    # by xcdat. We rename it to "bnds" to make it compatible.
-                    ds_ref = _rename_bounds_dim(ds_ref)
-
                     ds_ref_climo = ds_ref.temporal.climatology(target_var, "month")
                     da_ref = vars_funcs[(target_var,)](ds_ref_climo[target_var]).rename(
                         var
@@ -536,33 +531,3 @@ def _save_metrics_to_json(parameter: ARMDiagsParameter, metrics_dict: dict[str,
         json.dump(metrics_dict, outfile)
 
     logger.info(f"Metrics saved in: {abs_path}")
-
-
-def _rename_bounds_dim(ds: xr.Dataset) -> xr.Dataset:
-    """
-    Renames the bounds dimension "nv" to "bnds" in the given xarray.Dataset for
-    xCDAT compatibility.
-
-    This is a temporary workaround to ensure compatibility with xCDAT's bounds
-    handling. The bounds dimension "nv" is commonly used in datasets to
-    represent the number of vertices in a polygon, but xCDAT expects the
-    bounds dimension to be in `xcdat.bounds.VALID_BOUNDS_DIMS`. This function
-    renames "nv" to "bnds" to align with xCDAT's expectations.
-
-    Parameters
-    ----------
-    ds : xr.Dataset
-        The input xarray.Dataset which may contain a bounds dimension named "nv".
-
-    Returns
-    -------
-    xr.Dataset
-        A new xarray.Dataset with the "nv" dimension renamed to "bnds" if it
-        existed; otherwise, the original dataset copy.
-    """
-    ds_new = ds.copy()
-
-    if "nv" in ds_new.dims:
-        ds_new = ds_new.rename({"nv": "bnds"})
-
-    return ds_new
diff --git a/pyproject.toml b/pyproject.toml
index 0891a5210..11c3552c3 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -43,7 +43,7 @@ dependencies = [
   "pywavelets",
   "scipy",
   "xarray >=2024.03.0",
-  "xcdat >=0.7.3,<1.0.0",
+  "xcdat >=0.9.1,<1.0.0",
   "xesmf >=0.8.7",
   "xskillscore >=0.0.20",
 ]