An error occurred while loading the file. Please try again.
-
Thibault Hallouin authoredda78292b
#include <xtensor/xmath.hpp>
#include <xtensor/xview.hpp>
#include <xtensor/xsort.hpp>
#include "probabilist/evaluator.hpp"
namespace evalhyd
{
namespace probabilist
{
// Determine observed realisation of threshold(s) exceedance.
//
// \require q_obs:
// Streamflow observations.
// shape: (time,)
// \require q_thr:
// Streamflow exceedance threshold(s).
// shape: (thresholds,)
// \assign o_k:
// Event observed outcome.
// shape: (thresholds, time)
void Evaluator::calc_o_k()
{
// determine observed realisation of threshold(s) exceedance
o_k = q_obs >= xt::view(q_thr, xt::all(), xt::newaxis());
}
// Determine mean observed realisation of threshold(s) exceedance.
//
// \require o_k:
// Event observed outcome.
// shape: (thresholds, time)
// \require t_msk:
// Temporal subsets of the whole record.
// shape: (subsets, time)
// \assign bar_o:
// Mean event observed outcome.
// shape: (subsets, samples, thresholds)
void Evaluator::calc_bar_o()
{
// apply the mask
// (using NaN workaround until reducers work on masked_view)
auto o_k_masked = xt::where(
xt::view(t_msk, xt::all(), xt::newaxis(), xt::all()),
o_k, NAN
);
// compute variable one sample at a time
bar_o = xt::zeros<double>({n_msk, n_exp, n_thr});
for (int e = 0; e < n_exp; e++)
{
// apply the bootstrap sampling
auto o_k_masked_sampled =
xt::view(o_k_masked, xt::all(), xt::all(), b_exp[e]);
// compute mean "climatology" relative frequency of the event
// $\bar{o} = \frac{1}{n} \sum_{k=1}^{n} o_k$
xt::view(bar_o, xt::all(), e, xt::all()) =
xt::nanmean(o_k_masked_sampled, -1);
}
}
// Determine forecast probability of threshold(s) exceedance to occur.
//
// \require q_prd:
// Streamflow predictions.
// shape: (members, time)
// \require q_thr:
// Streamflow exceedance threshold(s).