ENH: add shading='nearest' and 'auto' to pcolormesh

[jklymak]

PR Summary

Redo of #9629. Lots of discussion there.

See built example here: https://29535-1385122-gh.circle-artifacts.com/0/home/circleci/project/doc/build/html/gallery/images_contours_and_fields/pcolormesh_grids.html#sphx-glr-gallery-images-contours-and-fields-pcolormesh-grids-py

Se new pcolormesh docs here:

https://29535-1385122-gh.circle-artifacts.com/0/home/circleci/project/doc/build/html/api/_as_gen/matplotlib.axes.Axes.pcolormesh.html#matplotlib.axes.Axes.pcolormesh

Suppose len(x) is N, len(y) is M, and size(C) is M, N. Previously, if shading='flat' then pcolor/mesh would silently drop the last row and column of C and x and y would specify the edges of the colorer quadrilaterals.

While compatible with Matlab, this is both arbitrary, and probably not what the user wanted.

This PR changes a few related things:

1. if shading='flat' and the dimensions of x, y and C match, a Deprecation Warning is raised the data is no longer silently dropped.
2. A new shading='nearest' is provided that only allows for x, y and C to match, and centers the colors on the grid points.
3. A new shading='auto' is provided that will choose 'flat' or 'nearest' depending on the size of the data.
4. A new rcParam['pcolor.shading'] is provided that defaults to 'flat' (for back compatibility).

Demo:

import numpy as np
import matplotlib.pyplot as plt

fig, axs = plt.subplots(3, 2, sharex=True, sharey=True)

nx = 4
ny = 5
x = np.arange(nx)
y = np.arange(ny)

X, Y = np.meshgrid(x, y)
np.random.seed(19680808)

C = np.random.rand(ny, nx)

# deprecation warning:
axs[0, 0].pcolormesh(x, y, C, vmin=0, vmax=1, shading='flat')
axs[0, 0].scatter(X[:], Y[:], c=C.flat, ec='w')
axs[0, 0].set_title('flat: Emits deprecation warning')

# no deprecation warning
axs[0, 1].pcolormesh(x, y, C[:-1, :-1], vmin=0, vmax=1, shading='flat')
axs[0, 1].scatter(X[:], Y[:], c=C.flat, ec='w')
axs[0, 1].set_title('flat: No warning')

# no deprecation warning
axs[1, 0].pcolormesh(x, y, C, shading='nearest', vmin=0, vmax=1)
axs[1, 0].scatter(X[:], Y[:], c=C.flat, ec='w')
axs[1, 0].set_title('nearest')

# no deprecation warning
axs[1, 1].pcolormesh(x, y, C, shading='gouraud', vmin=0, vmax=1)
axs[1, 1].scatter(X[:], Y[:], c=C.flat, ec='w')
axs[1, 1].set_title('gouraud')

axs[1, 1].set_xlim([-0.5, 3.5])
axs[1, 1].set_ylim([-0.5, 4.5])

# no deprecation warning:
axs[2, 0].pcolormesh(x, y, C, vmin=0, vmax=1, shading='auto')
axs[2, 0].scatter(X[:], Y[:], c=C.flat, ec='w')
axs[2, 0].set_title('auto: matching size')

# no deprecation warning
axs[2, 1].pcolormesh(x, y, C[:-1, :-1], vmin=0, vmax=1, shading='auto')
axs[2, 1].scatter(X[:], Y[:], c=C.flat, ec='w')
axs[2, 1].set_title('auto: c smaller')

The first call emits:

testpcolormesh.py:17: MatplotlibDeprecationWarning: shading='flat' when X and Y have 
the same dimensions as C is deprecated since 3.3.  Either specify the corners of 
the quadrilaterals with X and Y, or pass shading='auto', 'nearest' or 'gouraud'
or set rcParams['pcolor.shading']
  axs[0, 0].pcolormesh(x, y, C, vmin=0, vmax=1, shading='flat')

PR Checklist

• Has Pytest style unit tests
• Code is Flake 8 compliant
• New features are documented, with examples if plot related
• Documentation is sphinx and numpydoc compliant
• Added an entry to doc/users/next_whats_new/ if major new feature (follow instructions in README.rst there)
• Documented in doc/api/api_changes.rst if API changed in a backward-incompatible way

[jklymak]

Note this fails a bunch of tests because apparently we often call x, y, C with same size and the new Deprecation Warning causes the tests to fail.

[efiring]

It looks like you are missing the mixed cases like C.shape == (Ny+1, Nx). Did we handle them before?

[jklymak]

No we didn't handle that, and I think thats OK for simplicity's sake. If someone has a more complicated grid I feel they can wrangle it themselves....

[efiring]

Actually, we did handle it, because in the "flat" case we always did the indexing (your line 5661), which "worked" with the mixed case. In the code you have here, the mixed case will fail somewhere downstream instead of being caught with an exception, or simply handled as it was originally.

