updated map code to instantly refloat when refloat_halflife == 0.

Author: ChrisBarker-NOAA

py_gnome/gnome/maps/map.py

@@ -9,10 +9,8 @@
 As well as a few implementations of actual maps
 """
 
-# This is a re-write of the C++ raster map approach
-
-# Features:
-#  - Map now handles spillable area and map bounds as polygons
+# RasterMap is the primary implimentation:
+#  - Map handles spillable area and map bounds as polygons
 #  - raster is the same aspect ratio as the land
 #  - internally, raster is a numpy array
 #  - land raster is only as big as the land -- if the map bounds are bigger,
@@ -22,7 +20,6 @@
 #  - Perhaps we just use non-projected coordinates for the raster map?
 #    It makes for a little less computation at every step.
 
-
 import os
 import math
 import warnings
@@ -104,13 +101,13 @@ class GnomeMap(GnomeId):
     """
     The very simplest map for GNOME -- all water
     with only a bounding box for the map bounds.
-    """
+
+    Also the base class for all maps.
     # This also serves as a description of the interface and
     # base class for more complex maps
+    """
 
     _schema = GnomeMapSchema
-
-    refloat_halflife = None  # note -- no land, so never used
     _ref_as = 'map'
 
     def __init__(self,
@@ -154,6 +151,14 @@ def __init__(self,
         self.spillable_area = spillable_area
         self.land_polys = land_polys
 
+    @property
+    def refloat_halflife(self):
+        return self._refloat_halflife / self.seconds_in_hour
+
+    @refloat_halflife.setter
+    def refloat_halflife(self, value):
+        self._refloat_halflife = value * self.seconds_in_hour
+
     def __add__(self, other):
         # Just so pyghnome users will get a helpful message
         # if the do: model.map += a_map
@@ -868,13 +873,13 @@ def raster(self, arr):
         self._raster = np.ascontiguousarray(arr)
         self.build_coarser_rasters()
 
-    @property
-    def refloat_halflife(self):
-        return self._refloat_halflife / self.seconds_in_hour
+    # @property
+    # def refloat_halflife(self):
+    #     return self._refloat_halflife / self.seconds_in_hour
 
-    @refloat_halflife.setter
-    def refloat_halflife(self, value):
-        self._refloat_halflife = value * self.seconds_in_hour
+    # @refloat_halflife.setter
+    # def refloat_halflife(self, value):
+    #     self._refloat_halflife = value * self.seconds_in_hour
 
     @property
     def approximate_raster_interval(self):
@@ -1035,6 +1040,14 @@ def beach_elements(self, sc, model_time=None):
         last_water_positions[beached, :2] = \
             self.projection.to_lonlat(last_water_pos_pixel[beached, :2])
 
+        # zero half-life will not stick at all.
+        # so put the back to last water position right away.
+        if self.refloat_halflife == 0.0:
+            # set the beached elements to the last water position
+            next_pos[beached, :2] = last_water_positions[beached, :2]
+            # set them back to in_water
+            status_codes[beached] = oil_status.in_water
+
         self._set_off_map_status(sc)
 
         # # todo: need a prepare_for_model_run() so map adds these keys to
@@ -1130,16 +1143,18 @@ def to_pixel_array(self, coords):
         return self.projection.to_pixel(coords)
 
 
+# Fixme -- use utilities.convert_longitude
+
 def ShiftLon360(points):
     try:
-        points[points[:,0]<0,0] = points[:,0]+360
+        points[points[:,0] < 0, 0] = points[:, 0] + 360
     except ValueError:
         pass
     return points
 
 def ShiftLon180(points):
     try:
-        points[points[:,0]>180,0] = points[:,0]-360
+        points[points[:,0] > 180, 0] = points[:, 0] - 360
     except ValueError:
         pass
     return points
@@ -1157,6 +1172,7 @@ def __init__(self,
                  raster_size=4096 * 4096,
                  map_bounds=None,
                  spillable_area=None,
+                 land_polys=None,
                  shift_lons=0,
                  **kwargs):
         """
@@ -1177,7 +1193,9 @@ def __init__(self,
                           180, or 360 are valid inputs
         :type shiftLons: integer
 
-        Optional arguments (kwargs):
+        :param shift_lons=0: Whether to shift the longitude reference frame:
+                           Accepted values:
+                           0: do nothing. 180: shift to -180--180. 360: shift to 0--360
 
         :param refloat_halflife: the half-life (in hours) for the re-floating.
 
@@ -1186,6 +1204,8 @@ def __init__(self,
 
         :param spillable_area: The polygon bounding the spillable_area
 
+        Optional arguments (kwargs):
+
         :param id: unique ID of the object. Using UUID as a string.
                    This is only used when loading object from save file.
         :type id: string
@@ -1208,7 +1228,7 @@ def __init__(self,
         land_polys = PolygonSet()  # and lakes....
         spillable_area_bna = PolygonSet()
 
-        #add if based on input param
+        # add if based on input param
         tf = ShiftLon360 if shift_lons == 360 else ShiftLon180 if shift_lons == 180 else None
         if tf is not None:
             polygons.TransformData(tf)
@@ -1257,7 +1277,7 @@ def __init__(self,
         # get the raster as a numpy array:
         raster, projection = self.build_raster(land_polys, BB)
 
-        super(MapFromBNA, self).__init__(
+        super().__init__(
             raster=raster,
             projection=projection,
             map_bounds=map_bounds,
@@ -1299,7 +1319,6 @@ def build_raster(self, land_polys=None, BB=None):
 
         # draw the land to the background
         for poly in land_polys:
-            # fixme -- this should be something like "land"
             if poly.metadata[2] == '1':
                 canvas.draw_polygon(poly,
                                     line_color='land',
@@ -1391,6 +1410,8 @@ def to_geojson(self):
 class MapFromUGrid(RasterMap):
     """
     A raster land-water map, created from netcdf File of a UGrid
+
+    NOTE: not complete or well tested
     """
     _schema = MapFromUGridSchema
 

py_gnome/tests/unit_tests/test_maps/test_map.py

@@ -494,7 +494,7 @@ def test_refloat_some_onland(self):
 
 class Test_MapfromBNA:
 
-    print("instaniating map:", testbnamap)
+    print("instantiating map:", testbnamap)
     # NOTE: this is a pretty course map -- for testing
     bna_map = MapFromBNA(testbnamap, refloat_halflife=6, raster_size=1000)
 
@@ -910,6 +910,181 @@ def test_bna_no_map_bounds():
                                          (6., 10.),
                                          ])
 
+class TestZeroRefloat():
+    """
+    some tests for setting refloat to zero
+
+    Depending on the flag set -- it should
+    either refloat them all at the next timestep, or right away
+
+    """
+    # a very simple raster:
+    (w, h) = (20, 10)
+    raster = np.zeros((w, h), dtype=np.uint8)
+    # a single skinny vertical line:
+    raster[10, :] = 1
+
+    # old style no longer exists ...
+    # def test_zero_refloat_old_style(self):
+    #     """
+    #     test a few LEs
+    #     """
+    #     time_step = 900  # 15 minutes in seconds
+    #     gmap = RasterMap(refloat_halflife=0,
+    #                      raster=self.raster,
+    #                      map_bounds=((-50, -30), (-50, 30),
+    #                                  (50, 30), (50, -30)),
+    #                      projection=NoProjection())
+
+    #     gmap.instant_refloat=False
+    #     # all moving left to right:
+
+    #     spill = sample_sc_release(2)
+
+    #     spill['positions'] = np.array(((5.0, 5.0, 0.),
+    #                                    (5.0, 7.0, 0.),
+    #                                    ),
+    #                                  dtype=np.float64)
+    #     spill['next_positions'] = np.array(((15.0, 5.0, 0.),
+    #                                         (15.0, 7.0, 0.),
+    #                                         ),
+    #                                        dtype=np.float64)
+
+    #     gmap.beach_elements(spill)
+
+    #     pos = spill['positions']
+    #     npos = spill['next_positions']
+    #     lwp = spill['last_water_positions']
+    #     sc = spill['status_codes']
+
+    #     print("after beaching")
+    #     print("before refloating")
+
+    #     print(f"{pos=}")
+    #     print(f"{npos=}")
+    #     print(f"{lwp=}")
+    #     print(f"{sc=}")
+
+    #     # elements should be beached
+    #     assert np.all(sc == oil_status.on_land)
+    #     # last water position should be set to right before the land
+    #     assert np.all(lwp == [[9., 5., 0.],
+    #                          [9., 7., 0.]])
+    #     # next positions are on land
+    #     assert np.all(npos == [[10., 5., 0.],
+    #                           [10., 7., 0.]])
+
+    #     # positions are unchanged
+    #     assert np.all(pos == [[5., 5., 0.],
+    #                           [5., 7., 0.]])
+
+    #     gmap.refloat_elements(spill, time_step)
+
+    #     pos2 = spill['positions']
+    #     npos2 = spill['next_positions']
+    #     lwp2 = spill['last_water_positions']
+    #     sc2 = spill['status_codes']
+
+    #     print("after refloating")
+    #     print(f"{pos2=}")
+    #     print(f"{npos2=}")
+    #     print(f"{lwp2=}")
+    #     print(f"{sc2=}")
+
+    #     # after refloating:
+    #     # elements should not be beached
+    #     assert np.all(sc2 == oil_status.in_water)
+    #     # last water position should still be set to right before the land
+    #     assert np.all(lwp2 == [[9., 5., 0.],
+    #                            [9., 7., 0.]])
+    #     # next positions are irrelevent
+    #     # assert np.all(npos2 == [[10., 5., 0.],
+    #     #                       [10., 7., 0.]])
+
+    #     # positions are at last water position
+    #     assert np.all(pos2 == [[9., 5., 0.],
+    #                            [9., 7., 0.]])
+
+    def test_zero_refloat_instant_refloat(self):
+        """
+        test a few LEs
+        """
+        time_step = 900  # 15 minutes in seconds
+        gmap = RasterMap(refloat_halflife=0,
+                         raster=self.raster,
+                         map_bounds=((-50, -30), (-50, 30),
+                                     (50, 30), (50, -30)),
+                         projection=NoProjection())
+
+        gmap.instant_refloat=True
+
+        # all moving left to right:
+
+        spill = sample_sc_release(2)
+
+        spill['positions'] = np.array(((5.0, 5.0, 0.),
+                                       (5.0, 7.0, 0.),
+                                       ),
+                                     dtype=np.float64)
+        spill['next_positions'] = np.array(((15.0, 5.0, 0.),
+                                            (15.0, 7.0, 0.),
+                                            ),
+                                           dtype=np.float64)
+
+        gmap.beach_elements(spill)
+
+        pos = spill['positions']
+        npos = spill['next_positions']
+        lwp = spill['last_water_positions']
+        sc = spill['status_codes']
+
+        print("after beaching")
+        print("before refloating")
+
+        print(f"{pos=}")
+        print(f"{npos=}")
+        print(f"{lwp=}")
+        print(f"{sc=}")
+
+        # no elements should be beached
+        assert np.all(sc == oil_status.in_water)
+        # last water position should be set to right before the land
+        assert np.all(lwp == [[9., 5., 0.],
+                              [9., 7., 0.]])
+        # next positions are at last water position
+        assert np.all(npos == [[9., 5., 0.],
+                               [9., 7., 0.]])
+        # positions are unchanged
+        assert np.all(pos == [[5., 5., 0.],
+                              [5., 7., 0.]])
+
+        # refloat should not make a difference:
+        gmap.refloat_elements(spill, time_step)
+
+        pos2 = spill['positions']
+        npos2 = spill['next_positions']
+        lwp2 = spill['last_water_positions']
+        sc2 = spill['status_codes']
+
+        print("after refloating")
+        print(f"{pos2=}")
+        print(f"{npos2=}")
+        print(f"{lwp2=}")
+        print(f"{sc2=}")
+
+        # after refloating: -- no difference if none beached
+        # no elements should be beached
+        assert np.all(sc2 == oil_status.in_water)
+        # last water position should be set to right before the land
+        assert np.all(lwp2 == [[9., 5., 0.],
+                              [9., 7., 0.]])
+        # next positions are at last water position
+        assert np.all(npos2 == [[9., 5., 0.],
+                               [9., 7., 0.]])
+        # positions are unchanged
+        assert np.all(pos2 == [[5., 5., 0.],
+                              [5., 7., 0.]])
+
 
 class Test_lake():
     """