Hey @Try2Code ! Thanks again for the test data 😃 I played around a bit and indeed, the visualization using the PolyCollection is pretty much straightforward.

import xarray as xr
from matplotlib.collections import PolyCollection
import matplotlib.pyplot as plt
from matplotlib.colors import Normalize

fig, ax = plt.subplots()
ds = xr.open_dataset('massFlux_r2b6.nc')
xb = ds.elon_bnds.values
yb = ds.elat_bnds.values
coll = PolyCollection(
    np.dstack([xb, yb]), array=ds.mass_flux[0, 0].values,
    norm=Normalize(-1, 1))
ax.add_collection(coll)
ax.set_xlim(-np.pi, np.pi)
ax.set_ylim(-np.pi/2, np.pi/2)
plt.show()

The output image however does not look good because of the same reason that makes the triangular plotting not straight-forward. The problem here comes from the grid cells at the boundaries of our coordinate system, because they have vertices with x-coordinates close to -π as well as x-coordinates close to π.

So if we exclude these cells, the image looks much nicer

mask = ~(((xb < -1).any(axis=1) & (xb > 1).any(axis=1)) |
         ((yb < -1).any(axis=1) & (yb > 1).any(axis=1)))
coll2 = PolyCollection(
    np.dstack([xb[mask], yb[mask]]), array=ds.mass_flux[0, 0].values[mask],
    norm=Normalize(-1, 1))
ax.add_collection(coll2)
coll.set_visible(False)
plt.draw()

You can also use the projections from cartopy

import cartopy.crs as ccrs
fig, ax = plt.subplots(subplot_kw=dict(projection=ccrs.PlateCarree()))
ax.coastlines()
coll2 = PolyCollection(
    np.dstack([xb[mask], yb[mask]]) * 180/np.pi, array=ds.mass_flux[0, 0].values[mask],
    norm=Normalize(-1, 1), transform=ccrs.PlateCarree())
ax.add_collection(coll2)
ax.set_global()
plt.show()

and then you need a separate collection for the corners at the boundaries, e.g. via

xb2 = xb[~mask]
yb2 = yb[~mask]
xb2[xb2 > 1] -= 2*np.pi
boundaries_coll = PolyCollection(
    np.dstack([xb2, yb2]) * 180/np.pi, array=ds.mass_flux[0, 0].values[~mask],
    norm=Normalize(-1, 1), transform=ccrs.PlateCarree())
ax.add_collection(boundaries_coll)
plt.draw()

A closeup then looks like

This methodology would be straight-forward to be implemented in psyplot and it works for any type of grid, whether you 3, 4 or more vertices. However, as I assumed in the beginning, it is pretty slow. On my laptop, making the last figure, takes about 20 seconds.

So my question to you is: Do you think, it would be useful to have this implemented in psyplot? If yes, I would work on it during the next week.

plotting other non-triangular grids #6

Description

Metadata

Metadata

Assignees

Labels

Type

Projects

Milestone

Relationships

Development

Issue actions