comparison xarray_mapplot.py @ 4:9bbaab36a5d4 draft

"planemo upload for repository https://github.com/galaxyecology/tools-ecology/tree/master/tools/data_manipulation/xarray/ commit 2166974df82f97557b082a9e55135098e61640c4"
author ecology
date Thu, 20 Jan 2022 17:09:40 +0000
parents 663268794710
children
comparison
equal deleted inserted replaced
3:663268794710 4:9bbaab36a5d4
17 # [--xlim "x1,x2"] 17 # [--xlim "x1,x2"]
18 # [--ylim "y1,y2"] 18 # [--ylim "y1,y2"]
19 # [--range "valmin,valmax"] 19 # [--range "valmin,valmax"]
20 # [--threshold VAL] 20 # [--threshold VAL]
21 # [--label label-colorbar] 21 # [--label label-colorbar]
22 # [--config config-file]
22 # [--shift] 23 # [--shift]
23 # [-v] 24 # [-v]
24 # input varname 25 # input varname
25 # 26 #
26 # positional arguments: 27 # positional arguments:
44 # --xlim limited geographical area longitudes "x1,x2" 45 # --xlim limited geographical area longitudes "x1,x2"
45 # --ylim limited geographical area latitudes "y1,y2" 46 # --ylim limited geographical area latitudes "y1,y2"
46 # --range "valmin,valmax" for plotting 47 # --range "valmin,valmax" for plotting
47 # --threshold do not plot values below threshold 48 # --threshold do not plot values below threshold
48 # --label set a label for colormap 49 # --label set a label for colormap
50 # --config plotting parameters are passed via a config file
51 # (overwrite other plotting options)
49 # --shift shift longitudes if specified 52 # --shift shift longitudes if specified
50 # -v, --verbose switch on verbose mode 53 # -v, --verbose switch on verbose mode
51 # 54 #
52 55
53 import argparse 56 import argparse
66 69
67 import xarray as xr # noqa: E402 70 import xarray as xr # noqa: E402
68 71
69 72
70 class MapPlotXr (): 73 class MapPlotXr ():
71 def __init__(self, input, proj, varname, cmap, output, verbose=False, 74 def __init__(self, input, varname, output, verbose=False,
72 time=[], title="", latitude="latitude", 75 config_file="", proj="", shift=False):
73 longitude="longitude", land=0, ocean=0, 76
74 coastline=0, borders=0, xlim=[], ylim=[], 77 li = list(input.split(","))
75 threshold="", label="", shift=False, 78 if len(li) > 1:
76 range_values=[]): 79 self.input = li
77 self.input = input 80 else:
78 print("PROJ", proj) 81 self.input = input
79 if proj != "" and proj is not None: 82
80 self.proj = proj.replace('X', ':') 83 if proj != "" and proj is not None and Path(proj).exists():
81 else: 84 f = open(proj)
82 self.proj = proj 85 sdict = ''.join(
86 f.read().replace("\n", "").split('{')[1].split('}')[0]
87 )
88 self.proj = '{' + sdict.strip() + '}'
89 else:
90 self.proj = None
83 self.varname = varname 91 self.varname = varname
84 self.get_cmap(cmap)
85 self.time = time
86 self.latitude = latitude
87 self.longitude = longitude
88 self.land = land
89 self.ocean = ocean
90 self.coastline = coastline
91 self.borders = borders
92 self.xlim = xlim
93 self.ylim = ylim
94 self.range = range_values
95 self.threshold = threshold
96 self.shift = shift 92 self.shift = shift
97 self.xylim_supported = False 93 self.xylim_supported = False
98 self.colorbar = True 94 self.colorbar = True
99 self.title = title
100 if output is None: 95 if output is None:
101 self.output = Path(input).stem + '.png' 96 if type(self.input) is list:
97 self.output = Path(self.input[0]).stem + '.png'
98 else:
99 self.output = Path(self.input).stem + '.png'
102 else: 100 else:
103 self.output = output 101 self.output = output
104 self.verbose = verbose 102 self.verbose = verbose
105 self.dset = xr.open_dataset(self.input, use_cftime=True)
106
107 self.label = {} 103 self.label = {}
108 if label != "" and label is not None: 104 self.time = []
109 self.label['label'] = label 105 self.xlim = []
106 self.ylim = []
107 self.range = []
108 self.latitude = "latitude"
109 self.longitude = "longitude"
110 self.land = 0
111 self.ocean = 0
112 self.coastline = 0
113 self.borders = 0
114 self.cmap = "coolwarm"
115 self.threshold = ""
116 self.title = ""
117
118 if config_file != "" and config_file is not None:
119 with open(config_file) as f:
120 sdict = ''.join(
121 f.read().replace("\n", "").split('{')[1].split('}')[0]
122 )
123 tmp = ast.literal_eval('{' + sdict.strip() + '}')
124 for key in tmp:
125 if key == 'time':
126 time = tmp[key]
127 self.time = list(map(int, time.split(",")))
128 if key == 'cmap':
129 self.get_cmap(tmp[key])
130 if key == 'latitude':
131 self.latitude = tmp[key]
132 if key == 'longitude':
133 self.longitude = tmp[key]
134 if key == 'land':
135 self.land = float(tmp[key])
136 if key == 'ocean':
137 self.ocean = float(tmp[key])
138 if key == 'coastline':
139 self.coastline = float(tmp[key])
140 if key == 'borders':
141 self.borders = float(tmp[key])
142 if key == 'xlim':
143 xlim = tmp[key]
144 self.xlim = list(map(float, xlim.split(",")))
145 if key == 'ylim':
146 ylim = tmp[key]
147 self.ylim = list(map(float, ylim.split(",")))
148 if key == 'range':
149 range_values = tmp[key]
150 self.range = list(map(float, range_values.split(",")))
151 if key == 'threshold':
152 self.threshold = float(tmp[key])
153 if key == 'label':
154 self.label['label'] = tmp[key]
155 if key == 'title':
156 self.title = tmp[key]
157
158 if type(self.input) is list:
159 self.dset = xr.open_mfdataset(self.input, use_cftime=True)
160 else:
161 self.dset = xr.open_dataset(self.input, use_cftime=True)
162
110 if verbose: 163 if verbose:
111 print("input: ", self.input) 164 print("input: ", self.input)
112 print("proj: ", self.proj) 165 print("proj: ", self.proj)
113 print("varname: ", self.varname) 166 print("varname: ", self.varname)
114 print("time: ", self.time) 167 print("time: ", self.time)
135 def projection(self): 188 def projection(self):
136 if self.proj is None: 189 if self.proj is None:
137 return ccrs.PlateCarree() 190 return ccrs.PlateCarree()
138 191
139 proj_dict = ast.literal_eval(self.proj) 192 proj_dict = ast.literal_eval(self.proj)
140
141 user_proj = proj_dict.pop("proj") 193 user_proj = proj_dict.pop("proj")
142 if user_proj == 'PlateCarree': 194 if user_proj == 'PlateCarree':
143 self.xylim_supported = True 195 self.xylim_supported = True
144 return ccrs.PlateCarree(**proj_dict) 196 return ccrs.PlateCarree(**proj_dict)
145 elif user_proj == 'AlbersEqualArea': 197 elif user_proj == 'AlbersEqualArea':
316 parser = argparse.ArgumentParser() 368 parser = argparse.ArgumentParser()
317 parser.add_argument( 369 parser.add_argument(
318 'input', 370 'input',
319 help='input filename with geographical coordinates (netCDF format)' 371 help='input filename with geographical coordinates (netCDF format)'
320 ) 372 )
321
322 parser.add_argument( 373 parser.add_argument(
323 '--proj', 374 '--proj',
324 help='Specify the projection on which we draw' 375 help='Config file with the projection on which we draw'
325 ) 376 )
326 parser.add_argument( 377 parser.add_argument(
327 'varname', 378 'varname',
328 help='Specify which variable to plot (case sensitive)' 379 help='Specify which variable to plot (case sensitive)'
329 ) 380 )
330 parser.add_argument( 381 parser.add_argument(
331 '--cmap',
332 help='Specify which colormap to use for plotting'
333 )
334 parser.add_argument(
335 '--output', 382 '--output',
336 help='output filename to store resulting image (png format)' 383 help='output filename to store resulting image (png format)'
337 ) 384 )
338 parser.add_argument( 385 parser.add_argument(
339 '--time', 386 '--config',
340 help='list of times to plot for multiple plots' 387 help='pass plotting parameters via a config file'
341 )
342 parser.add_argument(
343 '--title',
344 help='plot title'
345 )
346 parser.add_argument(
347 '--latitude',
348 help='variable name for latitude'
349 )
350 parser.add_argument(
351 '--longitude',
352 help='variable name for longitude'
353 )
354 parser.add_argument(
355 '--land',
356 help='add land on plot with alpha value [0-1]'
357 )
358 parser.add_argument(
359 '--ocean',
360 help='add oceans on plot with alpha value [0-1]'
361 )
362 parser.add_argument(
363 '--coastline',
364 help='add coastline with alpha value [0-1]'
365 )
366 parser.add_argument(
367 '--borders',
368 help='add country borders with alpha value [0-1]'
369 )
370 parser.add_argument(
371 '--xlim',
372 help='limited geographical area longitudes "x1,x2"'
373 )
374 parser.add_argument(
375 '--ylim',
376 help='limited geographical area latitudes "y1,y2"'
377 )
378 parser.add_argument(
379 '--range',
380 help='min and max values for plotting "minval,maxval"'
381 )
382 parser.add_argument(
383 '--threshold',
384 help='do not plot values below threshold'
385 )
386 parser.add_argument(
387 '--label',
388 help='set a label for colorbar'
389 ) 388 )
390 parser.add_argument( 389 parser.add_argument(
391 '--shift', 390 '--shift',
392 help='shift longitudes if specified', 391 help='shift longitudes if specified',
393 action="store_true" 392 action="store_true"
396 "-v", "--verbose", 395 "-v", "--verbose",
397 help="switch on verbose mode", 396 help="switch on verbose mode",
398 action="store_true") 397 action="store_true")
399 args = parser.parse_args() 398 args = parser.parse_args()
400 399
401 if args.time is None: 400 dset = MapPlotXr(input=args.input, varname=args.varname,
402 time = [] 401 output=args.output, verbose=args.verbose,
403 else: 402 config_file=args.config, proj=args.proj,
404 time = list(map(int, args.time.split(","))) 403 shift=args.shift)
405 if args.xlim is None:
406 xlim = []
407 else:
408 xlim = list(map(float, args.xlim.split(",")))
409 if args.ylim is None:
410 ylim = []
411 else:
412 ylim = list(map(float, args.ylim.split(",")))
413 if args.range is None:
414 range_values = []
415 else:
416 range_values = list(map(float, args.range.split(",")))
417 if args.latitude is None:
418 latitude = "latitude"
419 else:
420 latitude = args.latitude
421 if args.longitude is None:
422 longitude = "longitude"
423 else:
424 longitude = args.longitude
425 if args.land is None:
426 land = 0
427 else:
428 land = float(args.land)
429 if args.ocean is None:
430 ocean = 0
431 else:
432 ocean = float(args.ocean)
433 if args.coastline is None:
434 coastline = 0
435 else:
436 coastline = float(args.coastline)
437 if args.borders is None:
438 borders = 0
439 else:
440 borders = float(args.borders)
441
442 dset = MapPlotXr(input=args.input, proj=args.proj, varname=args.varname,
443 cmap=args.cmap, output=args.output, verbose=args.verbose,
444 time=time, title=args.title,
445 latitude=latitude, longitude=longitude, land=land,
446 ocean=ocean, coastline=coastline, borders=borders,
447 xlim=xlim, ylim=ylim, threshold=args.threshold,
448 label=args.label, shift=args.shift,
449 range_values=range_values)
450 404
451 if dset.time == []: 405 if dset.time == []:
452 dset.plot() 406 dset.plot()
453 else: 407 else:
454 for t in dset.time: 408 for t in dset.time: