ENH: Dimension-dependent weighting factor flexibility in gridding (#1545

)

* ENH: Dimension-dependent weighting factor flexibility in gridding
(previously, weighting scaling was optional only for the z-dim)

* TST: modify tests to align with new mapping scaling factor flexibility

* FIX: dist2 bug fix + STY

* FIX: optional conversion of  h_factor and dist_factor to array (compatibility)

examples/mapping/plot_grid_single_sweep_ppi.py

@@ -81,7 +81,7 @@
 
 print(
     "===\n\nUsing 2D weighting by "
-    "setting h_factor and zdist_factor to 0.0, the max range looks OK now\n\n==="
+    "setting h_factor and dist_factor z component to 0.0, the max range looks OK now\n\n==="
 )
 fig2 = plt.figure(figsize=(12, 6), tight_layout=True)
 fig2.suptitle("Cartesian gridded")
@@ -92,8 +92,8 @@
         grid_shape=(1, 1601, 1601),
         grid_limits=((0, 10000.0), [-40000, 40000], [-40000, 40000]),
         fields=[processed_field],
-        h_factor=0.0,
-        zdist_factor=0.0,
+        h_factor=(0.0, 1.0, 1.0),
+        dist_factor=(0.0, 1.0, 1.0),
     )
     ax = fig2.add_subplot(1, 2, ax_ind + 1)
     ax.imshow(

pyart/map/_gate_to_grid_map.pyx

@@ -96,11 +96,12 @@ cdef class DistBeamRoI(RoIFunction):
     Radius of influence which expands with distance from multiple radars.
     """
 
-    cdef float h_factor, min_radius, beam_factor
+    cdef float min_radius, beam_factor
     cdef int num_offsets
     cdef float[:, :] offsets
+    cdef float[::1] h_factor
 
-    def __init__(self, h_factor, nb, bsp, min_radius, offsets):
+    def __init__(self, float[::1] h_factor, float nb, float bsp, float min_radius, offsets):
         """ initalize. """
         cdef int i
         self.h_factor = h_factor
@@ -133,8 +134,9 @@ cdef class DistBeamRoI(RoIFunction):
             z_offset = self.offsets[i, 0]
             y_offset = self.offsets[i, 1]
             x_offset = self.offsets[i, 2]
-            roi = (self.h_factor * ((z - z_offset) / 20.0) +
-                   sqrt((y - y_offset)**2 + (x - x_offset)**2) *
+            roi = (self.h_factor[0] * ((z - z_offset) / 20.0) +
+                   sqrt((self.h_factor[1] * (y - y_offset))**2 +
+                        (self.h_factor[2] * (x - x_offset))**2) *
                    self.beam_factor)
             if roi < self.min_radius:
                 roi = self.min_radius
@@ -233,7 +235,7 @@ cdef class GateToGridMapper:
             float[:, :, ::1] field_data,
             char[:, :, ::1] field_mask, char[:, ::1] excluded_gates,
             RoIFunction roi_func, int weighting_function,
-            float zdist_factor):
+            float[::1] dist_factor):
         """
         Map radar gates unto the regular grid.
 
@@ -265,10 +267,11 @@ cdef class GateToGridMapper:
             Function to use for weighting gates based upon distance.
             0 for Barnes, 1 for Cressman, 2 for Nearest and 3 for Barnes 2
             neighbor weighting.
-        zdist_factor: float
-            Scaling factor for squared z difference in distance calculation.
-            A value of 0.0 combined with an h_factor=0.0 (if calling
-            DistBeamRoI) or z_factor=0.0 (if calling DistRoI) results in
+        dist_factor: 3-element float32 array
+            Scaling factors for squared z,y,x difference in distance calculation.
+            For example:
+            A value of (0.0, 1.0, 1.0)  combined with an h_factor=(0.0, 1.0, 1.0)
+            (if calling DistBeamRoI) or z_factor=0.0 (if calling DistRoI) results in
             the exclusion of the z dimension in gridding weighting and could
             serve as a potential solution for gridding a single PPI sweep.
 
@@ -294,15 +297,15 @@ cdef class GateToGridMapper:
                 masks = field_mask[nray, ngate]
 
                 self.map_gate(x, y, z, roi, values, masks, weighting_function,
-                              zdist_factor)
+                              dist_factor)
 
     @cython.initializedcheck(False)
     @cython.cdivision(True)
     @cython.boundscheck(False)
     @cython.wraparound(False)
     cdef int map_gate(self, float x, float y, float z, float roi,
                       float[:] values, char[:] masks,
-                      int weighting_function, float zdist_factor):
+                      int weighting_function, float[:] dist_factor):
         """ Map a single gate to the grid. """
 
         cdef float xg, yg, zg, weight, roi2, dist2, min_dist2
@@ -348,7 +351,9 @@ cdef class GateToGridMapper:
                         xg = self.x_step * xi
                         yg = self.y_step * yi
                         zg = self.z_step * zi
-                        dist2 = ((xg - x)**2 + (yg - y)**2 + zdist_factor * (zg - z)**2)
+                        dist2 = (dist_factor[2] * (xg - x)**2 +
+                                 dist_factor[1] * (yg - y)**2 +
+                                 dist_factor[0] * (zg - z)**2)
                         if dist2 >= roi2:
                             continue
                         for i in range(self.nfields):
@@ -370,7 +375,9 @@ cdef class GateToGridMapper:
                         xg = self.x_step * xi
                         yg = self.y_step * yi
                         zg = self.z_step * zi
-                        dist2 = (xg-x)*(xg-x) + (yg-y)*(yg-y) + zdist_factor * (zg-z)*(zg-z)
+                        dist2 = (dist_factor[2] * (xg-x)*(xg-x) +
+                                 dist_factor[1] * (yg-y)*(yg-y) +
+                                 dist_factor[0] * (zg-z)*(zg-z))
 
                         if dist2 > roi2:
                             continue

pyart/map/gates_to_grid.py

@@ -36,10 +36,10 @@ def map_gates_to_grid(
     z_factor=0.05,
     xy_factor=0.02,
     min_radius=500.0,
-    h_factor=1.0,
+    h_factor=(1.0, 1.0, 1.0),
     nb=1.5,
     bsp=1.0,
-    zdist_factor=1.0,
+    dist_factor=(1.0, 1.0, 1.0),
     **kwargs
 ):
     """
@@ -105,6 +105,16 @@ def map_gates_to_grid(
         except TypeError:
             grid_origin_alt = np.mean(radars[0].altitude["data"])
 
+    # convert input h_factor and dist_factor from scalar, tuple, or list to array
+    if isinstance(h_factor, (tuple, list)):
+        h_factor = np.array(h_factor, dtype="float32")
+    elif isinstance(h_factor, float):
+        h_factor = np.full(3, h_factor, dtype="float32")
+    if isinstance(dist_factor, (tuple, list)):
+        dist_factor = np.array(dist_factor, dtype="float32")
+    elif isinstance(dist_factor, float):
+        dist_factor = np.full(3, dist_factor, dtype="float32")
+
     gatefilters = _parse_gatefilters(gatefilters, radars)
     cy_weighting_function = _detemine_cy_weighting_func(weighting_function)
     projparams = _find_projparams(grid_origin, radars, grid_projection)
@@ -175,7 +185,7 @@ def map_gates_to_grid(
             excluded_gates,
             roi_func,
             cy_weighting_function,
-            zdist_factor,
+            dist_factor,
         )
 
     # create and return the grid dictionary

tests/map/test_gates_to_grid.py

@@ -103,7 +103,7 @@ def test_map_to_grid_dist_beam_roi():
         fields=["reflectivity"],
         min_radius=30,
         bsp=0.0,
-        h_factor=0.0,
+        h_factor=(0.0, 1.0, 1.0),
     )
     center_slice = grids["reflectivity"][1, 4, :]
     assert_almost_equal(np.round(center_slice), EXPECTED_CENTER_SLICE)

tests/map/test_grid_mapper.py

@@ -103,7 +103,7 @@ def test_map_to_grid_dist_beam_roi():
         fields=["reflectivity"],
         min_radius=30,
         bsp=0.0,
-        h_factor=0.0,
+        h_factor=(0.0, 1.0, 1.0),
     )
     center_slice = grids["reflectivity"][1, 4, :]
     assert_almost_equal(center_slice, EXPECTED_CENTER_SLICE)