Skip to content

Commit

Permalink
Add itslive to prepro (#1539)
Browse files Browse the repository at this point in the history
---------

Co-authored-by: Lilian Schuster <37300567+lilianschuster@users.noreply.github.com>
fmaussion and lilianschuster authored Mar 9, 2023
1 parent 410bc78 commit ee6904a
Showing 16 changed files with 181 additions and 73 deletions.
4 changes: 2 additions & 2 deletions docs/api.rst
Original file line number Diff line number Diff line change
@@ -150,7 +150,7 @@ the majority of OGGM's tasks). They are parallelizable.
tasks.mass_conservation_inversion
tasks.filter_inversion_output
tasks.get_inversion_volume
tasks.compute_velocities
tasks.compute_inversion_velocities
tasks.distribute_thickness_per_altitude
tasks.distribute_thickness_interp
tasks.find_inversion_calving_from_any_mb
@@ -444,7 +444,7 @@ including the model flowlines. This is achieved by choosing preprocessing level
# The working directory is where OGGM will store the run's data
cfg.PATHS['working_dir'] = os.path.join(gettempdir(), 'Docs_GlacierDir2')
# The base url is where to find the pre-processed directories
base_url = 'https://cluster.klima.uni-bremen.de/~oggm/gdirs/oggm_v1.6/exps/L3-L5_files/W5E5_melt_calib'
base_url = 'https://cluster.klima.uni-bremen.de/~oggm/gdirs/oggm_v1.6/L3-L5_files/2023.1/elev_bands/W5E5'
gdirs = workflow.init_glacier_directories('RGI60-11.00897',
from_prepro_level=5,
prepro_base_url=base_url,
2 changes: 1 addition & 1 deletion docs/shop.rst
Original file line number Diff line number Diff line change
@@ -499,7 +499,7 @@ have access to the timeseries through the glacier directory:

.. ipython:: python
base_url = 'https://cluster.klima.uni-bremen.de/~oggm/gdirs/oggm_v1.6/exps/L3-L5_files/W5E5_melt_calib/'
base_url = 'https://cluster.klima.uni-bremen.de/~oggm/gdirs/oggm_v1.6/L3-L5_files/2023.1/elev_bands/W5E5'
gdir = workflow.init_glacier_directories('RGI60-11.00897',
from_prepro_level=3,
prepro_base_url=base_url,
2 changes: 1 addition & 1 deletion docs/structure.rst
Original file line number Diff line number Diff line change
@@ -80,7 +80,7 @@ downloading and extracting these data locally:
# Where to fetch the pre-processed directories - this can be changed
server_url = 'https://cluster.klima.uni-bremen.de/~oggm/gdirs/'
experiment_url = 'oggm_v1.6/exps/L3-L5_files/W5E5_melt_calib'
experiment_url = 'oggm_v1.6/L3-L5_files/2023.1/elev_bands/W5E5'
base_url = server_url + experiment_url
.. ipython:: python
2 changes: 1 addition & 1 deletion docs/whats-new.rst
Original file line number Diff line number Diff line change
@@ -528,7 +528,7 @@ Enhancements
- New function ``cfg.add_to_basenames`` now allows users to define their own
entries in glacier directories (:issue:`731`).
By `Fabien Maussion <https://github.com/fmaussion>`_.
- New function ``inversion.compute_velocities`` writes the section and
- New function ``inversion.compute_inversion_velocities`` writes the section and
surface veloicites in the inversion output (:issue:`876`).
By `Beatriz Recinos <https://github.com/bearecinos>`_.
- Added ASTER v3 as optional DEM. Requires credentials to
22 changes: 20 additions & 2 deletions oggm/cli/prepro_levels.py
Original file line number Diff line number Diff line change
@@ -80,8 +80,9 @@ def run_prepro_levels(rgi_version=None, rgi_reg=None, border=None,
elev_bands=False, centerlines=False,
override_params=None,
mb_calibration_strategy='informed_threestep',
add_consensus_thickness=False, add_millan_thickness=False,
add_millan_velocity=False, add_hugonnet_dhdt=False,
add_consensus_thickness=False, add_itslive_velocity=False,
add_millan_thickness=False, add_millan_velocity=False,
add_hugonnet_dhdt=False,
start_level=None, start_base_url=None, max_level=5,
logging_level='WORKFLOW', disable_dl_verify=False,
dynamic_spinup=False, err_dmdtda_scaling_factor=1,
@@ -133,6 +134,9 @@ def run_prepro_levels(rgi_version=None, rgi_reg=None, border=None,
add_consensus_thickness : bool
adds (reprojects) the consensus estimates thickness to the glacier
directories. With elev_bands=True, the data will also be binned.
add_itslive_velocity : bool
adds (reprojects) the ITS_LIVE velocity to the glacier
directories. With elev_bands=True, the data will also be binned.
add_millan_thickness : bool
adds (reprojects) the millan thickness to the glacier
directories. With elev_bands=True, the data will also be binned.
@@ -420,6 +424,10 @@ def _time_log():
from oggm.shop.bedtopo import add_consensus_thickness
workflow.execute_entity_task(add_consensus_thickness, gdirs)
bin_variables.append('consensus_ice_thickness')
if add_itslive_velocity:
from oggm.shop.its_live import velocity_to_gdir
workflow.execute_entity_task(velocity_to_gdir, gdirs)
bin_variables.append('itslive_v')
if add_millan_thickness:
from oggm.shop.millan22 import thickness_to_gdir
workflow.execute_entity_task(thickness_to_gdir, gdirs)
@@ -479,6 +487,10 @@ def _time_log():
opath = os.path.join(sum_dir, 'glacier_statistics_{}.csv'.format(rgi_reg))
utils.compile_glacier_statistics(gdirs, path=opath)

if add_itslive_velocity:
from oggm.shop.its_live import compile_itslive_statistics
opath = os.path.join(sum_dir, 'itslive_statistics_{}.csv'.format(rgi_reg))
compile_itslive_statistics(gdirs, path=opath)
if add_millan_thickness or add_millan_velocity:
from oggm.shop.millan22 import compile_millan_statistics
opath = os.path.join(sum_dir, 'millan_statistics_{}.csv'.format(rgi_reg))
@@ -793,6 +805,11 @@ def parse_args(args):
'estimates to the glacier directories. '
'With --elev-bands, the data will also be '
'binned.')
parser.add_argument('--add-itslive-velocity', nargs='?', const=True, default=False,
help='adds (reprojects) the ITS_LIVE velocity '
'estimates to the glacier directories. '
'With --elev-bands, the data will also be '
'binned.')
parser.add_argument('--add-millan-thickness', nargs='?', const=True, default=False,
help='adds (reprojects) the millan thickness '
'estimates to the glacier directories. '
@@ -889,6 +906,7 @@ def parse_args(args):
centerlines=args.centerlines,
add_consensus_thickness=args.add_consensus_thickness,
add_millan_thickness=args.add_millan_thickness,
add_itslive_velocity=args.add_itslive_velocity,
add_millan_velocity=args.add_millan_velocity,
add_hugonnet_dhdt=args.add_hugonnet_dhdt,
disable_dl_verify=args.disable_dl_verify,
8 changes: 3 additions & 5 deletions oggm/core/flowline.py
Original file line number Diff line number Diff line change
@@ -3254,9 +3254,8 @@ def run_random_climate(gdir, nyears=1000, y0=None, halfsize=15,
the glacier directory to process
nyears : int
length of the simulation
y0 : int, optional
central year of the random climate period. The default is to be
centred on t*.
y0 : int
central year of the random climate period. Has to be set!
halfsize : int, optional
the half-size of the time window (window size = 2 * halfsize + 1)
bias : float
@@ -3370,8 +3369,7 @@ def run_constant_climate(gdir, nyears=1000, y0=None, halfsize=15,
length of the simulation (default: as long as needed for reaching
equilibrium)
y0 : int
central year of the requested climate period. The default is to be
centred on t*.
central year of the random climate period. Has to be set!
halfsize : int, optional
the half-size of the time window (window size = 2 * halfsize + 1)
bias : float
46 changes: 31 additions & 15 deletions oggm/core/inversion.py
Original file line number Diff line number Diff line change
@@ -51,7 +51,7 @@


@entity_task(log, writes=['inversion_input'])
def prepare_for_inversion(gdir, add_debug_var=False,
def prepare_for_inversion(gdir,
invert_with_rectangular=True,
invert_all_rectangular=False,
invert_with_trapezoid=True,
@@ -454,9 +454,12 @@ def mass_conservation_inversion(gdir, glen_a=None, fs=None, write=True,
volume[i] = out_thick[i] * w[i] * cl['dx']

# Sanity check
if np.any(out_thick <= 0):
log.warning("Found zero or negative thickness: "
"this should not happen.")
if np.any(out_thick <= -1e-2):
log.warning(f"{gdir.rgi_id} Found negative thickness: "
f"n={(out_thick <= -1e-2).sum()}, "
f"v={np.min(out_thick)}.")

out_thick = utils.clip_min(out_thick, 0)

if write:
cl['is_trapezoid'] = is_trap
@@ -615,8 +618,17 @@ def get_inversion_volume(gdir):


@entity_task(log, writes=['inversion_output'])
def compute_velocities(gdir, glen_a=None, fs=None, filesuffix='',
with_sliding=False):
def compute_velocities(*args, **kwargs):
"""Deprecated. Use compute_inversion_velocities instead."""
warnings.warn("`compute_velocities` has been renamed to "
"`compute_inversion_velocities`. Prefer to use the new"
"name from now on.")
return compute_inversion_velocities(*args, **kwargs)


@entity_task(log, writes=['inversion_output'])
def compute_inversion_velocities(gdir, glen_a=None, fs=None, filesuffix='',
with_sliding=None):
"""Surface velocities along the flowlines from inverted ice thickness.
Computed following the methods described in Cuffey and Paterson (2010)
@@ -631,19 +643,20 @@ def compute_velocities(gdir, glen_a=None, fs=None, filesuffix='',
The output is written in 'inversion_output.pkl' in m yr-1
You'll need to call prepare_for_inversion with the `add_debug_var=True`
kwarg for this to work!
Parameters
----------
gdir : Glacier directory
gdir : :py:class:`oggm.GlacierDirectory`
the glacier directory to process
with_sliding : bool
default is False, if set to True we will compute using
the sliding component.
glen_a : Glen A, defaults to PARAMS
fs : sliding, defaults to PARAMS
if set to True, we will compute velocities the sliding component.
the default is to check if sliding was used for the inversion
(fs != 0)
glen_a : float
Glen's deformation parameter A. Defaults to PARAMS['inversion_glen_a']
fs : float
sliding paramter, defaults to PARAMS['inversion_fs']
filesuffix : str
add a suffix to the output file
add a suffix to the output file (optional)
"""

# Defaults
@@ -655,6 +668,9 @@ def compute_velocities(gdir, glen_a=None, fs=None, filesuffix='',
rho = cfg.PARAMS['ice_density']
glen_n = cfg.PARAMS['glen_n']

if with_sliding is None:
with_sliding = fs != 0

# Getting the data for the main flowline
cls = gdir.read_pickle('inversion_output')

13 changes: 2 additions & 11 deletions oggm/core/massbalance.py
Original file line number Diff line number Diff line change
@@ -1428,20 +1428,11 @@ def mb_calibration_from_geodetic_mb(gdir, *,
The data table can be obtained with utils.get_geodetic_mb_dataframe().
It is equivalent to the original data from Hugonnet, but has some outlier
values filtered. See `this notebook <>`_ for more details.
values filtered. See `this notebook` for more details.
The problem of calibrating many unknown parameters on geodetic data is
currently unsolved. This is OGGM's current take, based on trial and
error and based on ideas from the litterature.
Here are the new defaults for the W5E5 dataset:
- for precipitation: if PARAMS['use_winter_prcp_fac'] is True (the default)
and the baseline dataset is W5E5 (the default), then we use a precipitation
factor obtained by deriving a relation between the winter precipitation
and at a glacier and the necessary correction to match WGMS observations
(see Schuster et al., 2023 for details)
- for temperature: use_temp_bias_from_file
error and based on ideas from the literature.
Parameters
----------
8 changes: 3 additions & 5 deletions oggm/sandbox/edu.py
Original file line number Diff line number Diff line change
@@ -44,9 +44,8 @@ def __init__(
bias : float, optional
set to the alternative value of the annual bias [mm we yr-1]
you want to use (the default is to use the calibrated value)
y0 : int, optional, default: tstar
the year at the center of the period of interest. The default
is to use tstar as center.
y0 : int
the year at the center of the period of interest. Has to be set!
halfsize : int, optional
the half-size of the time window (window size = 2 * halfsize + 1)
filename : str, optional
@@ -176,8 +175,7 @@ def run_constant_climate_with_bias(
is to take the last year available
bias : float
bias of the mb model. Default is to use the calibrated one, which
is often a better idea. For t* experiments it can be useful to set it
to zero
is zero usually anyways.
kwargs : dict
kwargs to pass to the FluxBasedModel instance
"""
14 changes: 7 additions & 7 deletions oggm/shop/gcm_climate.py
Original file line number Diff line number Diff line change
@@ -196,7 +196,7 @@ def process_monthly_isimip_data(gdir, output_filesuffix='',
member : str
ensemble member gcm that you want to process
ssp : str
ssp scenario to process (only 'ssp126' or 'ssp585' are available)
ssp scenario to process (only 'ssp126', 'ssp370' or 'ssp585' are available)
year_range : tuple of str
the year range for which the anomalies are computed
(passed to process_gcm_gdata). Default for ISIMIP3b `('1979', '2014')
@@ -218,6 +218,7 @@ def process_monthly_isimip_data(gdir, output_filesuffix='',
# Glacier location
glon = gdir.cenlon
glat = gdir.cenlat

if testing:
gcm_server = 'https://cluster.klima.uni-bremen.de/~oggm/test_climate/'
else:
@@ -229,14 +230,13 @@ def process_monthly_isimip_data(gdir, output_filesuffix='',
fpath_spec = path + '{}_w5e5_'.format(member) + '{ssp}_{var}' + add
fpath_temp = fpath_spec.format(var='tasAdjust', ssp=ssp)
fpath_temp_h = fpath_spec.format(var='tasAdjust', ssp='historical')

fpath_precip = fpath_spec.format(var='prAdjust', ssp=ssp)
fpath_precip_h = fpath_spec.format(var='prAdjust', ssp='historical')
fpath_temp = utils.file_downloader(fpath_temp)
fpath_temp_h = utils.file_downloader(fpath_temp_h)

fpath_precip = utils.file_downloader(fpath_precip)
fpath_precip_h = utils.file_downloader(fpath_precip_h)
with utils.get_lock():
fpath_temp = utils.file_downloader(fpath_temp)
fpath_temp_h = utils.file_downloader(fpath_temp_h)
fpath_precip = utils.file_downloader(fpath_precip)
fpath_precip_h = utils.file_downloader(fpath_precip_h)

# Read the GCM files
with xr.open_dataset(fpath_temp_h, use_cftime=True) as tempds_hist, \
Loading
Oops, something went wrong.

0 comments on commit ee6904a

Please sign in to comment.