Run OGGM with GCM data
Contents
Run OGGM with GCM data#
In this example, we illustrate how to do a typical “projection run”, i.e. using GCM data (here CMIP5 and CMIP6). There are two important steps:
simulate all glaciers from their inventory date up to a “present day” geometry (New in version 1.4: this is already available in the pre-processed directories)
simulate their future evolution from the present day state
# Libs
import xarray as xr
import matplotlib.pyplot as plt
# Locals
import oggm.cfg as cfg
from oggm import utils, workflow, tasks
from oggm.core import gcm_climate
/usr/local/pyenv/versions/3.10.10/lib/python3.10/site-packages/oggm/core/gcm_climate.py:3: FutureWarning: The module `oggm.core.gcm_climate` has moved to oggm.shop.gcm_climate. This compatibility module will be removed in future OGGM versions
warnings.warn('The module `oggm.core.gcm_climate` has moved to '
Pre-processed directories#
Let’s do a run for two Himalayan glaciers: Ngojumba and Khumbu.
# Initialize OGGM and set up the default run parameters
cfg.initialize(logging_level='WARNING')
# Here we override some of the default parameters
# How many grid points around the glacier?
# Make it large if you expect your glaciers to grow large:
# here, 40 is more than enough
cfg.PARAMS['border'] = 40
# Local working directory (where OGGM will write its output)
cfg.PATHS['working_dir'] = utils.gettempdir('OGGM_gcm_run', reset=True)
# RGI glaciers: Ngojumba and Khumbu
rgi_ids = utils.get_rgi_glacier_entities(['RGI60-15.03473', 'RGI60-15.03733'])
# Go - get the pre-processed glacier directories
gdirs = workflow.init_glacier_directories(rgi_ids, from_prepro_level=5)
2023-03-07 12:46:23: oggm.cfg: Reading default parameters from the OGGM `params.cfg` configuration file.
2023-03-07 12:46:23: oggm.cfg: Multiprocessing switched OFF according to the parameter file.
2023-03-07 12:46:23: oggm.cfg: Multiprocessing: using all available processors (N=2)
2023-03-07 12:46:30: oggm.workflow: init_glacier_directories from prepro level 5 on 2 glaciers.
2023-03-07 12:46:30: oggm.workflow: Execute entity tasks [gdir_from_prepro] on 2 glaciers
---------------------------------------------------------------------------
TypeError Traceback (most recent call last)
Cell In[2], line 17
14 rgi_ids = utils.get_rgi_glacier_entities(['RGI60-15.03473', 'RGI60-15.03733'])
16 # Go - get the pre-processed glacier directories
---> 17 gdirs = workflow.init_glacier_directories(rgi_ids, from_prepro_level=5)
File /usr/local/pyenv/versions/3.10.10/lib/python3.10/site-packages/oggm/workflow.py:533, in init_glacier_directories(rgidf, reset, force, from_prepro_level, prepro_border, prepro_rgi_version, prepro_base_url, from_tar, delete_tar)
531 if cfg.PARAMS['dl_verify']:
532 utils.get_dl_verify_data('cluster.klima.uni-bremen.de')
--> 533 gdirs = execute_entity_task(gdir_from_prepro, entities,
534 from_prepro_level=from_prepro_level,
535 prepro_border=prepro_border,
536 prepro_rgi_version=prepro_rgi_version,
537 base_url=prepro_base_url)
538 else:
539 # We can set the intersects file automatically here
540 if (cfg.PARAMS['use_intersects'] and
541 len(cfg.PARAMS['intersects_gdf']) == 0 and
542 not from_tar):
File /usr/local/pyenv/versions/3.10.10/lib/python3.10/site-packages/oggm/workflow.py:191, in execute_entity_task(task, gdirs, **kwargs)
187 if ng > 3:
188 log.workflow('WARNING: you are trying to run an entity task on '
189 '%d glaciers with multiprocessing turned off. OGGM '
190 'will run faster with multiprocessing turned on.', ng)
--> 191 out = [pc(gdir) for gdir in gdirs]
193 return out
File /usr/local/pyenv/versions/3.10.10/lib/python3.10/site-packages/oggm/workflow.py:191, in <listcomp>(.0)
187 if ng > 3:
188 log.workflow('WARNING: you are trying to run an entity task on '
189 '%d glaciers with multiprocessing turned off. OGGM '
190 'will run faster with multiprocessing turned on.', ng)
--> 191 out = [pc(gdir) for gdir in gdirs]
193 return out
File /usr/local/pyenv/versions/3.10.10/lib/python3.10/site-packages/oggm/workflow.py:108, in _pickle_copier.__call__(self, arg)
106 for func in self.call_func:
107 func, kwargs = func
--> 108 res = self._call_internal(func, arg, kwargs)
109 return res
File /usr/local/pyenv/versions/3.10.10/lib/python3.10/site-packages/oggm/workflow.py:102, in _pickle_copier._call_internal(self, call_func, gdir, kwargs)
99 gdir, gdir_kwargs = gdir
100 kwargs.update(gdir_kwargs)
--> 102 return call_func(gdir, **kwargs)
File /usr/local/pyenv/versions/3.10.10/lib/python3.10/site-packages/oggm/workflow.py:251, in gdir_from_prepro(entity, from_prepro_level, prepro_border, prepro_rgi_version, base_url)
248 tar_base = utils.get_prepro_gdir(prepro_rgi_version, rid, prepro_border,
249 from_prepro_level, base_url=base_url)
250 from_tar = os.path.join(tar_base.replace('.tar', ''), rid + '.tar.gz')
--> 251 return oggm.GlacierDirectory(entity, from_tar=from_tar)
File /usr/local/pyenv/versions/3.10.10/lib/python3.10/site-packages/oggm/utils/_workflow.py:2529, in GlacierDirectory.__init__(self, rgi_entity, base_dir, reset, from_tar, delete_tar)
2525 raise RuntimeError('RGI Version not supported: '
2526 '{}'.format(self.rgi_version))
2528 # remove spurious characters and trailing blanks
-> 2529 self.name = filter_rgi_name(name)
2531 # region
2532 reg_names, subreg_names = parse_rgi_meta(version=self.rgi_version[0])
File /usr/local/pyenv/versions/3.10.10/lib/python3.10/site-packages/oggm/utils/_funcs.py:734, in filter_rgi_name(name)
728 def filter_rgi_name(name):
729 """Remove spurious characters and trailing blanks from RGI glacier name.
730
731 This seems to be unnecessary with RGI V6
732 """
--> 734 if name is None or len(name) == 0:
735 return ''
737 if name[-1] in ['À', 'È', 'è', '\x9c', '3', 'Ð', '°', '¾',
738 '\r', '\x93', '¤', '0', '`', '/', 'C', '@',
739 'Å', '\x06', '\x10', '^', 'å', ';']:
TypeError: object of type 'float' has no len()
New in version 1.4: the _historical runs#
As of v1.4, the level 5 files now come with a pre-computed model run from the RGI outline date to the last possible date given by the historical climate data. These files are stored in the directory with a _historical suffix. Let’s compile them for our two glaciers:
ds = utils.compile_run_output(gdirs, input_filesuffix='_historical')
ds.volume.plot(hue='rgi_id');
Each RGI glacier has an “inventory date”, the time at which the ouline is valid. It can change between glaciers:
gdirs[0].rgi_date, gdirs[1].rgi_date
ds.time[-1]
The command that was run on preprocessing was the run_from_climate_data which will run OGGM from this inventory date up to a desired date (here, the last year with valid historical climate data). So the data above is strictly equivalent to:
workflow.execute_entity_task(tasks.run_from_climate_data, gdirs,
output_filesuffix='_my_spinup', # to use the files as input later on
);
ds = utils.compile_run_output(gdirs, input_filesuffix='_my_spinup')
ds.volume.plot(hue='rgi_id');
One thing to remember here is that although we try hard to avoid spin-up issues, the glacier after the inversion is not in a perfect dynamical state. Some variable (in particular glacier length) might need some years to adjust:
ds.length.plot(hue='rgi_id');
Download and process GCM data from CMIP5#
A typical use case for OGGM will be to use climate model output (here CMIP5). We use some of the files we store online here, but you can use whichever you want. The above mentioned storage contains information about the data, how to cite them and tabular summaries of the available GCMs.
We download the files and bias correct them to the calibration data (here: CRU, see process_gcm_data which is called by process_cmip_data). This step is very important, as the model is very sensitive to temperature variability.
bp = 'https://cluster.klima.uni-bremen.de/~oggm/cmip5-ng/pr/pr_mon_CCSM4_{}_r1i1p1_g025.nc'
bt = 'https://cluster.klima.uni-bremen.de/~oggm/cmip5-ng/tas/tas_mon_CCSM4_{}_r1i1p1_g025.nc'
for rcp in ['rcp26', 'rcp45', 'rcp60', 'rcp85']:
# Download the files
ft = utils.file_downloader(bt.format(rcp))
fp = utils.file_downloader(bp.format(rcp))
# bias correct them
workflow.execute_entity_task(gcm_climate.process_cmip_data, gdirs,
filesuffix='_CMIP5_CCSM4_{}'.format(rcp), # recognize the climate file for later
fpath_temp=ft, # temperature projections
fpath_precip=fp, # precip projections
);
Projection runs#
We now run OGGM under various scenarios and start from the end year of the historical run:
for rcp in ['rcp26', 'rcp45', 'rcp60', 'rcp85']:
rid = '_CMIP5_CCSM4_{}'.format(rcp)
workflow.execute_entity_task(tasks.run_from_climate_data, gdirs, ys=2020,
climate_filename='gcm_data', # use gcm_data, not climate_historical
climate_input_filesuffix=rid, # use the chosen scenario
init_model_filesuffix='_historical', # this is important! Start from 2020 glacier
output_filesuffix=rid, # recognize the run for later
);
Plot model output#
f, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 4))
for rcp in ['rcp26', 'rcp45', 'rcp60', 'rcp85']:
rid = '_CMIP5_CCSM4_{}'.format(rcp)
ds = utils.compile_run_output(gdirs, input_filesuffix=rid)
ds.isel(rgi_id=0).volume.plot(ax=ax1, label=rcp);
ds.isel(rgi_id=1).volume.plot(ax=ax2, label=rcp);
plt.legend();
New: and now from CMIP6!#
Very similar steps, although some links have changed and RCPs are now SSPs! See this table for a summary of available GCMs.
bp = 'https://cluster.klima.uni-bremen.de/~oggm/cmip6/GCM/CESM2/CESM2_{}_r1i1p1f1_pr.nc'
bt = 'https://cluster.klima.uni-bremen.de/~oggm/cmip6/GCM/CESM2/CESM2_{}_r1i1p1f1_tas.nc'
for ssp in ['ssp126', 'ssp585']: # Removed 'ssp245', 'ssp370' because the files are large!
# Download the files
ft = utils.file_downloader(bt.format(ssp))
fp = utils.file_downloader(bp.format(ssp))
# bias correct them
workflow.execute_entity_task(gcm_climate.process_cmip_data, gdirs,
filesuffix='_CMIP6_CESM2_{}'.format(ssp), # recognize the climate file for later
fpath_temp=ft, # temperature projections
fpath_precip=fp, # precip projections
);
for ssp in ['ssp126', 'ssp585']:
rid = '_CMIP6_CESM2_{}'.format(ssp)
workflow.execute_entity_task(tasks.run_from_climate_data, gdirs, ys=2020,
climate_filename='gcm_data', # use gcm_data, not climate_historical
climate_input_filesuffix=rid, # use the chosen scenario
init_model_filesuffix='_historical', # this is important! Start from 2020 glacier
output_filesuffix=rid, # recognize the run for later
);
f, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 4))
for ssp in ['ssp126', 'ssp585']:
rid = '_CMIP6_CESM2_{}'.format(ssp)
ds = utils.compile_run_output(gdirs, input_filesuffix=rid)
ds.isel(rgi_id=0).volume.plot(ax=ax1, label=ssp);
ds.isel(rgi_id=1).volume.plot(ax=ax2, label=ssp);
plt.legend();
What’s next?#
see also the tutorial on Merge, analyse and visualize OGGM GCM runs
return to the OGGM documentation
back to the table of contents