Fast and accurate method to compute the following precipitation indices from climate data organised in yearly blocks of NetCDF files:
- temporal mean
- maximum yearly/seasonal daily or hourly precipitation
- wet day frequency
- percentiles
Daily and hourly temporal input frequencies can be processed and percentiles can be computed as all-day/hour or wet-day/hour according to Schär et al. (2016). Indices can be derived for entire years or seasons (JJA, SON, DJF, MAM).
To reduce memory consumption, yearly blocks are processed sequentially and the relevant data is kept in memory. To e.g. compute the 90% percentile, the largest 10% of data per grid cell is kept in memory and sequentially updated. After all yearly blocks have been processed, the percentile is computed via interpolation. This method is identical to applying the function numpy.percentile(). Part of the percentile calculation is parallelised with Numba – the performance of the script depends thus on the available number of CPUs (respectively threads).
The following Python packages are required to use this repository: NumPy, SciPy, netCDF4, Xarray, Dask, Numba. First, create a Conda environment with all required Python packages:
conda create -n prec_percentiles -c conda-forge numpy scipy matplotlib netcdf4 xarray dask numbaand activate this environment. In the script precipitation_indices.py, the path to the script auxiliary.py must then be adjusted. Afterwards, precipitation_indices.py can be used. The Python scripts in the folder comparison additionally require the packages Matplotlib, cartopy and cmcrameri for plotting.
Example data for the script can be downloaded from ESGF. Select the following options: Project → CORDEX, Domain → EUR-11, Experiment → evaluation, RCM Model → COSMO-crCLIM-v1-1, Time Frequency → 1hr, Variable → pr. Then download the first 10 files:
pr_EUR-11_ECMWF-ERAINT_evaluation_r1i1p1_CLMcom-ETH-COSMO-crCLIM-v1-1_v1_1hr_197901010030-197912312330.nc
..
pr_EUR-11_ECMWF-ERAINT_evaluation_r1i1p1_CLMcom-ETH-COSMO-crCLIM-v1-1_v1_1hr_198801010030-198812312330.nc
The script submit_daint_multithread.sh allows to submit precipitation_indices.py as a job on the CSCS cluster Daint.
The folder comparison contains Python scripts to compare obtained results with other methods:
-
precipitation_indices_CDO.sh: Compute the same precipitation indices for the example data with CDO. This script is approximately ~7 times slower than precipitation_indices.py (depending on the available number of CPUs/threads). With CDO, percentiles are estimated by sorting the data into bins (see CDO Manual) to reduce memory consumption. By default, 101 bins are applied. To obtain more accurate results, the bin size can be increase with e.g. export CDO_PCTL_NBINS=1001.
-
python_versus_CDO.py: Compare output from precipitation_indices.py with results from the above shell script applying CDO. The plot below shows the result for eight precipitation indices computed from a COSMO simulation with a 0.11° grid spacing and derived from hourly data for the period 1979 - 1988. The output from precipitation_indices.py is shown in the left and the absolute difference (CDO - Python) in the right column. The maximum absolute error occurring in the field is displayed in the upper right of the figure panels. The non-random deviations for temporal mean are likely caused by the slight difference in computing the average over the 10 years (in the Python script, the average is weighted by the number of days in a year). For the yearly maximum, the wet day frequency and the intensity, the deviations are very small (and random) compared to the actual values. For percentiles, the deviations are considerably larger and particularly occur in the southern part of the domain (at least for the percentiles in the range of 90% - 99%). These deviations originate from two two different computational methods: The Python scripts considers all individual values whereas in CDO, values are temporally aggregated into 1001 bins for computing percentiles.
- python_evaluate.py: Check that the fast script precipitation_indices.py yields identical results (apart from numerical errors) compared to a simpler (but less efficient) method in Python that computes indices from unpartitioned time series (e.g. by applying the function numpy.percentiles() to compute percentiles).
- Ban, N., Caillaud, C., Coppola, E. et al. (2021): The first multi-model ensemble of regional climate simulations at kilometer-scale resolution, part I: evaluation of precipitation. Clim Dyn 57, 275–302. https://doi.org/10.1007/s00382-021-05708-w
- Chinita, M. J., Richardson, .M, Teixeira, J. and Miranda, P. M. A. (2021): Global mean frequency increases of daily and sub-daily heavy precipitation in ERA 5. Environ Res Lett 16(7):74035. https://doi.org/10.1088/1748-9326/ac0caa
- Schär, C., Ban, N., Fischer, E. M. et al. (2016): Percentile indices for assessing changes in heavy precipitation events. Climatic Change 137, 201–216. https://doi.org/10.1007/s10584-016-1669-2