Skip to content

Commit

Permalink
Merge pull request #526 from mraspaud/fix-frequency
Browse files Browse the repository at this point in the history
  • Loading branch information
djhoese authored Aug 15, 2023
2 parents dafbb34 + bea2050 commit 04ea38d
Showing 1 changed file with 56 additions and 30 deletions.
86 changes: 56 additions & 30 deletions pyresample/geometry.py
Original file line number Diff line number Diff line change
Expand Up @@ -279,13 +279,16 @@ def get_boundary_lonlats(self):
return (SimpleBoundary(s1_lon.squeeze(), s2_lon.squeeze(), s3_lon.squeeze(), s4_lon.squeeze()),
SimpleBoundary(s1_lat.squeeze(), s2_lat.squeeze(), s3_lat.squeeze(), s4_lat.squeeze()))

def get_bbox_lonlats(self, frequency: Optional[int] = None, force_clockwise: bool = True) -> tuple:
def get_bbox_lonlats(self, vertices_per_side: Optional[int] = None, force_clockwise: bool = True,
frequency: Optional[int] = None) -> tuple:
"""Return the bounding box lons and lats.
Args:
frequency:
vertices_per_side:
The number of points to provide for each side. By default (None)
the full width and height will be provided.
frequency:
Deprecated, use vertices_per_side
force_clockwise:
Perform minimal checks and reordering of coordinates to ensure
that the returned coordinates follow a clockwise direction.
Expand All @@ -311,15 +314,19 @@ def get_bbox_lonlats(self, frequency: Optional[int] = None, force_clockwise: boo
pyresample (ex. :class:`pyresample.spherical.SphPolygon`).
"""
lons, lats = self._get_bbox_elements(self.get_lonlats, frequency)
if frequency is not None:
warnings.warn("The `frequency` argument is pending deprecation, use `vertices_per_side` instead",
PendingDeprecationWarning, stacklevel=2)
vertices_per_side = vertices_per_side or frequency
lons, lats = self._get_bbox_elements(self.get_lonlats, vertices_per_side)
if force_clockwise and not self._corner_is_clockwise(
lons[0][-2], lats[0][-2], lons[0][-1], lats[0][-1], lons[1][1], lats[1][1]):
# going counter-clockwise
lons, lats = self._reverse_boundaries(lons, lats)
return lons, lats

