Versions Compared

Key

  • This line was added.
  • This line was removed.
  • Formatting was changed.

...

Notice that the point at the top of South America still shows the local minimum that is evident in the cell-centered plots above, thus this method should still preserve gradients and oscillations in the data.

Note that the triangulation method can produce weird results (holes in the plot area) when the data spans a plot boundary. For example, as shown above, we need to “fix” the longitude coordinates to go from -180 to 180, or we get a big hole in the plot. This also happens with polar plots, whether we fix the longitudes first or not. A quick example to demonstrate:

Code Block
vname = 'PS'
data = ds[vname]
lon = ds['lon']
lat = ds['lat']

# Plot using triangulation in lat/lon coordinates
data_proj = crs.PlateCarree()  # data is always assumed to be in lat/lon
projections = (crs.PlateCarree(), crs.Orthographic(central_latitude=60))
for plot_proj in projections:
    figure = pyplot.figure()
    ax = figure.add_subplot(111, projection=plot_proj)
    ax.coastlines(linewidth=0.2)
    ax.set_global()
    # Note, we need fix_lon here using this method
    pl = ax.tripcolor(lon, lat, data, transform=data_proj)

results in the following two figures:

...

We can actually fix both of these examples by first projecting the coordinates explicitly, and then doing the triangulation, rather than doing the triangulation in lat/lon coordinates and then projecting:

Code Block
for plot_proj, proj_name in zip(plot_projs, proj_names):
    figure = pyplot.figure()
    ax = figure.add_subplot(111, projection=plot_proj)
    ax.set_global()
    ax.coastlines(linewidth=0.2)

    # Transform coordinates to target projection
    tcoords = plot_proj.transform_points(data_proj,lon.values,lat.values)
    
    # Subset data
    xi = tcoords[:,0] != numpy.inf
    tc = tcoords[xi,:]
    datai = data[:][xi]  # convert to numpy array, then subset
    
    # Plot in transformed coordinates
    pl = ax.tripcolor(tc[:,0],tc[:,1], datai, shading='flat') # 'gouraud')
    
    # Save figure
    figure.savefig(f'{proj_name}_map-triangulation.png', bbox_inches='tight')

...

The main downside of this method (triangulation) is that it becomes painfully slow for large grids. We may have a couple ways around this. First, we might expect that a significant amount of time is spent calculating the Delaunay triangulation. It turns out this is only partly true, as the following timing results show (using meshes ne4, ne30, ne120, ne240, ne512):

...