Commit 3b41a20e authored by POLCHER Jan's avatar POLCHER Jan 🚴🏾
Browse files

First set of adjusted parameter and new calculation of maximum acceptable time step.

parent d0a8f271
......@@ -25,6 +25,43 @@ import RPPtools as RPP
#
NCFillValue=1.0e20
#
# Function to compute quantiles with weights
#
def weighted_quantile(values, quantiles, sample_weight=None,
values_sorted=False, old_style=False):
""" Very close to numpy.percentile, but supports weights.
NOTE: quantiles should be in [0, 1]!
:param values: numpy.array with data
:param quantiles: array-like with many quantiles needed
:param sample_weight: array-like of the same length as `array`
:param values_sorted: bool, if True, then will avoid sorting of
initial array
:param old_style: if True, will correct output to be consistent
with numpy.percentile.
:return: numpy.array with computed quantiles.
"""
values = np.array(values)
quantiles = np.array(quantiles)
if sample_weight is None:
sample_weight = np.ones(len(values))
sample_weight = np.array(sample_weight)
assert np.all(quantiles >= 0) and np.all(quantiles <= 1), \
'quantiles should be in [0, 1]'
if not values_sorted:
sorter = np.argsort(values)
values = values[sorter]
sample_weight = sample_weight[sorter]
weighted_quantiles = np.cumsum(sample_weight) - 0.5 * sample_weight
if old_style:
# To be convenient with numpy.percentile
weighted_quantiles -= weighted_quantiles[0]
weighted_quantiles /= weighted_quantiles[-1]
else:
weighted_quantiles /= np.sum(sample_weight)
return np.interp(quantiles, weighted_quantiles, values)
#
# Add time constants and possibly other prameters to the graph file
#
def addparameters(outnf, version, part, tcst, vtyp, NCFillValue) :
......@@ -204,7 +241,7 @@ class HydroGraph :
ncvar.units = "count"
ncvar[:] = thist[:]
#
maxtstep = np.quantile(gvar[~np.isnan(gvar)], 0.1)
maxtstep = weighted_quantile(gvar[~np.isnan(gvar)]*tcst.stream_tcst, 0.50, sample_weight=garea[~np.isnan(gvar)])
mtstep = outnf.createVariable("MaxTimeStep", vtyp, ('dimpara'), fill_value=NCFillValue)
mtstep.title = "The maximum time-step possible given the distribution of HTU properties"
mtstep.units = "s"
......
......@@ -317,10 +317,10 @@ class HydroParameter :
self.swamp_cst = 0.2
elif hydrogrid.hd >= 0.016 :
# Case for MERIT
self.stream_tcst = 0.0014/1000*RPP.OneDay
self.fast_tcst = 0.07/1000*RPP.OneDay
self.slow_tcst = 0.98/1000*RPP.OneDay
self.flood_tcst = 4.0/1000*RPP.OneDay
self.stream_tcst = 2.6
self.fast_tcst = 550.0
self.slow_tcst = 7500.0
self.flood_tcst = 345.6
self.swamp_cst = 0.2
elif hydrogrid.hd >= 0.008 :
# Case for HydroSEHDS
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment