MAGIC ExtremeEvents namelist

[earnone]

Namelist to calculate ETCCDI climate extremes indices calling the climdex.pcic.ncdf R library based on temperature and precipitation model data.

old version: https://github.com/ESMValGroup/ESMValTool-private/tree/CRESCENDO_ExtremeEvents_revised

More recent work in issue #906

A temporary conversion to ESMValTool2 showed the new output from ESMVal preprocessor (former climo file) is not any more accepted by the climdex library. This was tested on data from several CMIP5 models. Calling the library by hand with a new ESMVal preprocessed file, e.g. _> sfile <- "pr_day_historical_EC-Earth_r8i1p11990-2000.nc" _> author.data <- list(institution="ISAC-CNR", institutionid="ISAC") _> create.indices.from.files(sfile, workdir, sfile[1], author.data,base.range=c(1990,2000))

gives the following error message: Error in as.POSIXlt.character(x, tz, format, ...) : character string is not in a standard unambiguous format In addition: Warning messages: 1: In min(args.flat, na.rm = na.rm) : no non-missing arguments to min; returning Inf 2: In max(args.flat, na.rm = na.rm) : no non-missing arguments to max; returning -Inf

It seems the problem depends on the bnds dimension added to the new ESMValTool2 preprocessed file.

[earnone]

@bjoernbroetz Could we revise here what happened to the former extremeEvents namelist/recipe? Was this eventually ported into v2.0? I recall we read output from this to run further analysis with our tool #324.

The porting was halted because of the climdex.pcic.ncdf library being unable to deal with files pre-processed by the v.2.0 preproc. I guess also the work performed by our Norwegian colleagues is dependent on this.

[bjoernbroetz]

Yes, I will check again into this and come back here with more information.

[bouweandela]

@bjoernbroetz @earnone What is the status of this?

[earnone]

The files preprocessed by v2 are not readable by the climdex.pcic.ncdf library at the core of this dlr namelist. So that halted the porting into v2. There was also an overlap with hyint, in the sense that hyint can read the results from this namelist and do further processing. If this is not going to be ported, we can run directly climdex adding a subroutine to hyint if we find a way to make climdex.pcic.ncdf reading the v2 data. @bjoernbroetz what are your plans?

[jhardenberg]

To be precise, the problem is the following: when the climdex R library attempts to process a file generated by the ESMValTool preprocessor it gets confused by the time_bnds, lon_bnds and lat_bnds and spits the following error message:

Error in nc_parse_group_structure(vars, verbose = verbose) : 
  Error, trying to add dim named time  to the root group, but there is already  a dim named  time in the root group, with different length.  Have you mistakenly resued  a dimension name for two different dims?  Have you correctly specifed  the groups of any dims with duplicated names?

[jhardenberg]

I just found a solution: to remove the bound information from the input files. If an input file is further preprocessed with nco, for example using the following script, the climdex library works as expected and does not complain:

#!/bin/bash
infile=$1
ncks -C -O -x -v lat_bnds,lon_bnds,time_bnds $infile temp$$.nc
ncatted -O -a bounds,time,d,, temp$$.nc
ncatted -O -a bounds,lat,d,, temp$$.nc
ncatted -O -a bounds,lon,d,, temp$$.nc
ncatted -O -a coordinates,pr,d,, temp$$.nc
mv temp$$.nc $infile

[maritsandstad]

Now I am very confused. Should I consult a different branch? Is someone else doing the climdex porting (that would be wonderful really:-))? Do you want it to start over from scratch? Have you guys been working on some things to port this and then I should work from there, or just stay away? Also should this be part of the MAGIC tools?

[earnone]

