Merge remote-tracking branch 'remotes/origin/master'

.gitattributes

@@ -1,2 +1,2 @@
 # Treat Jupyter notebook JSON files as binary
-*.ipynb binary
+# *.ipynb binary

.gitignore

@@ -52,3 +52,10 @@ stats.dat
 .hg
 *.tmp*
 tests/_LargeFileStash
+
+# Specific files
+run.py
+tests/GLWACSO/HSP2results/hspp007.hdf
+tests/GLWACSO/HSPFresults/hspf006.HBN
+tests/GLWACSO/HSP2results/hspp007.uci
+tests/test_report_conversion.html

HSP2/ADCALC.py

@@ -3,7 +3,7 @@
 License: LGPL2
 '''
 
-from numpy import zeros
+from numpy import zeros, array
 from numba import njit
 from HSP2.utilities import make_numba_dict
 
@@ -29,6 +29,7 @@ def adcalc(store, siminfo, uci, ts):
 	nexits = int(ui['NEXITS'])          # table type GEN-INFO
 	ui['simlen'] = siminfo['steps']
 	ui['delts'] = siminfo['delt'] * 60.0     # delts is the simulation interval in seconds
+	ui['uunits']  = siminfo['units']
 
 	# calculated timeseries for advect()
 	if 'SROVOL' not in ts:
@@ -67,13 +68,23 @@ def _adcalc_(ui, ts):
 	simlen = int(ui['simlen'])
 	nexits = int(ui['NEXITS'])
 	delts  = ui['delts']
+	uunits = ui['uunits']
+
+	# units conversion constants, 1 ACRE is 43560 sq ft. assumes input in acre-ft
+	VFACT = 43560.0
+	AFACT = 43560.0
+	if uunits == 2:
+		# si units conversion constants, 1 hectare is 10000 sq m, assumes area input in hectares, vol in Mm3
+		VFACT = 1.0e6
+		AFACT = 10000.0
+
 	# table ADCALC-DATA
 	if 'CRRAT' in ui:
 		crrat = ui['CRRAT']
 	else:
 		crrat = 1.5
 	if 'VOL' in ui:
-		vol = ui['VOL']
+		vol = ui['VOL'] * VFACT
 	else:
 		vol = 0.0
 
@@ -98,6 +109,8 @@ def _adcalc_(ui, ts):
 		for index in range(nexits):
 			ts['EOVOL' + str(index + 1)] = EOVOL[:, index]
 
+	ROS = ui['ROS']
+
 	# external time series
 	O = zeros((simlen, nexits))
 	if nexits > 1:
@@ -107,10 +120,10 @@ def _adcalc_(ui, ts):
 		O[:, 0] = ts['RO']
 
 	for loop in range(simlen):
-		vols = VOL[loop-1] * 43560  if loop > 0 else vol
+		vols = VOL[loop-1] * VFACT  if loop > 0 else vol
 
 		o  = O[loop]
-		os = O[loop-1] if loop > 0 else O[loop]
+		os = O[loop-1] if loop > 0 else array([ROS])
 		ro = 0.0
 		ros= 0.0
 		for index in range(nexits):

HSP2/COPY.py

@@ -0,0 +1,48 @@
+from HSP2.utilities import get_timeseries
+import pandas as pd
+from typing import List, Dict
+
+class Copy():
+    """
+    Partial implementation of the COPY module. 
+    In HSPF, COPY module supports ability able to 'copy' in this case output timeseries to 
+    locations specified in NETWORK block or to EXTERNAL SOURCES. 
+    This functionality is not currently implemented, presently only loading from EXT SOURCES
+    """
+
+    def __init__(self, store: pd.HDFStore, sim_info: Dict, ext_sources: List) -> None:
+        self._ts = {}
+        self._ts['MEAN'] = {}
+        self._ts['POINT'] = {}
+
+        ts = get_timeseries(store, ext_sources, sim_info)
+        for source in ext_sources:
+            themn = source.TMEMN
+            themsb = source.TMEMSB
+            self.set_ts(ts[f'{themn}{themsb}'], themn, themsb)
+
+    def set_ts(self, ts: pd.Series, themn: str, themsb: str) -> None:
+        """Set the provided timeseries to ts dictionary
+        
+        ts: pd.Series
+            pandas Series class instance corresponding to a timeseries
+        tmemn: str, {'MEAN', 'POINT'} 
+            Target member name, specifies if target timeseries is in mean-valued 
+            or point-valued dictionaries
+        tmemsb: str, 
+            Target member name subscripts, acts as key for mean-valued and point-valued dictionaries
+            Original HSPF restricts this to 0-20 but no limit enforced in HSP2
+        """
+        self._ts[themn][themsb] = ts
+
+    def get_ts(self, tmemn: str, tmemsb: str) -> pd.Series:
+        """Gets the specified timeseries from copy class instance based
+        
+        tmemn: str, {'MEAN', 'POINT'} 
+            Target member name, specifies if target timeseries is in mean-valued 
+            or point-valued dictionaries
+        tmemsb: str, 
+            Target member name subscripts, acts as key for mean-valued and point-valued dictionaries
+            Original HSPF restricts this to 0-20 but no limit enforced in HSP2
+        """
+        return self._ts[tmemn][tmemsb]

HSP2/GENER.py

@@ -0,0 +1,149 @@
+from numba import njit
+import numpy as np
+import pandas as pd
+from typing import Dict, List
+
+from HSP2.utilities import get_timeseries
+
+class Gener():
+    """
+    Partial implementation of the GENER module. 
+    Currently supports OPCODES 1-7, 9-23, & 25-26
+    """
+
+    def __init__(self, segment: str, copies: Dict, geners: Dict, ddlinks: Dict, ddgener: Dict) -> None:
+        self.ts_input_1 = pd.Series() # type: pd.Series
+        self. ts_input_2 = pd.Series() # type: pd.Series
+        self.ts_output = pd.Series()  # type: pd.Series
+        self.opcode = 0               # type: int
+        self.k = -1.0E30              # type: float 
+
+        self.opcode = ddgener['OPCODE'][segment]
+        if self.opcode in [9,10,11,24,25,26]:
+            self.k = ddgener['PARM'][segment]
+
+        for link in ddlinks[segment]:
+            ts = pd.Series()
+            if link.SVOL == 'COPY': 
+                copy = copies[link.SVOLNO]
+                ts = copy.get_ts(link.SMEMN,link.SMEMSB1)
+                if link.MFACTOR != 1: ts *= link.MFACTOR
+            elif link.SVOL == 'GENER':
+                gener = geners[link.SVOLNO]
+                ts = gener.get_ts()
+                if link.MFACTOR != 1: ts *= link.MFACTOR
+            else:
+                raise NotImplementedError(f"Invalid SVOL. GENER module does not currently support reading TimeSeries for '{link.SVOL}'")
+
+            if link.TGRPN == 'INPUT' and link.TMEMN == 'ONE':
+                self.ts_input_1 = ts            
+            elif link.TGRPN == 'INPUT' and link.TMEMN == 'TWO':
+                self.ts_input_2 = ts
+            else:
+                raise AttributeError(f"No attribute {link.TGRPN}{link.THEMN} to assign TimeSeries. Should be either 'INPUTONE' or 'INPUTWO'")
+
+        self._execute_gener()
+
+    def get_ts(self) -> pd.Series:
+        """
+        Returns the result TimeSeries generated from executing the operation specified by the OPCODE.
+        """
+        return self.ts_output
+
+    def _execute_gener(self) -> None:
+        gener_op = getattr(self, f'_opcode{self.opcode}')
+        ts_result = gener_op()
+        #May need additional logic here to set default of 1.0E30 to be consistent with FORTRAN code
+        self.ts_output = ts_result
+
+    def _opcode1(self) -> pd.Series:
+        return np.abs(self.ts_input_1)
+
+    def _opcode2(self) -> pd.Series:
+        return np.sqrt(self.ts_input_1)
+
+    def _opcode3(self) -> pd.Series: 
+        return np.trunc(self.ts_input_1)
+
+    def _opcode4(self) -> pd.Series:
+        return np.ceil(self.ts_input_1)
+
+    def _opcode5(self) -> pd.Series:
+        return np.floor(self.ts_input_1)
+
+    def _opcode6(self) -> pd.Series:
+        return np.log(self.ts_input_1)
+
+    def _opcode7(self) -> pd.Series:
+        return np.log10(self.ts_input_1)
+
+    def _opcode8(self) -> pd.Series:
+		#Not presently implemented, read UCI would need to modify to 
+        #process NTERMS and COEFFS sub blocks of GENER block
+        raise NotImplementedError("GENER OPCODE 8 is not currently supported")
+
+    def _opcode9(self) -> pd.Series:
+        return np.power(self.k , self.ts_input_1)
+
+    def _opcode10(self) -> pd.Series:
+        return np.power(self.ts_input_1, self.k)
+
+    def _opcode11(self) -> pd.Series:
+        return np.add(self.ts_input_1, self.k)
+
+    def _opcode12(self) -> pd.Series:
+        return np.sin(self.ts_input_1)
+
+    def _opcode13(self) -> pd.Series:
+        return np.cos(self.ts_input_1)
+
+    def _opcode14(self) -> pd.Series:
+        return np.tan(self.ts_input_1)
+
+    def _opcode15(self) -> pd.Series:
+        return np.cumsum(self.ts_input_1)
+
+    def _opcode16(self) -> pd.Series:
+        return np.add(self.ts_input_1, self.ts_input_2)
+
+    def _opcode17(self) -> pd.Series:
+        return np.subtract(self.ts_input_1, self.ts_input_2)
+
+    def _opcode18(self) -> pd.Series:
+        return np.multiply(self.ts_input_1, self.ts_input_2)
+
+    def _opcode19(self) -> pd.Series:
+        return np.divide(self.ts_input_1, self.ts_input_2)
+
+    def _opcode20(self) -> pd.Series:
+        return np.maximum(self.ts_input_1, self.ts_input_2)
+
+    def _opcode21(self) -> pd.Series:
+        return np.minimum(self.ts_input_1, self.ts_input_2)
+
+    def _opcode22(self) -> pd.Series:
+        return np.power(self.ts_input_1, self.ts_input_2)
+
+    def _opcode23(self) -> pd.Series:
+        ts_out = pd.Series(index=self.ts_input_1.index)
+        s = 0        
+        for idx in self.ts_input_1.index:
+            result = self.ts_input_1[idx] - self.ts_input_2[idx] - s
+            if result < 0:
+                s = s - self.ts_input_1[idx] + self.ts_input_2[idx]
+                result = 0
+            else:
+                s = 0
+            ts_out[idx] = result
+        return ts_out
+
+    def _opcode24(self) -> pd.Series:
+        #skip for now 
+        #would need to figure out timeseries length component
+        raise NotImplementedError("GENER OPCODE 24 is not currently supported")
+
+    def _opcode25(self) -> pd.Series:
+        return np.maximum(self.ts_input_1, self.k)
+
+    def _opcode26(self) -> pd.Series:
+        return np.minimum(self.ts_input_1, self.k)