Merge pull request #98 from HydrologicEngineeringCenter/cumulus_speedup

Grid write performance boost

Author: oskarhurst

pyproject.toml

@@ -1,6 +1,6 @@
 [tool.poetry]
 name = "hec-dss-python"
-version = "0.1.24"
+version = "0.1.25"
 description = "Python wrapper for the HEC-DSS file database C library."
 authors = ["Hydrologic Engineering Center"]
 license = "MIT"

setup.cfg

@@ -1,6 +1,6 @@
 [metadata]
 name = hecdss
-version = 0.1.24
+version = 0.1.25
 author = Hydrologic Engineering Center
 author_email =hec.dss@usace.army.mil
 description = Python wrapper for the HEC-DSS file database C library.

src/hecdss/gridded_data.py

@@ -102,24 +102,17 @@ def range_limit_table(self, minval, maxval, range_, bins, datasize, data):
         if bins > max_bins:
             bins = max_bins
 
-        self.rangeLimitTable = [0] * bins
-        self.numberEqualOrExceedingRangeLimit = [0] * bins
-
+        self.rangeLimitTable = np.empty(bins, dtype=float)
         self.rangeLimitTable[0] = NULL_INT
         self.rangeLimitTable[1] = minval
 
         step = range_ / bins
-
-        for i in range(2, bins):
-            self.rangeLimitTable[i] = minval + step * i
+        self.rangeLimitTable[2:] = minval + step * np.arange(2, bins)
 
         self.rangeLimitTable[bins - 1] = maxval
         # Exceedance
-        sorted_data = np.sort(data)
-        n = len(sorted_data)
-        for jdx in range(bins):
-            idx = np.searchsorted(sorted_data, self.rangeLimitTable[jdx], side="left")
-            self.numberEqualOrExceedingRangeLimit[jdx] = n - idx
+        mask = data[None, :] >= self.rangeLimitTable[:, None]
+        self.numberEqualOrExceedingRangeLimit = mask.sum(axis=1)
 
     def update_grid_info(self):
         """
@@ -134,8 +127,8 @@ def update_grid_info(self):
         self.meanDataValue = np.nanmean(self.data)
 
         self.data = np.nan_to_num(self.data, nan=NULL_INT)
-        self.numberOfRanges = math.floor(1 + 3.322 * math.log10(n) + 1)
-        flat_data = self.data.ravel()
+        self.numberOfRanges = math.floor(2 + 3.322 * math.log10(n))
+        flat_data = self.data.flatten()
         self.range_limit_table(self.minDataValue, self.maxDataValue, bin_range, self.numberOfRanges, n, flat_data)
 
     @staticmethod

src/hecdss/native.py

@@ -331,7 +331,6 @@ def hec_dss_gridStore(
             self,
             gd,
     ):
-
         self.dll.hec_dss_pdStore.restype = c_int
         self.dll.hec_dss_pdStore.argtypes = [
             ctypes.POINTER(ctypes.c_void_p),  # dss (dss file pointer)
@@ -395,8 +394,9 @@ def hec_dss_gridStore(
         c_rangeLimitTable = (c_float * len(gd.rangeLimitTable))(*gd.rangeLimitTable)
         c_numberEqualOrExceedingRangeLimit = (c_int * len(gd.numberEqualOrExceedingRangeLimit))(
             *gd.numberEqualOrExceedingRangeLimit)
-        flat_list = gd.data.flatten()
-        c_data = (c_float * len(flat_list))(*flat_list)
+
+        arr = gd.data.astype('float32', copy=False)
+        c_data = arr.ctypes.data_as(ctypes.POINTER(ctypes.c_float))
 
         return self.dll.hec_dss_gridStore(self.handle, c_pathname, c_gridType, c_dataType,
                                           c_lowerLeftCellX, c_lowerLeftCellY,