NCAR / wrf-python

A collection of diagnostic and interpolation routines for use with output from the Weather Research and Forecasting (WRF-ARW) Model.
https://wrf-python.readthedocs.io
Apache License 2.0
411 stars 155 forks source link

MercatorWithLatTS no longer compatible with cartopy 0.17 #77

Closed honnorat closed 6 years ago

honnorat commented 6 years ago

Using cartopy 0.17, it is no longer possible to use wrf.get_cartopy() on a Mercator grid with non-zero true latitude

>>> import wrf
>>> import netCDF4 as cdf
>>> wrf.get_cartopy(wrfin=cdf.Dataset("geo_em.d01.nc"))
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/home/honnorat/conda/lib/python3.6/site-packages/wrf/util.py", line 3497, in get_cartopy
    squeeze, cache)
  File "/home/honnorat/conda/lib/python3.6/site-packages/wrf/util.py", line 3362, in _get_proj_obj
    proj_obj = wrf_proj.cartopy()
  File "/home/honnorat/conda/lib/python3.6/site-packages/wrf/projection.py", line 407, in cartopy
    return self._cartopy()
  File "/home/honnorat/conda/lib/python3.6/site-packages/wrf/projection.py", line 705, in _cartopy
    globe = self._globe())
  File "/home/honnorat/conda/lib/python3.6/site-packages/wrf/projection.py", line 81, in __init__
    self._threshold = np.diff(self.x_limits)[0] / 720
  File "/home/honnorat/conda/lib/python3.6/site-packages/cartopy/crs.py", line 1071, in x_limits
    return self._x_limits
AttributeError: 'MercatorWithLatTS' object has no attribute '_x_limits'
$ conda list "(wrf-python|cartopy)"
# packages in environment at /home/honnorat/conda:
#
# Name                    Version                   Build  Channel
cartopy                   0.17.0           py36h81b52dc_0    conda-forge
wrf-python                1.2.0            py36h18b3941_1    conda-forge

The culprit is https://github.com/SciTools/cartopy/commit/9dce9084b67f30f94d8d878ff508a45f69fecd3b (included in v0.17), where @QuLogic changed internal variable cartopy.crs.Mercator._xlimits to cartopy.crs.Mercator._x_limits. Hence MercatorWithLatTS is broken since it modifies this internal state.

I do not see any clean solution to change MercatorWithLatTS while still being compatible with cartopy <= 0.16 and cartopy == 0.17.

QuLogic commented 6 years ago

Cartopy 0.17 support latitude_true_scale in Mercator, so something like:

if cartopy.__version__ >= '0.17':
   class MercatorWithLatTS(ccrs.Mercator):
        pass
else:
    class MercatorWithLatTS(ccrs.Mercator):
        # The old implementation

should work.

honnorat commented 6 years ago

Ok, thanks. PR #78 proposes a quick fix to make wrf-python 0.17 usable with cartopy 0.17.

But beyond that, what is the purpose of MercatorWithLatTS ? What does this class do that ccrs.Mercator doesn't ?

bladwig1 commented 6 years ago

Thanks for the bug report. This was originally written before ccrs.Mercator supported latitude_true_scale. However, there's an additional problem that it fixes. For one of my Mercator test files, the minimum and maximum xlimits get calculated to the same values, which causes a crash when plotting. I'll have to see if this is still an issue before merging the PR. If it is, I'll submit a fix for it and close the PR.

bladwig1 commented 6 years ago

The following example illustrates the problem, and it doesn't look like the new cartopy addresses it:

from cartopy.crs import Mercator, Globe
import numpy as np

stand_lon = np.float32(143.63)
truelat1 = np.float32(39.694)

m = Mercator(central_longitude=stand_lon, 
             latitude_true_scale=truelat1, 
             globe=Globe(ellipse=None, semimajor_axis=6370000, 
                         semiminor_axis=6370000.))

print(m._x_limits)

This prints:

(-15398519.86417731, -15398519.864185866)

These values are invalid for plotting.

For reference, if you create the default Mercator, you get:

(-20037508.342789244, 20037508.342789244)

This is why the MercatorWithLatTS has the sign flipping. It looks like this is still needed.

bladwig1 commented 6 years ago

I just looked at your PR, and it should be fine, so I'll go ahead and merge. Thanks!