Failed to fetch fork details. Try again later.
-
Guillaume Blanchy authoredce3939f7
Forked from
reversaal / OhmPi
Source project has a limited visibility.
#ifndef EVALHYD_DETERMINIST_EVALUATOR_HPP
#define EVALHYD_DETERMINIST_EVALUATOR_HPP
#include <vector>
#include <xtensor/xexpression.hpp>
#include <xtensor/xtensor.hpp>
#include "../maths.hpp"
namespace evalhyd
{
namespace determinist
{
class Evaluator
{
private:
// members for input data
const xt::xtensor<double, 2>& q_obs;
const xt::xtensor<double, 2>& q_prd;
xt::xtensor<bool, 3> t_msk;
const std::vector<xt::xkeep_slice<int>>& b_exp;
// members for dimensions
std::size_t n_tim;
std::size_t n_msk;
std::size_t n_srs;
std::size_t n_exp;
std::vector<std::size_t> inner_dims;
std::vector<std::size_t> inter_dims;
std::vector<std::size_t> mean_dims;
std::vector<std::size_t> final_dims;
// members for computational elements
xt::xtensor<double, 4> mean_obs;
xt::xtensor<double, 4> mean_prd;
xt::xtensor<double, 2> quad_err;
xt::xtensor<double, 4> quad_obs;
xt::xtensor<double, 4> quad_prd;
xt::xtensor<double, 3> r_pearson;
xt::xtensor<double, 3> alpha;
xt::xtensor<double, 3> bias;
public:
// constructor method
Evaluator(const xt::xtensor<double, 2>& obs,
const xt::xtensor<double, 2>& prd,
const xt::xtensor<bool, 2>& msk,
const std::vector<xt::xkeep_slice<int>>& exp) :
q_obs{obs}, q_prd{prd}, b_exp{exp}
{
// determine size for recurring dimensions
n_tim = q_prd.shape(1);
n_srs = q_prd.shape(0);
n_msk = msk.shape(0);
n_exp = b_exp.size();
// determine dimensions for inner components, intermediate
// elements, for time average elements, and for final metrics:
// -> inner components shape: (samples, subsets, series)
inner_dims = {n_msk, n_exp, n_srs};
// -> intermediate elements shape: (samples, subsets, series, time)
inter_dims = {n_msk, n_exp, n_srs, n_tim};
// -> time average elements shape: (samples, subsets, series, 1)
mean_dims = {n_msk, n_exp, n_srs, 1};
// -> final metrics shape: (series, subsets, samples)
final_dims = {n_srs, n_msk, n_exp};
// drop time steps where observations or predictions are NaN
auto obs_nan = xt::isnan(q_obs);
7172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140
auto prd_nan = xt::isnan(q_prd);
t_msk = xt::ones<bool>({n_msk, n_srs, n_tim});
xt::view(t_msk, xt::all()) =
xt::view(msk, xt::all(), xt::newaxis(), xt::all());
for (int m = 0; m < n_msk; m++) {
xt::view(t_msk, m) =
xt::where(obs_nan | prd_nan,
false, xt::view(t_msk, m));
}
};
// members for evaluation metrics
xt::xtensor<double, 3> RMSE;
xt::xtensor<double, 3> NSE;
xt::xtensor<double, 3> KGE;
xt::xtensor<double, 3> KGEPRIME;
// methods to compute elements
void calc_mean_obs();
void calc_mean_prd();
void calc_quad_err();
void calc_quad_obs();
void calc_quad_prd();
void calc_r_pearson();
void calc_alpha();
void calc_bias();
// methods to compute metrics
void calc_RMSE();
void calc_NSE();
void calc_KGE();
void calc_KGEPRIME();
};
// Compute the mean of the observations.
//
// \require q_obs:
// Streamflow observations.
// shape: (series, time)
// \assign mean_obs:
// Mean observed streamflow.
// shape: (subsets, samples, series, 1)
void Evaluator::calc_mean_obs()
{
mean_obs = xt::zeros<double>(mean_dims);
for (int m = 0; m < n_msk; m++) {
// apply the mask
// (using NaN workaround until reducers work on masked_view)
auto obs_masked = xt::where(xt::view(t_msk, m), q_obs, NAN);
// compute variable one sample at a time
for (int e = 0; e < n_exp; e++) {
// apply the bootstrap sampling
auto obs = xt::view(obs_masked, xt::all(), b_exp[e]);
xt::view(mean_obs, m, e) =
xt::nanmean(obs, -1, xt::keep_dims);
}
}
}
// Compute the mean of the predictions.
//
// \require q_prd:
// Streamflow predictions.
// shape: (series, time)
// \assign mean_prd:
// Mean predicted streamflow.
// shape: (subsets, samples, series, 1)
void Evaluator::calc_mean_prd()