Mollweide projection - completely corrupted contourf when there are gridpoints at the pole #2137
Description
Description
This is a similar issue to #1333 or #2016, but the effect is much more dramatic, not only limited to the poles.
It seems to happen when there are gridpoints exactly at the poles, not when they are further north than 89.14,
so #1334 might not fix it.
Therefore, I'm thinking that it might be a different problem.
Basically, plotting the exact same data, this is what I'm getting:
Setting the latitude grid to go from -89.9 to 89.9 degrees fixes the problem, but it still seems like a bug.
Code to reproduce
import numpy as np
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
def function(theta, phi):
return 0.5 * (1 + np.cos(theta) ** 2) * np.cos(2 * phi)
nlon = 240
nlat = 120
lon_degrees = np.linspace(0, 360.0, nlon)
lat_degrees = np.linspace(-90., 90., nlat)
lon = np.deg2rad(lon_degrees)
lat = np.deg2rad(lat_degrees)
lon2d, lat2d = np.meshgrid(lon, lat)
data = function(np.pi/2 - lat2d, lon2d)
fig = plt.figure(figsize=(10, 5))
ax1 = fig.add_subplot(2, 2, 1)
mappable = ax1.contourf(
np.rad2deg(lon2d),
np.rad2deg(lat2d),
data,
cmap=plt.get_cmap("seismic"),
levels=np.linspace(-1, 1, num=50),
)
plt.colorbar(mappable, ax=ax1)
ax1.set_title('Bare contourf')
proj = ccrs.Mollweide()
ax2 = fig.add_subplot(2, 2, 2, projection=proj)
mappable = ax2.contourf(
np.rad2deg(lon2d),
np.rad2deg(lat2d),
data,
cmap=plt.get_cmap("seismic"),
levels=np.linspace(-1, 1, num=50),
transform=ccrs.PlateCarree(),
)
ax2.coastlines()
ax2.set_global()
plt.colorbar(mappable, ax=ax2)
ax2.set_title('Mollweide')
proj = ccrs.PlateCarree()
ax3 = fig.add_subplot(2, 2, 3, projection=proj)
mappable = ax3.contourf(
np.rad2deg(lon2d),
np.rad2deg(lat2d),
data,
cmap=plt.get_cmap("seismic"),
levels=np.linspace(-1, 1, num=50),
transform=ccrs.PlateCarree(),
)
ax3.coastlines()
ax3.set_global()
plt.colorbar(mappable, ax=ax3)
ax3.set_title('PlateCarree')
proj = ccrs.Robinson()
ax4 = fig.add_subplot(2, 2, 4, projection=proj)
mappable = ax4.contourf(
np.rad2deg(lon2d),
np.rad2deg(lat2d),
data,
cmap=plt.get_cmap("seismic"),
levels=np.linspace(-1, 1, num=50),
transform=ccrs.PlateCarree(),
)
ax4.coastlines()
ax4.set_global()
plt.colorbar(mappable, ax=ax4)
ax4.set_title('Robinson')Versions
Cartopy version 0.21.0, matplotlib version 3.6.1.