def _get_bbox_elements(self, coord_fun, frequency: Optional[int] = None) -> tuple:
s1_slice, s2_slice, s3_slice, s4_slice = self._get_bbox_slices(frequency)
def _get_bbox_elements(self, coord_fun, vertices_per_side: Optional[int] = None) -> tuple:
s1_slice, s2_slice, s3_slice, s4_slice = self._get_bbox_slices(vertices_per_side)
s1_dim1, s1_dim2 = coord_fun(data_slice=s1_slice)
s2_dim1, s2_dim2 = coord_fun(data_slice=s2_slice)
s3_dim1, s3_dim2 = coord_fun(data_slice=s3_slice)
Expand Down Expand Up @@ -348,14 +355,14 @@ def _filter_bbox_nans(
new_dim2_sides.append(dim2_side[is_valid_mask])
return new_dim1_sides, new_dim2_sides

def _get_bbox_slices(self, frequency):
def _get_bbox_slices(self, vertices_per_side):
height, width = self.shape
if frequency is None:
if vertices_per_side is None:
row_num = height
col_num = width
else:
row_num = frequency
col_num = frequency
row_num = vertices_per_side
col_num = vertices_per_side
s1_slice = (0, np.linspace(0, width - 1, col_num, dtype=int))
s2_slice = (np.linspace(0, height - 1, row_num, dtype=int), -1)
s3_slice = (-1, np.linspace(width - 1, 0, col_num, dtype=int))
Expand Down Expand Up @@ -398,25 +405,34 @@ def _corner_is_clockwise(lon1, lat1, corner_lon, corner_lat, lon2, lat2):
is_clockwise = -np.pi < angle < 0
return is_clockwise

def get_edge_lonlats(self, frequency=None):
def get_edge_lonlats(self, vertices_per_side=None, frequency=None):
"""Get the concatenated boundary of the current swath."""
lons, lats = self.get_bbox_lonlats(frequency=frequency, force_clockwise=False)
if frequency is not None:
warnings.warn("The `frequency` argument is pending deprecation, use `vertices_per_side` instead",
PendingDeprecationWarning, stacklevel=2)
vertices_per_side = vertices_per_side or frequency
lons, lats = self.get_bbox_lonlats(vertices_per_side=vertices_per_side, force_clockwise=False)
blons = np.ma.concatenate(lons)
blats = np.ma.concatenate(lats)
return blons, blats

def get_edge_bbox_in_projection_coordinates(self, frequency: Optional[int] = None):
def get_edge_bbox_in_projection_coordinates(self, vertices_per_side: Optional[int] = None,
frequency: Optional[int] = None):
"""Return the bounding box in projection coordinates."""
x, y = self._get_bbox_elements(self.get_proj_coords, frequency)
if frequency is not None:
warnings.warn("The `frequency` argument is pending deprecation, use `vertices_per_side` instead",
PendingDeprecationWarning, stacklevel=2)
vertices_per_side = vertices_per_side or frequency
x, y = self._get_bbox_elements(self.get_proj_coords, vertices_per_side)
return np.hstack(x), np.hstack(y)

def boundary(self, frequency=None, force_clockwise=False):
def boundary(self, vertices_per_side=None, force_clockwise=False, frequency=None):
"""Retrieve the AreaBoundary object.
Parameters
----------
frequency:
The number of points to provide for each side. By default (None)
vertices_per_side:
(formerly `frequency`) The number of points to provide for each side. By default (None)
the full width and height will be provided.
force_clockwise:
Perform minimal checks and reordering of coordinates to ensure
Expand All @@ -428,7 +444,11 @@ def boundary(self, frequency=None, force_clockwise=False):
Default is False.
"""
from pyresample.boundary import AreaBoundary
lon_sides, lat_sides = self.get_bbox_lonlats(frequency=frequency,
if frequency is not None:
warnings.warn("The `frequency` argument is pending deprecation, use `vertices_per_side` instead",
PendingDeprecationWarning, stacklevel=2)
vertices_per_side = vertices_per_side or frequency
lon_sides, lat_sides = self.get_bbox_lonlats(vertices_per_side=vertices_per_side,
force_clockwise=force_clockwise)
return AreaBoundary.from_lonlat_sides(lon_sides, lat_sides)

Expand Down Expand Up @@ -1550,20 +1570,20 @@ def is_geostationary(self):
return False
return 'geostationary' in coord_operation.method_name.lower()

def _get_geo_boundary_sides(self, frequency=None):
def _get_geo_boundary_sides(self, vertices_per_side=None):
"""Retrieve the boundary sides list for geostationary projections."""
# Define default frequency
if frequency is None:
frequency = 50
if vertices_per_side is None:
vertices_per_side = 50
# Ensure at least 4 points are used
if frequency < 4:
frequency = 4
if vertices_per_side < 4:
vertices_per_side = 4
# Ensure an even number of vertices for side creation
if (frequency % 2) != 0:
frequency = frequency + 1
lons, lats = get_geostationary_bounding_box_in_lonlats(self, nb_points=frequency)
if (vertices_per_side % 2) != 0:
vertices_per_side = vertices_per_side + 1
lons, lats = get_geostationary_bounding_box_in_lonlats(self, nb_points=vertices_per_side)
# Retrieve dummy sides for GEO (side1 and side3 always of length 2)
side02_step = int(frequency / 2) - 1
side02_step = int(vertices_per_side / 2) - 1
lon_sides = [lons[slice(0, side02_step + 1)],
lons[slice(side02_step, side02_step + 1 + 1)],
lons[slice(side02_step + 1, side02_step * 2 + 1 + 1)],
Expand All @@ -1576,15 +1596,17 @@ def _get_geo_boundary_sides(self, frequency=None):
]
return lon_sides, lat_sides

def boundary(self, frequency=None, force_clockwise=False):
def boundary(self, *, vertices_per_side=None, force_clockwise=False, frequency=None):
"""Retrieve the AreaBoundary object.
Parameters
----------
frequency:
vertices_per_side:
The number of points to provide for each side. By default (None)
the full width and height will be provided, except for geostationary
projection where by default only 50 points are selected.
frequency:
Deprecated, use vertices_per_side
force_clockwise:
Perform minimal checks and reordering of coordinates to ensure
that the returned coordinates follow a clockwise direction.
Expand All @@ -1595,10 +1617,14 @@ def boundary(self, frequency=None, force_clockwise=False):
Default is False.
"""
from pyresample.boundary import AreaBoundary
if frequency is not None:
warnings.warn("The `frequency` argument is pending deprecation, use `vertices_per_side` instead",
PendingDeprecationWarning, stacklevel=2)
vertices_per_side = vertices_per_side or frequency
if self.is_geostationary:
lon_sides, lat_sides = self._get_geo_boundary_sides(frequency=frequency)
lon_sides, lat_sides = self._get_geo_boundary_sides(vertices_per_side=vertices_per_side)
else:
lon_sides, lat_sides = self.get_bbox_lonlats(frequency=frequency,
lon_sides, lat_sides = self.get_bbox_lonlats(vertices_per_side=vertices_per_side,
force_clockwise=force_clockwise)
boundary = AreaBoundary.from_lonlat_sides(lon_sides, lat_sides)
return boundary
Expand Down

0 comments on commit 04ea38d

Please sign in to comment.