@maritsandstad yes, we were looking at making the extremeEvents working in v2 since the hyint recipe uses climdex output files as input for further post-processing. The only issue that halted us was the fact that the climdex.pcic.ncdf library was not able to digest input files pre-processed by esmvaltool v2, and this was sorted by @jhardenberg with the above nco commands (removing bounds). I am afraid we had not realised you were working on the porting too. @jhardenberg can update shortly, he has in fact just pushed a first version of the extremeEvents ported into v2 following the standard we developed for all our r recipes. There are still quite a few refinements to implement, we can maybe do that together to speed things up?

[earnone]

See https://github.com/ESMValGroup/ESMValTool/tree/version2_ExtremeEvents_v2

[maritsandstad]

Ok, I'd be super happy for anyone else to help with the work! The original code not only included generation of ETCCDI index calculation, but also generation of some plots (Glecker and timeseries). I suppose that is still a feature that should be ported? Will you port that too, or just the index generation? Also what branch are you working in? Making use of your work would at least very much speed things up for me as esmvaltool is very much a side project for me...

[earnone]

Good, I think it is much faster if we finalise the porting and interact with you for what needed to avoid duplicates. From the R coding point of view you will easily understand the various changes, it just takes some time to get acquainted with v2. The plotting part was included in the porting. See the branch here above.

[jhardenberg]

Hi Marit, I also left a comment in #906 , let's continue here. For now https://github.com/ESMValGroup/ESMValTool/tree/version2_ExtremeEvents_v2 computes only the indices but only because I did not yet have time to test the plots. I will try and give a look later.

To document what is needed to test the diagnostic: some additional packages are needed in the environment. We will add those to environment.yml later, but for a quick test we can proceed as follows, installing them by hand:

conda install -c conda-forge r-snow
conda install -c conda-forge r-udunits2
conda install -c conda-forge r-proj4

and then from R doing

> install.packages("devtools")
> library(devtools)
>install_github('pacificclimate/climdex.pcic.ncdf', ref='release')

[maritsandstad]

Ok, thank you both, that sounds great! I will have a look and start testing soon.

[bouweandela]

We will add those to environment.yml later

Please add them to the file esmvaltool/install/R/r_requirements.txt instead, as this is the correct location for R dependencies.

[jhardenberg]

Please add them to the file esmvaltool/install/R/r_requirements.txt instead, as this is the correct location for R dependencies.

Not in this case: I actually tried that first, but at least when I tested in my environment, installing these directly from R would fail, while adding them from conda works .... but I will try again.

We actually have another important issue: the climdex.pcic.ncdf R library does not exist on CRAN, using install_github is the only official method of installation (see https://github.com/pacificclimate/climdex.pcic.ncdf ): install_github('pacificclimate/climdex.pcic.ncdf', ref='release')

So how do we proceed, should we modify setup.R so that it allows installing also packages from github directly (actually only this one ;) ... ) ?

[bouweandela]

I think we should contact the developers of this library and ask them nicely if they would be willing to release their software on CRAN. If necessary we can help them do that. Does anyone know them personally? @jhardenberg @earnone @maritsandstad @ESMValGroup/esmvaltool-coreteam

[jhardenberg]

Ok, the current version of branch version2_ExtremeEvents_v2 now also produces the timeseries and gleckler plots :)

@maritsandstad : there was probably a missing configuration file, since some parameters were undefined, I made a new internal config file "cfg_extreme.r" containing these with some educated guesses (for example parameters controlling the alignment of graphic elements such as gl_xscale_spacer_RMSE ). Do you have better values to suggest for these? (indeed the resulting Gleckler plot has some elements misaligned, I guess better values of these parameters would fix them). Some of the parameters in this file will later need to be moved to the recipe (I left them here for now for testing purpouses), such as the index selection lists timeseries_idx and gleckler_idx .... One which is mysterious for me is CMIP_name in make_timeseries_plot.r ... I understand this to be the name of the model ensemble (you plot the ensemble mean if I understood correctly). So what should this be?

Another question: the original code was searching among the datasets, separating those with project "OBS" from the others (these are treated separately in the plots). This does not necessarily work in ESMValtool2 since there could be different obs projects such as "OBS" and "obs4mips" .... Actually the current ESMValToolv2 standard is that in the recipe a "reference_model" can be declared and I am using that. For now I used that assuming it would be only one, but looking at the plotting routines I see that you intend these to be actually multiple datasets. Should we make reference_model a list and parse that?

Finally, for the record: another change I had to do is that in function createGrid I had to add a selectgrid,2 command to the cdo command string, since the netcdf files it operates on now actually contain two grids and the it the second one which we want.

[jhardenberg]

PS, additional comment on cdo: at a certain point we should convert the code to use the cdo wrapper using system2 which we devoloped and use in the quantilebias diagnostic.

[jhardenberg]

I think we should contact the developers of this library and ask them nicely if they would be willing to release their software on CRAN. If necessary we can help them do that. Does anyone know them personally? @jhardenberg @earnone @maritsandstad @ESMValGroup/esmvaltool-coreteam

I somehow doubt that this will be successful. This library is nowadays the official reference to compute the ETCCDI indices (there are some subtle details to be taken into account in their apparently simple definitions), so that everybody is using them. I am sure that this request has already been made....

Alternative: on their official webpage http://pacificclimate.github.io/climdex.pcic.ncdf/ they also offer a tarball to download. This is quite compact (64k), and maybe we could distribute it directly with ESMValtool (provided licensing issues allow that) ? In any case: @bouweandela why are you so firmly against allowing using install_github in this case?

[maritsandstad]

Know him personally would be an overstatement, but I have his contact information and have been in touch with him from time to time. My boss knows him though. If you would like, I can definitely send a request to him about releasing the code on CRAN. However, I also think it is not highly likely that he will concede to making CRAN-release.

All this developement sounds great!

Yes, the Comparison is meant to be done against several observational data (basically the whole diagnostic is supposed to reproduce the figures in figure 9.37 of the AR5 WG1 IPCC report). The Glecker plot should be modular so it can produce differently divided squares for different number of observations to be compared with, the timeseries plot should also accommodate differing numbers of observations to be plotted on top of the model ensemble spread. A list of reference models would hence be good. I'll look into the CMIP_name-issue shortly.

[bouweandela]

I have his contact information and have been in touch with him from time to time. My boss knows him though. If you would like, I can definitely send a request to him about releasing the code on CRAN. However, I also think it is not highly likely that he will concede to making CRAN-release.

@maritsandstad Yes, if you could ask him would be great. You can tell him that I would be willing to help out, if the reason that he does not want to do this is time constraints. Or do you think it is another reason?

In any case: @bouweandela why are you so firmly against allowing using install_github in this case?

@jhardenberg Because ESMValTool should be reliable and installing from anywhere else then the official repository for a language is not reliable.

If there will be no release on CRAN the alternative will be to copy the source code into the ESMValTool repository, as we've done for the CVPD package for example, but this is not preferred, because it will require manually updating that code occasionally.

[maritsandstad]

@bouweandela : I'll send him an email now. Time constraints may very well be the reason, but there may be other reasons too.

@jhardenberg : The CMIP_name is the legend to be used for the ensemble spread in the time series plot. So it would be whatever the user would like to dub his or her model ensemble. Maybe CMIP_name isn't the greatest denotation, MIP_name might be more apt, however, a description of this might be given in the diagnostic documentation, no?

[jhardenberg]

Ok, we can include MIP_name in the recipe then. What about the other parameters in cfg_extreme.r? Do you have better default values? (you will find the file in the version2_ExtremeEvents_v2 branch)

[maritsandstad]

@jhardenberg : I will have a look at the defaults this afternoon

[maritsandstad]

@jhardenberg , like the idiot I am, can I just ask how you have been running (which config file are you using)?

[maritsandstad]

