3. Diurnal Cycle
This notebook aims at inroducing new users on how to use the PCDMI diurnal cycle drivers.
This diagram shows how various drivers are linked together.
To run the following demo, it is expected that you have downloaded the sample data as demonstrated in the download notebook. For this tutorial we will be using 3 years worth of 3-hourly data resampled to a 5x5 grid.
The following cell loads the demo input and output directories set in the download notebook.
[1]:
from user_choices import demo_data_directory, demo_output_directory
Daily Means
Like all other drivers in the PCMDI Metrics Package, Diurnal Cycle relies on parameter input files.
Our first driver starts from CMIP5 data and computes the daily means.
[2]:
with open("basic_diurnal_compute_daily_mean.py") as f:
print(f.read())
Now to run this simply call the driver:
computeStdOfDailyMeans.py -p basic_diurnal_compute_daily_mean.py
[3]:
%%bash
computeStdOfDailyMeans.py -p basic_diurnal_compute_daily_mean.py
SKIPPING: []
FILES: ['demo_data/misc_demo_data/atm/3hr/pr/pr_3hr_IPSL-CM5A-LR_historical_r1i1p1_5x5_1997-1999.nc']
PARAMS: [<pcmdi_metrics.diurnal.common.INPUT object at 0x7ff368087a60>]
Data source: IPSL-CM5A-LR
Opening demo_data/misc_demo_data/atm/3hr/pr/pr_3hr_IPSL-CM5A-LR_historical_r1i1p1_5x5_1997-1999.nc ...
Data source: IPSL-CM5A-LR
Opening demo_data/misc_demo_data/atm/3hr/pr/pr_3hr_IPSL-CM5A-LR_historical_r1i1p1_5x5_1997-1999.nc ...
Year 1997:
Reading pr from demo_data/misc_demo_data/atm/3hr/pr/pr_3hr_IPSL-CM5A-LR_historical_r1i1p1_5x5_1997-1999.nc for time interval 1997-7-1 0:0:0.0 to 1997-7-31 23:59:0.0 ...
Year 1998:
Reading pr from demo_data/misc_demo_data/atm/3hr/pr/pr_3hr_IPSL-CM5A-LR_historical_r1i1p1_5x5_1997-1999.nc for time interval 1998-7-1 0:0:0.0 to 1998-7-31 23:59:0.0 ...
Year 1999:
Reading pr from demo_data/misc_demo_data/atm/3hr/pr/pr_3hr_IPSL-CM5A-LR_historical_r1i1p1_5x5_1997-1999.nc for time interval 1999-7-1 0:0:0.0 to 1999-7-31 23:59:0.0 ...
This generates a netcdf file pr_IPSL-CM5A-LR_Jul_1997-1999_std_of_dailymeans.nc which contains the daily standard deviation at each cell.
Looking at our diagram, the next driver to run is the one computing the mean of the standard deviation from daily means over a region of interest. First we open the parameter file.
[4]:
with open("basic_diurnal_std_daily_mean.py") as f:
print(f.read())
# output directory
results_dir = 'demo_output/diurnal/json'
# input directory which is actually the output of previous driver
modpath = 'demo_output/diurnal/nc'
# filenames template
filename_template = 'pr_%(model)_Jul_%(firstyear)-%(lastyear)_std_of_dailymeans.nc'
# model to use
model = 'IPSL-CM5A-LR'
experiment = 'historical'
realization = 'r1i1p1'
# Month to use
month = 7
# Period
firstyear = 1997 # included
lastyear = 1999 # included
# Latitudes/longitudes to use
lat1 = -50.
lat2 = 50.
lon1 = 0.
lon2 = 360.
# Name
region_name = "TRMM"
# Number of workers
num_workers = 4
Now to run this simply call the driver
std_of_dailymeans.py -p basic_diurnal_std_daily_mean.py
[5]:
%%bash
std_of_dailymeans.py -p basic_diurnal_std_daily_mean.py
TEMPLATE NAME: pr_IPSL-CM5A-LR_Jul_1997-1999_std_of_dailymeans.nc
Preparing to write output to JSON file ...
Initializing dictionary of statistical results ...
['demo_output/diurnal/nc/pr_IPSL-CM5A-LR_Jul_1997-1999_std_of_dailymeans.nc']
PARAMS: [<pcmdi_metrics.diurnal.common.INPUT object at 0x7fa1f8669dc0>]
Specifying latitude / longitude domain of interest ...
Reading demo_output/diurnal/nc/pr_IPSL-CM5A-LR_Jul_1997-1999_std_of_dailymeans.nc ...
Shape = (20, 72)
Finding RMS area-average ...
For IPSL-CM5A-LR in Jul, average variance of daily values = ( 3.92 mm/d)^2
Writing output to JSON file ... {'IPSL-CM5A-LR': {'TRMM': 3.9163177715285187}}
REG MASK: {}
done
INFO::2021-11-15 14:56::pcmdi_metrics:: Results saved to a json file: /Users/ordonez4/Documents/git/pcmdi_metrics/doc/jupyter/Demo/demo_output/diurnal/json/pr_Jul_1997_1999_std_of_dailymeans.json
This generates a json file: pr_Jul_1997_1999_std_of_dailymeans.json
OPTIONAL: You could also append a new region to this json file by overwritting some of our parameters from the command line. This example creates results for the tropics instead of the whole globe.
[6]:
%%bash
std_of_dailymeans.py -p basic_diurnal_std_daily_mean.py \
--region_name=TROPICS \
--lat1=-30. \
--lat2=30. \
--lon1=0. \
--lon2=360 \
--append
TEMPLATE NAME: pr_IPSL-CM5A-LR_Jul_1997-1999_std_of_dailymeans.nc
Preparing to write output to JSON file ...
['demo_output/diurnal/nc/pr_IPSL-CM5A-LR_Jul_1997-1999_std_of_dailymeans.nc']
PARAMS: [<pcmdi_metrics.diurnal.common.INPUT object at 0x7fe1386a3970>]
Specifying latitude / longitude domain of interest ...
Reading demo_output/diurnal/nc/pr_IPSL-CM5A-LR_Jul_1997-1999_std_of_dailymeans.nc ...
Shape = (12, 72)
Finding RMS area-average ...
For IPSL-CM5A-LR in Jul, average variance of daily values = ( 3.89 mm/d)^2
Writing output to JSON file ... {'IPSL-CM5A-LR': {'TRMM': 3.9163177715285187, 'TROPICS': 3.890735445846671}}
REG MASK: {'TRMM': {'id': 'TRMM', 'domain': {'TRMM': 'cdutil.region.domain(latitude=(-50.0, 50.0),longitude=(0.0, 360.0))'}}}
done
INFO::2021-11-15 14:56::pcmdi_metrics:: Results saved to a json file: /Users/ordonez4/Documents/git/pcmdi_metrics/doc/jupyter/Demo/demo_output/diurnal/json/pr_Jul_1997_1999_std_of_dailymeans.json
Diurnal Composite
Looking again at our diagram we can now start again from the original 3 hourly data, and run the composite script
[7]:
with open("basic_diurnal_composite.py") as f:
print(f.read())
# output directory
results_dir = 'demo_output/diurnal/nc'
# input directory
modpath = 'demo_data/misc_demo_data/atm/3hr/pr/'
# filenames template
filename_template = 'pr_3hr_%(model)_%(experiment)_%(realization)_5x5_1997-1999.nc'
# model to use
model = 'IPSL-CM5A-LR'
experiment = 'historical'
realization = 'r1i1p1'
# Month to use
month = 7
# Period
firstyear = 1997 # included
lastyear = 1999 # included
# Number of workers
num_workers = 4
[8]:
%%bash
compositeDiurnalStatistics.py -p basic_diurnal_composite.py
done
TEMPLATE: pr_3hr_IPSL-CM5A-LR_historical_r1i1p1_5x5_1997-1999.nc
FILES: ['demo_data/misc_demo_data/atm/3hr/pr/pr_3hr_IPSL-CM5A-LR_historical_r1i1p1_5x5_1997-1999.nc']
PARAMS: [<pcmdi_metrics.diurnal.common.INPUT object at 0x7faca8978d60>]
Data source: IPSL-CM5A-LR
Opening demo_data/misc_demo_data/atm/3hr/pr/pr_3hr_IPSL-CM5A-LR_historical_r1i1p1_5x5_1997-1999.nc ...
Year 1997:
Reading pr from demo_data/misc_demo_data/atm/3hr/pr/pr_3hr_IPSL-CM5A-LR_historical_r1i1p1_5x5_1997-1999.nc for time interval 1997-7-1 0:0:0.0 to 1997-7-31 23:59:0.0 ...
Shape: (248, 36, 72)
DATA FROM: 1997-7-1 1:30:0.0 to 1997-7-31 22:30:0.0
Shape = (8, 36, 72)
8 timepoints per day, 3 hr intervals between timepoints
Creating temporary storage and output fields ...
Computing Local Standard Times for GMT 1.50 ...
Computing Local Standard Times for GMT 4.50 ...
Computing Local Standard Times for GMT 7.50 ...
Computing Local Standard Times for GMT 10.50 ...
Computing Local Standard Times for GMT 13.50 ...
Computing Local Standard Times for GMT 16.50 ...
Computing Local Standard Times for GMT 19.50 ...
Computing Local Standard Times for GMT 22.50 ...
Choosing timepoints with GMT 1.50 ...
days per mo : 31
Choosing timepoints with GMT 4.50 ...
days per mo : 31
Choosing timepoints with GMT 7.50 ...
days per mo : 31
Choosing timepoints with GMT 10.50 ...
days per mo : 31
Choosing timepoints with GMT 13.50 ...
days per mo : 31
Choosing timepoints with GMT 16.50 ...
days per mo : 31
Choosing timepoints with GMT 19.50 ...
days per mo : 31
Choosing timepoints with GMT 22.50 ...
days per mo : 31
Year 1998:
Reading pr from demo_data/misc_demo_data/atm/3hr/pr/pr_3hr_IPSL-CM5A-LR_historical_r1i1p1_5x5_1997-1999.nc for time interval 1998-7-1 0:0:0.0 to 1998-7-31 23:59:0.0 ...
Shape: (248, 36, 72)
Choosing timepoints with GMT 1.50 ...
days per mo : 31
Choosing timepoints with GMT 4.50 ...
days per mo : 31
Choosing timepoints with GMT 7.50 ...
days per mo : 31
Choosing timepoints with GMT 10.50 ...
days per mo : 31
Choosing timepoints with GMT 13.50 ...
days per mo : 31
Choosing timepoints with GMT 16.50 ...
days per mo : 31
Choosing timepoints with GMT 19.50 ...
days per mo : 31
Choosing timepoints with GMT 22.50 ...
days per mo : 31
Year 1999:
Reading pr from demo_data/misc_demo_data/atm/3hr/pr/pr_3hr_IPSL-CM5A-LR_historical_r1i1p1_5x5_1997-1999.nc for time interval 1999-7-1 0:0:0.0 to 1999-7-31 23:59:0.0 ...
Shape: (248, 36, 72)
Choosing timepoints with GMT 1.50 ...
days per mo : 31
Choosing timepoints with GMT 4.50 ...
days per mo : 31
Choosing timepoints with GMT 7.50 ...
days per mo : 31
Choosing timepoints with GMT 10.50 ...
days per mo : 31
Choosing timepoints with GMT 13.50 ...
days per mo : 31
Choosing timepoints with GMT 16.50 ...
days per mo : 31
Choosing timepoints with GMT 19.50 ...
days per mo : 31
Choosing timepoints with GMT 22.50 ...
days per mo : 31
Computing mean and standard deviation over all GMT 1.50 timepoints ...
Computing mean and standard deviation over all GMT 4.50 timepoints ...
Computing mean and standard deviation over all GMT 7.50 timepoints ...
Computing mean and standard deviation over all GMT 10.50 timepoints ...
Computing mean and standard deviation over all GMT 13.50 timepoints ...
Computing mean and standard deviation over all GMT 16.50 timepoints ...
Computing mean and standard deviation over all GMT 19.50 timepoints ...
Computing mean and standard deviation over all GMT 22.50 timepoints ...
This produces 3 output files:
pr_IPSL-CM5A-LR_Jul_1997-1999_diurnal_avg.nc
pr_IPSL-CM5A-LR_Jul_1997-1999_diurnal_std.nc
pr_IPSL-CM5A-LR_LocalSolarTimes.nc
These contain respectively the 1997-1999 mean diurnal cycle for July, the standard deviation of these results across individual years, and the local solar time. Results for each of these are available for the entire domain.
We can now generate ASCII files for composite diurnal cycle (w/ error bars) at selected grid points using the fourierDiurnalGridpoints.py script.
[9]:
%%bash
fourierDiurnalGridpoints.py -p basic_diurnal_fourier.py
LSTFILES: ['demo_output/diurnal/nc/pr_IPSL-CM5A-LR_LocalSolarTimes.nc']
TMPL pr_IPSL-CM5A-LR_LocalSolarTimes.nc
Results sent to: /Users/ordonez4/Documents/git/pcmdi_metrics/doc/jupyter/Demo/demo_output/diurnal/ascii/pr_Jul_1997-1999_fourierDiurnalGridPoints.asc
This produces an ascii file: pr_Jul_1997-1999_fourierDiurnalGridPoints.asc
Standard Deviation of Hourly Values
Starting again from the composite results our diagram suggests we now compute the standard deviation of hourly values.
[10]:
with open("basic_diurnal_std_hourly_mean.py") as f:
print(f.read())
# output directory
results_dir = 'demo_output/diurnal/json'
# input directory which is actually the output of previous driver
modpath = 'demo_output/diurnal/nc'
# model to use
model = 'IPSL-CM5A-LR'
experiment = 'historical'
realization = 'r1i1p1'
# Month to use
month = 7
# Period
firstyear = 1997 # included
lastyear = 1999 # included
# Latitudes/longitudes to use
lat1 = -50.
lat2 = 50.
lon1 = 0.
lon2 = 360.
# Number of workers
num_workers = 4
[11]:
%%bash
std_of_hourlyvalues.py -p basic_diurnal_std_hourly_mean.py
TEMPLATE NAME: pr_IPSL-CM5A-LR_Jul_1997-1999_diurnal_std.nc
Specifying latitude / longitude domain of interest ...
Preparing to write output to JSON file ...
Initializing dictionary of statistical results ...
['demo_output/diurnal/nc/pr_IPSL-CM5A-LR_Jul_1997-1999_diurnal_std.nc']
PARAMS: [<pcmdi_metrics.diurnal.common.INPUT object at 0x7f89187349d0>]
Specifying latitude / longitude domain of interest ...
Reading demo_output/diurnal/nc/pr_IPSL-CM5A-LR_Jul_1997-1999_diurnal_std.nc ...
Shape = (8, 20, 72)
Finding RMS area-average ...
For IPSL-CM5A-LR in Jul, average variance of hourly values = ( 4.85 mm/d)^2
Writing output to JSON file ...
done
INFO::2021-11-15 14:56::pcmdi_metrics:: Results saved to a json file: /Users/ordonez4/Documents/git/pcmdi_metrics/doc/jupyter/Demo/demo_output/diurnal/json/pr_Jul_1997-1999_std_of_hourlymeans.json
This generated the following file: pr_Jul_1997-1999_std_of_hourlymeans.json
These results are used in Trenberth et al. (2017). They are a measure of the intermittancy of hourly values, which puts “error bars” on the mean diurnal cycle.
Going back to the results of the composite we can now run std_of_meandiurnalcycle.py which can use the same input parameter file as the daily mean computation.
[12]:
%%bash
std_of_meandiurnalcycle.py -p basic_diurnal_std_hourly_mean.py
TEMPLATE NAME: pr_IPSL-CM5A-LR_Jul_1997-1999_diurnal_avg.nc
Specifying latitude / longitude domain of interest ...
Preparing to write output to JSON file ...
Initializing dictionary of statistical results ...
['demo_output/diurnal/nc/pr_IPSL-CM5A-LR_Jul_1997-1999_diurnal_avg.nc']
PARAMS: [<pcmdi_metrics.diurnal.common.INPUT object at 0x7fd9e8929eb0>]
Specifying latitude / longitude domain of interest ...
Reading demo_output/diurnal/nc/pr_IPSL-CM5A-LR_Jul_1997-1999_diurnal_avg.nc ...
Shape = (8, 20, 72)
Finding standard deviation over first dimension (time of day) ...
Shape = (20, 72)
Finding r.m.s. average over 2nd-3rd dimensions (area) ...
For IPSL-CM5A-LR in Jul, average variance of hourly values = ( 2.15 mm/d)^2
Writing output to JSON file ...
KEYS AT END: ['DISCLAIMER', 'REFERENCE', 'RESULTS']
REG MASK: {}
done
INFO::2021-11-15 14:56::pcmdi_metrics:: Results saved to a json file: /Users/ordonez4/Documents/git/pcmdi_metrics/doc/jupyter/Demo/demo_output/diurnal/json/pr_Jul_1997-1999_std_of_meandiurnalcyc.json
This generates the following file: pr_Jul_1997-1999_std_of_meandiurnalcyc.json
Fourier Analysis
Again starting from the composite results let’s do the fourier analysis. This uses a new parameter file.
[13]:
with open("basic_diurnal_fourierAllGrid.py") as f:
print(f.read())
# output directory
results_dir = 'demo_output/diurnal/nc'
# input directory which is actually the output of previous driver
modpath = 'demo_output/diurnal/nc'
# model to use
model = 'IPSL-CM5A-LR'
experiment = 'historical'
realization = 'r1i1p1'
# Month to use
month = 7
# Period
firstyear = 1997 # included
lastyear = 1999 # included
# Number of workers
num_workers = 4
[14]:
%%bash
fourierDiurnalAllGrid.py -p basic_diurnal_fourierAllGrid.py
modpath demo_output/diurnal/nc
filename_template pr_%(model)_%(month)_%(firstyear)-%(lastyear)_diurnal_avg.nc
filename_template_LST pr_%(model)_LocalSolarTimes.nc
LSTFILES: ['demo_output/diurnal/nc/pr_IPSL-CM5A-LR_LocalSolarTimes.nc']
TMPL pr_IPSL-CM5A-LR_LocalSolarTimes.nc
Reading demo_output/diurnal/nc/pr_IPSL-CM5A-LR_LocalSolarTimes.nc ... pr_IPSL-CM5A-LR_LocalSolarTimes.nc
====================
IPSL-CM5A-LR
====================
Reading time series of mean diurnal cycle ...
Input shapes: (8, 36, 72) (8, 36, 72)
Getting latitude and longitude coordinates.
Taking fast Fourier transform of the mean diurnal cycle ...
Creating output arrays ...
Calling numpy FFT function ...
(8, 36, 72)
Converting from complex-valued FFT to real-valued amplitude and phase ...
Output:
cycmean (36, 72)
maxvalue (3, 36, 72)
tmax (3, 36, 72)
"Re-decorating" Fourier harmonics with grid info, etc., ...
... and writing to netCDF.
This generates 3 files:
pr_IPSL-CM5A-LR_Jul_1997-1999_tmean.nc
pr_IPSL-CM5A-LR_Jul_1997-1999_S.nc
pr_IPSL-CM5A-LR_Jul_1997-1999_tS.nc
We can now run the last script: savg_fourierWrappedInOut.py
[15]:
%%bash
savg_fourier.py -p basic_diurnal_std_hourly_mean.py
Specifying latitude / longitude domain of interest ...
Preparing to write output to JSON file ...
Initializing dictionary of statistical results ...
TEMPLATE: pr_IPSL-CM5A-LR_Jul_1997-1999_S.nc
['demo_output/diurnal/nc/pr_IPSL-CM5A-LR_Jul_1997-1999_S.nc']
Reading Amplitude from demo_output/diurnal/nc/pr_IPSL-CM5A-LR_Jul_1997-1999_S.nc ...
Reading Phase from demo_output/diurnal/nc/pr_IPSL-CM5A-LR_Jul_1997-1999_tS.nc ...
Reading sftlf from demo_output/diurnal/nc/cmip5.IPSL-CM5A-LR.historical.r0i0p0.fx.atm.fx.sftlf.*.latestX.xml ...
Failed reading sftlf from file (error was: list index out of range)
Creating one for you
Global mean land fraction = 0.252
Area-averaging globally, over land only, and over ocean only ...
Converting singleton transient variables to plain floating-point numbers ...
Jul 1-harmonic amplitude, phase = 0.583 mm/d, 9.681 hrsLST averaged globally
Jul 1-harmonic amplitude, phase = 2.243 mm/d, 11.588 hrsLST averaged over land
Jul 1-harmonic amplitude, phase = 0.380 mm/d, 4.854 hrsLST averaged over ocean
Area-averaging globally, over land only, and over ocean only ...
Converting singleton transient variables to plain floating-point numbers ...
Jul 2-harmonic amplitude, phase = 0.569 mm/d, 0.159 hrsLST averaged globally
Jul 2-harmonic amplitude, phase = 1.553 mm/d, 11.904 hrsLST averaged over land
Jul 2-harmonic amplitude, phase = 0.251 mm/d, 0.698 hrsLST averaged over ocean
Area-averaging globally, over land only, and over ocean only ...
Converting singleton transient variables to plain floating-point numbers ...
Jul 3-harmonic amplitude, phase = 0.223 mm/d, 3.454 hrsLST averaged globally
Jul 3-harmonic amplitude, phase = 0.797 mm/d, 3.681 hrsLST averaged over land
Jul 3-harmonic amplitude, phase = 0.058 mm/d, 2.230 hrsLST averaged over ocean
{'IPSL-CM5A-LR': {'TRMM': {'harmonic1': {'amp_avg_lnd': 2.242995862103426, 'pha_avg_lnd': 11.587881354946367, 'amp_avg_ocn': 0.38036934839820885, 'pha_avg_ocn': 4.8537920964277905}, 'harmonic2': {'amp_avg_lnd': 1.5533814759275162, 'pha_avg_lnd': 11.90414360066512, 'amp_avg_ocn': 0.2513985389784492, 'pha_avg_ocn': 0.6976628431643584}, 'harmonic3': {'amp_avg_lnd': 0.7967068716487153, 'pha_avg_lnd': 3.68122152040087, 'amp_avg_ocn': 0.05825477671034634, 'pha_avg_ocn': 2.230342743709479}}}}
done
INFO::2021-11-15 14:57::pcmdi_metrics:: Results saved to a json file: /Users/ordonez4/Documents/git/pcmdi_metrics/doc/jupyter/Demo/demo_output/diurnal/json/pr_Jul_1997-1999_savg_DiurnalFourier.json
This creates the following file:
pr_Jul_1997-1999_savg_DiurnalFourier.json
[ ]: