pydata
Maximum value not always included in pcolormesh output when levels are requested
What happened?
When plotting data with "pcolormesh" and requesting levels, if the calculated "vmax" if equal to the upper bound of the last bin, it is not included in the map, as if it was "over" or "missing".
What did you expect to happen?
I expected values == vmax to be painted with the color of the last bin, which is what contourf does.
Minimal Complete Verifiable Example
# /// script
# requires-python = ">=3.11"
# dependencies = [
# "xarray[complete]@git+https://github.com/pydata/xarray.git@main",
# ]
# ///
#
# This script automatically imports the development branch of xarray to check for issues.
# Please delete this header if you have _not_ tested this script with `uv run`!
import xarray as xr
xr.show_versions()
# your reproducer code ...
import xarray as xr
ds = xr.tutorial.open_dataset('air_temperature')
air = ds.air.isel(time=0).round(-1)
air.plot(levels=5)
Steps to reproduce
The example above gives the image below:
You can see that with the heavy discretization I imposed, the values at the bottom are equal to the calculated "vmax" and are not shown.
MVCE confirmation
- [x] Minimal example — the example is as focused as reasonably possible to demonstrate the underlying issue in xarray.
- [x] Complete example — the example is self-contained, including all data and the text of any traceback.
- [x] Verifiable example — the example copy & pastes into an IPython prompt or Binder notebook, returning the result.
- [x] New issue — a search of GitHub Issues suggests this is not a duplicate.
- [x] Recent environment — the issue occurs with the latest version of xarray and its dependencies.
Relevant log output
Anything else we need to know?
Matplotlib's from_levels_and_colors, it is said that :
The quantization levels used to construct the BoundaryNorm. Value v is quantized to level i if lev[i] <= v < lev[i+1].
So this explains why xarray's implementation gives no color. However, the result seems to be unexpected and incoherent with the output of contourf. I think it is logic to expect the maximum value to be included in the plot when simply passing levels as a scalar.
Environment
INSTALLED VERSIONS
commit: None python: 3.13.3 | packaged by conda-forge | (main, Apr 14 2025, 20:44:03) [GCC 13.3.0] python-bits: 64 OS: Linux OS-release: 6.17.4-200.fc42.x86_64 machine: x86_64 processor: byteorder: little LC_ALL: None LANG: fr_CA.UTF-8 LOCALE: ('fr_CA', 'UTF-8') libhdf5: 1.14.6 libnetcdf: 4.9.2
xarray: 2025.10.2.dev18+ge49cfc4f2 pandas: 2.2.3 numpy: 2.2.6 scipy: 1.16.1 netCDF4: 1.7.2 pydap: None h5netcdf: 1.6.1 h5py: 3.13.0 zarr: None cftime: 1.6.4 nc_time_axis: 1.4.1 iris: None bottleneck: 1.5.0 dask: 2025.4.1 distributed: 2025.4.1 matplotlib: 3.10.3 cartopy: None seaborn: None numbagg: 0.9.0 fsspec: 2025.3.2 cupy: None pint: 0.25 sparse: None flox: 0.10.5 numpy_groupies: 0.11.2 setuptools: 80.1.0 pip: 25.1.1 conda: None pytest: 8.3.5 mypy: 1.15.0 IPython: 9.2.0 sphinx: 8.2.3
This seems like a matplotlib feature/bug to me. xarrays sends the same cmap and norm to matplotlib as far as I can tell.
import matplotlib.pyplot as plt
import xarray as xr
plt.show()
xr.show_versions()
# your reproducer code ...
import xarray as xr
ds = xr.tutorial.open_dataset("air_temperature")
air = ds.air.isel(time=0).round(-1)
fig, ax = plt.subplots(1, 1)
air.plot.contourf(levels=5, ax=ax)
fig, ax = plt.subplots(1, 1)
air.plot.pcolormesh(levels=5, ax=ax)
import xarray.plot.utils
ckwargs = xarray.plot.utils._determine_cmap_params(
air.data,
vmin=None,
vmax=None,
cmap=None,
center=None,
robust=False,
extend=None,
levels=5,
filled=True,
norm=None,
_is_facetgrid=False,
)
ckwargs.pop("extend")
ckwargs.pop("levels")
fig, ax = plt.subplots(1, 1)
ax.pcolormesh(air.lon.data, air.lat.data, air.data, **ckwargs)
fig, ax = plt.subplots(1, 1)
ax.contourf(air.lon.data, air.lat.data, air.data, **ckwargs)
@Illviljan , do you think the fix in https://github.com/pydata/xarray/pull/10954 is fine, or should we get this fixed upstream?
