diff --git a/.github/workflows/build-pages.yml b/.github/workflows/build-pages.yml
index d39ed7a..ade88d5 100644
--- a/.github/workflows/build-pages.yml
+++ b/.github/workflows/build-pages.yml
@@ -23,6 +23,10 @@ jobs:
           ${PIP} install --upgrade pip setuptools wheel
           ${PIP} install -r requirements.txt
           ${PIP} uninstall -y progressbar2
+      - name: Set up Node.js
+        uses: actions/setup-node@v4
+        with:
+          node-version: "22"
       - name: Build Book
         run: |
           jupyter-book build .
diff --git a/oggm-edu/glacier_water_resources.ipynb b/oggm-edu/glacier_water_resources.ipynb
index afbe9e7..c834002 100644
--- a/oggm-edu/glacier_water_resources.ipynb
+++ b/oggm-edu/glacier_water_resources.ipynb
@@ -79,50 +79,50 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "import matplotlib.pyplot as plt\n",
     "import seaborn as sns\n",
     "sns.set_context('notebook')  # plot defaults"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "import xarray as xr\n",
     "import salem\n",
     "import numpy as np\n",
     "import pandas as pd"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "import oggm.cfg\n",
     "from oggm import utils, workflow, tasks, graphics\n",
     "from oggm.sandbox.edu import run_constant_climate_with_bias"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# OGGM options\n",
     "oggm.cfg.initialize(logging_level='WARNING')\n",
     "oggm.cfg.PATHS['working_dir'] = utils.gettempdir(dirname='WaterResources')\n",
     "oggm.cfg.PARAMS['min_ice_thick_for_length'] = 1  # a glacier is defined when ice is thicker than 1m\n",
     "oggm.cfg.PARAMS['store_model_geometry'] = True"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -150,13 +150,13 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# Hintereisferner\n",
     "rgi_id = 'RGI60-11.00897'"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -174,14 +174,14 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# We pick the elevation-bands glaciers because they run a bit faster - but they create more step changes in the area outputs\n",
     "base_url = 'https://cluster.klima.uni-bremen.de/~oggm/gdirs/oggm_v1.6/L3-L5_files/2025.6/elev_bands/W5E5/per_glacier_spinup'\n",
     "gdir = workflow.init_glacier_directories([rgi_id], from_prepro_level=5, prepro_border=160, prepro_base_url=base_url)[0]"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -201,7 +201,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# One interactive plot below requires Bokeh\n",
     "try:\n",
@@ -218,9 +220,7 @@
     "    # The rest of the notebook works without this dependency\n",
     "    out = None\n",
     "out"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -231,13 +231,13 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "fls = gdir.read_pickle('model_flowlines')\n",
     "graphics.plot_modeloutput_section(fls);"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -257,14 +257,14 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "years = np.arange(400)\n",
     "temp_bias_ts = pd.Series(years * 0. - 2, index=years)\n",
     "temp_bias_ts.plot(); plt.xlabel('Year'); plt.ylabel('Temperature bias (°C)');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -277,7 +277,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# file identifier where the model output is saved\n",
     "file_id = '_spinup'\n",
@@ -289,9 +291,7 @@
     "                     y0=2009, halfsize=10,  # Period which we will average and constantly repeat\n",
     "                     store_monthly_hydro=True,  # Monthly outputs provide additional information\n",
     "                     output_filesuffix=file_id);  # an identifier for the output file, to read it later"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -312,14 +312,14 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "with xr.open_dataset(gdir.get_filepath('model_diagnostics', filesuffix=file_id)) as ds:\n",
     "    # The last step of hydrological output is NaN (we can't compute it for this year)\n",
     "    ds = ds.isel(time=slice(0, -1)).load()"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -330,12 +330,12 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "ds"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -348,15 +348,15 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "fig, axs = plt.subplots(nrows=2, figsize=(10, 7), sharex=True)\n",
     "ds.volume_m3.plot(ax=axs[0])\n",
     "ds.length_m.plot(ax=axs[1])\n",
     "axs[0].set_xlabel(''); axs[0].set_title(f'{rgi_id}'); axs[1].set_xlabel('Years');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -383,7 +383,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "ny_s = 50  # Start\n",
     "ny_t = 150  # Trend\n",