[efiring]

I think you could fix it by changing "and" to "or" in your line 5654, and then out-denting line 5661.

[efiring]

Actually, 5654 could also change "and" to "or". Warn if either dimension isn't right.

[jklymak]

Yep that now emits the warning for:
ax.pcolormesh(x, y[:-1], C[:-1, :-1], vmin=0, vmax=1, shading='flat')

[efiring]

I've only taken a quick look so far, but I think this is on the right track. I'm OK with keeping the "shading" kwarg and adding the options as you have done. For the test failures, it sounds like the thing to do is modify the tests such that the test images don't change.

[jklymak]

The fact that so many tests are hitting the deprecation warning should tell us that even the folks who write examples and tests for Matplotlib make X, Y and C the same size.

[greglucas]

I think this is a great addition! I agree with your decision to update the shading keyword argument rather than adding another keyword for this use case, I think it makes logical sense.

I've just made some minor comments in the code below, but it is close to being ready in my opinion.

[jklymak]

The fact that so many tests are hitting the deprecation warning should tell us that even the folks who write examples and tests for Matplotlib make X, Y and C the same size.

Though I must admit I was responsible for about half of them 😉

[jklymak]

@efiring, and anyone else, this is ready for review I think, but happy to make more changes.

[anntzer]

size -> shape (for an array t, t.shape is its shape, t.size its number of elements, so it's really not the same thing) ("is size" -> "has shape"
(throughout)
I also prefer (M, N) to MxN but won't overly insist on it (but (M+1, N+1) reads nicer that M+1 x N+1, in particular, I think).

[jklymak]

I think I got all the shape/size changes...

[anntzer]

Symmetrically expanding the first and last pixels seem fine but needs to be documented below. Also I guess this will crash if X or Y has size 1? (Sure it's degenerate, but I guess you can either use a width of 1 (likewise needs doc), or error out more cleanly.)

[jklymak]

I guess we should error. Assuming width 1 seems undesirable....

[anntzer]

erroring is certainly fine with me

[efiring]

I think the general mpl strategy, though, is to try hard not to error out when data are simply missing:

In [2]: plt.plot([])
Out[2]: [<matplotlib.lines.Line2D at 0x119b13898>]

The idea is that if you have a long-running process making plots from a stream of data that might have gaps, it is better to show an empty window during a gap than to require that the user trap the gap condition and handle it manually.

[jklymak]

OK, fair enough. Instead of guessing a dX, I just make it [X, X] and nothing visible gets plotted. Anything else seems liable to just be wrong (dX = X/2, has problems if X=0, and just choosing dX=1/2 seems like it could mess with scaling of some plots inadvertently.

[anntzer]

There's something funny going on with sticky_edges and sharex/sharey; not sure whether it's due to this PR or something preexisting:

from pylab import *

Z = np.arange(15).reshape(3, 5)
x = np.arange(6)
y = np.arange(4)

fig, axs = plt.subplots(2, 2, sharex="col", sharey="col")
axs[0, 0].pcolormesh(x, y, Z, shading='auto')
axs[1, 1].pcolormesh(x, y, Z, shading='auto')

plt.show()

results in

i.e. sticky_edges are correctly applied on the top left ("master", from the PoV of sticky edges) but not on the bottom right ("follower").

[anntzer]

I guess the whole thing would read slightly nicer if you interpolated rows instead of columns (and then below did the double transpose in the ncols == Nx case and not in the nrows = Ny case):

if len(X) > 1:
    dX = np.diff(X, axis=0) / 2
    X = np.row_stack([X[0] - dX[0], X[:-1] + dX, X[-1] + dX[-1]])
else:
    X = np.row_stack([X, X])

(untested)

[jklymak]

That only reads nicer if you can remember the rule as to whether X[a] refers to that a-th row or column. I prefer the more explicit form.

[anntzer]

I think it's not too hard if you think of X as "list of lists" (then X[a] is just the a-th list, which is clearly the a-th row). Anyways, not insisting on it, this is fine too.

[jklymak]

which is clearly the a-th row

Why is that clear?

I think when possible, its helpful to be explicit that we are using a 2-D array versus a 1-D. If I read X[0]-dX[0] I immediately think X is 1-D.

[anntzer]

Because

[[1, 2, 3],
 [4, 5, 6]]

is just

[blah,  # row
 blah]  # row

? i.e. literally how numpy writes out the thing?
Anyways, that's a tangent, not really relevant here.

[jklymak]

Yeah, its probably just that I have to switch between python and Fortran too much to keep it straight.

[jklymak]

i.e. sticky_edges are correctly applied on the top left ("master", from the PoV of sticky edges) but not on the bottom right ("follower").

That is also present in master.

[anntzer]

Moved that to #16448.

[anntzer]

lgtm post-ci.

[anntzer]

Let's go :)

[jklymak]

Thanks for everyone's help and patience with this!