Merge pull request #512 from djhoese/bugfix-ignore-proj-warn-areastr

Ignore pyproj to_proj4 warning when converting an AreaDefinition to a string

pyresample/area_config.py

@@ -150,7 +150,7 @@ def _read_yaml_area_file_content(area_file_name):
     """Read one or more area files in to a single dict object."""
     from pyresample.utils import recursive_dict_update
 
-    if isinstance(area_file_name, (str, pathlib.Path)):
+    if isinstance(area_file_name, (str, pathlib.Path, io.IOBase)):
         area_file_name = [area_file_name]
 
     area_dict = {}

pyresample/future/geometry/area.py

@@ -30,6 +30,7 @@
     get_geostationary_angle_extent,
     get_geostationary_bounding_box_in_lonlats,
     get_geostationary_bounding_box_in_proj_coords,
+    ignore_pyproj_proj_warnings,
 )
 
 

pyresample/geo_filter.py

@@ -63,10 +63,10 @@ def get_valid_index(self, geometry_def):
         # Get projection coords
         proj_kwargs = {}
         if self.nprocs > 1:
-            proj = _spatial_mp.Proj_MP(**self.area_def.proj_dict)
+            proj = _spatial_mp.Proj_MP(self.area_def.crs)
             proj_kwargs["nprocs"] = self.nprocs
         else:
-            proj = Proj(**self.area_def.proj_dict)
+            proj = Proj(self.area_def.crs)
 
         x_coord, y_coord = proj(lons, lats, **proj_kwargs)
 

pyresample/geometry.py

@@ -18,6 +18,7 @@
 """Classes for geometry operations."""
 from __future__ import annotations
 
+import contextlib
 import hashlib
 import math
 import warnings
@@ -28,8 +29,9 @@
 from typing import Optional, Sequence, Union
 
 import numpy as np
+import pyproj
 import yaml
-from pyproj import Geod, Proj, transform
+from pyproj import Geod, Proj
 from pyproj.aoi import AreaOfUse
 
 from pyresample import CHUNK_SIZE
@@ -59,6 +61,25 @@
 HashType = hashlib._hashlib.HASH
 
 
+@contextlib.contextmanager
+def ignore_pyproj_proj_warnings():
+    """Wrap operations that we know will produce a PROJ.4 precision warning.
+
+    Only to be used internally to Pyresample when we have no other choice but
+    to use PROJ.4 strings/dicts. For example, serialization to YAML or other
+    human-readable formats or testing the methods that produce the PROJ.4
+    versions of the CRS.
+
+    """
+    with warnings.catch_warnings():
+        warnings.filterwarnings(
+            "ignore",
+            "You will likely lose important projection information",
+            UserWarning,
+        )
+        yield
+
+
 class DimensionError(ValueError):
     """Wrap ValueError."""
 
@@ -669,14 +690,13 @@ def geocentric_resolution(self, ellps='WGS84', radius=None, nadir_factor=2):
         if self.ndim == 1:
             raise RuntimeError("Can't confidently determine geocentric "
                                "resolution for 1D swath.")
-        from pyproj import transform
         rows = self.shape[0]
         start_row = rows // 2  # middle row
-        src = Proj('+proj=latlong +datum=WGS84')
+        src = CRS('+proj=latlong +datum=WGS84')
         if radius:
-            dst = Proj("+proj=cart +a={} +b={}".format(radius, radius))
+            dst = CRS("+proj=cart +a={} +b={}".format(radius, radius))
         else:
-            dst = Proj("+proj=cart +ellps={}".format(ellps))
+            dst = CRS("+proj=cart +ellps={}".format(ellps))
         # simply take the first two columns of the middle of the swath
         lons = self.lons[start_row: start_row + 1, :2]
         lats = self.lats[start_row: start_row + 1, :2]
@@ -692,7 +712,8 @@ def geocentric_resolution(self, ellps='WGS84', radius=None, nadir_factor=2):
         lats = lats.ravel()
         alt = np.zeros_like(lons)
 
-        xyz = np.stack(transform(src, dst, lons, lats, alt), axis=1)
+        transformer = pyproj.Transformer.from_crs(src, dst)
+        xyz = np.stack(transformer.transform(lons, lats, alt), axis=1)
         dist = np.linalg.norm(xyz[1] - xyz[0])
         dist = dist[np.isfinite(dist)]
         if not dist.size:
@@ -795,7 +816,8 @@ def copy(self):
     @staticmethod
     def _do_transform(src, dst, lons, lats, alt):
         """Run pyproj.transform and stack the results."""
-        x, y, z = transform(src, dst, lons, lats, alt)
+        transformer = pyproj.Transformer.from_crs(src, dst)
+        x, y, z = transformer.transform(lons, lats, alt)
         return np.dstack((x, y, z))
 
     def aggregate(self, **dims):
@@ -807,8 +829,8 @@ def aggregate(self, **dims):
         import dask.array as da
         import pyproj
 