@@ -391,13 +393,13 @@
     "temp_bias_ts = np.concatenate([np.full(ny_s, -2.), np.linspace(-2, 0.5, ny_t), np.full(ny_e, 0.5)])\n",
     "temp_bias_ts = pd.Series(temp_bias_ts, index=np.arange(ny_s + ny_t + ny_e))\n",
     "temp_bias_ts.plot(); plt.xlabel('Year'); plt.ylabel('Temperature bias (°C)');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# file identifier where the model output is saved\n",
     "file_id = '_lin_temp'\n",
@@ -410,9 +412,7 @@
     "                     store_monthly_hydro=True,  # Monthly outputs provide additional information\n",
     "                     init_model_filesuffix='_spinup',  # We want to start from the glacier in equilibrium we created earlier\n",
     "                     output_filesuffix=file_id);  # an identifier for the output file, to read it later"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -423,27 +423,27 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "with xr.open_dataset(gdir.get_filepath('model_diagnostics', filesuffix=file_id)) as ds:\n",
     "    # The last step of hydrological output is NaN (we can't compute it for this year)\n",
     "    ds = ds.isel(time=slice(0, -1)).load()"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "fig, axs = plt.subplots(nrows=3, figsize=(10, 7), sharex=True)\n",
     "ds.volume_m3.plot(ax=axs[0])\n",
     "ds.length_m.plot(ax=axs[1])\n",
     "temp_bias_ts.plot(ax=axs[2], c='C3')\n",
     "axs[0].set_xlabel(''); axs[0].set_title(f'{rgi_id}'); axs[1].set_xlabel(''); axs[2].set_xlabel('Years'); axs[2].set_ylabel('Temperature bias (°C)');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -470,13 +470,13 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "sel_vars = [v for v in ds.variables if 'month_2d' not in ds[v].dims]\n",
     "df_annual = ds[sel_vars].to_dataframe()"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -487,7 +487,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# Select only the runoff variables\n",
     "runoff_vars = ['melt_off_glacier', 'melt_on_glacier', 'liq_prcp_off_glacier', 'liq_prcp_on_glacier']\n",
@@ -498,9 +500,7 @@
     "fig, ax = plt.subplots(figsize=(10, 3.5), sharex=True)\n",
     "df_runoff.sum(axis=1).plot(ax=ax)\n",
     "plt.ylabel('Runoff (Mt)'); plt.xlabel('Years'); plt.title(f'Total annual runoff for {rgi_id}');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -516,13 +516,13 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "f, ax = plt.subplots(figsize=(10, 6))\n",
     "df_runoff.plot.area(ax=ax, color=sns.color_palette('rocket')); plt.xlabel('Years'); plt.ylabel('Runoff (Mt)'); plt.title(rgi_id);"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -544,7 +544,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# Solution of the last question above\n",
     "glacier_mass = ds.volume_m3.to_series() * oggm.cfg.PARAMS['ice_density'] * 1e-9  # In Megatonnes, Mt\n",
@@ -553,9 +555,7 @@
     "plt.axhline(y=0, color='k', ls=':')\n",
     "plt.ylabel('Annual glacier mass change (Mt yr$^{-1}$)')\n",
     "plt.xlabel('Years'); plt.title('Glacier contribution to annual runoff');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -573,15 +573,15 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# Select only the runoff variables and convert them to megatonnes (instead of kg)\n",
     "monthly_runoff = ds['melt_off_glacier_monthly'] + ds['melt_on_glacier_monthly'] + ds['liq_prcp_off_glacier_monthly'] + ds['liq_prcp_on_glacier_monthly']\n",
     "monthly_runoff = monthly_runoff.rolling(time=31, center=True, min_periods=1).mean() * 1e-9\n",
     "monthly_runoff.clip(0).plot(cmap='Blues', cbar_kwargs={'label': 'Runoff (Mt)'}); plt.xlabel('Months'); plt.ylabel('Years');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -594,18 +594,20 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# Compute total precipitation (Snow + Liquid)\n",
     "tot_precip = ds['liq_prcp_off_glacier_monthly'] + ds['liq_prcp_on_glacier_monthly'] + ds['snowfall_off_glacier_monthly'] + ds['snowfall_on_glacier_monthly']\n",
     "tot_precip *= 1e-9  # in Mt"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# Plot these data at year 0\n",
     "yr = 0\n",
