Work on figures

pism · Jan 11, 2025 · 3cbc60a · 3cbc60a
1 parent 0c266e6
commit 3cbc60a
Show file tree

Hide file tree

Showing 3 changed files with 614 additions and 18 deletions.
diff --git a/analysis/analyze_scalar.py b/analysis/analyze_scalar.py
@@ -134,7 +134,7 @@
         "--reference_date",
         help="""Reference date.""",
         type=str,
-        default="2020-01-1",
+        default="1986-01-1",
     )
     parser.add_argument(
         "--n_jobs",
@@ -275,7 +275,7 @@
         obs_mankoff_basins = set(observed_mankoff_ds.basin.values)
         obs_grace_basins = set(observed_grace_ds.basin.values)
 
-        simulated = filtered_ds
+        simulated = simulated_retreat_ds
 
         sim_basins = set(simulated.basin.values)
         sim_grace = set(simulated.basin.values)
@@ -305,12 +305,21 @@
             {"time": resampling_frequency}
         ).mean()
 
-        obs_mean_vars_mankoff: list[str] = ["grounding_line_flux", "mass_balance"]
+        obs_mean_vars_mankoff: list[str] = [
+            "grounding_line_flux",
+            "mass_balance",
+            "surface_mass_balance",
+        ]
         obs_std_vars_mankoff: list[str] = [
             "grounding_line_flux_uncertainty",
             "mass_balance_uncertainty",
+            "surface_mass_balance_uncertainty",
+        ]
+        sim_vars_mankoff: list[str] = [
+            "grounding_line_flux",
+            "mass_balance",
+            "surface_mass_balance",
         ]
-        sim_vars_mankoff: list[str] = ["grounding_line_flux", "mass_balance"]
 
         sim_plot_vars = (
             [ragis_config["Cumulative Variables"]["cumulative_mass_balance"]]
@@ -332,17 +341,18 @@
             fudge_factor=mankoff_fudge_factor,
             params=params,
         )
-
         for filter_var in obs_mean_vars_mankoff:
             plot_basins(
-                observed_mankoff_basins_resampled_ds,
-                simulated_prior_mankoff[sim_plot_vars],
+                observed_mankoff_basins_resampled_ds.load(),
+                simulated_prior_mankoff[sim_plot_vars].load(),
                 simulated_posterior_mankoff.sel({"filtered_by": filter_var})[
                     sim_plot_vars
-                ],
+                ].load(),
                 filter_var=filter_var,
                 filter_range=filter_range,
+                figsize=(4.2, 3.2),
                 fig_dir=fig_dir,
+                fontsize=6,
                 fudge_factor=mankoff_fudge_factor,
                 reference_date=reference_date,
                 config=config,
@@ -369,15 +379,16 @@
 
         for filter_var in obs_mean_vars_grace:
             plot_basins(
-                observed_grace_basins_resampled_ds,
-                simulated_prior_grace[sim_plot_vars],
+                observed_grace_basins_resampled_ds.load(),
+                simulated_prior_grace[sim_plot_vars].load(),
                 simulated_posterior_grace.sel({"filtered_by": filter_var})[
                     sim_plot_vars
-                ],
+                ].load(),
                 filter_var=filter_var,
                 filter_range=filter_range,
                 fudge_factor=grace_fudge_factor,
                 fig_dir=fig_dir,
+                figsize=(4.2, 3.2),
                 config=config,
             )
 

diff --git a/calibration/adjust_tillwat.py b/calibration/adjust_tillwat.py
@@ -23,6 +23,7 @@
 from argparse import ArgumentDefaultsHelpFormatter, ArgumentParser
 from pathlib import Path
 
+import pint_xarray  # pylint: disable=unused-import
 import xarray as xr
 from dask.diagnostics import ProgressBar
 
@@ -60,14 +61,28 @@
     speed_threshold = options.speed_threshold
 
     pism_ds = xr.open_dataset(infile)
+
+    pism_config = pism_ds.pism_config
+    dx = pism_config.attrs["grid.dx"]
+    dy = pism_config.attrs["grid.dy"]
+    area = xr.DataArray(dx * dy).pint.quantify("m^2")
+    area_km2 = area.pint.to("km^2")
+
+    tillwat_mask = (pism_ds["tillwat"] > 0.0).where(pism_ds["mask"] == 2)
     speed_ds = xr.open_dataset(speed_file)
     speed_da = speed_ds[speed_var].interp_like(pism_ds["tillwat"])
     tillwat_da = pism_ds["tillwat"]
     pism_ds["tillwat"] = tillwat_da.where(
         speed_da <= speed_threshold, other=tillwat_max
     )
+    tillwat_mask_updated = (pism_ds["tillwat"] > 0.0).where(pism_ds["mask"] == 2)
 
+    twa_o = tillwat_mask.sum() * area_km2
+    twa_u = tillwat_mask_updated.sum() * area_km2
+    print(f"Original temperate area {twa_o.values}")
+    print(f"New temperate area {twa_u.values}")
+    print(f"Difference {(twa_u.values / twa_o.values) * 100 - 100 }%")
     comp = {"zlib": True, "complevel": 2, "_FillValue": None}
     encoding = {var: comp for var in pism_ds.data_vars}
     with ProgressBar():
-        pism_ds.to_netcdf(outfile, encoding=encoding)
+        pism_ds.pint.dequantify().to_netcdf(outfile, encoding=encoding)