-        geocent = pyproj.Proj(proj='geocent')
-        latlong = pyproj.Proj(proj='latlong')
+        geocent = pyproj.CRS(proj='geocent')
+        latlong = pyproj.CRS(proj='latlong')
         res = da.map_blocks(self._do_transform, latlong, geocent,
                             self.lons.data, self.lats.data,
                             da.zeros_like(self.lons.data), new_axis=[2],
@@ -1641,7 +1663,7 @@ def from_epsg(cls, code, resolution):
         up = max(up1, up2)
         low = min(low1, low2)
         area_extent = (left, low, right, up)
-        return create_area_def(crs.name, crs.to_dict(), area_extent=area_extent, resolution=resolution)
+        return create_area_def(crs.name, crs, area_extent=area_extent, resolution=resolution)
 
     @classmethod
     def from_extent(cls, area_id, projection, shape, area_extent, units=None, **kwargs):
@@ -1889,7 +1911,8 @@ def proj_str(self):
     def __str__(self):
         """Return string representation of the AreaDefinition."""
         # We need a sorted dictionary for a unique hash of str(self)
-        proj_dict = self.proj_dict
+        with ignore_pyproj_proj_warnings():
+            proj_dict = self.proj_dict
         proj_param_str = ', '.join(["'%s': '%s'" % (str(k), str(proj_dict[k])) for k in sorted(proj_dict.keys())])
         proj_str = '{' + proj_param_str + '}'
         if not self.proj_id:
@@ -1953,7 +1976,8 @@ def dump(self, filename=None):
         if self.crs.to_epsg() is not None:
             proj_dict = {'EPSG': self.crs.to_epsg()}
         else:
-            proj_dict = self.crs.to_dict()
+            with ignore_pyproj_proj_warnings():
+                proj_dict = self.crs.to_dict()
 
         res = OrderedDict(description=self.description,
                           projection=OrderedDict(proj_dict),
@@ -2198,7 +2222,7 @@ def colrow2lonlat(self, cols, rows):
         Both scalars and arrays are supported. To be used with scarse
         data points instead of slices (see get_lonlats).
         """
-        p = Proj(self.proj_str)
+        p = Proj(self.crs)
         x = self.projection_x_coords
         y = self.projection_y_coords
         return p(x[cols], y[rows], inverse=True)

pyresample/gradient/__init__.py

@@ -70,7 +70,8 @@ def create_gradient_search_resampler(source_geo_def, target_geo_def):
 @da.as_gufunc(signature='(),()->(),()')
 def transform(x_coords, y_coords, src_prj=None, dst_prj=None):
     """Calculate projection coordinates."""
-    return pyproj.transform(src_prj, dst_prj, x_coords, y_coords)
+    transformer = pyproj.Transformer.from_crs(src_prj, dst_prj)
+    return transformer.transform(x_coords, y_coords)
 
 
 class StackingGradientSearchResampler(BaseResampler):
@@ -104,32 +105,32 @@ def _get_projection_coordinates(self, datachunks):
             try:
                 self.src_x, self.src_y = self.source_geo_def.get_proj_coords(
                     chunks=datachunks)
-                src_prj = pyproj.Proj(**self.source_geo_def.proj_dict)
+                src_crs = self.source_geo_def.crs
                 self.use_input_coords = True
             except AttributeError:
                 self.src_x, self.src_y = self.source_geo_def.get_lonlats(
                     chunks=datachunks)
-                src_prj = pyproj.Proj("+proj=longlat")
+                src_crs = pyproj.CRS.from_string("+proj=longlat")
                 self.use_input_coords = False
             try:
                 self.dst_x, self.dst_y = self.target_geo_def.get_proj_coords(
                     chunks=CHUNK_SIZE)
-                dst_prj = pyproj.Proj(**self.target_geo_def.proj_dict)
+                dst_crs = self.target_geo_def.crs
             except AttributeError:
                 if self.use_input_coords is False:
                     raise NotImplementedError('Cannot resample lon/lat to lon/lat with gradient search.')
                 self.dst_x, self.dst_y = self.target_geo_def.get_lonlats(
                     chunks=CHUNK_SIZE)
-                dst_prj = pyproj.Proj("+proj=longlat")
+                dst_crs = pyproj.CRS.from_string("+proj=longlat")
             if self.use_input_coords:
                 self.dst_x, self.dst_y = transform(
                     self.dst_x, self.dst_y,
-                    src_prj=dst_prj, dst_prj=src_prj)
+                    src_prj=dst_crs, dst_prj=src_crs)
                 self.prj = pyproj.Proj(**self.source_geo_def.proj_dict)
             else:
                 self.src_x, self.src_y = transform(
                     self.src_x, self.src_y,
-                    src_prj=src_prj, dst_prj=dst_prj)
+                    src_prj=src_crs, dst_prj=dst_crs)
                 self.prj = pyproj.Proj(**self.target_geo_def.proj_dict)
 
     def _get_src_poly(self, src_y_start, src_y_end, src_x_start, src_x_end):