@@ -619,9 +621,7 @@
     "ax.fill_between(r.month_2d, r, p, where=(r > p), facecolor='C3', interpolate=True, alpha=0.5)\n",
     "plt.ylabel('Mt yr$^{-1}$'); plt.legend(loc='best')\n",
     "plt.xlabel('Month'); plt.title(f'Total monthly runoff and precipitation at year {yr}');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -641,7 +641,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "f, ax = plt.subplots(figsize=(10, 6))\n",
     "cmap = sns.color_palette('magma', 3)\n",
@@ -649,9 +651,7 @@
     "    monthly_runoff.sel(time=yr).plot(ax=ax, color=cmap[i], label=f'Year {yr}')\n",
     "plt.ylabel('Mt yr$^{-1}$'); plt.legend(loc='best')\n",
     "plt.xlabel('Month'); plt.title('Total monthly runoff change with time');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -692,7 +692,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# Make a dataframe out of the xarray dataset for simplicity\n",
     "df_monthly_runoff = pd.DataFrame(monthly_runoff.data, index=ds.time, columns=monthly_runoff.month_2d)\n",
@@ -708,9 +710,7 @@
     "f, ax = plt.subplots(figsize=(10, 6))\n",
     "df_monthly_runoff[['Spring (Apr-Jun)', 'Summer (Jul-Sep)', 'Winter (Oct-Mar)']].plot.area(ax=ax, color=sns.color_palette('rocket'))\n",
     "plt.xlabel('Years'); plt.ylabel('Runoff (Mt)'); plt.title('Runoff by season');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -721,7 +721,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "f, (ax1, ax2) = plt.subplots(2, 1, gridspec_kw={'height_ratios': [2, 1]}, figsize=(10, 8), sharex=True)\n",
     "\n",
@@ -735,9 +737,7 @@
     "(temp_bias_ts + 2).plot(ax=ax2, color='C2', label='Temperature', linewidth=3)\n",
     "ax2.legend(loc='lower right'); ax2.set_xlabel('Years'); ax2.set_ylabel('$\\Delta T (°C)$')\n",
     "sns.despine();"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
diff --git a/oggm-edu/glacier_water_resources_projections.ipynb b/oggm-edu/glacier_water_resources_projections.ipynb
index 9940ed6..e9f1ebb 100644
--- a/oggm-edu/glacier_water_resources_projections.ipynb
+++ b/oggm-edu/glacier_water_resources_projections.ipynb
@@ -50,49 +50,49 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "import matplotlib.pyplot as plt\n",
     "import seaborn as sns\n",
     "sns.set_context('notebook')  # plot defaults"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "import xarray as xr\n",
     "import salem\n",
     "import numpy as np\n",
     "import pandas as pd"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "import oggm.cfg\n",
     "from oggm import utils, workflow, tasks, graphics"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# OGGM options\n",
     "oggm.cfg.initialize(logging_level='WARNING')\n",
     "oggm.cfg.PATHS['working_dir'] = utils.gettempdir(dirname='WaterResources-Proj')\n",
     "oggm.cfg.PARAMS['min_ice_thick_for_length'] = 1  # a glacier is when ice thicker than 1m\n",
     "oggm.cfg.PARAMS['store_model_geometry'] = True"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -120,13 +120,13 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# Hintereisferner\n",
     "rgi_id = 'RGI60-11.00897'"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -144,14 +144,14 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# We pick the elevation-bands glaciers because they run a bit faster - but they create more step changes in the area outputs\n",
     "base_url = 'https://cluster.klima.uni-bremen.de/~oggm/gdirs/oggm_v1.6/L3-L5_files/2025.6/elev_bands/W5E5/per_glacier_spinup'\n",
     "gdir = workflow.init_glacier_directories([rgi_id], from_prepro_level=5, prepro_border=160, prepro_base_url=base_url)[0]"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -171,7 +171,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# One interactive plot below requires Bokeh\n",
     "try:\n",
@@ -188,20 +190,18 @@
     "    # The rest of the notebook works without this dependency\n",
     "    out = None\n",
     "out"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {
     "tags": []
    },
