Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

pcolormesh does not work with south pole ortho projection #409

Open
avipersin opened this issue Jun 5, 2018 · 1 comment
Open

pcolormesh does not work with south pole ortho projection #409

avipersin opened this issue Jun 5, 2018 · 1 comment

Comments

@avipersin
Copy link

avipersin commented Jun 5, 2018

The north pole works fine when setting lat_0=90 but when setting lat_0=-90 it doesn't work. The data shown in the southern projection is nonsense.

Related to #347 and #350 but the workaround described there for the "stereo" projection does not work for the "ortho" projection.

Code to replicate:

from mpl_toolkits.basemap import Basemap
import numpy as np
import matplotlib.pyplot as plt


lons = np.arange(-179.95, 180.01, 0.1)
lats = np.arange(-89.95, 90.01, 0.1)

# m = Basemap(projection='ortho', lon_0=-270, lat_0=-90., resolution='c')
m = Basemap(projection='ortho', lon_0=0, lat_0=90., resolution='c')
data = np.random.random((1800, 3600))
lons, lats = np.meshgrid(lons, lats)
x, y = m(lons, lats)

# data = np.ma.masked_where((x < m.xmin) | (x > m.xmax) | (y < m.ymin) | (y > m.ymax), data)

im = m.pcolormesh(x, y, data, latlon=False, cmap='RdBu')

m.colorbar(im)
m.drawcoastlines()
m.drawcountries()
plt.show()
@molinav
Copy link
Member

molinav commented Jan 6, 2024

For the record, if I run your two examples for North Pole and South Pole using the latest basemap 1.3.9 and (infamous) Python 2.7 (and therefore matplotlib 2.2.5), it seems to work fine:
north_pole
south_pole

However, if I use basemap 1.3.9 and Python 3.11 (together with matplotlib 3.7.4), then I get an error and the plots are not drawn:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
Cell In[1], line 17
     13 x, y = m(lons, lats)
     15 # data = np.ma.masked_where((x < m.xmin) | (x > m.xmax) | (y < m.ymin) | (y > m.ymax), data)
---> 17 im = m.pcolormesh(x, y, data, latlon=False, cmap='RdBu')
     19 m.colorbar(im)
     20 m.drawcoastlines()

File C:\ProgramData\pyenv-win\pyenv-win\versions\py311\Lib\site-packages\mpl_toolkits\basemap\__init__.py:543, in _transform.<locals>.with_transform(self, x, y, data, *args, **kwargs)
    541     # convert lat/lon coords to map projection coords.
    542     x, y = self(x,y)
--> 543 return plotfunc(self,x,y,data,*args,**kwargs)

File C:\ProgramData\pyenv-win\pyenv-win\versions\py311\Lib\site-packages\mpl_toolkits\basemap\__init__.py:3451, in Basemap.pcolormesh(self, x, y, data, **kwargs)   
   3449 self._save_use_hold(ax, kwargs)
   3450 try:
-> 3451     ret =  ax.pcolormesh(x,y,data,**kwargs)
   3452 finally:
   3453     self._restore_hold(ax)

File C:\ProgramData\pyenv-win\pyenv-win\versions\py311\Lib\site-packages\matplotlib\__init__.py:1459, in _preprocess_data.<locals>.inner(ax, data, *args, **kwargs) 
   1456 @functools.wraps(func)
   1457 def inner(ax, *args, data=None, **kwargs):
   1458     if data is None:
-> 1459         return func(ax, *map(sanitize_sequence, args), **kwargs)
   1461     bound = new_sig.bind(ax, *args, **kwargs)
   1462     auto_label = (bound.arguments.get(label_namer)
   1463                   or bound.kwargs.get(label_namer))

File C:\ProgramData\pyenv-win\pyenv-win\versions\py311\Lib\site-packages\matplotlib\axes\_axes.py:6218, in Axes.pcolormesh(self, alpha, norm, cmap, vmin, vmax, shading, antialiased, *args, **kwargs)
   6215 shading = shading.lower()
   6216 kwargs.setdefault('edgecolors', 'none')
-> 6218 X, Y, C, shading = self._pcolorargs('pcolormesh', *args,
   6219                                     shading=shading, kwargs=kwargs)
   6220 coords = np.stack([X, Y], axis=-1)
   6221 # convert to one dimensional array, except for 3D RGB(A) arrays

File C:\ProgramData\pyenv-win\pyenv-win\versions\py311\Lib\site-packages\matplotlib\axes\_axes.py:5715, in Axes._pcolorargs(self, funcname, shading, *args, **kwargs)
   5713 if funcname == 'pcolormesh':
   5714     if np.ma.is_masked(X) or np.ma.is_masked(Y):
-> 5715         raise ValueError(
   5716             'x and y arguments to pcolormesh cannot have '
   5717             'non-finite values or be of type '
   5718             'numpy.ma.core.MaskedArray with masked values')
   5719     # safe_masked_invalid() returns an ndarray for dtypes other
   5720     # than floating point.
   5721     if isinstance(X, np.ma.core.MaskedArray):

ValueError: x and y arguments to pcolormesh cannot have non-finite values or be of type numpy.ma.core.MaskedArray with masked values

This seems to appear because at some poing basemap is converting the input arrays into masked arrays, which is not supported by pcolormesh. If we do the input preprocessing in advance with the code snippet below:

outside = (x < m.xmin) | (x > m.xmax) | (y < m.ymin) | (y > m.ymax)
x = np.where(outside, 1E20, x)
y = np.where(outside, 1E20, y)

then your code snippet is also working with newer Python and matplotlib versions, and I think I can now see the bug you were describing:
north_pole_mpl374
south_pole_mpl374

Since both wheels for basemap 1.3.9 on Python 2.7 and Python 3.11 bundle the same GEOS library (3.6.5), I would guess that the source of change is coming from another place (pyproj or matplotlib upgrades). I remember that some upgrades on pyproj have caused unexpected behaviour in the past (see #443 and #463). I would need to investigate this more thoroughly.

@molinav molinav changed the title pcolormesh does not work with south pole ortho projeciton pcolormesh does not work with south pole ortho projection Jan 6, 2024
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants