Repository: climate Updated Branches: refs/heads/master 4cf79f3d0 -> 69793310f
CLIMATE-934 Fixed error with old Bounds constructor. Fixed error with matplotlib and date starting at epoch. Added Apache license. Minor Pylint. Project: http://git-wip-us.apache.org/repos/asf/climate/repo Commit: http://git-wip-us.apache.org/repos/asf/climate/commit/69793310 Tree: http://git-wip-us.apache.org/repos/asf/climate/tree/69793310 Diff: http://git-wip-us.apache.org/repos/asf/climate/diff/69793310 Branch: refs/heads/master Commit: 69793310f3ba76fa9165e20a7fd8c8275ade7555 Parents: 4cf79f3 Author: Michael Anderson <michaelanderson@Michaels-iMac.local> Authored: Wed Nov 22 13:27:00 2017 -0500 Committer: Michael Anderson <michaelanderson@Michaels-iMac.local> Committed: Wed Nov 22 13:27:00 2017 -0500 ---------------------------------------------------------------------- examples/time_series_with_regions.py | 135 +++++++++++++++++------------- 1 file changed, 79 insertions(+), 56 deletions(-) ---------------------------------------------------------------------- http://git-wip-us.apache.org/repos/asf/climate/blob/69793310/examples/time_series_with_regions.py ---------------------------------------------------------------------- diff --git a/examples/time_series_with_regions.py b/examples/time_series_with_regions.py index 3bb133c..168da90 100644 --- a/examples/time_series_with_regions.py +++ b/examples/time_series_with_regions.py @@ -1,19 +1,41 @@ +# Licensed to the Apache Software Foundation (ASF) under one +# or more contributor license agreements. See the NOTICE file +# distributed with this work for additional information +# regarding copyright ownership. The ASF licenses this file +# to you under the Apache License, Version 2.0 (the +# "License"); you may not use this file except in compliance +# with the License. You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, +# software distributed under the License is distributed on an +# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +# KIND, either express or implied. See the License for the +# specific language governing permissions and limitations +# under the License. + +''' + Download three netCDF files, process the files to be the same shape, + divide the data into 13 subregions and plot a monthly time series for each sub region. +''' + +import sys +import datetime +from os import path +from calendar import monthrange +import ssl + +import numpy as np + # Apache OCW lib immports -from ocw.dataset import Dataset, Bounds +from ocw.dataset import Bounds import ocw.data_source.local as local import ocw.data_source.rcmed as rcmed import ocw.dataset_processor as dsp -import ocw.evaluation as evaluation -import ocw.metrics as metrics import ocw.plotter as plotter import ocw.utils as utils -import datetime -import numpy as np -import numpy.ma as ma -from os import path -import sys - if sys.version_info[0] >= 3: from urllib.request import urlretrieve else: @@ -21,17 +43,17 @@ else: # But note that this might need an update when Python 4 # might be around one day from urllib import urlretrieve -import ssl + if hasattr(ssl, '_create_unverified_context'): ssl._create_default_https_context = ssl._create_unverified_context - # File URL leader -FILE_LEADER = "http://zipper.jpl.nasa.gov/dist/" +FILE_LEADER = 'http://zipper.jpl.nasa.gov/dist/' + # Three Local Model Files -FILE_1 = "AFRICA_KNMI-RACMO2.2b_CTL_ERAINT_MM_50km_1989-2008_pr.nc" -FILE_2 = "AFRICA_ICTP-REGCM3_CTL_ERAINT_MM_50km-rg_1989-2008_pr.nc" -FILE_3 = "AFRICA_UCT-PRECIS_CTL_ERAINT_MM_50km_1989-2008_pr.nc" +FILE_1 = 'AFRICA_KNMI-RACMO2.2b_CTL_ERAINT_MM_50km_1989-2008_pr.nc' +FILE_2 = 'AFRICA_ICTP-REGCM3_CTL_ERAINT_MM_50km-rg_1989-2008_pr.nc' +FILE_3 = 'AFRICA_UCT-PRECIS_CTL_ERAINT_MM_50km_1989-2008_pr.nc' LAT_MIN = -45.0 LAT_MAX = 42.24 @@ -56,55 +78,56 @@ region_counter = 0 # Download necessary NetCDF file if not present if not path.exists(FILE_1): - print("Downloading %s" % (FILE_LEADER + FILE_1)) + print('Downloading %s' % (FILE_LEADER + FILE_1)) urlretrieve(FILE_LEADER + FILE_1, FILE_1) if not path.exists(FILE_2): - print("Downloading %s" % (FILE_LEADER + FILE_2)) + print('Downloading %s' % (FILE_LEADER + FILE_2)) urlretrieve(FILE_LEADER + FILE_2, FILE_2) if not path.exists(FILE_3): - print("Downloading %s" % (FILE_LEADER + FILE_3)) + print('Downloading %s' % (FILE_LEADER + FILE_3)) urlretrieve(FILE_LEADER + FILE_3, FILE_3) -""" Step 1: Load Local NetCDF File into OCW Dataset Objects and store in list""" -target_datasets.append(local.load_file(FILE_1, varName, name="KNMI")) -target_datasets.append(local.load_file(FILE_2, varName, name="REGCM")) -target_datasets.append(local.load_file(FILE_3, varName, name="UCT")) - +# Step 1: Load Local NetCDF File into OCW Dataset Objects and store in list +target_datasets.append(local.load_file(FILE_1, varName, name='KNMI')) +target_datasets.append(local.load_file(FILE_2, varName, name='REGCM')) +target_datasets.append(local.load_file(FILE_3, varName, name='UCT')) -""" Step 2: Fetch an OCW Dataset Object from the data_source.rcmed module """ -print("Working with the rcmed interface to get CRU3.1 Daily Precipitation") +# Step 2: Fetch an OCW Dataset Object from the data_source.rcmed module +print('Working with the rcmed interface to get CRU3.1 Daily Precipitation') # the dataset_id and the parameter id were determined from # https://rcmes.jpl.nasa.gov/content/data-rcmes-database CRU31 = rcmed.parameter_dataset( 10, 37, LAT_MIN, LAT_MAX, LON_MIN, LON_MAX, START, END) - -""" Step 3: Processing datasets so they are the same shape ... """ -print("Processing datasets so they are the same shape") +# Step 3: Processing datasets so they are the same shape +print('Processing datasets so they are the same shape') CRU31 = dsp.water_flux_unit_conversion(CRU31) CRU31 = dsp.normalize_dataset_datetimes(CRU31, 'monthly') for member, each_target_dataset in enumerate(target_datasets): target_datasets[member] = dsp.subset(target_datasets[member], EVAL_BOUNDS) - target_datasets[member] = dsp.water_flux_unit_conversion(target_datasets[ - member]) + target_datasets[member] = dsp.water_flux_unit_conversion(target_datasets[member]) target_datasets[member] = dsp.normalize_dataset_datetimes( target_datasets[member], 'monthly') -print("... spatial regridding") +print('... spatial regridding') new_lats = np.arange(LAT_MIN, LAT_MAX, gridLatStep) new_lons = np.arange(LON_MIN, LON_MAX, gridLonStep) CRU31 = dsp.spatial_regrid(CRU31, new_lats, new_lons) - for member, each_target_dataset in enumerate(target_datasets): target_datasets[member] = dsp.spatial_regrid( target_datasets[member], new_lats, new_lons) -# find climatology monthly for obs and models +# Find climatology monthly for obs and models. CRU31.values, CRU31.times = utils.calc_climatology_monthly(CRU31) +# Shift the day of the month to the end of the month as matplotlib does not handle +# the xticks elegantly when the first date is the epoch and tries to determine +# the start of the xticks based on a value < 1. +for index, item in enumerate(CRU31.times): + CRU31.times[index] = datetime.date(item.year, item.month, monthrange(item.year, item.month)[1]) for member, each_target_dataset in enumerate(target_datasets): target_datasets[member].values, target_datasets[ @@ -112,41 +135,41 @@ for member, each_target_dataset in enumerate(target_datasets): # make the model ensemble target_datasets_ensemble = dsp.ensemble(target_datasets) -target_datasets_ensemble.name = "ENS" +target_datasets_ensemble.name = 'ENS' # append to the target_datasets for final analysis target_datasets.append(target_datasets_ensemble) -""" Step 4: Subregion stuff """ +# Step 4: Subregion stuff list_of_regions = [ - Bounds(-10.0, 0.0, 29.0, 36.5), - Bounds(0.0, 10.0, 29.0, 37.5), - Bounds(10.0, 20.0, 25.0, 32.5), - Bounds(20.0, 33.0, 25.0, 32.5), - Bounds(-19.3, -10.2, 12.0, 20.0), - Bounds(15.0, 30.0, 15.0, 25.0), - Bounds(-10.0, 10.0, 7.3, 15.0), - Bounds(-10.9, 10.0, 5.0, 7.3), - Bounds(33.9, 40.0, 6.9, 15.0), - Bounds(10.0, 25.0, 0.0, 10.0), - Bounds(10.0, 25.0, -10.0, 0.0), - Bounds(30.0, 40.0, -15.0, 0.0), - Bounds(33.0, 40.0, 25.0, 35.0)] - -region_list = [["R" + str(i + 1)] for i in xrange(13)] + Bounds(lat_min=-10.0, lat_max=0.0, lon_min=29.0, lon_max=36.5), + Bounds(lat_min=0.0, lat_max=10.0, lon_min=29.0, lon_max=37.5), + Bounds(lat_min=10.0, lat_max=20.0, lon_min=25.0, lon_max=32.5), + Bounds(lat_min=20.0, lat_max=33.0, lon_min=25.0, lon_max=32.5), + Bounds(lat_min=-19.3, lat_max=-10.2, lon_min=12.0, lon_max=20.0), + Bounds(lat_min=15.0, lat_max=30.0, lon_min=15.0, lon_max=25.0), + Bounds(lat_min=-10.0, lat_max=10.0, lon_min=7.3, lon_max=15.0), + Bounds(lat_min=-10.9, lat_max=10.0, lon_min=5.0, lon_max=7.3), + Bounds(lat_min=33.9, lat_max=40.0, lon_min=6.9, lon_max=15.0), + Bounds(lat_min=10.0, lat_max=25.0, lon_min=0.0, lon_max=10.0), + Bounds(lat_min=10.0, lat_max=25.0, lon_min=-10.0, lon_max=0.0), + Bounds(lat_min=30.0, lat_max=40.0, lon_min=-15.0, lon_max=0.0), + Bounds(lat_min=33.0, lat_max=40.0, lon_min=25.0, lon_max=35.0)] + +region_list = [['R' + str(i + 1)] for i in xrange(13)] for regions in region_list: firstTime = True - subset_name = regions[0] + "_CRU31" - # labels.append(subset_name) #for legend, uncomment this line + subset_name = regions[0] + '_CRU31' + labels.append(subset_name) subset = dsp.subset(CRU31, list_of_regions[region_counter], subset_name) tSeries = utils.calc_time_series(subset) results.append(tSeries) tSeries = [] firstTime = False for member, each_target_dataset in enumerate(target_datasets): - subset_name = regions[0] + "_" + target_datasets[member].name - # labels.append(subset_name) #for legend, uncomment this line + subset_name = regions[0] + '_' + target_datasets[member].name + labels.append(subset_name) subset = dsp.subset(target_datasets[member], list_of_regions[region_counter], subset_name) @@ -154,8 +177,8 @@ for regions in region_list: results.append(tSeries) tSeries = [] - plotter.draw_time_series(np.array(results), CRU31.times, labels, regions[ - 0], ptitle=regions[0], fmt='png') + plotter.draw_time_series(np.array(results), CRU31.times, labels, regions[0], + label_month=True, ptitle=regions[0], fmt='png') results = [] tSeries = [] labels = []