@bouweandela : I'm in contact with the climdex.pcic person now. Apparantly they had tried to put the climdex.pcic.ncdf package on CRAN previously, but got rejected from CRAN because the package was not of "general use". Could you possibly offer some help on that count?

[jhardenberg]

@jhardenberg , like the idiot I am, can I just ask how you have been running (which config file are you using)?

Hi, I am not sure that I understood your question: At the moment I am using my default version2_development conda environment, modified with some additional packages as described above https://github.com/ESMValGroup/ESMValTool/issues/270#issuecomment-488625895 in order to allow installing devtools and the climdex library. Is this what you meant? If you mean instead the missing parameters, these are now listed in esmvaltool/diag_scripts/extreme_events/cfg_extreme.r .....

[jhardenberg]

I updated the branch: now the diagnostic does not use direct system calls anymore and uses the wrapper cdo() function as we do in the quantilebias diagnostic (@maritsandstad there was a long discussion in the PR for that metric on why this would be better ....)

[jhardenberg]

Ok, I fixed the default offsets and scales incfg_extfremes.r a bit. Now both the timeseries and Gleckler plot look good ;) (tested only in png) Multiple obs datasets are now supported providing a list of "reference_dataset". @maritsandstad if you could check if all is ok, we could then probably move to the next step and make this into a PR ?

[jhardenberg]

Hello @maritsandstad , we noticed one thing which may need to be corrected, but would like your opinion first. As it is now, the diagnostic assumes that the climdex base.range (the reference period in which reference percentiles etc are computed) goes from the maximum of the model start years to the minimum of the model end years (i.e. the largest common period). So essentially the metric is computed over the same period which is used as reference period. We (I and Enrico) believe this to be wrong. In fact also fig 9.37 of the IPCC report was computed over different period (1981-2000) using as a base period 1961-1990. So essentially base.period should be an option in the configuration file ....

PS: I added provenance, testing now before pushing the change.

[maritsandstad]

Yes, you are correct that the base period probably shouldn't be the whole period. However, setting a default can be somewhat difficult... In case of the two plots the diagnostic now produces. The most important thing is to use the same base periods and periods of interest (those two would normally be different) for each of the models and reanalysis being compared. Ideally a base period should be of 30 years or so. Good values would be the AR5 base period of 1961-1990, or an extended version of the new present day AR6 period of 1995-2014 to 1985-2014. However, users may use the tool with inputs that do not contain any of these year ranges (at least not for all the models). They may even input data that contain fewer than 30 years of data.

The question is how to handle these issues. It could also very well be that they want to study a period of interest that is much shorter than the entire timespan of the data. One option could be to choose some defaults that we find reasonable, and that make it clear that for different uses these options must be given in the user-config file?

Also, my dumb question was how you run. I load the environment in the branch and used the esmvaltool command with config-user-extremetest.yml configfile and the recipe_extreme_events.yml recipe. I guess you have a different config-file, maybe?

[jhardenberg]

The question is how to handle these issues. It could also very well be that they want to study a period of interest that is much shorter than the entire timespan of the data. One option could be to choose some defaults that we find reasonable, and that make it clear that for different uses these options must be given in the user-config file?

You are right that making this selection compulsory in the recipe is really needed. Actually the user will need to specify two things: the base period and the period of interest. The model data provided will need to cover both in full. Ok, I will try and add this quickly.

Also, my dumb question was how you run. I load the environment in the branch and used the esmvaltool command with config-user-extremetest.yml configfile and the recipe_extreme_events.yml recipe. I guess you have a different config-file, maybe?

Well, I have my own config-user.yml pointing to the original model data on our machine, but that should not be a problem. The environment I use is the current default esmvaltool2 environment, I have not tested yet with the environment in the branch ( i.e. I did not do "conda env update" ) + I installed by hand the additional conda packages needed as described in my post above. I was planning to fix the environment of the branch later today. Des your esmvaltool not run because of the environment or because of the config file ?

[maritsandstad]