+   "outputs": [],
    "source": [
     "type(out)"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -213,30 +213,32 @@
   {
    "cell_type": "markdown",
    "metadata": {},
-   "source": "The glacier directories we are using come with a run called a \"spinup\", which is a simulation running from 1980 to 2020 while approximately matching glacier volume and area at the glacier outline date:"
+   "source": [
+    "The glacier directories we are using come with a run called a \"spinup\", which is a simulation running from 1980 to 2020 while approximately matching glacier volume and area at the glacier outline date:"
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "with xr.open_dataset(gdir.get_filepath('model_diagnostics', filesuffix='_spinup_historical')) as ds:\n",
     "    ds = ds.sel(time=slice(2000, 2020)).load()"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "fig, axs = plt.subplots(nrows=2, figsize=(10, 7), sharex=True)\n",
     "ds.volume_m3.plot(ax=axs[0])\n",
     "ds.area_m2.plot(ax=axs[1])\n",
     "axs[1].scatter(gdir.rgi_date+1, gdir.rgi_area_m2, label=f'Glacier area in {gdir.rgi_date+1}')\n",
     "axs[0].set_xlabel(''); axs[0].set_title(f'{rgi_id}'); axs[1].set_xlabel('Years'); plt.legend();"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -263,7 +265,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# file identifier where the model output is saved\n",
     "file_id = '_ct'\n",
@@ -276,9 +280,7 @@
     "                     init_model_filesuffix='_spinup_historical',  # Start from the glacier in 2020 after spinup\n",
     "                     store_monthly_hydro=True,  # Monthly outputs provide additional information\n",
     "                     output_filesuffix=file_id);  # an identifier for the output file to read it later"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -289,14 +291,14 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "with xr.open_dataset(gdir.get_filepath('model_diagnostics', filesuffix=file_id)) as ds:\n",
     "    # The last step of hydrological output is NaN (we can't compute it for this year)\n",
     "    ds = ds.isel(time=slice(0, -1)).load()"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -307,12 +309,12 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "ds"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -325,15 +327,15 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "fig, axs = plt.subplots(nrows=2, figsize=(10, 7), sharex=True)\n",
     "ds.volume_m3.plot(ax=axs[0])\n",
     "ds.area_m2.plot(ax=axs[1])\n",
     "axs[0].set_xlabel(''); axs[0].set_title(f'{rgi_id}'); axs[1].set_xlabel('Years');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -360,22 +362,26 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "sel_vars = [v for v in ds.variables if 'month_2d' not in ds[v].dims]\n",
     "df_annual = ds[sel_vars].to_dataframe()"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
    "metadata": {},
-   "source": "Then we can select the hydrological variables and sum them to get the total annual runoff:"
+   "source": [
+    "Then we can select the hydrological variables and sum them to get the total annual runoff:"
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# Select only the runoff variables\n",
     "runoff_vars = ['melt_off_glacier', 'melt_on_glacier', 'liq_prcp_off_glacier', 'liq_prcp_on_glacier']\n",
@@ -384,9 +390,7 @@
     "fig, ax = plt.subplots(figsize=(10, 3.5), sharex=True)\n",
     "df_runoff.sum(axis=1).plot(ax=ax)\n",
     "plt.ylabel('Mt'); plt.xlabel('Years'); plt.title(f'Total annual runoff for {rgi_id}');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -402,13 +406,13 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "f, ax = plt.subplots(figsize=(10, 6))\n",
     "df_runoff.plot.area(ax=ax, color=sns.color_palette(\"rocket\")); plt.xlabel('Years'); plt.ylabel('Runoff (Mt)'); plt.title(rgi_id);"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -435,7 +439,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "glacier_mass = ds.volume_m3.to_series() * oggm.cfg.PARAMS['ice_density'] * 1e-9  # In Megatonnes, Mt\n",
     "\n",
@@ -443,9 +449,7 @@
     "plt.axhline(y=0, color='k', ls=':')\n",
     "plt.ylabel('Annual glacier mass change (Mt yr$^{-1}$)')\n",
     "plt.xlabel('Years'); plt.title('Glacier contribution to annual runoff');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -471,16 +475,16 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# Select only the runoff variables and convert them to megatonnes (instead of kg)\n",
     "monthly_runoff = (ds['melt_off_glacier_monthly'] + ds['melt_on_glacier_monthly'] +\n",
     "                  ds['liq_prcp_off_glacier_monthly'] + ds['liq_prcp_on_glacier_monthly'])\n",
     "monthly_runoff *= 1e-9\n",
     "monthly_runoff.clip(0).plot(cmap='Blues', cbar_kwargs={'label': 'Mt'}); plt.xlabel('Months'); plt.ylabel('Years');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -493,15 +497,15 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "monthly_runoff.sel(time=[0, 30, 99]).plot(hue='time')\n",
     "plt.title('Annual cycle')\n",
     "plt.xlabel('Month')\n",
     "plt.ylabel('Runoff (Mt)');"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -523,22 +527,22 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# We keep Hintereisferner from earlier, but also add a new glacier\n",
     "rgi_ids = [rgi_id, 'RGI60-14.23809']"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "gdirs = workflow.init_glacier_directories(rgi_ids, from_prepro_level=5, prepro_border=160, prepro_base_url=base_url)"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -549,12 +553,12 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "gdirs"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -565,7 +569,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# One interactive plot below requires Bokeh\n",
     "try:\n",
@@ -577,9 +583,7 @@
     "    # The rest of the notebook works without this dependency\n",
     "    out = None\n",
     "out"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -599,7 +603,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# you can choose one of these 5 different GCMs:\n",
     "# 'gfdl-esm4_r1i1p1f1', 'mpi-esm1-2-hr_r1i1p1f1', 'mri-esm2-0_r1i1p1f1' (\"low sensitivity\" models, within typical ranges from AR6)\n",
@@ -617,9 +623,7 @@
     "                                 # recognize the climate file for later\n",
     "                                 output_filesuffix=f'_ISIMIP3b_{member}_{ssp}'\n",
     "                                 );"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -639,7 +643,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "for ssp in ['ssp126', 'ssp370', 'ssp585']:\n",
     "    rid = f'_ISIMIP3b_{member}_{ssp}'\n",
@@ -655,13 +661,13 @@
     "                             output_filesuffix=rid,\n",
     "                             # add monthly diagnostics\n",
     "                             store_monthly_hydro=True);"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# Create the figure\n",
     "f, ax = plt.subplots(figsize=(18, 7), sharex=True)\n",
@@ -687,9 +693,7 @@
     "    df_roll.plot(ax=ax, label=ssp, color=sns.color_palette(\"rocket\")[i])\n",
     "\n",
     "ax.set_ylabel('Annual runoff (Mt)'); ax.set_xlabel('Year'); plt.title(gdirs[0].rgi_id); plt.legend();"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
@@ -700,7 +704,9 @@
   },
   {
    "cell_type": "code",
+   "execution_count": null,
    "metadata": {},
+   "outputs": [],
    "source": [
     "# Create the figure\n",
     "f, ax = plt.subplots(figsize=(18, 7), sharex=True)\n",
@@ -727,14 +733,14 @@
     "    df_roll.plot(ax=ax, label=ssp, color=sns.color_palette(\"rocket\")[i])\n",
     "\n",
     "ax.set_ylabel('Annual runoff (Mt)'); ax.set_xlabel('Year'); plt.title(gdirs[1].rgi_id); plt.legend();"
-   ],
-   "outputs": [],
-   "execution_count": null
+   ]
   },
   {
    "cell_type": "markdown",
    "metadata": {},
-   "source": "Unlike for Hintereisferner, these projections indicate that the annual runoff will increase in all the scenarios for the first half of the century. The higher SSP scenarios can reach **peak water** later in the century, since the excess melt can continue to increase. For the lower SSP scenarios on the other hand, the glacier might be approaching a new equilibrium, which reduces the runoff earlier in the century ([Rounce et al., 2020](https://www.frontiersin.org/articles/10.3389/feart.2019.00331/full)). After **peak water** is reached (SSP126: ~2050, SSP585: ~2060 in these projections), the annual runoff begins to decrease. This decrease occurs because the shrinking glacier is no longer able to support the high levels of melt."
+   "source": [
+    "Unlike for Hintereisferner, these projections indicate that the annual runoff will increase in all the scenarios for the first half of the century. The higher SSP scenarios can reach **peak water** later in the century, since the excess melt can continue to increase. For the lower SSP scenarios on the other hand, the glacier might be approaching a new equilibrium, which reduces the runoff earlier in the century ([Rounce et al., 2020](https://www.frontiersin.org/articles/10.3389/feart.2019.00331/full)). After **peak water** is reached (SSP126: ~2050, SSP585: ~2060 in these projections), the annual runoff begins to decrease. This decrease occurs because the shrinking glacier is no longer able to support the high levels of melt."
+   ]
   },
   {
    "cell_type": "markdown",
diff --git a/requirements.txt b/requirements.txt
index 3df5565..3f5baca 100644
--- a/requirements.txt
+++ b/requirements.txt
@@ -1,4 +1,4 @@
-jupyter-book
+jupyter-book<2
 oggm
 bokeh
 panel