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

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....

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.

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

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

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

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) toMxN but won't overly insist on it (but(M+1, N+1) reads nicer thatM+1 x N+1, in particular, I think).

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

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.)

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

erroring is certainly fine with me

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.Line2Dat0x119b13898>]

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.

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.

I guess the whole thing would read slightly nicer if you interpolated rows instead of columns (and then below did the double transpose in thencols == Nx case and not in thenrows = 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)

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

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.

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 readX[0]-dX[0] I immediately thinkX is 1-D.

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.

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