You are right that making this selection compulsory in the recipe is really needed. Actually the user will need to specify two things: the base period and the period of interest. The model data provided will need to cover both in full. Ok, I will try and add this quickly.

Yes, that would be great if you could add that!

Well, I have my own config-user.yml pointing to the original model data on our machine, but that should not be a problem. The environment I use is the current default esmvaltool2 environment, I have not tested yet with the environment in the branch ( i.e. I did not do "conda env update" ) + I installed by hand the additional conda packages needed as described in my post above. I was planning to fix the environment of the branch later today. Des your esmvaltool not run because of the environment or because of the config file ?

config file, I think, but I'll try again now.

[jhardenberg]

Hi @maritsandstad : currently the plots are created regridding all datasets over a common grid. Currently simply the first file found in the working directory is used (in function createGrid). This is a bit random. It would be better to ask the user to specify the dataset to use as reference and as default maybe use the first reference dataset, no?

[maritsandstad]

@jhardenberg : Yes, I agree that would be a better choice.

[mattiarighi]

It would be better to ask the user to specify the dataset to use as reference and as default maybe use the first reference dataset, no?

I would also suggest to make reference_dataset a mandatory setting in this case. Random (i.e. not reproducible) results should be avoided if possible.

[maritsandstad]

Also I updated the environment to test and it seems the environment setup does not work well with the reference_dataset list: ValueError: Error validating data /net/pdo/div/pdo/extreme/masan/V2/ESMValTool/esmvaltool/recipes/recipe_extreme_events.yml with schema /uio/kant/div-cicero-u1/masan/.conda/envs/esmvaltool/lib/python3.6/site-packages/ESMValTool-2.0a1-py3.6.egg/esmvaltool/recipe_schema.yml diagnostics.ExtremeEvents.variables.pr.reference_dataset: '['ACCESS1-0', 'ACCESS1-3', 'IPSL-CM5A-MR']' is not a str.

[bouweandela]

Yes, reference_dataset when defined in the variable can only be a single dataset. Maybe you could use reference_datasets? Note that in that case we do not check if the data is available, though we could in the future. You could also define the reference_datasets in the diagnostic script section, if they are the same for all variables anyway.

[jhardenberg]

@maritsandstad did you do a pull to the latest version of the branch before testing? That had already been fixed.

Actually I took the liberty to make one small change in esmvaltool/recipe_schema.yml writing reference_dataset: any(str(), list(str()), required=False) instead of reference_dataset: str(required=False)

@bouweandela please notice that this change maintains existing functionality (it does not break other recipes/diagnostics) and allows passing multiple reference datasets . ....

[jhardenberg]

Fixed the envirionment, now after doing conda update and Rscript setup.R you only still need to do install_github('pacificclimate/climdex.pcic.ncdf', ref='release') by hand from R for the time being.

[bouweandela]

Actually I took the liberty to make one small change in esmvaltool/recipe_schema.yml writing reference_dataset: any(str(), list(str()), required=False) instead of reference_dataset: str(required=False)

Unfortunately I don't think that is a very good idea: the reference_dataset value can e.g. be used as a target grid for the regridding function. However, if you specify multiple datasets with differing grids, what dataset should we regrid to? By changing the schema, the user will no longer get an easy to understand error message, but probably the program will crash later on if you try to use the reference_dataset as a regridding target. It would be better to follow my previous suggestion https://github.com/ESMValGroup/ESMValTool/issues/270#issuecomment-489589202.

[jhardenberg]

I do not agree, since I would find it very confusing for users to have both a reference_datasetand a reference_datasets tag. As for the regridding target this could be simply the first reference dataset given (this is actually the solution which I will implement). In any case no big deal, I will follow your suggestion, also because it is indeed better to have one single tag in the diagnostic script section rather than have to repeat it for each variable ....

[jhardenberg]

