diff --git a/include/evalhyd/detail/determinist/efficiencies.hpp b/include/evalhyd/detail/determinist/efficiencies.hpp
index df13462708dfcb79a3c8bbc4e924da9d787e1c76..034c84fdd30347cc390721aa57aaed48d430d029 100644
--- a/include/evalhyd/detail/determinist/efficiencies.hpp
+++ b/include/evalhyd/detail/determinist/efficiencies.hpp
@@ -17,22 +17,16 @@ namespace evalhyd
{
/// Compute the Pearson correlation coefficient.
///
- /// \param q_obs
- /// Streamflow observations.
- /// shape: (series, time)
- /// \param q_prd
- /// Streamflow predictions.
- /// shape: (series, time)
- /// \param mean_obs
- /// Mean observed streamflow.
- /// shape: (subsets, samples, series, 1)
- /// \param mean_prd
- /// Mean predicted streamflow.
- /// shape: (subsets, samples, series, 1)
- /// \param quad_obs
+ /// \param err_obs
+ /// Errors between observations and mean observation.
+ /// shape: (subsets, samples, series, time)
+ /// \param err_prd
+ /// Errors between predictions and mean prediction.
+ /// shape: (subsets, samples, series, time)
+ /// \param quad_err_obs
/// Quadratic errors between observations and mean observation.
/// shape: (subsets, samples, series, time)
- /// \param quad_prd
+ /// \param quad_err_prd
/// Quadratic errors between predictions and mean prediction.
/// shape: (subsets, samples, series, time)
/// \param t_msk
@@ -50,14 +44,11 @@ namespace evalhyd
/// \return
/// Pearson correlation coefficients.
/// shape: (subsets, samples, series)
- template <class XD2>
inline xt::xtensor<double, 3> calc_r_pearson(
- const XD2& q_obs,
- const XD2& q_prd,
- const xt::xtensor<double, 4>& mean_obs,
- const xt::xtensor<double, 4>& mean_prd,
- const xt::xtensor<double, 4>& quad_obs,
- const xt::xtensor<double, 4>& quad_prd,
+ const xt::xtensor<double, 4>& err_obs,
+ const xt::xtensor<double, 4>& err_prd,
+ const xt::xtensor<double, 4>& quad_err_obs,
+ const xt::xtensor<double, 4>& quad_err_prd,
const xt::xtensor<bool, 3>& t_msk,
const std::vector<xt::xkeep_slice<int>>& b_exp,
std::size_t n_srs,
@@ -72,26 +63,15 @@ namespace evalhyd
for (std::size_t m = 0; m < n_msk; m++)
{
- // apply the mask
- // (using NaN workaround until reducers work on masked_view)
- auto prd_masked = xt::where(xt::view(t_msk, xt::all(), m),
- q_prd, NAN);
- auto obs_masked = xt::where(xt::view(t_msk, xt::all(), m),
- q_obs, NAN);
-
// compute variable one sample at a time
for (std::size_t e = 0; e < n_exp; e++)
{
- auto prd = xt::view(prd_masked, xt::all(), b_exp[e]);
- auto obs = xt::view(obs_masked, xt::all(), b_exp[e]);
- auto r_num = xt::nansum(
- (prd - xt::view(mean_prd, m, e))
- * (obs - xt::view(mean_obs, m, e)),
- -1
- );
+ auto prd = xt::view(err_prd, m, e, xt::all(), b_exp[e]);
+ auto obs = xt::view(err_obs, m, e, xt::all(), b_exp[e]);
+ auto r_num = xt::nansum(prd * obs, -1);
- auto prd2 = xt::view(quad_prd, m, e, xt::all(), b_exp[e]);
- auto obs2 = xt::view(quad_obs, m, e, xt::all(), b_exp[e]);
+ auto prd2 = xt::view(quad_err_prd, m, e, xt::all(), b_exp[e]);
+ auto obs2 = xt::view(quad_err_obs, m, e, xt::all(), b_exp[e]);
auto r_den = xt::sqrt(
xt::nansum(prd2, -1) * xt::nansum(obs2, -1)
);
@@ -269,7 +249,7 @@ namespace evalhyd
/// \param quad_err
/// Quadratic errors between observations and predictions.
/// shape: (series, time)
- /// \param quad_obs
+ /// \param quad_err_obs
/// Quadratic errors between observations and mean observation.
/// shape: (subsets, samples, series, time)
/// \param t_msk
@@ -289,7 +269,7 @@ namespace evalhyd
/// shape: (series, subsets, samples)
inline xt::xtensor<double, 3> calc_NSE(
const xt::xtensor<double, 2>& quad_err,
- const xt::xtensor<double, 4>& quad_obs,
+ const xt::xtensor<double, 4>& quad_err_obs,
const xt::xtensor<bool, 3>& t_msk,
const std::vector<xt::xkeep_slice<int>>& b_exp,
std::size_t n_srs,
@@ -314,7 +294,7 @@ namespace evalhyd
auto err2 = xt::view(quad_err_masked, xt::all(), b_exp[e]);
xt::xtensor<double, 1> f_num =
xt::nansum(err2, -1);
- auto obs2 = xt::view(quad_obs, m, e, xt::all(), b_exp[e]);
+ auto obs2 = xt::view(quad_err_obs, m, e, xt::all(), b_exp[e]);
xt::xtensor<double, 1> f_den =
xt::nansum(obs2, -1);
diff --git a/include/evalhyd/detail/determinist/errors.hpp b/include/evalhyd/detail/determinist/errors.hpp
index f6e12354f9d45d05ebf34a2bef08793750fd50e1..3b2a66b50c7c469e6d491008628a0298420a0f2a 100644
--- a/include/evalhyd/detail/determinist/errors.hpp
+++ b/include/evalhyd/detail/determinist/errors.hpp
@@ -171,7 +171,7 @@ namespace evalhyd
return xt::square(err);
}
- /// Compute the quadratic error between observations and mean observation.
+ /// Compute the error between observations and mean observation.
///
/// \param q_obs
/// Streamflow observations.
@@ -191,10 +191,10 @@ namespace evalhyd
/// \param n_exp
/// Number of bootstrap samples.
/// \return
- /// Quadratic errors between observations and mean observation.
+ /// Errors between observations and mean observation.
/// shape: (subsets, samples, series, time)
template <class XD2>
- inline xt::xtensor<double, 4> calc_quad_obs(
+ inline xt::xtensor<double, 4> calc_err_obs(
const XD2& q_obs,
const xt::xtensor<double, 4>& mean_obs,
const xt::xtensor<bool, 3>& t_msk,
@@ -204,7 +204,7 @@ namespace evalhyd
std::size_t n_exp
)
{
- xt::xtensor<double, 4> quad_obs =
+ xt::xtensor<double, 4> err_obs =
xt::zeros<double>({n_msk, n_exp, n_srs, n_tim});
for (std::size_t m = 0; m < n_msk; m++)
@@ -216,16 +216,31 @@ namespace evalhyd
for (std::size_t e = 0; e < n_exp; e++)
{
- xt::view(quad_obs, m, e) = xt::square(
+ xt::view(err_obs, m, e) = (
obs_masked - xt::view(mean_obs, m, e)
);
}
}
- return quad_obs;
+ return err_obs;
}
- /// Compute the quadratic error between predictions and mean prediction.
+ /// Compute the quadratic error between observations and mean observation.
+ ///
+ /// \param err_obs
+ /// Errors between observations and mean observation.
+ /// shape: (subsets, samples, series, time)
+ /// \return
+ /// Quadratic errors between observations and mean observation.
+ /// shape: (subsets, samples, series, time)
+ inline xt::xtensor<double, 4> calc_quad_err_obs(
+ const xt::xtensor<double, 4>& err_obs
+ )
+ {
+ return xt::square(err_obs);
+ }
+
+ /// Compute the error between predictions and mean prediction.
///
/// \param q_prd
/// Streamflow predictions.
@@ -245,10 +260,10 @@ namespace evalhyd
/// \param n_exp
/// Number of bootstrap samples.
/// \return
- /// Quadratic errors between predictions and mean prediction.
+ /// Errors between predictions and mean prediction.
/// shape: (subsets, samples, series, time)
template <class XD2>
- inline xt::xtensor<double, 4> calc_quad_prd(
+ inline xt::xtensor<double, 4> calc_err_prd(
const XD2& q_prd,
const xt::xtensor<double, 4>& mean_prd,
const xt::xtensor<bool, 3>& t_msk,
@@ -258,7 +273,7 @@ namespace evalhyd
std::size_t n_exp
)
{
- xt::xtensor<double, 4> quad_prd =
+ xt::xtensor<double, 4> quad_err_prd =
xt::zeros<double>({n_msk, n_exp, n_srs, n_tim});
for (std::size_t m = 0; m < n_msk; m++)
@@ -270,13 +285,28 @@ namespace evalhyd
for (std::size_t e = 0; e < n_exp; e++)
{
- xt::view(quad_prd, m, e) = xt::square(
+ xt::view(quad_err_prd, m, e) = (
prd_masked - xt::view(mean_prd, m, e)
);
}
}
- return quad_prd;
+ return quad_err_prd;
+ }
+
+ /// Compute the quadratic error between predictions and mean prediction.
+ ///
+ /// \param err_obs
+ /// Errors between predictions and mean prediction.
+ /// shape: (subsets, samples, series, time)
+ /// \return
+ /// Quadratic errors between predictions and mean prediction.
+ /// shape: (subsets, samples, series, time)
+ inline xt::xtensor<double, 4> calc_quad_err_prd(
+ const xt::xtensor<double, 4>& err_prd
+ )
+ {
+ return xt::square(err_prd);
}
}
diff --git a/include/evalhyd/detail/determinist/evaluator.hpp b/include/evalhyd/detail/determinist/evaluator.hpp
index a5d5a3e0abc5d13334fc430cb2f373f4937c9228..8276179cc3508b48ff73ce15cb86e88dadb07b2a 100644
--- a/include/evalhyd/detail/determinist/evaluator.hpp
+++ b/include/evalhyd/detail/determinist/evaluator.hpp
@@ -42,8 +42,10 @@ namespace evalhyd
xtl::xoptional<xt::xtensor<double, 2>, bool> err;
xtl::xoptional<xt::xtensor<double, 2>, bool> abs_err;
xtl::xoptional<xt::xtensor<double, 2>, bool> quad_err;
- xtl::xoptional<xt::xtensor<double, 4>, bool> quad_obs;
- xtl::xoptional<xt::xtensor<double, 4>, bool> quad_prd;
+ xtl::xoptional<xt::xtensor<double, 4>, bool> err_obs;
+ xtl::xoptional<xt::xtensor<double, 4>, bool> quad_err_obs;
+ xtl::xoptional<xt::xtensor<double, 4>, bool> err_prd;
+ xtl::xoptional<xt::xtensor<double, 4>, bool> quad_err_prd;
xtl::xoptional<xt::xtensor<double, 3>, bool> r_pearson;
xtl::xoptional<xt::xtensor<double, 3>, bool> alpha;
xtl::xoptional<xt::xtensor<double, 3>, bool> gamma;
@@ -116,28 +118,50 @@ namespace evalhyd
return quad_err.value();
};
- xt::xtensor<double, 4> get_quad_obs()
+ xt::xtensor<double, 4> get_err_obs()
{
- if (!quad_obs.has_value())
+ if (!err_obs.has_value())
{
- quad_obs = elements::calc_quad_obs(
+ err_obs = elements::calc_err_obs(
q_obs, get_mean_obs(), t_msk,
n_srs, n_tim, n_msk, n_exp
);
}
- return quad_obs.value();
+ return err_obs.value();
};
- xt::xtensor<double, 4> get_quad_prd()
+ xt::xtensor<double, 4> get_quad_err_obs()
{
- if (!quad_prd.has_value())
+ if (!quad_err_obs.has_value())
{
- quad_prd = elements::calc_quad_prd(
+ quad_err_obs = elements::calc_quad_err_obs(
+ get_err_obs()
+ );
+ }
+ return quad_err_obs.value();
+ };
+
+ xt::xtensor<double, 4> get_err_prd()
+ {
+ if (!err_prd.has_value())
+ {
+ err_prd = elements::calc_err_prd(
q_prd, get_mean_prd(), t_msk,
n_srs, n_tim, n_msk, n_exp
);
}
- return quad_prd.value();
+ return err_prd.value();
+ };
+
+ xt::xtensor<double, 4> get_quad_err_prd()
+ {
+ if (!quad_err_prd.has_value())
+ {
+ quad_err_prd = elements::calc_quad_err_prd(
+ get_err_prd()
+ );
+ }
+ return quad_err_prd.value();
};
xt::xtensor<double, 3> get_r_pearson()
@@ -145,8 +169,9 @@ namespace evalhyd
if (!r_pearson.has_value())
{
r_pearson = elements::calc_r_pearson(
- q_obs, q_prd, get_mean_obs(), get_mean_prd(),
- get_quad_obs(), get_quad_prd(), t_msk, b_exp,
+ get_err_obs(), get_err_prd(),
+ get_quad_err_obs(), get_quad_err_prd(),
+ t_msk, b_exp,
n_srs, n_msk, n_exp
);
}
@@ -272,7 +297,7 @@ namespace evalhyd
if (!NSE.has_value())
{
NSE = metrics::calc_NSE(
- get_quad_err(), get_quad_obs(), t_msk, b_exp,
+ get_quad_err(), get_quad_err_obs(), t_msk, b_exp,
n_srs, n_msk, n_exp
);
}