diff --git a/CHANGELOG b/CHANGELOG index 9a5b110fd0620a70b6a309c91481f6262fa34268..96decfe671451df3d271b5818cc4da4f53c669ed 100644 --- a/CHANGELOG +++ b/CHANGELOG @@ -539,3 +539,13 @@ v2.13 Updates for latest python versions. Improved error checks. + + +v2.14 +----- + +Plot pixels als polygons in case track information is not present. + src/cso/cso_catalogue.py + src/cso/cso_plot.py + + diff --git a/doc/source/history.rst b/doc/source/history.rst index 9294260129788f8d23c55d2c3887f368bd8c697d..4e999294353c1c9326a02f764d1219a2a2bba61f 100644 --- a/doc/source/history.rst +++ b/doc/source/history.rst @@ -121,15 +121,10 @@ A summary of the versions and changes. * | *v2.13* | Minor updates to trap errors and support latest Python versions. - - -To be included -============== - -Features available in branches, that should be included in master version. - -* Create map plots as a list of pixels, in case track information is not available. *(E.D.)* +* | *v2.14* + | Plot pixels als polygons in case track information is not present. + To be done ========== @@ -162,7 +157,7 @@ A wishlist of developments. How to create a version tag =========================== -Example for tag "v2.13". +Example for tag "v2.14". 1. Update file formatting following conventions. @@ -182,17 +177,17 @@ Example for tag "v2.13". c. Change the release tag that will also be imported to the documenation in the `pyproject.toml <../../../pyproject.toml>`_ file:: - release = "2.13" + release = "2.14" d. Commit and push these changes:: - git commit -m 'Defined v2.13' . + git commit -m 'Defined v2.14' . git push 4. Tag and push the code to the git repository:: - git tag -a 'v2.13' -m 'Tagged v2.13' - git push origin v2.13 + git tag -a 'v2.14' -m 'Tagged v2.14' + git push origin v2.14 diff --git a/pyproject.toml b/pyproject.toml index 5c89a7cbdf1967724798337cfdb28439e8b66970..d52ded81ee906d352e34f72452be2ab2f3918c33 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -4,7 +4,7 @@ build-backend = "setuptools.build_meta" [project] name = "cso" -version = "2.13" +version = "2.14" authors = [ { name = "Arjo Segers" }, { name = "Lewis Blake" }, diff --git a/src/cso/cso_catalogue.py b/src/cso/cso_catalogue.py index 0643db50f15d3c23dfb3f646d744c5eb7f91df53..916a56520b5f2f3f67c3bd1786353afa69cd7653 100644 --- a/src/cso/cso_catalogue.py +++ b/src/cso/cso_catalogue.py @@ -35,6 +35,11 @@ # Extended unit conversions. # Fixed reading of settings for map styles. # +# 2025-09, Arjo Segers +# Plot pixels als polygons in case track information is not present. +# Extended unit conversions. +# Fixed string formatting. +# ######################################################################## @@ -156,6 +161,30 @@ class CSO_CatalogueBase(utopya.UtopyaRc): # apply factor: values = values * 1e6 + # vertical column density: + elif conversion == "molec/cm2 -> umol/m2": + # apply factor: + # molec/cm2 / (molec/mol) umol/mol cm2/m2 + values = values / cso_constants.Avog * 1e6 * 1e4 + + # vertical column density: + elif conversion == "molec/cm2 -> 1e15 mlc/cm2": + # apply factor: + # molec/cm2 + values = values / 1e15 + + # vertical column density: + elif conversion == "(molec/cm2)**2 -> 1e15 mlc/cm2": + # apply factor: + # molec/cm2 + values = numpy.sqrt(values) / 1e15 + + # vertical column density, squared + elif conversion == "(molec/cm2)**2 -> umol/m2": + # apply factor: + # molec/cm2 / (molec/mol) umol/mol cm2/m2 + values = numpy.sqrt(values) / cso_constants.Avog * 1e6 * 1e4 + # vertical column density: elif conversion == "mol m-2 -> mmol/m2": # apply factor: @@ -383,7 +412,11 @@ class CSO_Catalogue(CSO_CatalogueBase): filename = os.path.join(inputdir, fname) # orbit key: - orbit = xlst.GetValue(fname, "orbit") + if "orbit" in xlst: + orbit = xlst.GetValue(fname, "orbit") + else: + orbit = "" + #endif ## testing ... # if '1200' not in filename : continue @@ -444,23 +477,32 @@ class CSO_Catalogue(CSO_CatalogueBase): # settings for this variable: vkey = "var.%s" % varname # plot type: - ptype = self.GetSetting("var.{varname}.type", default="map") + ptype = self.GetSetting(f"var.{varname}.type", default="map") # input variable: - vname = self.GetSetting("var.{varname}.source", default=varname) + vname = self.GetSetting(f"var.{varname}.source", default=varname) # target units: - vunits = self.GetSetting("var.{varname}.units", default="None") + vunits = self.GetSetting(f"var.{varname}.units", default="None") # long name used in labels: - long_name = self.GetSetting("var.{varname}.long_name", default=varname) + long_name = self.GetSetting(f"var.{varname}.long_name", default=varname) + + # track defined in files? + on_track = "track_longitude" in orb.ds # switch: if ptype == "map": # style: - vmin = eval(self.GetSetting("var.{varname}.vmin", default="None")) - vmax = eval(self.GetSetting("var.{varname}.vmax", default="None")) - colors = eval(self.GetSetting("var.{varname}.colors", default="None")) + vmin = eval(self.GetSetting(f"var.{varname}.vmin", default="None")) + vmax = eval(self.GetSetting(f"var.{varname}.vmax", default="None")) + colors = eval(self.GetSetting(f"var.{varname}.colors", default="None")) # extract corner grids and values: - xx, yy, values, attrs = orb.GetTrack(vname) + if on_track: + xx, yy, values, attrs = orb.GetTrack(sstate.ds[vname]) + else: + xx, attrs = orb.GetData( "longitude_bounds" ) + yy, attrs = orb.GetData( "latitude_bounds" ) + values, attrs = orb.GetData( vname ) + #endif # adhoc, some chocho files have nans in lon/lat bound arrays, # these are still present because the scan line has some valid values ... @@ -512,19 +554,35 @@ class CSO_Catalogue(CSO_CatalogueBase): bmp.update(bmp_kwargs) # create map figure: - fig = cso_plot.QuickPat( - values, - x=xx, - y=yy, - vmin=vmin, - vmax=vmax, - cmap=dict(colors=colors), # color_bad=color_nan), - bmp=bmp, - domain=domain, - cbar=dict(label=label), - title=title, - figsize=figsize, - ) + if on_track: + fig = cso_plot.QuickPat( + values, + x=xx, + y=yy, + vmin=vmin, + vmax=vmax, + cmap=dict(colors=colors), # color_bad=color_nan), + bmp=bmp, + domain=domain, + cbar=dict(label=label), + title=title, + figsize=figsize, + ) + else: + fig = cso_plot.QuickPoly( + values, + x=xx, + y=yy, + vmin=vmin, + vmax=vmax, + cmap=dict(colors=colors, color_bad=color_nan), + bmp=bmp, + domain=domain, + cbar=dict(label=label), + title=title, + figsize=figsize, + ) + #endif # save: fig.Export(figfile) @@ -812,25 +870,35 @@ class CSO_SimCatalogue(CSO_CatalogueBase): else: vfile, vname = "data", vsource # endif + # track defined in files? + on_track = "track_longitude" in sdata.ds # extract corner grids and values: if vfile == "state": # check ... if vname not in sstate.ds.keys(): - logging.error( - "variable `%s` not found in state from file: %s" - % (vname, statefile) - ) + logging.error(f"variable '{vname}' not found in state from file: {statefile}") raise Exception # endif # extract corner grids and values: - xx, yy, values, attrs = sdata.GetTrack(sstate.ds[vname]) + if on_track: + xx, yy, values, attrs = sdata.GetTrack(sstate.ds[vname]) + else: + xx, attrs = sdata.GetData( "longitude_bounds" ) + yy, attrs = sdata.GetData( "latitude_bounds" ) + values, attrs = sstate.GetData( vname ) + #endif elif vfile == "data": # extract corner grids and values from data file: - xx, yy, values, attrs = sdata.GetTrack(vname) + if on_track: + xx, yy, values, attrs = sdata.GetTrack(vname) + else: + xx, attrs = sdata.GetData( "longitude_bounds" ) + yy, attrs = sdata.GetData( "latitude_bounds" ) + values, attrs = sdata.GetData( vname ) + #endtry else: logging.error( - 'unsupported type keyword "%s" in source description "%s" for variable "%s"' - % (vfile, vsource, varname) + f"unsupported type keyword '{vfile}' in source description '{vsource}' for variable '{varname}'" ) raise Exception # endif @@ -849,9 +917,7 @@ class CSO_SimCatalogue(CSO_CatalogueBase): if (len(shp) == 3) and (shp[2] == 1): values = values[:, :, 0] else: - logging.error( - 'could not squeeze "%s" %s to 2D track' % (vname, str(shp)) - ) + logging.error(f"could not squeeze '{vname}' {shp} to 2D track") raise Exception # endif # endif @@ -866,20 +932,38 @@ class CSO_SimCatalogue(CSO_CatalogueBase): bmp = dict(countries=True) bmp.update(bmp_kwargs) + # info ... + logging.info(indent + " plot pattern ...") # create map figure: - fig = cso_plot.QuickPat( - values, - x=xx, - y=yy, - vmin=vmin, - vmax=vmax, - cmap=dict(colors=colors, color_bad=color_nan), - bmp=bmp, - domain=domain, - cbar=dict(label=label), - title=title, - figsize=figsize, - ) + if on_track: + fig = cso_plot.QuickPat( + values, + x=xx, + y=yy, + vmin=vmin, + vmax=vmax, + cmap=dict(colors=colors, color_bad=color_nan), + bmp=bmp, + domain=domain, + cbar=dict(label=label), + title=title, + figsize=figsize, + ) + else: + fig = cso_plot.QuickPoly( + values, + x=xx, + y=yy, + vmin=vmin, + vmax=vmax, + cmap=dict(colors=colors, color_bad=color_nan), + bmp=bmp, + domain=domain, + cbar=dict(label=label), + title=title, + figsize=figsize, + ) + #endif # save: fig.Export(figfile) @@ -902,7 +986,7 @@ class CSO_SimCatalogue(CSO_CatalogueBase): # break else: - logging.info(indent + " no state file: %s" % statefile) + logging.info(f"{indent} no state file: {statefile}") # endif # state file present # next: diff --git a/src/cso/cso_gridded.py b/src/cso/cso_gridded.py index 5774a2c0ab1a14489044a9d12a006572a52e9c13..cd54bdbc3d8f4fd207097d89b364da44cffe49ec 100644 --- a/src/cso/cso_gridded.py +++ b/src/cso/cso_gridded.py @@ -40,6 +40,9 @@ # 2025-09, Arjo Segers # Updated calculation of temporal means for files with time records. # +# 2026-01, Arjo Segers +# Fixed timerange setup. +# ######################################################################## @@ -593,14 +596,12 @@ class CSO_GriddedAverage(utopya.UtopyaRc): elif tstep in ["day", "daily"]: # start and end of days: rtt1 = pandas.date_range(start=outfile_t1, end=outfile_t2, freq="D")[:-1] - rtt2 = pandas.date_range( - start=rtt1[0] + pandas.Timedelta(days=1), periods=len(rtt1) - ) + rtt2 = pandas.date_range(start=rtt1[0] + pandas.Timedelta(days=1), periods=len(rtt1)) # elif tstep in ["hour", "hourly"]: # hours, end is now the same as start ... rtt1 = pandas.date_range(start=outfile_t1, end=outfile_t2, freq="h")[:-1] - rtt2 = tt1 + rtt2 = rtt1 else: logging.error(f"unsupported output file frequency {tstep}") raise Exception diff --git a/src/cso/cso_plot.py b/src/cso/cso_plot.py index 7227b145659e7339335f4ebe0a506d2bf869b8cd..91aecc73b609137eb01a9c7dde47a73606f05b06 100644 --- a/src/cso/cso_plot.py +++ b/src/cso/cso_plot.py @@ -38,6 +38,9 @@ # Improved supported for defined number of color entries. # Support masked arrays for selecting no-data values. # +# 2025-09, Arjo Segers +# Implemented 'QuickPoly' method to plot list of polygons. +# # # NOTES @@ -1195,6 +1198,7 @@ class ColorbarFigure(Figure): # * def Polygons(self, xx, yy, values, vmin=None, vmax=None, vbounds=None, **kwargs): + """ Add colored polygons. @@ -1211,81 +1215,36 @@ class ColorbarFigure(Figure): """ - # obtain mappable for colorbar, - # easiest was just to mark the centers with a tiny dot; - # this will define the "cmap" and "cbar" attributes: - p = self.Scatter( - xx.mean(axis=1), - yy.mean(axis=1), - c=values, - vmin=vmin, - vmax=vmax, - vbounds=vbounds, - marker=".", - markersize=0, - linewidths=0, - ) - - # map values to colors: - colors = self.cmap(self.cnorm(values)) - # loop: - for ip in range(len(values)): - # add patch: - self.ax.fill(xx[ip, :], yy[ip, :], color=colors[ip]) - # endfor - - # enddef Polygons - - # * - - def PolygonsMesh(self, xx, yy, values, vmin=None, vmax=None, vbounds=None, **kwargs): - """ - Add colored polygons. - - Arguments: - - * ``xx``, ``yy`` : 2D corner fields, shapes (ny+1,nx+1) - - Optional arguments: - - * ``vmin``,``vmax`` : data range for coloring, min and max values by default, - or the settings for PColor if that was called before - * ``vbounds`` : value boundaries (alternative for vmin,vmax) - - """ + # modules: + import numpy + import matplotlib.pyplot as plt + import matplotlib.cm as cm + import matplotlib.collections - # to obtain mappable for colorbar, - # easiest was just to mark the centers with a tiny dot; - # this will define the "cmap" and "cbar" attributes: - p = self.Scatter( - (xx[:-1, :-1] + xx[1:, :-1] + xx[:-1, 1:] + xx[1:, 1:]) / 4.0, - (yy[:-1, :-1] + yy[1:, :-1] + yy[:-1, 1:] + yy[1:, 1:]) / 4.0, - c=values, - vmin=vmin, - vmax=vmax, - vbounds=vbounds, - marker=".", - linewidths=0, - ) + # collect vertices into array with shape (npoly,ncorner,2) + verts = numpy.array( [xx.T,yy.T] ).T + # set mapping from data to color if not done yet: + if self.cnorm is None: + # set norm function: + self.Set_ColorNorm(values, vmin=vmin, vmax=vmax, vbounds=vbounds) + # endif + # map values to colors: colors = self.cmap(self.cnorm(values)) - # shape: - ny, nx = values.shape - # loop: - for j in range(ny): - for i in range(nx): - # add patch: - self.ax.fill( - [xx[j, i], xx[j, i + 1], xx[j + 1, i + 1], xx[j + 1, i]], - [yy[j, i], yy[j, i + 1], yy[j + 1, i + 1], yy[j + 1, i]], - color=colors[j, i], - ) - # endfor # i - # endfor # j + # define polygons: + collection = matplotlib.collections.PolyCollection( verts, + facecolors=colors, edgecolors="face", + cmap=self.cmap, norm=self.cnorm ) + + # add to figure: + p = self.ax.add_collection( collection ) - # enddef PolygonsMesh + # add colorbar if not done yet: + self.Add_ColorBar(p) + + # enddef Polygons # * @@ -1734,6 +1693,137 @@ def QuickPat( # * +def QuickPoly( + cc, + x=None, + y=None, + domain=None, + title=None, + vmin=None, + vmax=None, + dmax=None, + vbounds=None, + bmp=None, + cbar=dict(), + border=None, + fig=None, + **kwargs, +): + """ + + Plot series of polygons. + + Arguments passed to :py:meth:`ColorbarFigure.Polygons` method: + + * ``x``, ``y`` : polygon coordinates with shape ``(npoly,ncorner)``; + * ``domain=(left,right,bottom,top)`` : coordinate bounds; if ``domain`` is ``None`` or + some elements are ``None``, values are based on the ``x`` or ``y`` coordinates if defined. + * ``cc`` : data values, shape (ny,nx) + * ``vmin``,``vmax`` : data range for coloring + * ``vbounds`` : value bounds for colorbar (defines interval with same color) + + Other arguments are passed to 'ColorbarFigure' : + + * ``figsize = (width,height)`` + * ``domain=(left,right,bottom,top)`` : axes extent + * ``title`` : text above axes + * ``cmap = dict( 'colors'=['blue','white','red'], ncolor=9 )`` + * ``cbar = dict( 'ntick'=9, extend='both', label='entity [units]' )``, + or ``False`` to skip colorbar because it is plotted elsewhere. + * ``bmp`` : map projection, see :py:meth:`.Figure.AddAxes`:: + + bmp=dict( countries=True, raster=[30.0,30.0] ) + + Return value: + + * ``fig`` : Figure instance + + If an already existing Figure instance is passed as the ``fig`` argument, + the pattern will be added to the plot. + + """ + + # modules: + import numpy + + # set domain if not present yet: + if domain is None: + domain = [x.min(), x.max(), y.min(), y.max()] + else: + if domain[0] is None: + domain[0] = xx.min() + if domain[1] is None: + domain[1] = xx.max() + if domain[2] is None: + domain[2] = yy.min() + if domain[3] is None: + domain[3] = yy.max() + # endif + + # replace or add domain in bmp argument: + if type(bmp) == dict: + bmp["domain"] = domain + + # symetric range? + if dmax is not None: + # range: + vmin = -dmax + vmax = dmax + # colors: + if "cmap" in kwargs.keys(): + if "colors" not in kwargs["cmap"].keys(): + kwargs["cmap"]["colors"] = "pwb" + # endif + else: + kwargs["cmap"] = {"colors": "pwb"} + # endif + # endif + + # value range: + if hasattr(cc, "values"): + cmin = numpy.nanmin(cc.values) + cmax = numpy.nanmax(cc.values) + else: + cmin = numpy.nanmin(cc) + cmax = numpy.nanmax(cc) + # endif + + # auto extend? + if "extend" not in cbar.keys(): + # extend below minimum? + extend_min = (vmin is not None) and (cmin < vmin) + # extend above maximum? + extend_max = (vmax is not None) and (cmax > vmax) + # set keyword: + if extend_min and (not extend_max): + cbar["extend"] = "min" + elif extend_min and extend_max: + cbar["extend"] = "both" + elif (not extend_min) and extend_max: + cbar["extend"] = "max" + else: + cbar["extend"] = "neither" + # endif + # endif + + # setup figure: + if fig is None: + fig = ColorbarFigure(cbar=cbar, bmp=bmp, domain=domain, title=title, **kwargs) + # endif + + # add polygons: + fig.Polygons(x, y, cc, vmin=vmin, vmax=vmax, vbounds=vbounds) + + # ok + return fig + + +# enddef QuickPoly + + +# * + + def QuickDens(x, y, bins=10, normed=False, **kwargs): """ Density plot. diff --git a/src/utopya/utopya_index.py b/src/utopya/utopya_index.py index 9759db29b999cbf9e86431b4f252ed3d7b2e438a..9dd2284ce4b66f16fa4cb30f00029e5d84167375 100644 --- a/src/utopya/utopya_index.py +++ b/src/utopya/utopya_index.py @@ -758,8 +758,6 @@ class IndexPart(utopya_rc.UtopyaRc): if len(line) > 0: line = line + sep # add content: - # print( 'xxx line = ', line) - # print( '..x level_value= ', level_value ) line = line + level_value # endfor