• now usingreference_datasets in the diagnostic section to specify the 'obs' files
• added the two tags base_range and analysis_range in the recipe to specify the base.range for climdex and the range to be used for plotting the timeseries and computing the Gleckler plot. The model data will need to include both periods (this is checked)

[jhardenberg]

Now a tag regrid_dataset can be used in the recipe to specify which dataset to use for regridding all datasets. If not specified the first reference_dataset is used.

[jhardenberg]

To start testing the branch I started a PR (please ignore #1051 (closed) , I was a bit clumsy in creating the PR)

[bouweandela]

To be precise, the problem is the following: when the climdex R library attempts to process a file generated by the ESMValTool preprocessor it gets confused by the time_bnds, lon_bnds and lat_bnds and spits the following error message: ...

I just found a solution: to remove the bound information from the input files. If an input file is further preprocessed with nco, for example using the following script, the climdex library works as expected and does not complain: ...

The output of the ESMValTool is written by iris, so it can be expected to comply with the cf-conventions. Instead of using external tools to modify the data so the climdex R library understands it, I think it would be better to file an issue (and possibly a pull request) with the climdex library to make it understand the bounds. After all, such a library can be expected to be able to read data that complies with the cf-conventions.

[maritsandstad]

The output of the ESMValTool is written by iris, so it can be expected to comply with the cf-conventions. Instead of using external tools to modify the data so the climdex R library understands it, I think it would be better to file an issue (and possibly a pull request) with the climdex library to make it understand the bounds. After all, such a library can be expected to be able to read data that complies with the cf-conventions.

You can always try, but I would bet any substantial sum that they will make the necessary changes. However, might be worth a try...

[bouweandela]

I have created an issue to ask the authors of the climdex.pcic.ncdf package about publication on CRAN: https://github.com/pacificclimate/climdex.pcic.ncdf/issues/16.

[bouweandela]

@jhardenberg Could you please open an issue about this https://github.com/ESMValGroup/ESMValTool/issues/270#issuecomment-487994554 in the relevant repository?

[bouweandela]

The recipe/diagnostic script has now been merged. I noticed however, that provenance is not correctly recorded for a number of files. After running the default recipe, the following files do not have a reference to any preprocessor output files, so very probably have a wrong provenance record (i.e. missing/wrong ancestors):

./plots/extreme_events/main/Gleckler_CMIP_24-models_11-idx_4-obs_1981-2000_provenance.xml
./work/extreme_events/main/timeseries/r95pETCCDI_yr_25enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/sdiiETCCDI_yr_75enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/tnxETCCDI_yr_ensmean_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/rx5dayETCCDI_yr_ensmean_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/txxETCCDI_yr_25enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/rx5dayETCCDI_yr_25enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/rx5dayETCCDI_yr_75enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/sdiiETCCDI_yr_ensmean_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/cddETCCDI_yr_25enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/fdETCCDI_yr_25enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/txxETCCDI_yr_75enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/cddETCCDI_yr_75enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/tnnETCCDI_yr_25enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/txnETCCDI_yr_25enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/sdiiETCCDI_yr_25enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/txnETCCDI_yr_75enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/rx1dayETCCDI_yr_25enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/r95pETCCDI_yr_75enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/tnnETCCDI_yr_75enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/rx1dayETCCDI_yr_ensmean_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/cddETCCDI_yr_ensmean_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/txnETCCDI_yr_ensmean_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/fdETCCDI_yr_ensmean_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/trETCCDI_yr_ensmean_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/tnnETCCDI_yr_ensmean_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/fdETCCDI_yr_75enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/trETCCDI_yr_75enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/trETCCDI_yr_25enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/r95pETCCDI_yr_ensmean_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/rx1dayETCCDI_yr_75enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/tnxETCCDI_yr_25enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/tnxETCCDI_yr_75enspctl_for_timeseries_provenance.xml
./work/extreme_events/main/timeseries/txxETCCDI_yr_ensmean_for_timeseries_provenance.xml