diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000000000000000000000000000000000000..27bf1b458cba0aa4c27d2ab6ec971b627beaa862 --- /dev/null +++ b/.gitignore @@ -0,0 +1,9 @@ +*.pyc +.vscode +*.png +test-reports +mapgenerator.egg-info +.vscode +.pytest_cache +.eggs +.coverage \ No newline at end of file diff --git a/Makefile b/Makefile index ace5605976e35c521ef5c02bb92a7f5bd0abd0bd..23879d80ceef4e466f7cee0ea96e4ee76b3e64f0 100644 --- a/Makefile +++ b/Makefile @@ -2,6 +2,6 @@ init: pip install -r requirements.txt test: - py.test tests + pytest tests .PHONY: init test diff --git a/bin/map_main.py b/bin/map_main.py index 4a052f7f585f77c1f5e81fc1e7a1574e41c46bb2..b174da71dc2f741c583a96b84c2eaee37585491a 100755 --- a/bin/map_main.py +++ b/bin/map_main.py @@ -17,17 +17,22 @@ # You should have received a copy of the GNU General Public License # along with MapGenerator. If not, see . -from __future__ import print_function -from mapgenerator.plotting.config import ArgumentParser - """ Main module for MapGenerator. Only contains an interface class to all functionality implemented on MapGenerator. """ + +from __future__ import print_function + +import logging +import matplotlib +mpl_logger = logging.getLogger('matplotlib') +mpl_logger.setLevel('CRITICAL') + +from mapgenerator.plotting.config import ArgumentParser from mapgenerator.plotting import plotmap import sys -import logging from datetime import datetime, timedelta -logging.basicConfig(level=logging.DEBUG) + log = logging.getLogger(__name__) @@ -44,7 +49,7 @@ class MapGenerator: """ try: args = self.parser.parse_args() - log.info(args) + #log.info(args) req = vars(args) # print help if no args if req.values()==[None for i in range(len(req.values()))]: @@ -56,18 +61,27 @@ class MapGenerator: # req.pop('srcvars') # pass only valid values and cast boolean strings to boolean res = {k:eval(v) if v in ('True','False') else v for k,v in req.items() if v} + + loglevel = 'loglevel' in res and \ + getattr(logging, res['loglevel']) or \ + getattr(logging, 'WARNING') + + mpl_logger.setLevel(loglevel) + log.setLevel(loglevel) + log.info(res) + # various plot types - if res.has_key('pltype'): + if 'pltype' in res: if res['pltype'] == 'cross': - plotmap.PlotCross().plot(**res) + plotmap.PlotCross(loglevel).plot(**res) elif res['pltype'] == 'timeseries': - plotmap.PlotSeries().plot(**res) + plotmap.PlotSeries(loglevel).plot(**res) elif res['pltype'] == 'map': - plotmap.PlotMap().plot(**res) + plotmap.PlotMap(loglevel).plot(**res) else: res['pltype'] = 'map' - plotmap.PlotMap().plot(**res) + plotmap.PlotMap(loglevel).plot(**res) except Exception as e: log.error('Unhandled exception on MapGenerator: %s' % e, exc_info=True) diff --git a/mapgenerator/evaluation/DAevaluator.py b/mapgenerator/evaluation/DAevaluator.py index 9f71d27b9737117e9f1f9c0aa49d28b77d501be8..2f6fe95d368e1189d966aec8f8d44fdf60d35472 100644 --- a/mapgenerator/evaluation/DAevaluator.py +++ b/mapgenerator/evaluation/DAevaluator.py @@ -559,25 +559,25 @@ if __name__ == "__main__": nargs = len(sys.argv) if nargs != 2: - print "Error: maximum of 2 args" + print("Error: maximum of 2 args") sys.exit(1) def exec_all(): #print "EXEC_ALL" try: testDA() - except Exception, e: - print "Error (DA):", e + except Exception as e: + print("Error (DA):", e) try: testPollen() - except Exception, e: - print "Error (Pollen):", e + except Exception as e: + print("Error (Pollen):", e) if nargs == 1: exec_all() else: attr = sys.argv[1] - print "ATTR:", attr + print("ATTR:", attr) current_module = sys.modules[__name__] # print attr, current_module.__name__ # for i in dir(current_module): diff --git a/mapgenerator/plotting/config.py b/mapgenerator/plotting/config.py index 9a35ee9656a343954c1b01eaa61b4acfadc7b022..261bb8afbdfafc56207b860cf29ded006a2f33aa 100644 --- a/mapgenerator/plotting/config.py +++ b/mapgenerator/plotting/config.py @@ -6,7 +6,7 @@ import os import configargparse as argparse -import ConfigParser +import configparser from optparse import OptionParser import mapgenerator @@ -32,7 +32,7 @@ class ArgumentParser(object): """ try: - self.parser = argparse.ArgumentParser(description='Main parser for MapGenerator.') + self.parser = argparse.ArgumentParser(description='') self.parser.add_argument('-V', '--version', action='version', version=mapgenerator.__version__, help="returns MapGenerator version number and exit") @@ -53,9 +53,15 @@ class ArgumentParser(object): help="main variable") self.parser.add_argument("--plottype", dest="pltype", - choices=['map','cross','timeseries'], + choices=['map', 'cross', 'timeseries'], help="type of plots to generate") # + self.parser.add_argument("--loglevel", + dest="loglevel", + choices=['DEBUG', 'INFO', 'WARNING', 'ERROR', + 'CRITICAL'], + help="logging levelbug (default=WARNING)") + self.parser.add_argument("--indir", dest="indir", help="directory where is/are stored the input file(s)") @@ -65,9 +71,6 @@ class ArgumentParser(object): self.parser.add_argument("--file_format", dest="filefmt", help="format of the output file(s) (default=png)") -# self.parser.add_argument("--debug", -# dest="debug", -# help="debug (default=False)") self.parser.add_argument("--timesteps", dest="timesteps", nargs="+", @@ -102,6 +105,9 @@ class ArgumentParser(object): self.parser.add_argument("--colorbar", dest="colorbar", help="toggle colorbar (default=True)") + self.parser.add_argument("--formats", + dest="formats", + help="colorbar labels format (default=auto)") self.parser.add_argument("--extend", dest="extend", choices=['both','max','min'], @@ -182,7 +188,7 @@ class ArgumentParser(object): help="dimension name and operation in the form DIMENSION,OPERATION (mean, max, min, ...) as documented in the numpy reference guide. If OPERATION is a number it will extract the correspondent element") self.parser.add_argument("--alpha", dest="alpha", - help="alpha value for colormap netween 0 and 1 (default=1") + help="alpha value for colormap netween 0 and 1 (default=1)") ############ self.parser.add_argument("--img_template", dest="img_template", @@ -278,7 +284,7 @@ class ArgumentParser(object): def readConf(section=None, fpath=None): """ Read configuration """ - config = ConfigParser.RawConfigParser() + config = configparser.RawConfigParser() config.read(fpath) if section == None: return config.sections() diff --git a/mapgenerator/plotting/definitions.py b/mapgenerator/plotting/definitions.py index b0a1e40268190bb62e3a1f96a63c20f826d5f17c..d391cb319b2ef69da5125e680a1122b192190c66 100644 --- a/mapgenerator/plotting/definitions.py +++ b/mapgenerator/plotting/definitions.py @@ -9,6 +9,7 @@ from netCDF4 import Dataset as nc #from grads import GaLab import pandas as pd import numpy as np +import numpy.ma as ma from scipy.interpolate import griddata from datetime import datetime from datetime import timedelta @@ -18,10 +19,35 @@ import copy import logging -logging.basicConfig(level=logging.DEBUG) +#logging.basicConfig(level=logging.WARNING) log = logging.getLogger(__name__) +def is_grid_regular(lon, lat): + """ Check spacing if uniform. Returns True if grid is regular """ + + #print("Checking if lon is strictly_increasing ...") + if not np.array([strictly_increasing(lon[x]) \ + for x in np.arange(lon.shape[0])]).all(): + return False + + #print("Checking if lat is strictly_increasing ...") + if not np.array([strictly_increasing(lat[:,y]) \ + for y in np.arange(lat.shape[1])]).all(): + return False + + #print("Checking lon spacing ...") + is_x_reg = np.array([(np.spacing(lon[x]) == np.spacing(lon[x+1])).all() \ + for x in np.arange(lon.shape[0]-1)]).all() + + print(is_x_reg) + if not is_x_reg: return is_x_reg + + #print("Checking lat spacing ...") + return np.array([(np.spacing(lat[:,y]) == np.spacing(lat[:,y+1])).all() \ + for y in np.arange(lat.shape[1]-1)]).all() + + def strictly_increasing(L): """ Check if vector values are strictly increasing """ return all(x= cur_coord[0]) & \ (var1[:] <= cur_coord[-1])) - #log.info("idxs %s" % (str(idxs))) + log.info("idxs %s %s %s %s" % (str(idxs), var1.shape, len(idxs), idxs[0].shape)) return var1[idxs], idxs[0] log.info("Not subsetting") return var[:], None @@ -92,7 +123,8 @@ def parse_parameters_list(plist): except: rg = tmp steps = None - start, end = [float(i.replace('m','-')) for i in rg.split('-')] + start, end = [float(i.replace('m','-')) if i.find('.')>=0 else + int(i.replace('m','-')) for i in rg.split('-')] if not steps: steps = (end-start)/10 else: @@ -150,13 +182,14 @@ def parse_path(dir, f): class MapData(object): """Wrapper class to handle all the data needed to draw a map""" - def __init__(self, **kwargs): + def __init__(self, loglevel='WARNING', **kwargs): + log.setLevel(loglevel) log.info("mapData initializer") self.reset() self.maxReset() # Process optional arguments for dtype in ('map_data','wind_data','contour_data','scatter_data'): - if kwargs.has_key(dtype): + if dtype in kwargs: setattr(self, dtype, kwargs[dtype]) def reset(self): @@ -286,9 +319,9 @@ class MapDataHandler(object): self.ncdf[0].append(fin) for v in fin.variables: - if v in ('lon','longitude','nav_lon'): + if v in ('longitude','lon','nav_lon'): x, x_idxs = extract_coords(self, fin.variables[v]) - elif v in ('lat','latitude','nav_lat'): + elif v in ('latitude','lat','nav_lat'): y, y_idxs = extract_coords(self, fin.variables[v], coord='lats') self.grid = (x, y) @@ -402,6 +435,7 @@ class MapDataHandler(object): values = srcvar[tstep,:,:][self.grid_idxs[1],:][:,self.grid_idxs[0]] else: values = srcvar[tstep,:,:] + log.info("1. SRCVAR: %s - GRID_IDXS: %s" % (str(values.shape), str(self.grid_idxs))) else: # otherwise search for a variable of the dataset into the # expression given @@ -417,6 +451,7 @@ class MapDataHandler(object): newvar = srcvar[tstep,:,:] values = eval(expvar) break + log.info("2. SRCVAR: %s - GRID_IDXS: %s" % (str(values.shape), str(self.grid_idxs))) if mapData: #exp_res = self.ncdf.exp(EXPVAR) mapData.append(values) @@ -520,6 +555,8 @@ class DataFrameHandler(object): if type(self.reader) == dict: self.reader = pd.DataFrame(self.reader) + self.reader = self.reader.dropna() + if datecol not in self.reader.columns: return dict(self.reader) @@ -569,7 +606,7 @@ class DataTransform(object): class MapDrawOptions(object): """ """ - def __init__(self, **kwargs): + def __init__(self, loglevel='WARNING', **kwargs): ## General Options self.area_thresh = None # Default area threshold for maps self.resolution = None # Default resolution is intermediate @@ -614,6 +651,8 @@ class MapDrawOptions(object): self.shapef_width_init = 0.75 # Initial line width for shapefiles self.shapef_width_step = 0.25 # Line width difference between two consecutive shapefiles + log.setLevel(loglevel) + def hasShapeFiles(self): return (len(self.shapefiles) > 0) @@ -660,7 +699,7 @@ class MapGenerator(object): class MapCross(MapGenerator): """ Define common attributes to maps and cross sections """ - def __init__(self, **kwargs): + def __init__(self, loglevel='WARNING', **kwargs): """ Initialize class with MapGenerator attributes plus some others """ MapGenerator.__init__(self, **kwargs) self.bounds = kwargs.get('bounds', None) @@ -668,7 +707,7 @@ class MapCross(MapGenerator): self.colors = kwargs.get('colors', 'jet') # self.over = kwargs.get('over', None) # self.under = kwargs.get('under', None) - self.bad = kwargs.get('bad', 'gray') + self.bad = kwargs.get('bad', None) self.lat = kwargs.get('lat', []) self.lon = kwargs.get('lon', []) self.wind = kwargs.get('wind', []) @@ -679,6 +718,7 @@ class MapCross(MapGenerator): self.contours_int = kwargs.get('contours_int', 1) self.smooth = kwargs.get('smooth', False) self.colorbar = kwargs.get('colorbar', True) + self.formats = kwargs.get('formats', None) self.noruntime = kwargs.get('noruntime', False) self.timesteps = kwargs.get('timesteps', '0') self.limits = kwargs.get('limits', None) diff --git a/mapgenerator/plotting/plotmap.py b/mapgenerator/plotting/plotmap.py index 235f1e9ec32018a2ea9c8b3ad976020b456a7d30..a4ab4c29c9fcd498775a0564f371673770ba3ca9 100644 --- a/mapgenerator/plotting/plotmap.py +++ b/mapgenerator/plotting/plotmap.py @@ -6,7 +6,7 @@ import matplotlib as mpl -#mpl.use('Agg') +mpl.use('Qt5Agg') import matplotlib.pyplot as plt from mpl_toolkits.basemap import Basemap @@ -17,7 +17,6 @@ from datetime import datetime #from datetime import timedelta #import time import os.path -import commands import subprocess #import sys from collections import Iterable @@ -25,48 +24,48 @@ from PIL import Image #from mapgenerator.plotting.mapgenerator import MapCross #from mapgenerator.plotting.mapgenerator import MapDrawOptions #from mapgenerator.plotting.mapgenerator import MapDataHandler -from definitions import MapData -from definitions import MapCross -from definitions import MapDrawOptions -from definitions import MapDataHandler -from definitions import DataFrameHandler -from definitions import strictly_increasing -from definitions import do_interpolation -from definitions import parse_parameter -from definitions import parse_parameters -from definitions import parse_parameters_list -from definitions import parse_path -import urllib2 +from .definitions import MapData +from .definitions import MapCross +from .definitions import MapDrawOptions +from .definitions import MapDataHandler +from .definitions import DataFrameHandler +from .definitions import is_grid_regular +from .definitions import do_interpolation +from .definitions import parse_parameter +from .definitions import parse_parameters +from .definitions import parse_parameters_list +from .definitions import parse_path +from urllib.request import urlopen import copy from glob import glob import logging -logging.basicConfig(level=logging.DEBUG) +#logging.basicConfig(level=logging.WARNING) log = logging.getLogger(__name__) class PlotSeries(): """ Main class for plotting time series """ - def __init__(self, **kwargs): + def __init__(self, loglevel='WARNING', **kwargs): pass class PlotCross(): """ Main class for plotting cross sections """ - def __init__(self, **kwargs): + def __init__(self, loglevel='WARNING', **kwargs): pass class PlotMap(MapCross, MapDrawOptions): """ Main class for plotting maps """ - def __init__(self, **kwargs): + def __init__(self, loglevel='WARNING', **kwargs): """ Initialize class with attributes """ - MapCross.__init__(self, **kwargs) - MapDrawOptions.__init__(self, **kwargs) + MapCross.__init__(self, loglevel, **kwargs) + MapDrawOptions.__init__(self, loglevel, **kwargs) self.nocontourf = False self.maxdata = False self.maxtitle = None @@ -77,7 +76,7 @@ class PlotMap(MapCross, MapDrawOptions): self.kmz = False self.norm = None self.cmap = None - self.formats = None + #self.formats = None self.map = None #Map #self.first_image = True # Is is the first image of the batch? self.mapName = None # name of the figure to map animation colors @@ -86,6 +85,7 @@ class PlotMap(MapCross, MapDrawOptions): # 'countries': # {'linewidth': 0.3, 'color': 'grey', }, # } + log.setLevel(loglevel) def __setattr__(self, key, value): #log.info("SETATTR: {} - {}".format(key, value)) @@ -175,24 +175,24 @@ class PlotMap(MapCross, MapDrawOptions): # a = plt.axes([nap, 0, .14, 1])#, frameon=False) # plt.axis('off') mpl.rcParams['axes.linewidth'] = 0.1 + mpl.rcParams['axes.formatter.useoffset'] = False if self.ticks: #log.info("ticks %s" % str(self.ticks)) self.ticks = parse_parameters_list(self.ticks) elif self.bounds: self.ticks = self.bounds + else: + try: + self.ticks = mco.levels + except: + pass + log.info("***** {} *****".format(self.formats)) if self.map: cb = self.map.colorbar(mco, ticks=self.ticks, format=self.formats, pad="2%", size="3%", extend=self.extend, drawedges=drawedges) else: cb = plt.colorbar(mco, ticks=self.ticks, format=self.formats, pad=.06, extend=self.extend, drawedges=drawedges) - if self.ticks: - cb.set_ticklabels(self.ticks) - else: - try: - cb.set_ticklabels(mco.levels) - except: - pass cb.ax.tick_params(labelsize=float(self.coordsopts[1])) for lin in cb.ax.yaxis.get_ticklines(): lin.set_visible(False) @@ -204,7 +204,7 @@ class PlotMap(MapCross, MapDrawOptions): shrink = 0.75 upper = 0.863 - 0.45 * (0.8 - shrink) # formula seems to work ok lower = 0.137 + 0.45 * (0.8 - shrink) - span = upper - lower + span = upper - lower arrows = self.limits if self.bounds: @@ -248,8 +248,8 @@ class PlotMap(MapCross, MapDrawOptions): } try: - pTitle = unicode(title.decode('utf-8')) % titleDic #("%sh" % sTime, valid) - except Exception, e: + pTitle = str(title) % titleDic #("%sh" % sTime, valid) + except Exception as e: print("Title error: %s", str(e)) pTitle = title % titleDic #("%sh" % sTime, valid) @@ -266,7 +266,7 @@ class PlotMap(MapCross, MapDrawOptions): now_t = datetime.now() diff_t = now_t - last_t if (msg != None): - print 'TIME:', msg, ' done in ', diff_t.seconds, ' s' + print('TIME:', msg, ' done in ', diff_t.seconds, ' s') last_t = now_t def runCommand(self, commStr, fatal=True): @@ -323,7 +323,6 @@ class PlotMap(MapCross, MapDrawOptions): # else: # strs.append(str(b)) # -# self.formats = mpl.ticker.FixedFormatter(strs) # Set the colormap self.setColorMap() if not self.resolution: @@ -334,6 +333,13 @@ class PlotMap(MapCross, MapDrawOptions): def genImageMap(self, figName, grid, data, imgTitle, cur_scatter_data=None, **kwargs): """ Generate image map """ + + # overwrite option + if os.path.exists(figName + ".png") and not self.overwrite: + print(figName, " already exists.") + plt.clf() + return figName + map_data = data.getMapData() # FIXME scatter_data = data.getScatterData() or cur_scatter_data @@ -370,29 +376,25 @@ class PlotMap(MapCross, MapDrawOptions): glon, glat = np.meshgrid(glon, glat) #print(glon, glat) - #log.info("GLON: %s, GLAT: %s" % (str(glon.shape), str(glat.shape))) + log.info("1. GLON: %s, GLAT: %s" % (str(glon.shape), str(glat.shape))) # if curvilinear grid do interpolation, return already gridded coords - if not strictly_increasing(glon[0]) or not strictly_increasing(glat[:,0]): + if not is_grid_regular(glon, glat): map_data[0], glon, glat = do_interpolation(map_data[0], glon, glat) + log.info("2. GLON: %s, GLAT: %s" % (str(glon.shape), str(glat.shape))) + if not self.nomap: x, y = self.map(*(glon, glat)) else: x, y = glon, glat -# print "X", x.shape -# print "Y", y.shape + + log.info("3. GLON: %s, GLAT: %s" % (str(x.shape), str(y.shape))) # print "map_data", map_data[0].shape #self.printTime("meshgrid") #log.info("X: %s, Y: %s" % (str(x.shape), str(y.shape))) - # overwrite option - if os.path.exists(figName + ".png") and not self.overwrite: - print figName, " already exists." - plt.clf() - return figName - # # nomap option # if self.nomap and os.path.exists("%s-%s/%s.png" % (runDate, varName, figName)) and not self.overwrite: # print figName, " already exists." @@ -421,28 +423,61 @@ class PlotMap(MapCross, MapDrawOptions): self.mgplot.add_axes([xstart_norm, ystart_norm, xsize_norm, ysize_norm]) #print "Switching AXES for MAP:", ax_map - if self.alpha!=None and self.bounds!=None: + if self.alpha != None and self.bounds != None: if self.extend in ('min', 'both'): - map_data[0] = np.ma.masked_where(map_data[0]=self.bounds[0])&(map_data[0]= self.bounds[0]) & (map_data[0] < self.bounds[1]), map_data[0]) mco = None + log.info(type(map_data[0])) + map_data[0] = np.ma.filled(map_data[0], np.nan) + log.info(type(map_data[0])) + if not self.nocontourf and not self.nomap and self.smooth: - #log.info("X: %s, Y: %s, DATA: %s" % (x.shape, y.shape, map_data[0].shape)) - mco = self.map.contourf(x, y, map_data[0], cmap=self.cmap, norm=self.norm, levels=self.bounds, extend=self.extend, horizontalalignment='center', alpha=self.alpha, antialiased=True) + log.info("X: %s, Y: %s, DATA: %s" % (x.shape, y.shape, map_data[0])) + mco = self.map.contourf(x, y, + map_data[0], + cmap=self.cmap, + norm=self.norm, + levels=self.bounds, + extend=self.extend, + horizontalalignment='center', + alpha=self.alpha, + antialiased=True) elif not self.nocontourf and not self.nomap and not self.smooth: - mco = self.map.pcolormesh(x, y, map_data[0], cmap=self.cmap, norm=self.norm, alpha=self.alpha, antialiased=True) + log.info("X: %s, Y: %s, DATA: %s" % (x.shape, y.shape, map_data[0])) + mco = self.map.pcolormesh(x, y, + map_data[0], + cmap=self.cmap, + norm=self.norm, + alpha=self.alpha, + antialiased=True) elif not self.nocontourf and self.nomap and self.smooth: - mco = plt.contourf(x, y, map_data[0], cmap=self.cmap, norm=self.norm, levels=self.bounds, extend=self.extend, horizontalalignment='center', alpha=self.alpha, antialiased=True) + mco = plt.contourf(x, y, + map_data[0], + cmap=self.cmap, + norm=self.norm, + levels=self.bounds, + extend=self.extend, + horizontalalignment='center', + alpha=self.alpha, + antialiased=True) + plt.axis([self.lon[0], self.lon[-1], self.lat[0], self.lat[-1]]) elif not self.nocontourf and self.nomap and not self.smooth: - mco = plt.pcolormesh(x, y, map_data[0], cmap=self.cmap, norm=self.norm, alpha=self.alpha, antialiased=True) + mco = plt.pcolormesh(x, y, + map_data[0], + cmap=self.cmap, + norm=self.norm, + alpha=self.alpha, + antialiased=True) + plt.axis([self.lon[0], self.lon[-1], self.lat[0], self.lat[-1]]) # Draw shape files @@ -459,18 +494,54 @@ map_data[0]) if scatter_data is not None: log.info("Plotting scatter data: {} of keys {}".format(str(scatter_data), str(scatter_data.keys()))) if mco and self.bounds: - self.map.scatter(scatter_data['lon'].tolist(), scatter_data['lat'].tolist(), s=scatter_data['size'].tolist(), c=scatter_data['color'].tolist(), marker='o', linewidth=0.3, vmin=self.bounds[0], vmax=self.bounds[-1], cmap=self.cmap, zorder=10) + self.map.scatter(scatter_data['lon'].tolist(), + scatter_data['lat'].tolist(), + s=scatter_data['size'].tolist(), + c=scatter_data['color'].tolist(), + marker='o', + linewidth=0.3, + vmin=self.bounds[0], + vmax=self.bounds[-1], + cmap=self.cmap, + norm=self.norm, + zorder=10) elif mco: - self.map.scatter(scatter_data['lon'].tolist(), scatter_data['lat'].tolist(), s=scatter_data['size'].tolist(), c=scatter_data['color'].tolist(), marker='o', linewidth=0.3, cmap=self.cmap, zorder=10) + self.map.scatter(scatter_data['lon'].tolist(), + scatter_data['lat'].tolist(), + s=scatter_data['size'].tolist(), + c=scatter_data['color'].tolist(), + marker='o', + linewidth=0.3, + cmap=self.cmap, + norm=self.norm, + zorder=10) elif self.bounds: - mco = self.map.scatter(scatter_data['lon'].tolist(), scatter_data['lat'].tolist(), s=scatter_data['size'].tolist(), c=scatter_data['color'].tolist(), marker='o', linewidth=0.3, vmin=self.bounds[0], vmax=self.bounds[-1], cmap=self.cmap, zorder=10) + mco = self.map.scatter(scatter_data['lon'].tolist(), + scatter_data['lat'].tolist(), + s=scatter_data['size'].tolist(), + c=scatter_data['color'].tolist(), + marker='o', + linewidth=0.3, + vmin=self.bounds[0], + vmax=self.bounds[-1], + cmap=self.cmap, + norm=self.norm, + zorder=10) else: - mco = self.map.scatter(scatter_data['lon'].tolist(), scatter_data['lat'].tolist(), s=scatter_data['size'].tolist(), c=scatter_data['color'].tolist(), marker='o', linewidth=0.3, cmap=self.cmap, zorder=10) + mco = self.map.scatter(scatter_data['lon'].tolist(), + scatter_data['lat'].tolist(), + s=scatter_data['size'].tolist(), + c=scatter_data['color'].tolist(), + marker='o', + linewidth=0.3, + cmap=self.cmap, + norm=self.norm, + zorder=10) if data.hasWindData(): winds = data.getWindData() X,Y,U,V = delete_masked_points(x.ravel(), y.ravel(), winds['u'].ravel(), winds['v'].ravel()) - print "Wind scale is", self.wind_units, self.wind_scale + #print("Wind scale is", self.wind_units, self.wind_scale) if winds.has_key('barbs'): self.map.barbs(X, Y, U, V, #units=self.wind_units, #headlength=self.wind_head_length, @@ -519,26 +590,40 @@ map_data[0]) for ex in exclude: lvls = filter (lambda ex: ex != 0, lvls) if (len(lvls) > 0) and map_data != None: - mco = self.map.contourf(x, y, map_data[0], cmap=self.cmap, norm=self.norm, levels=self.bounds, extend=self.extend, horizontalalignment='center', alpha=self.alpha) + mco = self.map.contourf(x, y, + map_data[0], + cmap=self.cmap, + norm=self.norm, + levels=self.bounds, + extend=self.extend, + horizontalalignment='center', + alpha=self.alpha) #print "-----------------", map_data, cdata, "--------------------" if map_data != None and (map_data[0] == cdata).all(): - print ":::::::::::::::::::::::: SAME !!! :::::::::::::::::::" - CS = plt.contour(x, y, cdata, levels=self.bounds, - colors=self.contours_color, - linewidths=self.contours_linewidth, alpha=self.alpha) + print(":::::::::::::::::::::::: SAME !!! :::::::::::::::::::") + CS = plt.contour(x, y, + cdata, + levels=self.bounds, + colors=self.contours_color, + linewidths=self.contours_linewidth, + alpha=self.alpha) else: - print ":::::::::::::::::::::::: DIFFERENT !!! :::::::::::::::::::" - print "MIN:", cMin, "MAX:", cMax - CS = plt.contour(x, y, cdata, levels=lvls, - colors=self.contours_color, - linewidths=self.contours_linewidth, alpha=self.alpha) + print(":::::::::::::::::::::::: DIFFERENT !!! :::::::::::::::::::") + print("MIN:", cMin, "MAX:", cMax) + CS = plt.contour(x, y, + cdata, + levels=lvls, + colors=self.contours_color, + linewidths=self.contours_linewidth, + alpha=self.alpha) if self.contours_label: - self.mgplot.clabel(CS, inline=1, - fontsize=self.contours_lbl_fontsize, - #backgroundcolor='r', - fmt=self.contours_lbl_format) + self.mgplot.clabel(CS, + inline=1, + fontsize=self.contours_lbl_fontsize, + #backgroundcolor='r', + fmt=self.contours_lbl_format) # else: # print "Not drawing contours since levels has ", len(lvls), " elements" @@ -597,7 +682,7 @@ size=%s&\ dpi=%d&\ format=png24&\ f=image" % (self.lon[0], self.lat[0], self.lon[-1], self.lat[-1], size, self.dpi) - self.map.imshow(imread(urllib2.urlopen(basemap_url)), origin='upper') + self.map.imshow(imread(urlopen(basemap_url)), origin='upper') #logo if self.logo: @@ -615,10 +700,11 @@ f=image" % (self.lon[0], self.lat[0], self.lon[-1], self.lat[-1], size, self.dpi self.mgplot.figimage(nim, self.logo[1], self.logo[2], zorder=10) # SAVEFIG - log.info("printing %s" % figName) + fullname = "{}.{}".format(figName, self.filefmt) + log.warning("printing {}".format(fullname)) if self.kml or self.kmz: #self.mgplot.axes(frameon=0) - self.mgplot.savefig("%s.%s" % (figName, self.filefmt), + self.mgplot.savefig(fullname, bbox_inches='tight', frameon=0, pad_inches=0, dpi=self.dpi, transparent=True) @@ -641,7 +727,8 @@ f=image" % (self.lon[0], self.lat[0], self.lon[-1], self.lat[-1], size, self.dpi fontsize=self.fontsize, zorder=0, ) - self.mgplot.savefig("%s.%s" % (figName, self.filefmt), bbox_inches='tight', pad_inches=.2, dpi=self.dpi) + self.mgplot.savefig(fullname, bbox_inches='tight', pad_inches=.2, dpi=self.dpi) + #self.mgplot.show() #plt.tight_layout() #self.printTime("saving") @@ -692,7 +779,7 @@ f=image" % (self.lon[0], self.lat[0], self.lon[-1], self.lat[-1], size, self.dpi #comm = "montage %s -tile 2x -geometry 800x600+1+1 %s.gif" % (imgs, newMapNames[idx]) comm = "montage %s -tile 2x -geometry +1+1 %s.gif" % (imgs, newimg) #print "COMPARE COMMAND:", comm - st, out = commands.getstatusoutput(comm) + st, out = subprocess.getstatusoutput(comm) idx += 1 #montage -tile 2x -geometry 800x600+1+1 @@ -717,9 +804,9 @@ f=image" % (self.lon[0], self.lat[0], self.lon[-1], self.lat[-1], size, self.dpi if tpl: nimg0 = tpl % {'date':date, 'step':nimg0} nimg1 = tpl % {'date':date, 'step':nimg1} - st, out = commands.getstatusoutput("convert -delay 75 -loop 0 %s.gif %s.gif" % (nimg0, nimg1)) + st, out = subprocess.getstatusoutput("convert -delay 75 -loop 0 %s.gif %s.gif" % (nimg0, nimg1)) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) def genKML(self, figNames, srcvars, dims, online): @@ -826,8 +913,7 @@ f=image" % (self.lon[0], self.lat[0], self.lon[-1], self.lat[-1], size, self.dpi os.remove("%s.kml" % kmlName) - def aplot(self, x, y, map_data=None, wind_data=None, contour_data=None, -scatter_data=None, **kwargs): + def aplot(self, x, y, map_data=None, wind_data=None, contour_data=None, scatter_data=None, **kwargs): """ Plot one or more maps directly from numerical arrays """ # store options to be restored at the end @@ -841,7 +927,7 @@ scatter_data=None, **kwargs): # parse wind arguments if (self.windopts and (len(self.windopts) == 4)): - self.wind_scale = int(self.windopts[3]); + self.wind_scale = int(self.windopts[3]) log.info("SHP {}".format(str(self.shapefiles))) if self.shapefiles: self.shapefiles = [f.replace(".shp","") for f in glob(self.shapefiles)] @@ -899,6 +985,16 @@ scatter_data=None, **kwargs): # restore options without local function parameters vars(self).update(localvars) + def plot_cube(self, cube, **kwargs): + if 'title' not in kwargs: + kwargs['title'] = '{0.long_name} ({0.units})'.format(cube) + + self.aplot( + cube.coord('longitude').points, + cube.coord('latitude').points, + map_data=cube.data, + **kwargs + ) def plot(self, **kwargs): """ Plot one or more maps from netCDF file(s) """ @@ -915,12 +1011,12 @@ scatter_data=None, **kwargs): if isinstance(self.srcvars, str): self.srcvars = [self.srcvars] elif not isinstance(self.srcvars, list): - print('Error') + log.error('Error') if isinstance(self.srcfiles, str): self.srcfiles = [self.srcfiles] elif not isinstance(self.srcfiles, list): - print('Error') + log.error('Error') # parse wind arguments if (self.windopts and (len(self.windopts) == 4)): @@ -1050,7 +1146,7 @@ scatter_data=None, **kwargs): pTitle ) else: - print "Missing MAXTITLE, thus not generating MAX images!" + print("Missing MAXTITLE, thus not generating MAX images!") if self.kml or self.kmz: plt.clf() @@ -1082,7 +1178,7 @@ scatter_data=None, **kwargs): for lin in cb.ax.yaxis.get_ticklines(): lin.set_visible(False) - tit = unicode(pTitle.decode('utf-8')) + tit = str(pTitle) idx1 = tit.find('(') idx2 = tit.find(')') varUnit = tit[idx1:idx2+1] @@ -1111,12 +1207,12 @@ scatter_data=None, **kwargs): #generate KMZ - Offline if self.kmz: - print "Generating KMZ ..." + print("Generating KMZ ...") self.genKML(figNames, self.srcvars, dims, online=False) #generate KML - Online if self.kml: - print "Generating KML ..." + print("Generating KML ...") self.genKML(figNames, self.srcvars, dims) #generate animation. @@ -1156,7 +1252,7 @@ scatter_data=None, **kwargs): def loadConf(self, section=None, fpath=None, reset=False): """ Load existing configurations from file. """ - from config import readConf + from .config import readConf if fpath == None: fpath = parse_path(self.config_dir, self.config_file) @@ -1179,7 +1275,7 @@ scatter_data=None, **kwargs): def writeConf(self, section, fpath=None): """ Write configurations on file. """ - from config import writeConf + from .config import writeConf if fpath == None: fpath = parse_path(self.config_dir, self.config_file) diff --git a/mapgenerator/utils/cross.py b/mapgenerator/utils/cross.py index cb9560d85ef08ef6990f8beaf9219c8decbcd6bf..46a12e9d544073fab931abcd2e4eaf8cb37cfa53 100644 --- a/mapgenerator/utils/cross.py +++ b/mapgenerator/utils/cross.py @@ -25,7 +25,7 @@ UNDER = '#ffffff' FREQ = 3 # Timestep interval to plot INT = 2 -# Timesteps to plot (if < total timesteps of data) +# Timesteps to plot (if < total timesteps of data) TLIMIT = 24 # Lon/Lat line style @@ -121,20 +121,20 @@ if __name__ == "__main__": #LON:2, LAT:3, STNAME:5 if DEBUG: - print "Loading stations ..." + print("Loading stations ...") STLIST = np.loadtxt(TXNAME, delimiter=';', skiprows=1, dtype='str') if DEBUG: - print "Loading NetCDF ..." + print("Loading NetCDF ...") f = nc(NCNAME, 'r') if DEBUG: - print "Loading Coords ..." + print("Loading Coords ...") lons = f.variables['lon'][:].astype(DTYPE) lats = f.variables['lat'][:].astype(DTYPE) if DEBUG: - print "FindCoordIndexes ..." + print("FindCoordIndexes ...") LON0_INDEX = findCoordIndex(lons, LON0) LON1_INDEX = findCoordIndex(lons, LON1) LAT0_INDEX = findCoordIndex(lats, LAT0) @@ -145,86 +145,86 @@ if __name__ == "__main__": #level = level[LEV0_INDEX:LEV1_INDEX] if DEBUG: - print "Loading Main Var ..." + print("Loading Main Var ...") main_var = f.variables[VARNAME][:,LEV0_INDEX:LEV1_INDEX,LAT0_INDEX:LAT1_INDEX,LON0_INDEX:LON1_INDEX].astype(DTYPE) main_var = (np.ma.masked_where(main_var == MVAL, main_var)*MULTVAR).astype(DTYPE) if DEBUG: - print "Main Var Memory Usage:", main_var.nbytes + print("Main Var Memory Usage:", main_var.nbytes) if DEBUG: - print "Loading Z ..." + print("Loading Z ...") z = f.variables['fis'][:,LAT0_INDEX:LAT1_INDEX,LON0_INDEX:LON1_INDEX].astype(DTYPE) z = (np.ma.masked_where(z == MVAL, z)/9800).astype(DTYPE) if DEBUG: - print "Z Memory Usage:", z.nbytes + print("Z Memory Usage:", z.nbytes) if DEBUG: - print "Loading Time ..." + print("Loading Time ...") v_time = f.variables['time'] time = v_time[:] units = v_time.units if DEBUG: - print "Loading HSL ..." + print("Loading HSL ...") hlevel = f.variables['hsl'][:,LEV0_INDEX:LEV1_INDEX,LAT0_INDEX:LAT1_INDEX,LON0_INDEX:LON1_INDEX].astype(DTYPE) hlevel = (np.ma.masked_where(hlevel == MVAL, hlevel)/1000).astype(DTYPE) if DEBUG: - print "HSL Memory Usage:", hlevel.nbytes + print("HSL Memory Usage:", hlevel.nbytes) if DEBUG: - print "Loading TSL ..." + print("Loading TSL ...") tlevel = f.variables['tsl'][:,LEV0_INDEX:LEV1_INDEX,LAT0_INDEX:LAT1_INDEX,LON0_INDEX:LON1_INDEX].astype(DTYPE) tlevel = np.ma.masked_where(tlevel == MVAL, tlevel).astype(DTYPE) if DEBUG: - print "TSL Memory Usage:", tlevel.nbytes + print("TSL Memory Usage:", tlevel.nbytes) if DEBUG: - print "Loading PRES ..." + print("Loading PRES ...") plevel = (f.variables['pres'][:]*100).astype(DTYPE) if DEBUG: - print "Creating ones ..." + print("Creating ones ...") level = (ones(hlevel.shape)*MVAL).astype(DTYPE) th = (ones(hlevel.shape)*MVAL).astype(DTYPE) if DEBUG: - print "Populating tmplevel and tmpth ..." + print("Populating tmplevel and tmpth ...") for tt in range(0,len(time)): for l in range(0,hlevel.shape[1]): level[tt,l,:,:] = hlevel[tt,l,:,:] th[tt,l,:,:] = tlevel[tt,l,:,:]*((100000/plevel[l])**0.286) #taking reference pressure at 1000kPa if DEBUG: - print "Deleting tmp arrays ..." + print("Deleting tmp arrays ...") del plevel del hlevel del tlevel if DEBUG: - print "Creating level and th ..." + print("Creating level and th ...") level = np.ma.masked_where(level == MVAL, level).astype(DTYPE) th = np.ma.masked_where(th == MVAL, th).astype(DTYPE) if DEBUG: - print "++++++++++++++++++++++++++++++++++++++++++++++" - print "LEVEL", level.shape - print "++++++++++++++++++++++++++++++++++++++++++++++" - print "TH", th.shape - print "++++++++++++++++++++++++++++++++++++++++++++++" + print("++++++++++++++++++++++++++++++++++++++++++++++") + print("LEVEL", level.shape) + print("++++++++++++++++++++++++++++++++++++++++++++++") + print("TH", th.shape) + print("++++++++++++++++++++++++++++++++++++++++++++++") if DEBUG: - print "Loading USL ..." + print("Loading USL ...") udat = f.variables['usl_h'][:,LEV0_INDEX:LEV1_INDEX,LAT0_INDEX:LAT1_INDEX,LON0_INDEX:LON1_INDEX].astype(DTYPE) udat = (np.ma.masked_where(udat == MVAL, udat)*2).astype(DTYPE) if DEBUG: - print "UDAT Memory Usage:", udat.nbytes + print("UDAT Memory Usage:", udat.nbytes) if DEBUG: - print "Loading VSL ..." + print("Loading VSL ...") vdat = f.variables['vsl_h'][:,LEV0_INDEX:LEV1_INDEX,LAT0_INDEX:LAT1_INDEX,LON0_INDEX:LON1_INDEX].astype(DTYPE) vdat = (np.ma.masked_where(vdat == MVAL, vdat)*2).astype(DTYPE) if DEBUG: - print "VDAT Memory Usage:", vdat.nbytes + print("VDAT Memory Usage:", vdat.nbytes) cmap = mpl.colors.ListedColormap(COLORS) cmap.set_over(OVER) @@ -243,25 +243,25 @@ if __name__ == "__main__": STNAME = line[5] TIME = TIME_INIT - print "Station", STNAME, "LON", LON, "LAT", LAT + print("Station", STNAME, "LON", LON, "LAT", LAT) try: LON_INDEX = findCoordIndex(lons, LON) LAT_INDEX = findCoordIndex(lats, LAT) except: - print "Coords not found. Continue" + print("Coords not found. Continue") continue - print "TIME", range(int(time[0]),int(time[-1])+INT*FREQ, INT*FREQ) + print("TIME", range(int(time[0]),int(time[-1])+INT*FREQ, INT*FREQ)) for t in range(0,len(time),INT): plt.clf() step = TDGNUM % (t*FREQ) - print "TIME", t, "STEP", step + print("TIME", t, "STEP", step) MONTH = TIME.strftime("%b").upper() - # data plot (t, var, idx, typ + # data plot (t, var, idx, typ #dplot(t, lats, LON_INDEX, 'lat') #dplot(t, lons, LAT_INDEX, 'lon') @@ -279,7 +279,7 @@ if __name__ == "__main__": continue udat1[i,j] = np.ma.masked vdat1 = np.ma.masked_array(vdat[t,:,:,LON_INDEX]) - print "DATA", data1.shape, "LEV", level1.shape, "LATS", lats1.shape, "z1", z1.shape, "th1", th1.shape + print("DATA", data1.shape, "LEV", level1.shape, "LATS", lats1.shape, "z1", z1.shape, "th1", th1.shape) xtmp = createTicks(lats, 10) xts = createCLabel(xtmp, 'lat') @@ -338,11 +338,11 @@ if __name__ == "__main__": plt.subplots_adjust(bottom=.05, hspace=.3, left=.1, right=1) - print "Saving images", step, STNAME, "..." #date.crsc.station.XX.gif + print("Saving images", step, STNAME, "..." #date.crsc.station.XX.gif) #savefig(POUT + "%s%02d%02d12.crsc.%s.%s.png" % (TIME.year, TIME.month, TIME.day, STNAME, step)) - NFILE1 = POUT + "crsc.%s.%s.png" % (STNAME, step) - NFILE2 = POUT + "%s%02d%02d12.crsc.%s.%s.png" % (TIME_INIT.year, TIME_INIT.month, TIME_INIT.day, STNAME, step) + NFILE1 = "%scrsc.%s.%s.png" % POUT, STNAME, step + NFILE2 = "%s%s%02d%02d12.crsc.%s.%s.png" % POUT, TIME_INIT.year, TIME_INIT.month, TIME_INIT.day, STNAME, step plt.savefig(NFILE2) shutil.copy(NFILE2, NFILE1) diff --git a/mapgenerator/utils/iplot.py b/mapgenerator/utils/iplot.py index 1120022be7b31edabca2dc03aac4476f58b15270..bc5424ccef2e272eded67b5b2ac7a1044f812e51 100644 --- a/mapgenerator/utils/iplot.py +++ b/mapgenerator/utils/iplot.py @@ -1,11 +1,12 @@ -import plotly.plotly as py -import plotly.tools as tls +import chart_studio.plotly as py from plotly.graph_objs import Figure, Scatter, Contour, Data, Line +import numpy as np + trace1 = Contour( z=air, x=lon, - y=lat, + y=at, colorscale="RdBu", zauto=False, # custom contour levels zmin=-5, # first contour level @@ -27,7 +28,7 @@ def polygons_to_traces(poly_paths, N_poly): pos arg 1. (poly_paths): paths to polygons pos arg 2. (N_poly): number of polygon to convert ''' - traces = [] # init. plotting list + traces = [] # init. plotting list for i_poly in range(N_poly): poly_path = poly_paths[i_poly] diff --git a/mapgenerator/utils/map_generator.py b/mapgenerator/utils/map_generator.py index e6d5aa2e098ec04b19133faa5d6314cad7158871..863a9b22f4b1912d8aa2189e1f280e703a37ca04 100644 --- a/mapgenerator/utils/map_generator.py +++ b/mapgenerator/utils/map_generator.py @@ -65,10 +65,10 @@ class MapGenerator: """ create color map """ self.cmap = mpl.colors.ListedColormap(self.colors) if self.over: - print "Setting over to ", self.over + print("Setting over to ", self.over) self.cmap.set_over(color=self.over) if self.under: - print "Setting under to ", self.under + print("Setting under to ", self.under) self.cmap.set_under(color=self.under) self.cmap.set_bad('red') @@ -98,8 +98,8 @@ class MapGenerator: try: pTitle=unicode(title.decode('utf-8')) % titleDic #("%sh" % sTime, valid) - except Exception, e: - print e + except Exception as e: + print(e) pTitle=title % titleDic #("%sh" % sTime, valid) return pTitle @@ -109,27 +109,28 @@ class MapGenerator: try: start_t except NameError: - global start_t, last_t + global start_t + global last_t start_t = datetime.now() last_t = start_t now_t = datetime.now() diff_t = now_t - last_t if (msg != None): - print 'TIME:', msg, ' done in ', diff_t.seconds, ' s' + print('TIME:', msg, ' done in ', diff_t.seconds, ' s') last_t = now_t def runCommand(self, commStr, fatal=True): st, out = commands.getstatusoutput(commStr) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) if (fatal == True): sys.exit(1) def genAnim(self, indir, inpattern, outdir, outfile): ''' Generate an animation starting from a set of images, specified by the inpattern parameter ''' - print "Animation: ", indir, inpattern, outdir, outfile + print("Animation: ", indir, inpattern, outdir, outfile) # Create animation #self.runCommand("convert -delay 75 -loop 0 -crop 800x600+0+0 -layers Optimize %s/%s %s/%s" % (indir, inpattern, outdir, outfile)) self.runCommand("convert -delay 75 -loop 0 -layers Optimize %s/%s %s/%s" % (indir, inpattern, outdir, outfile)) @@ -293,9 +294,9 @@ class MapGenerator: x, y = self.map(*meshgrid(g.lon,g.lat)) else: x, y = meshgrid(g.lon,g.lat) - print "X", len(x) - print "Y", len(y) - print "mapData", mapData.shape + print("X", len(x)) + print("Y", len(y)) + print("mapData", mapData.shape) self.printTime("meshgrid") @@ -305,7 +306,7 @@ class MapGenerator: # overwrite option if os.path.exists(figName + ".png") and not OVERWRITE: - print figName, " already exists." + print(figName, " already exists.") clf() return figName @@ -313,7 +314,7 @@ class MapGenerator: if self.DRAW_OPTS.nodraw and \ os.path.exists("%s-%s/%s.png" % (runDate, varName, figName)) and \ not self.DRAW_OPTS.overwrite: - print figName, " already exists." + print(figName, " already exists.") clf() return figName @@ -322,7 +323,7 @@ class MapGenerator: # Draw filled contour if not self.DRAW_OPTS.keep_aspect: - print "Not keeping original aspect" + print("Not keeping original aspect") ysize_norm = self.DRAW_OPTS.ysize*0.8 xsize_norm = self.DRAW_OPTS.xsize*0.8 @@ -337,9 +338,9 @@ class MapGenerator: ystart_norm = (0.8 - ysize_norm) / 2 + ystart_base_norm xstart_norm = (0.8 - xsize_norm) / 2 + xstart_base_norm ax_map = axes([xstart_norm, ystart_norm, xsize_norm, ysize_norm]) - print "Switching AXES for MAP:", ax_map + print("Switching AXES for MAP:", ax_map) - print "******** INTERPOLATION **********:", self.DRAW_OPTS.nointerp + print("******** INTERPOLATION **********:", self.DRAW_OPTS.nointerp) if not self.DRAW_OPTS.nocontourf and not self.DRAW_OPTS.nomap and not self.DRAW_OPTS.nointerp: mco = self.map.contourf(x, y, mapData, cmap=self.cmap, norm=self.norm, levels=self.bounds, extend='both', horizontalalignment='center') @@ -349,28 +350,28 @@ class MapGenerator: mco = contourf(x, y, mapData, cmap=self.cmap, norm=self.norm, levels=self.bounds, extend='both', horizontalalignment='center') - # Draw shape files + # Draw shape files if self.DRAW_OPTS.hasShapeFiles(): line_w = self.DRAW_OPTS.shapef_width_init for shapef in self.DRAW_OPTS.shapefiles: - print "Processing shape file:", shapef, " with line width:", line_w + print("Processing shape file:", shapef, " with line width:", line_w) self.map.readshapefile(shapef, "%s" % os.path.basename(shapef), linewidth=line_w) line_w = max(self.DRAW_OPTS.shapef_width_step, line_w - self.DRAW_OPTS.shapef_width_step) self.printTime("contourf") - ## FIXME Modify to use scatter inside DATA + ## FIXME Modify to use scatter inside DATA ## if DATA.hasScatterData()... etc... if kwargs.has_key('scatter'): SCATTERD = kwargs['scatter'] if (SCATTERD != None): - print "Plotting scatter data", SCATTERD[0], SCATTERD[1] + print("Plotting scatter data", SCATTERD[0], SCATTERD[1]) self.map.scatter(SCATTERD[0], SCATTERD[1], s=SCATTERD[2], c=SCATTERD[3], marker='o', linewidth=0.5) if DATA.hasWindData(): WINDS = DATA.getWindData() X,Y,U,V = delete_masked_points(x.ravel(), y.ravel(), WINDS['u'].ravel(), WINDS['v'].ravel()) - print "Wind scale is", self.DRAW_OPTS.wind_units, self.DRAW_OPTS.wind_scale + print("Wind scale is", self.DRAW_OPTS.wind_units, self.DRAW_OPTS.wind_scale) if WINDS.has_key('barbs'): B = self.map.barbs(X, Y, U, V, #units=self.DRAW_OPTS.wind_units, #headlength=self.DRAW_OPTS.wind_head_length, @@ -420,15 +421,15 @@ class MapGenerator: if(len(lvls) > 0): mco = self.map.contourf(x, y, mapData, cmap=self.cmap, norm=self.norm, levels=self.bounds, extend='both', horizontalalignment='center') - print "-----------------", mapData, CDATA, "--------------------" + print("-----------------", mapData, CDATA, "--------------------") if (mapData == CDATA).all(): - print ":::::::::::::::::::::::: SAME !!! :::::::::::::::::::" + print(":::::::::::::::::::::::: SAME !!! :::::::::::::::::::") CS = contour(x, y, CDATA, levels=self.bounds, colors=self.DRAW_OPTS.contours_color, linewidths=self.DRAW_OPTS.contours_linewidth) else: - print ":::::::::::::::::::::::: DIFFERENT !!! :::::::::::::::::::" - print "MIN:", cMin, "MAX:", cMax + print(":::::::::::::::::::::::: DIFFERENT !!! :::::::::::::::::::") + print("MIN:", cMin, "MAX:", cMax) CS = contour(x, y, CDATA, levels=lvls, colors=self.DRAW_OPTS.contours_color, linewidths=self.DRAW_OPTS.contours_linewidth) @@ -518,7 +519,7 @@ class MapGenerator: self.first_image = False #convert from png to gif - print "printing ", figName + print("printing ", figName) # self.runCommand("convert -quality 80 +dither -remap %s %s.png %s.anim.png" % (self.mapName, figName, figName)) self.runCommand("convert +dither %s.png %s.gif" % (figName, figName)) self.runCommand("rm %s.png" % figName) @@ -582,5 +583,4 @@ class MapGenerator: nimg1 = tpl % {'date':date, 'step':nimg1} st, out = commands.getstatusoutput("convert -delay 75 -loop 0 %s.gif %s.gif" % (nimg0, nimg1)) if st != 0: - print "Error: %s" % str(out) - + print("Error: %s" % str(out)) diff --git a/mapgenerator/utils/old_map_utils.py b/mapgenerator/utils/old_map_utils.py index ee4c04c047d57697171d1022b63cd03bc393786d..cd6ddb69bad823cfa7b01302afb84f71c41cc5da 100644 --- a/mapgenerator/utils/old_map_utils.py +++ b/mapgenerator/utils/old_map_utils.py @@ -40,7 +40,7 @@ class MapGenerator: self.map = None #Map self.lat = None #Latitudes self.lon = None #Longitudes - self.res = None #Resolution + self.res = None #Resolution def setNorm(self): @@ -184,25 +184,26 @@ class MapGenerator: def setTitle(self, title): """ """ pass - + def printTime(self, msg=None): - # initial time + # initial time try: start_t except NameError: - global start_t, last_t + global start_t + global last_t start_t = datetime.now() last_t = start_t - + now_t = datetime.now() diff_t = now_t - last_t if (msg != None): - print 'TIME:', msg, ' done in ', diff_t.seconds, ' s' + print('TIME:', msg, ' done in ', diff_t.seconds, ' s') last_t = now_t def initMap(self, LAT, LON, RESOLUTION): """ initialize map """ - + self.lat = LAT self.lon = LON self.res = RESOLUTION @@ -233,23 +234,23 @@ class MapGenerator: g = DATA.grid x, y = self.map(*meshgrid(g.lon,g.lat)) - + self.printTime("meshgrid") - print "SHAPES: ", shape(DATA.grid.lon), shape(DATA.grid.lat), shape(x), shape(y) + print("SHAPES: ", shape(DATA.grid.lon), shape(DATA.grid.lat), shape(x), shape(y)) cm = mpl.colors.ListedColormap(COLORS) n = mpl.colors.BoundaryNorm(self.bounds, len(COLORS)) cm.set_under(color=UNDER) cm.set_over(color=OVER) - + self.map.contourf(x, y, DATA, cmap=cm, norm=n, levels=self.bounds, extend='both') self.printTime("contourf") if kwargs.has_key('winds'): WINDS = kwargs['winds'] - print WINDS + print(WINDS) if (WINDS != None): X,Y,U,V = delete_masked_points(x.ravel(), y.ravel(), WINDS['u'].ravel(), WINDS['v'].ravel()) Q = self.map.quiver(X, Y, U, V, units='width', headlength=2, width=0.0015, scale=400) @@ -306,16 +307,16 @@ class MapGenerator: arrow(0.83, pos + add, 0.15, 0, length_includes_head=True, head_length=0.15, head_width=0.025, fill=True, color='k') self.printTime("limits") - print 'FIGNAME: ', figName, ' ', type(figName) + print('FIGNAME: ', figName, ' ', type(figName)) savefig(figName + ".png") self.printTime("saving") clf() #convert from png to gif and delete png - print "printing ", figName, " ..." + print("printing ", figName, " ...") st, out = commands.getstatusoutput("convert %s.png %s.gif && rm %s.png" % (figName,figName,figName)) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) sys.exit(1) self.printTime("convert") @@ -358,15 +359,15 @@ class MapGenerator: newDate1 = time.strftime('%Y%m%d', time.strptime(runDate1,"%d%b%Y")) if (runDate0 != runDate1): - print "Error: different date, comparison impossible." - print "MOD0", mod0, " - MOD1", mod1 - print "DAT0", runDate0, " - DAT1", runDate1 + print("Error: different date, comparison impossible.") + print("MOD0", mod0, " - MOD1", mod1) + print("DAT0", runDate0, " - DAT1", runDate1) continue if (var0 != var1): - print "Error: different variable, comparison impossible." - print "MOD0", mod0, " - MOD1", mod1 - print "VAR0", var0, " - VAR1", var1 + print("Error: different variable, comparison impossible.") + print("MOD0", mod0, " - MOD1", mod1) + print("VAR0", var0, " - VAR1", var1) continue #print "****", runHour0, hOffset0, "****", runHour1, hOffset1, "****" @@ -378,7 +379,7 @@ class MapGenerator: (app, item0, item1, newMapName)) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) else: #print "==>", newMapName, "<==" retNames.append(newMapName) @@ -426,13 +427,13 @@ class MapGenerator: #rename from 00-06-... to 00-01 ... st, out = commands.getstatusoutput("ls %s*" % newMapNameTpl) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) else: names = out.split('\n') #print names i = 0 for num in range(0, len(names)): - print num, i + print(num, i) snum = "%02d" % i #print "mv %s %s.gif" % (names[num], newMapNameTpl.replace("??", snum)) st2, out2 = commands.getstatusoutput("mv -f %s %s.gif" % (names[num], newMapNameTpl.replace("??", snum))) @@ -440,13 +441,12 @@ class MapGenerator: break i += int(freq) if st2 != 0: - print "Error: %s" % str(out2) - + print("Error: %s" % str(out2)) newMapNameTpl2 = newMapNameTpl.replace("??", "loop") st, out = commands.getstatusoutput("convert -delay 75 -loop 0 %s.gif %s.gif" % \ (newMapNameTpl, newMapNameTpl2)) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) def run(): @@ -461,8 +461,8 @@ def run(): (options, args) = parser.parse_args() - print options - print args + print(options) + print(args) if len(args) < 1: parser.error("wrong number of arguments") @@ -488,12 +488,12 @@ def run(): os.chdir(curdir) ga.cmd("reinit") - print '@@@@@', fName, '; ', cFile + print('@@@@@', fName, '; ', cFile) OPTS = readConf(fName, cFile) - print OPTS + print(OPTS) if not OPTS: - print "Section '", fName, "' doesn't exist." + print("Section '", fName, "' doesn't exist.") continue INDIR = options.indir or OPTS['indir'] @@ -519,7 +519,7 @@ def run(): MAX = options.max or OPTS['max'] ANIM = options.anim or OPTS['anim'] - print "OPTIONS:", """ + print("OPTIONS:", """ INDIR = %s OUTDIR = %s SOURCES = %s @@ -543,13 +543,13 @@ def run(): MAX = %s """ % (INDIR, OUTDIR, SOURCES, SRCVAR, SRCGAP, BOUNDS, BOUNDARIES, COLORS, TOTAL, INTERVAL, FREQ, TITLE, LAT, LON, OVER, UNDER, LIMITS, WIND, -WINDOPTS, RESOLUTION, MAX) +WINDOPTS, RESOLUTION, MAX)) - # Open SOURCES + # Open SOURCES for f in SOURCES: ga.open(INDIR + '/' + f) - # Open WINDS if specified + # Open WINDS if specified if WIND: ga.open(INDIR + '/' + WIND) os.chdir(OUTDIR) @@ -567,27 +567,27 @@ WINDOPTS, RESOLUTION, MAX) valid_tmp = dims.time[1] valid = "%sz %s %s %s" % (valid_tmp[:-10], valid_tmp[-9:-7], valid_tmp[-7:-4], valid_tmp[-4:]) - + mg.initMap(LAT,LON,RESOLUTION) - + for nTime in range(0,TOTAL,INTERVAL): DATA = None #GET DATA from SOURCES - print "SOURCES: ", SOURCES, type(SOURCES), len(SOURCES) + print("SOURCES: ", SOURCES, type(SOURCES), len(SOURCES)) for src_n in range(0, len(SOURCES)): - print ">>>>> CURRENT SRC " , src_n, ':', SOURCES[src_n] - ga.cmd("set dfile %d" % (src_n + 1)) + print(">>>>> CURRENT SRC " , src_n, ':', SOURCES[src_n]) + ga.cmd("set dfile %d" % (src_n + 1)) ga.cmd("set lon %s %s" % (str(LON[0]), str(LON[1]))) ga.cmd("set lat %s %s" % (str(LAT[0]), str(LAT[1]))) ga.cmd("set t %s" % (nTime + int(SRCGAP[src_n]))) - if DATA != None: - print "DATALON: ", shape(DATA.data), " SHAPE LON: ", shape(DATA.grid.lon), " SHAPE LAT: ", shape(DATA.grid.lat) + if DATA != None: + print("DATALON: ", shape(DATA.data), " SHAPE LON: ", shape(DATA.grid.lon), " SHAPE LAT: ", shape(DATA.grid.lat)) lons, lats = meshgrid(DATA.grid.lon, DATA.grid.lat) - print "SHAPE lats/lons:", shape(lats), shape(lons), "lenghts: ", len(DATA.grid.lat), len(DATA.grid.lon) + print("SHAPE lats/lons:", shape(lats), shape(lons), "lenghts: ", len(DATA.grid.lat), len(DATA.grid.lon)) VVV, AAA = ga.interp(SRCVAR[src_n], lons, lats) #print "SHAPE: ", shape(VVV), type(VVV), shape(DATA.data), type(DATA.data) #print "DATA: ", VVV.data[0:10] - print shape(VVV), shape(AAA), type(VVV) + print(shape(VVV), shape(AAA), type(VVV)) DATA += VVV.reshape(shape(DATA)) else: #VVV = ga.exp(SRCVAR[src_n]) @@ -596,11 +596,11 @@ WINDOPTS, RESOLUTION, MAX) DATA = ga.exp(SRCVAR[src_n]) #print "DATA %d: " % (src_n), DATA.grid.lon[0:10], type(DATA), shape(DATA) - print "DATA: ", type(DATA), shape(DATA) - print "DATALON: ", shape(DATA.data), " SHAPE LON: ", shape(DATA.grid.lon), " SHAPE LAT: ", shape(DATA.grid.lat) - - - #Get wind data + print("DATA: ", type(DATA), shape(DATA)) + print("DATALON: ", shape(DATA.data), " SHAPE LON: ", shape(DATA.grid.lon), " SHAPE LAT: ", shape(DATA.grid.lat)) + + + #Get wind data WINDS=None if WIND: WINDS = {'u': None, 'v': None} @@ -615,11 +615,11 @@ WINDOPTS, RESOLUTION, MAX) sTime = "%02d" % pTime figName = "%s-%s-%s" % (fName, sTime, run_tmp) - print 'title: ' , TITLE, ' ', type(TITLE) + print('title: ' , TITLE, ' ', type(TITLE)) try: pTitle=unicode(TITLE.decode('utf-8')) % (sTime, valid) - except Exception, e: - print e + except Exception as e: + print(e) pTitle=TITLE % (sTime, valid) figNames.append(mg.genImageMap(ga, figName, valid, DATA, pTitle, COLORS, @@ -627,10 +627,10 @@ WINDOPTS, RESOLUTION, MAX) ## Calculate matrix of maxs if mx != None: mx = maximum(DATA, mx) - print "M [0, 0]: ", mx.data[0][0], DATA.data[0][0] + print("M [0, 0]: ", mx.data[0][0], DATA.data[0][0]) else: mx = DATA - print "E [0, 0]: ", mx.data[0][0] + print("E [0, 0]: ", mx.data[0][0]) # Create Max at requires intervals if MAX and nTime in MAX: @@ -647,7 +647,7 @@ WINDOPTS, RESOLUTION, MAX) if ANIM: st, out = commands.getstatusoutput("convert -delay 75 -loop 0 %s-%s-*-%s.gif %s-%s-%s-loop.gif" % (figName, varName, run_tmp, runDate, figName, varName)) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) sys.exit(1) mapNames.append(figNames) diff --git a/mapgenerator/utils/old_sds_map_generator.py b/mapgenerator/utils/old_sds_map_generator.py index 1278ca945c1266a1dab97cc16440f283c339e24d..e94baacb36dee8ef052f693f5e80e07a53423841 100755 --- a/mapgenerator/utils/old_sds_map_generator.py +++ b/mapgenerator/utils/old_sds_map_generator.py @@ -40,7 +40,7 @@ class MapGenerator: self.map = None #Map self.lat = None #Latitudes self.lon = None #Longitudes - self.res = None #Resolution + self.res = None #Resolution def setNorm(self): @@ -184,25 +184,25 @@ class MapGenerator: def setTitle(self, title): """ """ pass - + def printTime(self, msg=None): - # initial time + # initial time try: start_t except NameError: global start_t, last_t start_t = datetime.now() last_t = start_t - + now_t = datetime.now() diff_t = now_t - last_t if (msg != None): - print 'TIME:', msg, ' done in ', diff_t.seconds, ' s' + print('TIME:', msg, ' done in ', diff_t.seconds, ' s') last_t = now_t def initMap(self, LAT, LON, RESOLUTION): """ initialize map """ - + self.lat = LAT self.lon = LON self.res = RESOLUTION @@ -227,7 +227,7 @@ class MapGenerator: dims = ga.query("dims", Quiet=True) - print "nTime", nTime + print("nTime", nTime) if nTime == 0: global run_tmp global run @@ -238,7 +238,7 @@ class MapGenerator: else: runt = run_tmp[:-10] run = "%sh %s %s %s" % (runt, run_tmp[-9:-7], run_tmp[-7:-4], run_tmp[-4:]) - print "run", run + print("run", run) valid_tmp = dims.time[1] valid = "%sh %s %s %s" % (valid_tmp[:-10], valid_tmp[-9:-7], valid_tmp[-7:-4], valid_tmp[-4:]) @@ -288,12 +288,12 @@ class MapGenerator: figName = "%s-%s-%s-%s" % (runDate, modName, varName, sTime) if os.path.exists(figName + ".gif") and not OVERWRITE: - print figName, " already exists." + print(figName, " already exists.") clf() return figName if NODRAW and os.path.exists("%s-%s/%s.png" % (runDate, varName, figName)) and not OVERWRITE: - print figName, " already exists." + print(figName, " already exists.") clf() return figName @@ -322,12 +322,12 @@ class MapGenerator: ) # Z1 = dot(ga.exp(varName), multiply) -# g = Z1.grid +# g = Z1.grid # x, y = ga.map(*meshgrid(g.lon,g.lat)) # # ## Calculate matrix of maxs # if max != None: -# max = maximum(Z1, max) +# max = maximum(Z1, max) # print "[0, 0]: ", max.data[0][0], Z1.data[0][0] # else: # max = Z1 @@ -396,12 +396,12 @@ class MapGenerator: # a uint8 array between 0-255 nim = np.array(im).astype(np.float) / 255 - # With newer (1.0) versions of matplotlib, you can + # With newer (1.0) versions of matplotlib, you can # use the "zorder" kwarg to make the image overlay # the plot, rather than hide behind it... (e.g. zorder=10) figimage(nim, self.logo[1], self.logo[2], zorder=10) - print "printing ", figName, " ..." + print("printing ", figName, " ...") if NODRAW: savefig(figName + ".png", bbox_inches='tight', pad_inches=0, dpi=DPI, transparent=True) else: @@ -419,7 +419,7 @@ class MapGenerator: st, out = commands.getstatusoutput("convert %s.png %s.gif && rm %s.png" % (figName,figName,figName)) st, out = commands.getstatusoutput("convert %s.png %s.gif && rm %s.png" % (figName2,figName2,figName2)) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) sys.exit(1) clf() @@ -464,15 +464,15 @@ class MapGenerator: newDate1 = time.strftime('%Y%m%d', time.strptime(runDate1,"%d%b%Y")) if (runDate0 != runDate1): - print "Error: different date, comparison impossible." - print "MOD0", mod0, " - MOD1", mod1 - print "DAT0", runDate0, " - DAT1", runDate1 + print("Error: different date, comparison impossible.") + print("MOD0", mod0, " - MOD1", mod1) + print("DAT0", runDate0, " - DAT1", runDate1) continue if (var0 != var1): - print "Error: different variable, comparison impossible." - print "MOD0", mod0, " - MOD1", mod1 - print "VAR0", var0, " - VAR1", var1 + print("Error: different variable, comparison impossible.") + print("MOD0", mod0, " - MOD1", mod1) + print("VAR0", var0, " - VAR1", var1) continue #print "****", runHour0, hOffset0, "****", runHour1, hOffset1, "****" @@ -546,13 +546,13 @@ class MapGenerator: #rename from 00-06-... to 00-01 ... st, out = commands.getstatusoutput("ls %s*" % newMapNameTpl) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) else: names = out.split('\n') #print names i = 0 for num in range(0, len(names)): - print num, i + print(num, i) snum = "%02d" % i #print "mv %s %s.gif" % (names[num], newMapNameTpl.replace("??", snum)) st2, out2 = commands.getstatusoutput("mv -f %s %s.gif" % (names[num], newMapNameTpl.replace("??", snum))) @@ -560,13 +560,13 @@ class MapGenerator: break i += int(freq) if st2 != 0: - print "Error: %s" % str(out2) + print("Error: %s" % str(out2)) newMapNameTpl2 = newMapNameTpl.replace("??", "loop") st, out = commands.getstatusoutput("convert -delay 75 -loop 0 %s.gif %s.gif" % \ (newMapNameTpl, newMapNameTpl2)) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) def run(): @@ -581,7 +581,7 @@ def run(): Usage: %s """ % sys.argv[0] if len(sys.argv) not in (3,4): - print error_msg + print(error_msg) sys.exit(1) global run @@ -613,13 +613,13 @@ def run(): #os.chdir(curdir) if not os.path.exists(INDIR + '/' + fName): - print "File '", fName, "' doesn't exist." + print("File '", fName, "' doesn't exist.") continue ga.cmd("reinit") ga.open(fName) - ext = fName.split('.')[-1:][0] modName = '_'.join(fName.split('_')[1:]).replace('.' + ext,'') + ext = fName.split('.')[-1:][0] modName = '_'.join(fName.split('_')[1:]).replace('.' + ext,'') runDate = fName.split('_')[0] mg.rundate = runDate @@ -738,11 +738,11 @@ def run(): fig.savefig("%s-%s/%s-colorbar.png" % (runDate, varName, runDate), bbox_inches='tight', pad_inches=0, dpi=DPI, transparent=True) #generate KMZ - Offline - print "Generating KMZ ..." + print("Generating KMZ ...") mg.genKML(LON, LAT, figNames, runDate, online=False) #generate KML - Online - print "Generating KML ..." + print("Generating KML ...") mg.genKML(LON, LAT, figNames, runDate) #generate animation. TODO: use imagemagick libs @@ -750,7 +750,7 @@ def run(): st, out = commands.getstatusoutput("convert -delay 75 -loop 0 %s-%s-*-%s.gif %s-%s-%s-loop.gif" % (modName, varName, run_tmp, runDate, modName, varName)) st, out = commands.getstatusoutput("convert -delay 75 -loop 0 %s-%s-%s-??.gif %s-%s-%s--loop-.gif" % (runDate, modName, varName, runDate, modName, varName)) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) sys.exit(1) mapNames.append(figNames) diff --git a/mapgenerator/utils/tseries.py b/mapgenerator/utils/tseries.py index b75a285e245bf62309465c6602d323916dd8ec07..30d04f5508ac6e1acff0daa392e9838092a795d1 100644 --- a/mapgenerator/utils/tseries.py +++ b/mapgenerator/utils/tseries.py @@ -4,25 +4,25 @@ ## author: Francesco Benincasa # # -#import matplotlib -#matplotlib.use('Agg') +# import matplotlib +# matplotlib.use('Agg') # -#import matplotlib.pyplot as plt -#from pylab import plot -#from pylab import show -#import numpy as np -#from xml.dom.minidom import parseString -#from datetime import datetime -#from datetime import date -#from datetime import timedelta -#from calendar import monthrange -#import os.path +# import matplotlib.pyplot as plt +# from pylab import plot +# from pylab import show +# import numpy as np +# from xml.dom.minidom import parseString +# from datetime import datetime +# from datetime import date +# from datetime import timedelta +# from calendar import monthrange +# import os.path # -#from random import randrange +# from random import randrange # -#global FIGURE +# global FIGURE # -#def do_avg(l, k): +# def do_avg(l, k): # #print "LIST", l # #print "K", k # @@ -62,40 +62,40 @@ # return avg_p, avg_n # # -#class DataXMLParser: +# class DataXMLParser: # """ """ -# + # def __init__(self): # """ """ # pass -# + # def getCountryCode(self, dom): # """ """ # xmlTags = dom.getElementsByTagName('country_isocode') # elem = xmlTags[0] # return elem.firstChild.data -# + # def getStation(self, dom, code): # """ """ # station_name = '' # xmlTags = dom.getElementsByTagName('station') -# + # for elem in xmlTags: # for child in elem.childNodes: -## if child.tagName == 'name': -## station_name = child.firstChild.data +# # if child.tagName == 'name': +# # station_name = child.firstChild.data # if child.tagName == 'code' \ # and child.firstChild.data != code: # return None # return elem -# + # def getComponent(self, elem, name): # """ """ # meas = elem.getElementsByTagName('measurement') # result = {} # dtime = '' # val = '' -# + # for items in meas: # flag = False # for child in items.childNodes: @@ -109,70 +109,70 @@ # val = child.firstChild.data # if dtime and val: # result[dtime] = val -# + # return result -# + # def retrieveXMLData(self, filename, code): # f = open(filename) # data = f.read() # f.close() -# + # dom = parseString(data) -# + # datap = DataXMLParser() # elem = datap.getStation(dom, code) # pm10 = datap.getComponent(elem, 'PM10') # pm2_5 = datap.getComponent(elem, 'PM2_5') # cc = datap.getCountryCode(dom) -# + # #st_name = "%s, %s" % (st_name, cc) -# + # return pm10, pm2_5 -# -# -#class TextParser2: + + +# class TextParser2: # """ """ -# + # def __init__(self): # """ """ # pass -# + # def retrieveTXTData(self, filename, code, model=None, start_hour=0): # """ """ -# + # nomodels = ("AERONET", "AE1", "AE2", "VISIBILITY", "HAZE") -# + # f = open(filename) # data = f.read() # f.close() -# + # d1 = {} # lines = data.split('\n') -# + # for line in lines: -# + # if not line: # continue -# + # #data = line.split('\t') # #model station date hour value # data = line.split() -# + # mod = data[0] # if model != mod: # continue -# + # cod = data[1] # if code != cod: # continue -# + # fndat = os.path.basename(filename).split('.')[0] # if data[2] < fndat: # continue -# + # dat = datetime.strptime(data[2], "%Y%m%d").strftime("%Y-%m-%d") # hour = int(float(data[3])) -# + # # AERONET, AE1, AE2 # mdate = "%s %02d:00" % (dat, hour) # if model in nomodels: @@ -184,11 +184,11 @@ # else: # d1[mdate] = float(data[4]) # continue -# + # hour = start_hour + hour # if hour < START_H and model not in nomodels: # continue -# + # edelta = timedelta(days=1) # edat = (datetime.strptime(data[2], "%Y%m%d") + edelta).strftime("%Y-%m-%d") # if hour > 36: @@ -197,61 +197,61 @@ # delta = timedelta(days=hour/NHOURS) # dat = (datetime.strptime(data[2], "%Y%m%d") + delta).strftime("%Y-%m-%d") # hour = hour % NHOURS -# + # mdate = "%s %02d:00" % (dat, hour) # d1[mdate] = float(data[4]) -# + # #print "D1", d1 # return d1 -# -# -#class TextParser: + + +# class TextParser: # """ """ -# + # def __init__(self): # """ """ # pass -# + # def retrieveTXTData(self, filename, code): # """ """ # f = open(filename) # data = f.read() # f.close() -# + # d1 = {} # d2 = {} # lines = data.split('\n') -# + # for line in lines: -# + # if not line: # continue -# + # #data = line.split('\t') # data = line.split() # cod = data[0] -# + # if code != cod: # continue -# + # dat = datetime.strptime(data[1], "%Y%m%d").strftime("%Y-%m-%d") # mdate = "%s %02d:00" % (dat, int(data[2])) # pm10 = data[3] # pm2_5 = data[4] -# + # d1[mdate] = pm10 # d2[mdate] = pm2_5 -# + # return d1, d2 -# -# -#class PlotGenerator: -# + + +# class PlotGenerator: + # def __init__(self): # """ """ # plt.clf() # plt.cla() -# + # def setVar(self, xvals, var, mon, force_continue=False): # """ """ # hrs = var.keys() @@ -268,10 +268,10 @@ # st = ((MR+day)*NH)+(hour/INTERV) # #print "DAT", dat, "ST", st, "V", var[dat] # steps[st] = var[dat] -# + # res = [] # xvtmp = [] -# + # for k in range(NH,(MR+1)*NH+4): #xvals: # val = NULLVAL # if k in steps.keys(): @@ -283,42 +283,42 @@ # else: # if float(val) < 0: # val = NULLVAL -# + # if type(val) is list: # xvtmp.extend([k for i in range(0, len(val))]) # res.extend(val) # else: # xvtmp.append(k) # res.append(val) -# + # res = np.array(res) -# + # return xvtmp, res+.005 #_masked -# -# + + # def generatePlot(self, varN, title, fname, month, MR): # """ """ -# + # f = plt.gcf() # DefaultSize = f.get_size_inches() -# + # XTICKS_END = (MR+1)*NH + 4 -# -## XTICKS_LABELS = np.arange(1,(MR+1)) + +# # XTICKS_LABELS = np.arange(1,(MR+1)) # XTICKS_LABEL_START = XOFF / XTICKS_INT # XTICKS_LABEL_END = XTICKS_END / XTICKS_INT + 1 -# + # XTICKS_DATE_START = month + "01" # XTICKS_DATE_START_FORMAT = "%Y%m%d" # XTICKS_DATE_NUM = (XTICKS_END - XTICKS_START) / NH # XTICKS_DATE_INT = 1 -# + # dates = [datetime.strptime(XTICKS_DATE_START, XTICKS_DATE_START_FORMAT) + timedelta(days=i) for i in range(XTICKS_DATE_NUM)] -# + # XTL = [datetime.strftime(d, "%d") for d in dates] # XTICKS_LABELS = [(int(d) % 2) and d or '' for d in XTL] # XTICKS_LABELS.append('01') -# + # for attrs in varN: # x = attrs['xvals'] # y = attrs['var'] @@ -337,9 +337,9 @@ # except Exception, e: # print "ERROR", e # continue -# + # ax = plt.axis([XMIN,XTICKS_END,YMIN,YMAX]) -# + # ll = plt.legend( # loc=LEGEND, # fancybox=True, @@ -350,27 +350,27 @@ # lx = plt.xlabel(XLABEL % {'month': MONTH, 'year': YEAR}) # tx = plt.xticks(np.arange(XTICKS_START, XTICKS_END, XTICKS_INT), # XTICKS_LABELS) -# + # ty = plt.yticks(YTICKS, YTICKS_LABELS) -# + # if EVAL_TYPE == 'VISIBILITY': # ly = plt.ylabel(YLABEL2) # plt.yscale('log') # else: # ly = plt.ylabel(YLABEL1) -# + # tl = plt.title(title) # gr = plt.grid(GRID) -# + # f.set_size_inches((DefaultSize[0]*1.5, DefaultSize[1]*1.5)) # f.savefig(fname, dpi=72) # #if fname.find(FIGURE) >= 0: # if fname.find("Zinder") >= 0: # f.savefig("latest.png", bbox_inches='tight', pad_inches=.2, dpi=42) # plt.close('all') -# -# -#EVAL1_STATIONS = { + + +# EVAL1_STATIONS = { # "Barcelona": u"Barcelona (Spain)", # "Capo_Verde": u"Capo_Verde (Cabo Verde)", # "Dakar": u"Dakar (Senegal)", @@ -413,185 +413,185 @@ # "KAUST_Campus": u"KAUST_Campus (Saudi Arabia)", # "CUT-TEPAK": u"CUT-TEPAK (Cyprus)", # "Kuwait_University": u"Kuwait_University (Kuwait)", -#} -# -##name : (var name, start hour, color, type, , ) -# -#EVAL1_MODELS = { -#r'01-AE AERONET $>$ 0.6': ('AE1', 0, 'w', 'o', False, 1.0, 6, None), -#r'02-AE AERONET $\leq$ 0.6': ('AE2', 0, 'k', 'o', False, 1.0, 6, 'gray'), -#r'03-AOD$_{550}$ AERONET': ('AERONET', 0, 'y', '^', False, 1.8, 8, None), -##r'04-DOD$_{550}$ NMMB/BSC-Dust': ('NMMB-BSC_OPER', 12, 'r', '-', True, 1.8, 8, None), -#r'05-DOD$_{550}$ NMMB/BSC-Dust': ('SDSWAS_NMMB-BSC-v2_OPER', 12, 'r', '-', True, 2.5, 8, None), -##r'04-DOD$_{550}$ BSC_DREAM8b': ('BSC_DREAM8b', 12, 'r', '-', True, 1.8, 8, None), -##r'05-DOD$_{550}$ MACC-ECMWF': ('MACC-ECMWF', 0, 'b', '-', True, 1.8, 8, None), -##r'06-DOD$_{550}$ DREAM8-NMME-MACC': ('DREAM8-MACC', 12, 'm', '-', True, 1.8, 8, None), -##r'07-DOD$_{550}$ CHIMERE': ('CHIMERE', 0, 'c', '-', True, 1.8, 8, None), -##r'08-DOD$_{550}$ NMMB/BSC-Dust': ('NMMB-BSC', 12, 'g', '-', True, 1.8, 8, None), -##r'09-DOD$_{550}$ U.K. MetOffice MetUM': ('_UKMET.', 0, 'DarkOrange', '-', True, 1.8, 8, None), -##r'10-DOD$_{550}$ NASA GEOS-5': ('NASA-GE', 0, 'Brown', '-', True, 1.8, 8, None), -##r'11-DOD$_{550}$ NCEP NGAC': ('NCEP-NG', 0, 'Chartreuse', '-', True, 1.8, 8, None), -##r'12-DOD$_{550}$ MEDIAN': ('_MEDIAN', 0, 'k', '--', True, 2.5, 8, None), -#} -# -#EVAL2_STATIONS = { -#"DAUI": u"IN SALAH NORTH (Algeria)", -#"DAUZ": u"IN AMENAS/ZARZAI (Algeria)", -#"DATM": u"BORDJ MOKHTAR (Algeria)", -#"GQNN": u"NOUAKCHOTT (Mauritania)", -#"GQPP": u"NOUADHIBOU (Mauritania)", -#"LCRA": u"AKROTIRI (Cyprus)", -#"GGOV": u"BISSAU (Guinea-Bissau)", -#"HLLT": u"TRIPOLI INTL ARP (Libya)", -#"HLLB": u"BENINA/BENGHAZI (Libya)", -#"GMML": u"LAAYOUNE/HASSAN (Morocco)", -#"GMFK": u"ER-RACHIDIA (Morocco)", -#"DTTR": u"EL BORMA (Tunisia)", -#"GOTT": u"TAMBACOUNDA (Senegal)", -#"GOGS": u"CAPE SKIRING (Senegal)", -#"LTCJ": u"BATMAN (Turkey)", -#"OJAM": u"AMMAN/KING ABDUL (Jordan)", -#"HDAM": u"AMBOULI (Djibouti)", -#"HECA": u"CAIRO INTL AIRPORT (Egypt)", -#"HESH": u"ST CATHERINE (Egypt)", -#"OOMS": u"SEEB INTL (Oman)", -#"DNKN": u"KANO/MALLAM (Nigeria)", -#"OIKQ": u"GHESHM/DAYRESTAN (Iran)", -#"OIAW": u"AHWAZ (Iran)", -#"OIBB": u"BUSHEHR (Iran)", -#"FTTJ": u"NDJAMENA (Chad)", -#"OERF": u"RAFHA (Saudi Arabia)", -#"OERK": u"RIYADH (Saudi Arabia)", -#"OEAH": u"AL AHSA (Saudi Arabia)", -#"OEGS": u"GASSIM (Saudi Arabia)", -#"OENG": u"NEJRAN (Saudi Arabia)", -#"KQTZ": u"BAGHDAD (Irak)", -#"OMFJ": u"FUJAIRAH (U. Arab Emirates)", -#"OYHD": u"HODEIDAH (Yemen)", -#"OTBD": u"DOHA INTL AIRPOR (Qatar)", -#"OKBK": u"KUWAIT INTL (Kuwait)", -#"OBBI": u"BAHRAIN INTL ARP (Bahrain)", -#"GABS": u"BAMAKO/SENOU (Mali)", -#"DRRN": u"NIAMEY (Niger)", -#"DRZA": u"AGADEZ SOUTH (Niger)", -#} -# -#EVAL2_MODELS = { -#r'13-SCONC ($\mu$g/m$^3$) SAND/DUST': ('VISIBILITY', 0, 'r', '^', False, 1.8, 10, None), -#r'14-SCONC ($\mu$g/m$^3$) HAZE': ('HAZE', 0, 'r', '*', False, 1.8, 10, None), -#r'04-SCONC ($\mu$g/m$^3$) BSC_DREAM8b': ('BSC_DREAM8b', 12, 'r', '-', True, 1.8, 8, None), -#r'05-SCONC ($\mu$g/m$^3$) MACC-ECMWF': ('MACC-ECMWF', 0, 'b', '-', True, 1.8, 8, None), -#r'06-SCONC ($\mu$g/m$^3$) DREAM8-NMME-MACC': ('DREAM8-MACC', 12, 'm', '-', True, 1.8, 8, None), -#r'07-SCONC ($\mu$g/m$^3$) CHIMERE': ('CHIMERE', 0, 'c', '-', True, 1.8, 8, None), -#r'08-SCONC ($\mu$g/m$^3$) NMMB/BSC-Dust': ('NMMB-BSC', 12, 'g', '-', True, 1.8, 8, None), -#r'09-SCONC ($\mu$g/m$^3$) U.K. MetOffice MetUM': ('_UKMET.', 0, 'DarkOrange', '-', True, 1.8, 8, None), -#r'10-SCONC ($\mu$g/m$^3$) NASA GEOS-5': ('NASA-GE', 0, 'Brown', '-', True, 1.8, 8, None), -#r'11-SCONC ($\mu$g/m$^3$) NCEP NGAC': ('NCEP-NG', 0, 'Chartreuse', '-', True, 1.8, 8, None), -#r'12-SCONC ($\mu$g/m$^3$) MEDIAN': ('_MEDIAN', 0, 'k', '--', True, 2.5, 8, None), -#} -# -#EVAL_TYPE = '' -# -#NULLVAL = np.nan #None -# -#MONTH = YEAR = "" -# -#XLABEL = r"days" #r"%(month)s %(year)s" -#YLABEL1 = r"Aerosol Optical Depth (AOD$_{550}$), Dust Optical Depth (DOD$_{550}$)" -#YLABEL2 = r"Dust Surface Concentration - SCONC ($\mu$g/m$^3$)" -# -#LEGEND = 'upper center' #'upper right' -#GRID = False -# -#NUMDAYS = 31 -#DAYSTART = 1 -# -#NHOURS = 24 -#INTERV = 3 -#NH = NHOURS/INTERV #8 -# -#XOFF = NH * DAYSTART -#XMAX = NH * NUMDAYS + XOFF #768 - NH*32 -#XMIN = 0 + XOFF -# -#YOFF = 0 -#YMAX = 2.0 + YOFF #00 -#YMIN = 0 + YOFF -# -#XTICKS_START = XMIN -#XTICKS_END = XMAX -#XTICKS_INT = NH -#XTICKS_LABEL_START = XOFF / XTICKS_INT -#XTICKS_LABEL_END = XMAX / XTICKS_INT + 1 -#XTICKS_LABEL_INT = 1 -# -#XTICKS_DATE_START = "20110501" -#XTICKS_DATE_START_FORMAT = "%Y%m%d" -#XTICKS_DATE_NUM = (XTICKS_END - XTICKS_START) / NH -#XTICKS_DATE_INT = 1 -# -#dates = [datetime.strptime(XTICKS_DATE_START, +# } + +# #name : (var name, start hour, color, type, , ) + +# EVAL1_MODELS = { +# r'01-AE AERONET $>$ 0.6': ('AE1', 0, 'w', 'o', False, 1.0, 6, None), +# r'02-AE AERONET $\leq$ 0.6': ('AE2', 0, 'k', 'o', False, 1.0, 6, 'gray'), +# r'03-AOD$_{550}$ AERONET': ('AERONET', 0, 'y', '^', False, 1.8, 8, None), +# #r'04-DOD$_{550}$ NMMB/BSC-Dust': ('NMMB-BSC_OPER', 12, 'r', '-', True, 1.8, 8, None), +# r'05-DOD$_{550}$ NMMB/BSC-Dust': ('SDSWAS_NMMB-BSC-v2_OPER', 12, 'r', '-', True, 2.5, 8, None), +# #r'04-DOD$_{550}$ BSC_DREAM8b': ('BSC_DREAM8b', 12, 'r', '-', True, 1.8, 8, None), +# #r'05-DOD$_{550}$ MACC-ECMWF': ('MACC-ECMWF', 0, 'b', '-', True, 1.8, 8, None), +# #r'06-DOD$_{550}$ DREAM8-NMME-MACC': ('DREAM8-MACC', 12, 'm', '-', True, 1.8, 8, None), +# #r'07-DOD$_{550}$ CHIMERE': ('CHIMERE', 0, 'c', '-', True, 1.8, 8, None), +# #r'08-DOD$_{550}$ NMMB/BSC-Dust': ('NMMB-BSC', 12, 'g', '-', True, 1.8, 8, None), +# #r'09-DOD$_{550}$ U.K. MetOffice MetUM': ('_UKMET.', 0, 'DarkOrange', '-', True, 1.8, 8, None), +# #r'10-DOD$_{550}$ NASA GEOS-5': ('NASA-GE', 0, 'Brown', '-', True, 1.8, 8, None), +# #r'11-DOD$_{550}$ NCEP NGAC': ('NCEP-NG', 0, 'Chartreuse', '-', True, 1.8, 8, None), +# #r'12-DOD$_{550}$ MEDIAN': ('_MEDIAN', 0, 'k', '--', True, 2.5, 8, None), +# } + +# EVAL2_STATIONS = { +# "DAUI": u"IN SALAH NORTH (Algeria)", +# "DAUZ": u"IN AMENAS/ZARZAI (Algeria)", +# "DATM": u"BORDJ MOKHTAR (Algeria)", +# "GQNN": u"NOUAKCHOTT (Mauritania)", +# "GQPP": u"NOUADHIBOU (Mauritania)", +# "LCRA": u"AKROTIRI (Cyprus)", +# "GGOV": u"BISSAU (Guinea-Bissau)", +# "HLLT": u"TRIPOLI INTL ARP (Libya)", +# "HLLB": u"BENINA/BENGHAZI (Libya)", +# "GMML": u"LAAYOUNE/HASSAN (Morocco)", +# "GMFK": u"ER-RACHIDIA (Morocco)", +# "DTTR": u"EL BORMA (Tunisia)", +# "GOTT": u"TAMBACOUNDA (Senegal)", +# "GOGS": u"CAPE SKIRING (Senegal)", +# "LTCJ": u"BATMAN (Turkey)", +# "OJAM": u"AMMAN/KING ABDUL (Jordan)", +# "HDAM": u"AMBOULI (Djibouti)", +# "HECA": u"CAIRO INTL AIRPORT (Egypt)", +# "HESH": u"ST CATHERINE (Egypt)", +# "OOMS": u"SEEB INTL (Oman)", +# "DNKN": u"KANO/MALLAM (Nigeria)", +# "OIKQ": u"GHESHM/DAYRESTAN (Iran)", +# "OIAW": u"AHWAZ (Iran)", +# "OIBB": u"BUSHEHR (Iran)", +# "FTTJ": u"NDJAMENA (Chad)", +# "OERF": u"RAFHA (Saudi Arabia)", +# "OERK": u"RIYADH (Saudi Arabia)", +# "OEAH": u"AL AHSA (Saudi Arabia)", +# "OEGS": u"GASSIM (Saudi Arabia)", +# "OENG": u"NEJRAN (Saudi Arabia)", +# "KQTZ": u"BAGHDAD (Irak)", +# "OMFJ": u"FUJAIRAH (U. Arab Emirates)", +# "OYHD": u"HODEIDAH (Yemen)", +# "OTBD": u"DOHA INTL AIRPOR (Qatar)", +# "OKBK": u"KUWAIT INTL (Kuwait)", +# "OBBI": u"BAHRAIN INTL ARP (Bahrain)", +# "GABS": u"BAMAKO/SENOU (Mali)", +# "DRRN": u"NIAMEY (Niger)", +# "DRZA": u"AGADEZ SOUTH (Niger)", +# } + +# EVAL2_MODELS = { +# r'13-SCONC ($\mu$g/m$^3$) SAND/DUST': ('VISIBILITY', 0, 'r', '^', False, 1.8, 10, None), +# r'14-SCONC ($\mu$g/m$^3$) HAZE': ('HAZE', 0, 'r', '*', False, 1.8, 10, None), +# r'04-SCONC ($\mu$g/m$^3$) BSC_DREAM8b': ('BSC_DREAM8b', 12, 'r', '-', True, 1.8, 8, None), +# r'05-SCONC ($\mu$g/m$^3$) MACC-ECMWF': ('MACC-ECMWF', 0, 'b', '-', True, 1.8, 8, None), +# r'06-SCONC ($\mu$g/m$^3$) DREAM8-NMME-MACC': ('DREAM8-MACC', 12, 'm', '-', True, 1.8, 8, None), +# r'07-SCONC ($\mu$g/m$^3$) CHIMERE': ('CHIMERE', 0, 'c', '-', True, 1.8, 8, None), +# r'08-SCONC ($\mu$g/m$^3$) NMMB/BSC-Dust': ('NMMB-BSC', 12, 'g', '-', True, 1.8, 8, None), +# r'09-SCONC ($\mu$g/m$^3$) U.K. MetOffice MetUM': ('_UKMET.', 0, 'DarkOrange', '-', True, 1.8, 8, None), +# r'10-SCONC ($\mu$g/m$^3$) NASA GEOS-5': ('NASA-GE', 0, 'Brown', '-', True, 1.8, 8, None), +# r'11-SCONC ($\mu$g/m$^3$) NCEP NGAC': ('NCEP-NG', 0, 'Chartreuse', '-', True, 1.8, 8, None), +# r'12-SCONC ($\mu$g/m$^3$) MEDIAN': ('_MEDIAN', 0, 'k', '--', True, 2.5, 8, None), +# } + +# EVAL_TYPE = '' + +# NULLVAL = np.nan #None + +# MONTH = YEAR = "" + +# XLABEL = r"days" #r"%(month)s %(year)s" +# YLABEL1 = r"Aerosol Optical Depth (AOD$_{550}$), Dust Optical Depth (DOD$_{550}$)" +# YLABEL2 = r"Dust Surface Concentration - SCONC ($\mu$g/m$^3$)" + +# LEGEND = 'upper center' #'upper right' +# GRID = False + +# NUMDAYS = 31 +# DAYSTART = 1 + +# NHOURS = 24 +# INTERV = 3 +# NH = NHOURS/INTERV #8 + +# XOFF = NH * DAYSTART +# XMAX = NH * NUMDAYS + XOFF #768 - NH*32 +# XMIN = 0 + XOFF + +# YOFF = 0 +# YMAX = 2.0 + YOFF #00 +# YMIN = 0 + YOFF + +# XTICKS_START = XMIN +# XTICKS_END = XMAX +# XTICKS_INT = NH +# XTICKS_LABEL_START = XOFF / XTICKS_INT +# XTICKS_LABEL_END = XMAX / XTICKS_INT + 1 +# XTICKS_LABEL_INT = 1 + +# XTICKS_DATE_START = "20110501" +# XTICKS_DATE_START_FORMAT = "%Y%m%d" +# XTICKS_DATE_NUM = (XTICKS_END - XTICKS_START) / NH +# XTICKS_DATE_INT = 1 + +# dates = [datetime.strptime(XTICKS_DATE_START, # XTICKS_DATE_START_FORMAT) + timedelta(days=i) for i in range(XTICKS_DATE_NUM)] -# -#XTL = [datetime.strftime(d, "%d") for d in dates] -#XTICKS_LABELS = [(int(d) % 2) and d or '' for d in XTL] -# -#YTICKS_START = YMIN + 0.2 -#YTICKS_END = YMAX + 0.1 -#YTICKS_INT = 0.2 #25 -#YTICKS_LABEL_START = YTICKS_START -#YTICKS_LABEL_END = YTICKS_END -#YTICKS_LABEL_INT = YTICKS_INT -# -#YTICKS = () #np.arange(YTICKS_START, YTICKS_END, YTICKS_INT) -#YTICKS_LABELS = () #np.arange(YTICKS_LABEL_START, YTICKS_LABEL_END, YTICKS_LABEL_INT) -# -#START_H = 15 -# -# -#if __name__ == "__main__": -# + +# XTL = [datetime.strftime(d, "%d") for d in dates] +# XTICKS_LABELS = [(int(d) % 2) and d or '' for d in XTL] + +# YTICKS_START = YMIN + 0.2 +# YTICKS_END = YMAX + 0.1 +# YTICKS_INT = 0.2 #25 +# YTICKS_LABEL_START = YTICKS_START +# YTICKS_LABEL_END = YTICKS_END +# YTICKS_LABEL_INT = YTICKS_INT + +# YTICKS = () #np.arange(YTICKS_START, YTICKS_END, YTICKS_INT) +# YTICKS_LABELS = () #np.arange(YTICKS_LABEL_START, YTICKS_LABEL_END, YTICKS_LABEL_INT) + +# START_H = 15 + + +# if __name__ == "__main__": + # import glob # import sys -# + # msg_usage = """ -#Error: %s -#Usage: %s -#- : where find files to parse (absolute path) -#- : a date without day in format %%Y%%m -#- [AERONET|VISIBILITY]: evaluation type -#""" -# +# Error: %s +# Usage: %s +# - : where find files to parse (absolute path) +# - : a date without day in format %%Y%%m +# - [AERONET|VISIBILITY]: evaluation type +# """ + # if len(sys.argv) != 4: # print msg_usage % ("Bad number argument", sys.argv[0]) # sys.exit(1) -# + # path = sys.argv[1] # month = sys.argv[2] # eval_type = sys.argv[3] # if eval_type not in ('AERONET', 'VISIBILITY'): # print msg_usage % ("Wrong evaluation type: " + eval_type, sys.argv[0]) # sys.exit(1) -# + # EVAL_TYPE = eval_type -# + # try: # dt = datetime.strptime(month, "%Y%m") # except Exception, e: # print msg_usage % (str(e), sys.argv[0]) # sys.exit(1) -# + # MONTH = dt.strftime("%B") # YEAR = dt.strftime("%Y") -# + # filenames = glob.glob(path + month + '*obs') -# + # MR = monthrange(dt.year, dt.month)[1] # num = (MR+1)*NH #days in one month x hours in a day # xvals = np.arange(XTICKS_START, num) # p = PlotGenerator() # tp = TextParser2() -# + # if EVAL_TYPE == 'AERONET': # STATIONS = EVAL1_STATIONS # MODELS = EVAL1_MODELS @@ -602,12 +602,12 @@ # MODELS = EVAL2_MODELS # YTICKS = (0, 5, 20, 50, 200, 500, 2000, 5000, 20000, 100000) # YTICKS_LABELS = ('0', '5', '20', '50', '2000', '5000', '20000', '100000') -# + # #print STATIONS, MODELS # #print YTICKS, YTICKS_LABELS -# + # FIGURE = STATIONS.keys()[randrange(len(STATIONS.keys()))] -# + # for code in STATIONS.keys(): # print "STAT", code # aod = {} @@ -620,7 +620,7 @@ # aod[model] = {} # res = tp.retrieveTXTData(filename, code, modname, start) # aod[model].update(res) -# + # st_name = STATIONS[code] # varz = [] # mods = aod.keys() @@ -649,10 +649,10 @@ # 'mec': mec, # } # varz.append(item) -# + # title = unicode( "%s - %s %s\n" % (st_name, MONTH, YEAR)) #, date.strftime(dt, "%B %Y"))) # figname = "%s-%s.png" % (month, code) -# + # p.generatePlot( # varz, # title, @@ -660,23 +660,23 @@ # month, # MR, # ) -# -# -#### AIR QUALITY -# -#import matplotlib -#matplotlib.use('Agg') -# -#import matplotlib.pyplot as plt -#import numpy as np -#from xml.dom.minidom import parseString -#from datetime import datetime -#from datetime import date -#from calendar import monthrange -#import os.path -# -# -#STATIONS = { + + +# ### AIR QUALITY + +# import matplotlib +# matplotlib.use('Agg') + +# import matplotlib.pyplot as plt +# import numpy as np +# from xml.dom.minidom import parseString +# from datetime import datetime +# from datetime import date +# from calendar import monthrange +# import os.path + + +# STATIONS = { # "CY0002R": u"Cyprus: Ayia Marina", # "HU0002R": u"Hungary: k-Puszta", # "SI0033A": u"Slovenia: Murska Sobota", @@ -698,52 +698,52 @@ # "FR23124": u"La Tardière: France", # "FR35006": u"Tulle - Hugo: France", # "FR02005": u"Martigues l'Ile: France", -#} -# -#LABEL1 = '' -#LABEL2 = '' -# -#XLABEL = '' -#YLABEL = '' -# -#M = 50 -#OFFS = 10 -# -# -#class DataXMLParser: +# } + +# LABEL1 = '' +# LABEL2 = '' + +# XLABEL = '' +# YLABEL = '' + +# M = 50 +# OFFS = 10 + + +# class DataXMLParser: # """ """ -# + # def __init__(self): # """ """ # pass -# + # def getCountryCode(self, dom): # """ """ # xmlTags = dom.getElementsByTagName('country_isocode') # elem = xmlTags[0] # return elem.firstChild.data -# + # def getStation(self, dom, code): # """ """ # station_name = '' # xmlTags = dom.getElementsByTagName('station') -# + # for elem in xmlTags: # for child in elem.childNodes: -## if child.tagName == 'name': -## station_name = child.firstChild.data +# # if child.tagName == 'name': +# # station_name = child.firstChild.data # if child.tagName == 'code' \ # and child.firstChild.data != code: # return None # return elem -# + # def getComponent(self, elem, name): # """ """ # meas = elem.getElementsByTagName('measurement') # result = {} # dtime = '' # val = '' -# + # for items in meas: # flag = False # for child in items.childNodes: @@ -757,74 +757,74 @@ # val = child.firstChild.data # if dtime and val: # result[dtime] = val -# + # return result -# + # def retrieveXMLData(self, filename, code): # f = open(filename) # data = f.read() # f.close() -# + # dom = parseString(data) -# + # datap = DataXMLParser() # elem = datap.getStation(dom, code) # pm10 = datap.getComponent(elem, 'PM10') # pm2_5 = datap.getComponent(elem, 'PM2_5') # cc = datap.getCountryCode(dom) -# + # #st_name = "%s, %s" % (st_name, cc) -# + # return pm10, pm2_5 -# -# -#class TextParser: + + +# class TextParser: # """ """ -# + # def __init__(self): # """ """ # pass -# + # def retrieveTXTData(self, filename, code): # """ """ # f = open(filename) # data = f.read() # f.close() -# + # d1 = {} # d2 = {} # lines = data.split('\n') -# + # for line in lines: -# + # if not line: # continue -# + # #data = line.split('\t') # data = line.split() # cod = data[0] -# + # if code != cod: # continue -# + # dat = datetime.strptime(data[1], "%Y%m%d").strftime("%Y-%m-%d") # mdate = "%s %02d:00" % (dat, int(data[2])) # pm10 = data[3] # pm2_5 = data[4] -# + # d1[mdate] = pm10 # d2[mdate] = pm2_5 -# + # return d1, d2 -# -# -#class PlotGenerator: -# + + +# class PlotGenerator: + # def __init__(self): # """ """ # plt.clf() # plt.cla() -# + # def setVar(self, var, mon, MR): # """ """ # hrs = var.keys() @@ -840,7 +840,7 @@ # hour = dt.hour # if mon == dt.strftime("%Y%m"): # steps[(day*24)+hour] = var[dat] -# + # res = [] # for k in range(24,(MR+1)*24): # val = np.nan @@ -850,12 +850,12 @@ # val = np.nan # #if val != 0: print val, k # res.append(val) -# + # return res #np.array(res) -# + # def generatePlot(self, xvals, varN, title, fname, MR): # """ """ -# + # mxs = [] # for attrs in varN: # #print xvals, attrs['var'], type(attrs['var']), max(attrs['var']) @@ -863,10 +863,10 @@ # values = [float(i) for i in attrs['var']] # #if title.find('Cairo') >= 0: print values, "MAX:", np.max(values) # mxs.append(np.nanmax(values)) -# + # f = plt.gcf() # DefaultSize = f.get_size_inches() -# + # mx = np.nanmax(mxs) # #print title, mx # if np.isnan(mx): mx = 0 @@ -879,7 +879,7 @@ # ARGY = int((MAXY/OFFS)-(MAXY/OFFS)%INTY) # if ARGY == 0: ARGY = INTY # MAXX = (MR+1)*24 #768 -# + # ax = plt.axis([24,MAXX,0,MAXY]) # ll = plt.legend(loc='upper right') # lx = plt.xlabel(XLABEL) @@ -888,47 +888,47 @@ # ty = plt.yticks(np.arange(0,MAXY+1,ARGY)) # tl = plt.title(title) # gr = plt.grid(True) -# + # f.set_size_inches((DefaultSize[0]*1.5, DefaultSize[1]*1.5)) # f.savefig(fname, dpi=72) # plt.close('all') -# -# -#if __name__ == "__main__": -# + + +# if __name__ == "__main__": + # import glob # import sys -# + # msg_usage = """ -#Error: %s -#Usage: %s -#- : where find files to parse -#- : a date without day in format %%Y%%m -#""" -# +# Error: %s +# Usage: %s +# - : where find files to parse +# - : a date without day in format %%Y%%m +# """ + # if len(sys.argv) != 3: # print msg_usage % ("bad number argument", sys.argv[0]) # sys.exit(1) -# + # path = sys.argv[1] # month = sys.argv[2] -# + # try: # dt = datetime.strptime(month, "%Y%m") # except Exception, e: # print msg_usage % (str(e), sys.argv[0]) # sys.exit(1) -# + # filenames = glob.glob(path + '*') -# + # MR = monthrange(dt.year, dt.month)[1] # num = (MR+1)*24 #days in one month x hours in a day # xvals = np.arange(24,num) # p = PlotGenerator() # tp = TextParser() -# + # for code in STATIONS.keys(): -# + # if code == "IT0009R": # LABEL1 = 'Total' # LABEL2 = 'Coarse' @@ -943,7 +943,7 @@ # YLABEL = r"mass concentration $(\mu g/m^3)$" # M = 50 # OFFS = 10 -# + # pm10 = {} # pm2_5 = {} # for filename in filenames: @@ -953,14 +953,14 @@ # d1, d2 = tp.retrieveTXTData(filename, code) # pm10.update(d1) # pm2_5.update(d2) -# + # st_name = STATIONS[code] # var1 = p.setVar(pm10, month, MR) # var2 = p.setVar(pm2_5, month, MR) -# + # title = unicode( "%s\n%s\n" % (st_name, date.strftime(dt, "%B %Y"))) # figname = "%s-%s.png" % (month, code) -# + # p.generatePlot( # xvals, # ( @@ -993,10 +993,8 @@ logging.basicConfig(level=logging.DEBUG) class Evaluator(object): - def __init__(self, str_date, end_date, obs_path, mod_path, obs_tmpl, mod_tmpl, dir_outp="./out/", rej_path=None, rej_tmpl=None, en_debug=False): """Initialize attributes""" - self.str_date = str_date # start date YYYYMMDD self.end_date = end_date # end date YYYYMMDD self.obs_path = obs_path # observation path @@ -1019,10 +1017,11 @@ class Evaluator(object): self.plt_xlim = () # plot x bounds self.plt_ylim = () # plot y bounds self.plt_titl = '' # plot title - self.plt_type = 'xy' # types: - # 'xy': model in x and observation in y - # 'ts': time series with line - # 'tss': time series with scatter + self.plt_type = 'xy' + # types: + # 'xy': model in x and observation in y + # 'ts': time series with line + # 'tss': time series with scatter self.tmp_tmpl = '/tmp/%s_obs_file.txt' self.plt_imag = "image.png" self.fil_stat = "stats.txt" @@ -1032,22 +1031,18 @@ class Evaluator(object): if not os.path.exists(dir_outp): os.makedirs(dir_outp) - def getDateDict(self, dd): - sdict = { - 'date' : dd.strftime(self.date_fmt), - 'year' : dd.strftime('%Y'), + 'date': dd.strftime(self.date_fmt), + 'year': dd.strftime('%Y'), 'month': dd.strftime('%m'), - 'day' : dd.strftime('%d'), + 'day': dd.strftime('%d'), } return sdict - def checkFileMatch(self, obs_tmp, mod_tmp): """Choose only observation files matching with model date-timestep""" - # create dictionaries with date+timestep as keys otmp = {basename(i)[:10]: i for i in obs_tmp} mtmp = {basename(j)[:10]: j for j in mod_tmp} @@ -1057,25 +1052,24 @@ class Evaluator(object): return [otmp[i] for i in intr], [mtmp[j] for j in intr] - def readObsFile(self, obs_file): """Read observation file""" + obs = pd.read_table( + obs_file, + delimiter=' ', + skiprows=1, + # skipfooter=4, + usecols=(0, 1, 2, 3), + names=('datetime', 'lon', 'lat', self.obs_name), + warn_bad_lines=False, + error_bad_lines=False, + ) - obs = pd.read_table(obs_file, - delimiter=' ', - skiprows=1, - #skipfooter=4, - usecols=(0,1,2,3), - names=('datetime', 'lon', 'lat', self.obs_name), - warn_bad_lines=False, - error_bad_lines=False, - ) - - if self.mylogger: self.mylogger.debug("PARTIAL OBS:\n%s" % (obs.head())) + if self.mylogger: + self.mylogger.debug("PARTIAL OBS:\n%s", obs.head()) # extract the date+tstep obs_date = basename(obs_file)[:10] - # extract valid values obs = obs.loc[np.where(obs['datetime'].values == int(obs_date))] @@ -1091,17 +1085,18 @@ class Evaluator(object): # steps rstep = { - '06' : -3, - '12' : -2, - '18' : -1, - '00' : 0, # next day, so must go to previous + '06': -3, + '12': -2, + '18': -1, + '00': 0, # next day, so must go to previous } ostep = obs_date[-2:] - cstep = rstep[ostep] # current step + cstep = rstep[ostep] # current step if cstep == 0: # previous day - odate = datetime.strptime(obs_date[:-2], self.date_fmt)-timedelta(days=1) + odate = datetime.strptime(obs_date[:-2], self.date_fmt) + odate -= timedelta(days=1) else: # current day odate = datetime.strptime(obs_date[:-2], self.date_fmt) @@ -1116,20 +1111,29 @@ class Evaluator(object): for rej in rej_tmp: if self.mylogger: self.mylogger.debug("REJ FILE: %s" % (rej)) tmp = pd.read_table( - rej, - sep='\s*', - engine='python', - usecols=(0,1,2,3), - names=('datetime', 'lon', 'lat', self.obs_name), - ) - - if self.mylogger: self.mylogger.debug("CSTEP: %s REJ VALS: %s" % (type(cstep), tmp['datetime'].dtype)) + rej, + sep='\s*', + engine='python', + usecols=(0, 1, 2, 3), + names=('datetime', 'lon', 'lat', self.obs_name), + ) + + if self.mylogger: + self.mylogger.debug( + "CSTEP: %s REJ VALS: %s", type(cstep), + tmp['datetime'].dtype + ) # extract valid values tmp = tmp[tmp['datetime']==cstep] obs_rej_tmp.append(tmp) - obs_rej = pd.concat(obs_rej_tmp, names=('datetime', 'lon', 'lat', self.obs_name)) - if self.mylogger: self.mylogger.debug("SHAPE: %s\tOBS TO BE REJECTED 1:\n%s" % (obs_rej.shape, obs_rej)) + obs_rej = pd.concat( + obs_rej_tmp, names=('datetime', 'lon', 'lat', self.obs_name) + ) + if self.mylogger: + self.mylogger.debug( + "SHAPE: %s\tOBS TO BE REJECTED 1:\n%s", obs_rej.shape, obs_rej + ) # convert reject observation to the same format obs_rej['datetime'] = obs_rej['datetime'].replace(cstep, obs_date) @@ -1138,11 +1142,16 @@ class Evaluator(object): obs_rej['lat'] = (((180.0/(91-1)) * (91+4-obs_rej['lat']))-90).round(2) obs_rej[self.obs_name] = obs_rej[self.obs_name].round(2) - if self.mylogger: self.mylogger.debug("SHAPE: %s\tOBS TO BE REJECTED 2:\n%s" % (obs_rej.shape, obs_rej)) - if self.mylogger: self.mylogger.debug("OBS SHAPE BEFORE:\n%s" % (obs.shape,)) + if self.mylogger: + self.mylogger.debug( + "SHAPE: %s\tOBS TO BE REJECTED 2:\n%s", obs_rej.shape, obs_rej + ) + if self.mylogger: + self.mylogger.debug("OBS SHAPE BEFORE:\n%s", obs.shape) # remove from obs value present in obs_rej, rounded at the 2nd decimal obs = obs[~obs.round({'lon': 2, 'lat': 2}).isin(obs_rej)].dropna() - if self.mylogger: self.mylogger.debug("OBS SHAPE AFTER:\n%s" % (obs.shape,)) + if self.mylogger: + self.mylogger.debug("OBS SHAPE AFTER:\n%s", obs.shape) return obs @@ -1257,17 +1266,18 @@ class Evaluator(object): # call the checkFileMatch method if self.mylogger: - self.mylogger.debug("BEFORE\n\tOBSTMP: %s\n\tMODTMP: %s" %\ - (obs_tmp, mod_tmp)) + self.mylogger.debug( + "BEFORE\n\tOBSTMP: %s\n\tMODTMP: %s", obs_tmp, mod_tmp + ) obs_tmp, mod_tmp = self.checkFileMatch(obs_tmp, mod_tmp) if self.mylogger: - self.mylogger.debug("AFTER\n\tOBSTMP: %s\n\tMODTMP: %s" %\ - (obs_tmp, mod_tmp)) + self.mylogger.debug( + "AFTER\n\tOBSTMP: %s\n\tMODTMP: %s",obs_tmp, mod_tmp + ) self.obs_list.extend(obs_tmp) self.mod_list.extend(mod_tmp) - def createDateList(self): """Build the list of dates to search for""" @@ -1276,33 +1286,41 @@ class Evaluator(object): return [d.to_datetime() for d in pd.date_range(str_date, end_date)] - def selectRegion(self, lon, lat, mod): """Select a defined region""" if self.lon_bnds: - lon_idx = np.where((lon >= self.lon_bnds[0]) & (lon <= self.lon_bnds[1])) + lon_idx = np.where( + (lon >= self.lon_bnds[0]) & (lon <= self.lon_bnds[1]) + ) else: lon_idx = (np.arange(lon.size),) if self.lat_bnds: - lat_idx = np.where((lat >= self.lat_bnds[0]) & (lat <= self.lat_bnds[1])) + lat_idx = np.where( + (lat >= self.lat_bnds[0]) & (lat <= self.lat_bnds[1]) + ) else: lat_idx = (np.arange(lat.size),) - if self.mylogger: self.mylogger.debug("LON_IDX: %s LAT_IDX: %s" % (lon_idx[0].size, lat_idx[0].size)) + if self.mylogger: + self.mylogger.debug( + "LON_IDX: %s LAT_IDX: %s", lon_idx[0].size, lat_idx[0].size + ) # in case of observation lon and lat vector must have the same size if len(mod.shape) == 1: idx = np.intersect1d(lon_idx, lat_idx) return lon[idx], lat[idx], mod[idx] - return lon[lon_idx], lat[lat_idx], mod.squeeze()[lat_idx[0],:][:,lon_idx[0]].T - + return ( + lon[lon_idx], + lat[lat_idx], + mod.squeeze()[lat_idx[0], :][:,lon_idx[0]].T + ) def genStats(self, data): """Generate statistics""" - model = data[self.mod_name] obsrv = data[self.obs_name] @@ -1319,10 +1337,8 @@ class Evaluator(object): with open(outfile, 'w') as f: f.write(outputs) - def genPlots(self, data): """Generate plots""" - model = data[self.mod_name] obsrv = data[self.obs_name] datet = data['datetime'] @@ -1337,7 +1353,7 @@ class Evaluator(object): fmt = '%Y%m%d%H' data['datetime'] = pd.to_datetime(datet, format=fmt) data = data.set_index(pd.DatetimeIndex(datet)) - #if self.mylogger: self.mylogger.debug("DATA:\n%s" % (data)) + # if self.mylogger: self.mylogger.debug("DATA:\n%s" % (data)) if self.plt_type == 'tss': data.plot(style='.', ms=20) else: @@ -1352,17 +1368,13 @@ class Evaluator(object): plt.savefig(outimag) plt.clf() - def cleanTmp(self): """Clean temporary files""" - for f in glob(self.tmp_tmpl % '*'): os.remove(f) - def runEval(self): """Main method to run the evaluation""" - # remove tmp filelist self.cleanTmp() @@ -1399,115 +1411,144 @@ class Evaluator(object): # select region for model nlon, nlat, nmod = self.selectRegion(flon, flat, fmod) - if self.mylogger: self.mylogger.debug("NLON: %s NLAT: %s NMOD: %s" % (nlon.shape, nlat.shape, nmod.shape)) + if self.mylogger: + self.mylogger.debug( + "NLON: %s NLAT: %s NMOD: %s", nlon.shape, nlat.shape, + nmod.shape + ) # select region for obs xobs, yobs, nobs = self.selectRegion(obs['lon'].values, obs['lat'].values, obs[self.obs_name].values) - if self.mylogger: self.mylogger.debug("XOBS: %s YOBS: %s NOBS: %s" % (xobs.shape, yobs.shape, nobs.shape)) + if self.mylogger: + self.mylogger.debug( + "XOBS: %s YOBS: %s NOBS: %s", + xobs.shape, yobs.shape, nobs.shape + ) # interpolate netcdf values to observation coords interp_func = RectBivariateSpline(nlon, nlat, nmod) interp_mod = interp_func(xobs, yobs, grid=False) - if self.mylogger: self.mylogger.debug("INTERP_MOD: %s" % interp_mod.shape) - + if self.mylogger: + self.mylogger.debug("INTERP_MOD: %s", interp_mod.shape) f.close() # create dataframe with a column for observations and another one # with interpolated model values d = { - 'datetime' : pd.Series(obs['datetime']), - self.obs_name : pd.Series(nobs), - self.mod_name : pd.Series(interp_mod), - } + 'datetime': pd.Series(obs['datetime']), + self.obs_name: pd.Series(nobs), + self.mod_name: pd.Series(interp_mod), + } # append to the dataframe list frame_list.append(pd.DataFrame(d)) if frame_list: - if self.mylogger: self.mylogger.debug("FRAME LIST LEN: %s" % len(frame_list)) + if self.mylogger: + self.mylogger.debug("FRAME LIST LEN: %s" % len(frame_list)) # concatenate all dataframes into one data = pd.concat(frame_list) - if self.mylogger: self.mylogger.debug("DATA:\n%s" % data) + if self.mylogger: + self.mylogger.debug("DATA:\n%s", data) # generate statistics - self.genStats(data) #model, obsrv) + self.genStats(data) # model, obsrv) # generate plots - self.genPlots(data) #model, obsrv) + self.genPlots(data) # model, obsrv) else: - if self.mylogger: self.mylogger.debug("Empty file list. Exit.") + if self.mylogger: + self.mylogger.debug("Empty file list. Exit.") # remove tmp filelist self.cleanTmp() -def testDA(): - #/scratch/Earth/editomas/observations/modis_nrl_oceanland_aquaterra_AE_AI_count_selected/200705/2007051106_oceanland_aquaterra_obsnew.txt - #/scratch/Earth/editomas/observations/modis_nrl_oceanland_aquaterra_AE_AI_count_selected/200705/2007051112_oceanland_aquaterra_obsnew.txt - #/scratch/Earth/editomas/observations/modis_nrl_oceanland_aquaterra_AE_AI_count_selected/200705/2007051118_oceanland_aquaterra_obsnew.txt - #/scratch/Earth/editomas/observations/modis_nrl_oceanland_aquaterra_AE_AI_count_selected/200705/2007051200_oceanland_aquaterra_obsnew.txt - #/scratch/Earth/editomas/observations/modis_DB_land_aqua_AE_AI_Count_selected/200705/2007051112_deepblue_aqua_obsnew.txt - #/scratch/Earth/editomas/data/nmmb-bsc-ctm-v2.0_DA-2007-NRLDB-bincorr1-uss-1aer/ENS/departures/2007051106-2007051200/obs_dep_rejected.txt - - STR_DATE = "20070611" - END_DATE = "20070612" - - NUM = "84" - - OBS_PATH = ("/scratch/Earth/editomas/observations/modis_nrl_oceanland_aquaterra_AE_AI_count_selected/", - "/scratch/Earth/editomas/observations/modis_DB_land_aqua_AE_AI_Count_selected/",) - OBS_TMPL = ("%(year)s%(month)s/%(date)s*_oceanland_aquaterra_obsnew.txt", - "%(year)s%(month)s/%(date)s*_deepblue_aqua_obsnew.txt", ) - - REJ_PATH = None #("/scratch/Earth/editomas/data/nmmb-bsc-ctm-v2.0_DA-2007-NRLDB-bincorr1-uss-1aer/ENS/departures/",) - REJ_TMPL = None #("%(date)s*/obs_dep_rejected.txt",) - - #MOD_PATH = "/scratch/Earth/editomas/data/nmmb-bsc-ctm-v2.0_FC%s_CONTROL-2007-LR-5dd-IC-AN/ENS/nc_analysis_validation/" % NUM - MOD_PATH = "/scratch/Earth/editomas/data/nmmb-bsc-ctm-v2.0_FG-2007-NRLDB-bincorr1-uss-1aer/ENS/nc_analysis_validation/" - #FIX_TMPL = "_FC%(num)s_FC%(num)s_CONTROL-2007-LR-5dd-IC-AN.nc" % { 'num' : NUM } - FIX_TMPL = "_FG-2007-NRLDB-bincorr1-uss-1aer.nc" - MOD_TMPL = "%(date)s*" + FIX_TMPL - - DIR_OUTP = "./OUT/" - - DEBUG = True +def test_da(): + # /scratch/Earth/editomas/observations/modis_nrl_oceanland_aquaterra_AE_AI_count_selected/200705/2007051106_oceanland_aquaterra_obsnew.txt + # /scratch/Earth/editomas/observations/modis_nrl_oceanland_aquaterra_AE_AI_count_selected/200705/2007051112_oceanland_aquaterra_obsnew.txt + # /scratch/Earth/editomas/observations/modis_nrl_oceanland_aquaterra_AE_AI_count_selected/200705/2007051118_oceanland_aquaterra_obsnew.txt + # /scratch/Earth/editomas/observations/modis_nrl_oceanland_aquaterra_AE_AI_count_selected/200705/2007051200_oceanland_aquaterra_obsnew.txt + # /scratch/Earth/editomas/observations/modis_DB_land_aqua_AE_AI_Count_selected/200705/2007051112_deepblue_aqua_obsnew.txt + # /scratch/Earth/editomas/data/nmmb-bsc-ctm-v2.0_DA-2007-NRLDB-bincorr1-uss-1aer/ENS/departures/2007051106-2007051200/obs_dep_rejected.txt + + str_date = "20070611" + end_date = "20070612" + + # num = "84" + + obs_path = ( + "/scratch/Earth/editomas/observations/" + "modis_nrl_oceanland_aquaterra_AE_AI_count_selected/", + "/scratch/Earth/editomas/observations/" + "modis_DB_land_aqua_AE_AI_Count_selected/", + ) + obs_tmpl = ( + "%(year)s%(month)s/%(date)s*_oceanland_aquaterra_obsnew.txt", + "%(year)s%(month)s/%(date)s*_deepblue_aqua_obsnew.txt", + ) + + rej_path = None + # ("/scratch/Earth/editomas/data/" + # "nmmb-bsc-ctm-v2.0_DA-2007-NRLDB-bincorr1-uss-1aer/ENS/departures/",) + rej_tmpl = None # ("%(date)s*/obs_dep_rejected.txt",) + + # MOD_PATH = "/scratch/Earth/editomas/data/" + # "nmmb-bsc-ctm-v2.0_FC%s_CONTROL-2007-LR-5dd-IC-AN/ENS/" + # "nc_analysis_validation/" % NUM + mod_path = "/scratch/Earth/editomas/data/" \ + "nmmb-bsc-ctm-v2.0_FG-2007-NRLDB-bincorr1-uss-1aer/ENS/" \ + "nc_analysis_validation/" + # FIX_TMPL = "_FC%(num)s_FC%(num)s_CONTROL-2007-LR-5dd-IC-AN.nc" % \ + # { 'num' : NUM } + fix_tmpl = "_FG-2007-NRLDB-bincorr1-uss-1aer.nc" + mod_tmpl = "%(date)s*" + fix_tmpl + + dir_outp = "./OUT/" + + debug = True # create Class - ev = Evaluator(STR_DATE, END_DATE, OBS_PATH, MOD_PATH, OBS_TMPL, MOD_TMPL, DIR_OUTP, REJ_PATH, REJ_TMPL, DEBUG) + ev = Evaluator( + str_date, end_date, obs_path, mod_path, obs_tmpl, mod_tmpl, dir_outp, + rej_path, rej_tmpl, debug + ) # setting attributes ev.plt_xlim = (0, 4) ev.plt_ylim = (0, 4) - #ev.lon_bnds = (-25, 60) - #ev.lat_bnds = (0, 65) - #ev.mod_name = "" - #ev.obs_name = "" - ev.plt_titl = "Glob FG 2007 (REJ)" #% NUM - ev.plt_imag = "glob_fg2007_image_rej.png" #% NUM - ev.fil_stat = "glob_fg2007_stats_rej.txt" #% NUM + # ev.lon_bnds = (-25, 60) + # ev.lat_bnds = (0, 65) + # ev.mod_name = "" + # ev.obs_name = "" + ev.plt_titl = "Glob FG 2007 (REJ)" # % NUM + ev.plt_imag = "glob_fg2007_image_rej.png" # % NUM + ev.fil_stat = "glob_fg2007_stats_rej.txt" # % NUM # run ev.runEval() -def testPollen(): - - STR_DATE = "20130328" - END_DATE = "20130331" +def test_pollen(): + str_date = "20130328" + end_date = "20130331" - OBS_PATH = "/home/Earth/fbeninca/Programs/ACtools/evaluation/tests/polen/" - OBS_TMPL = "Pollen_barcelona_2013-14_Pinus_Platanus_massa.csv" + obs_path = "/home/Earth/fbeninca/Programs/ACtools/evaluation/tests/polen/" + obs_tmpl = "Pollen_barcelona_2013-14_Pinus_Platanus_massa.csv" - MOD_PATH = "/home/Earth/fbeninca/Programs/ACtools/evaluation/tests/polen/" - MOD_TMPL = "%(date)s00/%(date)s00_NMMB-BSC-CTM_regular.nc" + mod_path = "/home/Earth/fbeninca/Programs/ACtools/evaluation/tests/polen/" + mod_tmpl = "%(date)s00/%(date)s00_NMMB-BSC-CTM_regular.nc" - DIR_OUTP = "./OUT/" + dir_outp = "./OUT/" - DEBUG = True + debug = True # create Class - ev = Evaluator(STR_DATE, END_DATE, OBS_PATH, MOD_PATH, OBS_TMPL, MOD_TMPL, DIR_OUTP, en_debug=DEBUG) + ev = Evaluator( + str_date, end_date, obs_path, mod_path, obs_tmpl, mod_tmpl, dir_outp, + en_debug=debug + ) # setting attributes ev.obs_type = 'POLLEN' @@ -1519,9 +1560,9 @@ def testPollen(): ev.lat_bnds = (40, 42) ev.plt_xlim = (-0.0000005, 0.000003) ev.plt_ylim = (0, 0.000009) - ev.plt_titl = u"Pollen March 2013 (g/m$^3$)" #% NUM - ev.plt_imag = "pollen_201303_scatter.png" #% NUM - ev.fil_stat = "pollen_201303.txt" #% NUM + ev.plt_titl = u"Pollen March 2013 (g/m$^3$)" # % NUM + ev.plt_imag = "pollen_201303_scatter.png" # % NUM + ev.fil_stat = "pollen_201303.txt" # % NUM # run ev.runEval() @@ -1533,29 +1574,28 @@ if __name__ == "__main__": nargs = len(sys.argv) if nargs != 2: - print "Error: maximum of 2 args" + print("Error: maximum of 2 args") sys.exit(1) def exec_all(): - #print "EXEC_ALL" + # print "EXEC_ALL" try: - testDA() - except Exception, e: - print "Error (DA):", e + test_da() + except Exception as e: + print("Error (DA):", e) try: - testPollen() - except Exception, e: - print "Error (Pollen):", e + test_pollen() + except Exception as e: + print("Error (Pollen):", e) if nargs == 1: exec_all() else: attr = sys.argv[1] - print "ATTR:", attr + print("ATTR:", attr) current_module = sys.modules[__name__] # print attr, current_module.__name__ # for i in dir(current_module): # print i - hasattr(current_module, attr) and getattr(current_module, attr)() #or exec_all() - - + hasattr(current_module, attr) and getattr(current_module, attr)() + # or exec_all() diff --git a/setup.cfg b/setup.cfg new file mode 100644 index 0000000000000000000000000000000000000000..defd6ac9e7a519ae371274fb758db574f6adc314 --- /dev/null +++ b/setup.cfg @@ -0,0 +1,36 @@ +[aliases] +test=pytest + +[build_sphinx] +source-dir = docs/ +build-dir = docs/build +all_files = 1 +builder = html + +[tool:pytest] +addopts = + --flake8 + --doctest-modules + --ignore=mapgenerator/utils + --ignore=mapgenerator/test/unit/test_http.py + --ignore=tests + --cov=mapgenerator + --cov-report=term + --cov-report=xml:test-reports/coverage.xml + --cov-report=html:test-reports/coverage_html + --junit-xml=test-reports/report.xml + --html=test-reports/report.html +env = + MPLBACKEND = Agg +flake8-ignore = + docs/conf.py ALL + utils/iplot.py ALL + utils/cross.py ALL + tests/context ALL +log_level = WARNING + +[coverage:run] +parallel = true + +[pydocstyle] +convention = numpy diff --git a/setup.py b/setup.py index 90142ef881e4d70610423e1d7449e32c4cf6d235..81b93a9f8f5dd3f3ee98ad7688680ab4734723fa 100644 --- a/setup.py +++ b/setup.py @@ -4,7 +4,7 @@ # Copyright 2018 Barcelona Supercomputing Center - Centro Nacional de # Supercomputación (BSC-CNS) -# This file is part of MapGenerator. +# This file is part of MapGenerator. # MapGenerator is free software: you can redistribute it and/or modify # it under the terms of the GNU General Public License as published by @@ -24,6 +24,21 @@ from setuptools import setup from mapgenerator import __version__ +REQUIREMENTS = { + 'test': [ + 'pytest>=3.9', + 'pytest-cov', + 'pytest-env', + 'pytest-flake8', + 'pytest-html', + 'pytest-metadata>=1.5.1', + ], + 'setup': [ + 'pytest-runner', + 'setuptools_scm', + ], +} + setup( # Application name: name="mapgenerator", @@ -37,21 +52,23 @@ setup( # Packages #package_dir={'mapgenerator': 'mapgenerator/mapgenerator'}, - packages=['mapgenerator','mapgenerator.plotting'], + packages=['mapgenerator', 'mapgenerator.plotting'], # Include additional files into the package - #include_package_data=True, - scripts=['bin/mg',], + # include_package_data=True, + scripts=['bin/mg', ], # Details url="https://pypi.org/project/MapGenerator", - keywords=['earth sciences','weather','climate','air quality','2D maps'], + keywords=['earth sciences', 'weather', + 'climate', 'air quality', '2D maps'], description="Map Generator is a tool that provides easy to use 2D plotting functions for Earth sciences datasets.", -# long_description=open("README.rst").read(), -# long_description_content_type='text/x-rst', + # long_description=open("README.rst").read(), + # long_description_content_type='text/x-rst', # Dependent packages (distributions) + setup_requires=REQUIREMENTS['setup'], install_requires=[ "matplotlib", "pandas", @@ -60,6 +77,10 @@ setup( "netCDF4", "ConfigArgParse", "lxml", + "plotly", + "chart_studio", + "config", ], + tests_require=REQUIREMENTS['test'], ) diff --git a/tests/test_map_generator.py b/tests/test_map_generator.py index 395d945c9275fe82ed36f9cab53168f07bb44aa7..bb613261c97b5bf5f8b1d9762fc445e450007109 100644 --- a/tests/test_map_generator.py +++ b/tests/test_map_generator.py @@ -3,14 +3,14 @@ # # author: Francesco Benincasa +import subprocess +import time +from config import * +from pylab import * +from grads import * import matplotlib as mpl mpl.use('Agg') -from grads import * -from pylab import * - -import time, commands -from config import * VERSION = "0.1beta" @@ -28,13 +28,13 @@ def setColorMap(): cmap.set_bad('red') return cmap - -def setColorBar(drawedges=False,cax=None): + +def setColorBar(drawedges=False, cax=None): """ create color bar """ bounds = BOUNDS boundaries = BOUNDARIES - #print bounds, boundaries + # print bounds, boundaries extend = 'neither' if type(boundaries) is list: boundaries = [boundaries[0]] + bounds + [boundaries[1]] @@ -47,19 +47,19 @@ def setColorBar(drawedges=False,cax=None): boundaries = bounds + [boundaries] extend = 'max' - #print boundaries, extend + # print boundaries, extend norm = mpl.colors.BoundaryNorm(bounds, len(bounds)-1, clip=True) cb = mpl.colorbar.ColorbarBase(drawedges=drawedges, - ax=cax, - cmap=cmap, - norm=norm, - # to use 'extend', you must - # specify two extra boundaries: - boundaries=boundaries, - extend=extend, - ticks=bounds, # optional - spacing='uniform', - orientation='vertical') + ax=cax, + cmap=cmap, + norm=norm, + # to use 'extend', you must + # specify two extra boundaries: + boundaries=boundaries, + extend=extend, + ticks=bounds, # optional + spacing='uniform', + orientation='vertical') return cb @@ -71,33 +71,37 @@ def genImageMap(ga, modName, varName, multiply, nTime, imgTitle): ga.cmd("set lon %s %s" % (str(LON[0]), str(LON[1]))) ga.cmd("set lat %s %s" % (str(LAT[0]), str(LAT[1]))) ga.cmd("set t %s" % (nTime+1)) - + dims = ga.query("dims", Quiet=True) - print "nTime", nTime + print("nTime", nTime) if nTime == 0: global run_tmp global run - run_tmp = dims.time[0] #12Z30NOV2010 - run = "%sh %s %s %s" % (run_tmp[:2], run_tmp[3:5], run_tmp[5:8], run_tmp[8:]) - print "run", run + run_tmp = dims.time[0] # 12Z30NOV2010 + run = "%sh %s %s %s" % ( + run_tmp[:2], run_tmp[3:5], run_tmp[5:8], run_tmp[8:]) + print("run", run) valid_tmp = dims.time[1] - valid = "%sh %s %s %s" % (valid_tmp[:2], valid_tmp[3:5], valid_tmp[5:8], valid_tmp[8:]) + valid = "%sh %s %s %s" % ( + valid_tmp[:2], valid_tmp[3:5], valid_tmp[5:8], valid_tmp[8:]) ga.basemap(resolution='l') ga.transf = True ga.contourf(varName, mul=multiply, N=len(BOUNDS)-1, V=BOUNDS, - cbar=setColorBar, clines=False, colors=COLORS) #, cbar=None) + cbar=setColorBar, clines=False, colors=COLORS) # , cbar=None) - #add states and parallels + # add states and parallels coord_int = 10 - #coords normalization + # coords normalization lat_offset = abs(LAT[0]) % coord_int lon_offset = abs(LON[0]) % coord_int ga.map.drawcountries(linewidth=0.7) - ga.map.drawparallels(arange(LAT[0]+lat_offset, LAT[1], coord_int),labels=[1,0,0,0],linewidth=0.2) - ga.map.drawmeridians(arange(LON[0]+lon_offset, LON[1], coord_int),labels=[0,0,0,1],linewidth=0.2) + ga.map.drawparallels(arange( + LAT[0]+lat_offset, LAT[1], coord_int), labels=[1, 0, 0, 0], linewidth=0.2) + ga.map.drawmeridians(arange( + LON[0]+lon_offset, LON[1], coord_int), labels=[0, 0, 0, 1], linewidth=0.2) nTime = nTime*FREQ @@ -111,12 +115,14 @@ def genImageMap(ga, modName, varName, multiply, nTime, imgTitle): if ACMAD == '01': test_comp = "1-10 %s %s" % (month, year) - figName = "%s.%s.an%s.m%s.d01.10" % (varName.replace('_DUST','').lower(),modName,run_tmp[-4:],time.strftime('%m', time.strptime(run_tmp[-7:-4], '%b'))) + figName = "%s.%s.an%s.m%s.d01.10" % (varName.replace('_DUST', '').lower( + ), modName, run_tmp[-4:], time.strftime('%m', time.strptime(run_tmp[-7:-4], '%b'))) elif ACMAD == '02': test_comp = "11-20 %s %s" % (month, year) - figName = "%s.%s.an%s.m%s.d11.20" % (varName.replace('_DUST','').lower(),modName,run_tmp[-4:],time.strftime('%m', time.strptime(run_tmp[-7:-4], '%b'))) + figName = "%s.%s.an%s.m%s.d11.20" % (varName.replace('_DUST', '').lower( + ), modName, run_tmp[-4:], time.strftime('%m', time.strptime(run_tmp[-7:-4], '%b'))) elif ACMAD == '03': - if month in ('NOV','APR','JUN','SEP'): + if month in ('NOV', 'APR', 'JUN', 'SEP'): days = '21-30' elif month == 'FEB': if int(year) % 4 != 0: @@ -126,29 +132,31 @@ def genImageMap(ga, modName, varName, multiply, nTime, imgTitle): else: days = '21-31' test_comp = "%s %s %s" % (days, month, year) - figName = "%s.%s.an%s.m%s.d21.30" % (varName.replace('_DUST','').lower(),modName,run_tmp[-4:],time.strftime('%m', time.strptime(run_tmp[-7:-4], '%b'))) + figName = "%s.%s.an%s.m%s.d21.30" % (varName.replace('_DUST', '').lower( + ), modName, run_tmp[-4:], time.strftime('%m', time.strptime(run_tmp[-7:-4], '%b'))) elif ACMAD == '04': test_comp = "%s %s" % (month, year) - figName = "%s.%s.an%s.m%s" % (varName.replace('_DUST','').lower(),modName,run_tmp[-4:],time.strftime('%m', time.strptime(run_tmp[-7:-4], '%b'))) + figName = "%s.%s.an%s.m%s" % (varName.replace('_DUST', '').lower( + ), modName, run_tmp[-4:], time.strftime('%m', time.strptime(run_tmp[-7:-4], '%b'))) else: test_comp = (run, valid, 'H+'+sTime) - figName = "%s-%s-%s-%s" % (modName,varName,sTime,run_tmp) - + figName = "%s-%s-%s-%s" % (modName, varName, sTime, run_tmp) try: title(unicode(imgTitle.decode('utf-8')) % test_comp) except: title(imgTitle % test_comp) - print "printing ", figName, " ..." + print("printing ", figName, " ...") savefig(figName + ".png") clf() if TOGIF: - #convert from png to gif and delete png - st, out = commands.getstatusoutput("convert %s.png %s.gif && rm %s.png" % (figName,figName,figName)) + # convert from png to gif and delete png + st, out = subprocess.getstatusoutput( + "convert %s.png %s.gif && rm %s.png" % (figName, figName, figName)) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) sys.exit(1) return figName @@ -158,7 +166,7 @@ def compareMaps(mapNames): """ generate Maps Comparison """ for item0 in mapNames[0]: - #print "loop 1" + # print "loop 1" idx0 = str(mapNames[0].index(item0)) if len(idx0) == 1: idx0 = '0' + idx0 @@ -168,10 +176,10 @@ def compareMaps(mapNames): var0 = item0.split('-')[-3:-2][0] mod0 = '-'.join(item0.split('-')[:-3]) - newDate0 = time.strftime('%Y%m%d', time.strptime(runDate0,"%d%b%Y")) + newDate0 = time.strftime('%Y%m%d', time.strptime(runDate0, "%d%b%Y")) for item1 in mapNames[1]: - #print "loop 2" + # print "loop 2" idx1 = str(mapNames[1].index(item1)) if len(idx1) == 1: idx1 = '0' + idx1 @@ -181,39 +189,42 @@ def compareMaps(mapNames): var1 = item1.split('-')[-3:-2][0] mod1 = '-'.join(item1.split('-')[:-3]) - newDate1 = time.strftime('%Y%m%d', time.strptime(runDate1,"%d%b%Y")) + newDate1 = time.strftime( + '%Y%m%d', time.strptime(runDate1, "%d%b%Y")) if (runDate0 != runDate1): - print "Error: different date, comparison impossible." + print("Error: different date, comparison impossible.") continue if (var0 != var1): - print "Error: different variable, comparison impossible." + print("Error: different variable, comparison impossible.") continue if (int(runHour0)+int(hOffset0)) == (int(runHour1)+int(hOffset1)): - st, out = commands.getstatusoutput("convert +append %s.gif %s.gif %s-%s-%s-compared-map-%s-%s.gif" %\ - (item0, item1, newDate0, mod0, mod1, var0, idx0)) + st, out = subprocess.getstatusoutput("convert +append %s.gif %s.gif %s-%s-%s-compared-map-%s-%s.gif" % + (item0, item1, newDate0, mod0, mod1, var0, idx0)) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) else: break - st, out = commands.getstatusoutput("convert -delay 75 -loop 0 %s-%s-%s-compared-map-%s-??.gif %s-%s-%s-compared-map-%s-loop.gif" % \ - (newDate0, mod0, mod1, var0, newDate0, mod0, mod1, var0,)) + st, out = subprocess.getstatusoutput("convert -delay 75 -loop 0 %s-%s-%s-compared-map-%s-??.gif %s-%s-%s-compared-map-%s-loop.gif" % + (newDate0, mod0, mod1, var0, newDate0, mod0, mod1, var0,)) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) if __name__ == "__main__": - import os, sys, string + import os + import sys + import string error_msg = """\ Usage: %s """ % sys.argv[0] - if len(sys.argv) not in (3,4): - print error_msg + if len(sys.argv) not in (3, 4): + print(error_msg) sys.exit(1) global cmap @@ -223,7 +234,7 @@ if __name__ == "__main__": clf() - run = '' + run = '' run_tmp = '' mapNames = [] cmap = setColorMap() @@ -236,21 +247,22 @@ if __name__ == "__main__": curdir = os.getcwd() - ga = GaLab(Bin='grads',Window=False,Echo=False) + ga = GaLab(Bin='grads', Window=False, Echo=False) for fName in fNames: os.chdir(curdir) ga.cmd("reinit") - ga.open(INDIR + '/'+ fName) + ga.open(INDIR + '/' + fName) os.chdir(OUTDIR) - modName = '_'.join(fName.split('_')[1:]).replace('.nc','') + modName = '_'.join(fName.split('_')[1:]).replace('.nc', '') runDate = fName.split('_')[0] for sec in OPTS.keys(): data = sec.split('-') - secMod = '-'.join(data[:-1]) #unify all elements without the last one - secVar = data[-1:][0] #I suppose the last element is the variable + # unify all elements without the last one + secMod = '-'.join(data[:-1]) + secVar = data[-1:][0] # I suppose the last element is the variable if secMod == modName and secVar == varName: imgTitle = OPTS[sec]['title'] multiply = OPTS[sec]['mul'] @@ -264,25 +276,27 @@ if __name__ == "__main__": break figNames = [] - for nTime in range(0,TOTAL,INTERVAL): - figNames.append(genImageMap(ga, modName, varName, multiply, nTime, imgTitle)) + for nTime in range(0, TOTAL, INTERVAL): + figNames.append(genImageMap( + ga, modName, varName, multiply, nTime, imgTitle)) - #generate animation. TODO: use imagemagick libs + # generate animation. TODO: use imagemagick libs if anim: - st, out = commands.getstatusoutput("convert -delay 75 -loop 0 %s-%s-*-%s.gif %s-%s-%s-loop.gif" % (modName, varName, run_tmp, runDate, modName, varName)) + st, out = subprocess.getstatusoutput("convert -delay 75 -loop 0 %s-%s-*-%s.gif %s-%s-%s-loop.gif" % ( + modName, varName, run_tmp, runDate, modName, varName)) if st != 0: - print "Error: %s" % str(out) + print("Error: %s" % str(out)) sys.exit(1) mapNames.append(figNames) - #comparing different (2) models - #MODNAME-VAR-ISTANT-HHZDATE - #3H_MACC-ECMWF-SCONC_DUST-12-00Z18JAN2011.gif + # comparing different (2) models + # MODNAME-VAR-ISTANT-HHZDATE + # 3H_MACC-ECMWF-SCONC_DUST-12-00Z18JAN2011.gif num = len(mapNames) if num > 2: - print "Error: comparison only between 2 files." + print("Error: comparison only between 2 files.") sys.exit(1) if num > 1: compareMaps(mapNames)