From 90924065cc6365f65e3b33214884eb352a01e913 Mon Sep 17 00:00:00 2001
From: Thibault Hallouin <thibault.hallouin@inrae.fr>
Date: Mon, 1 Aug 2022 17:02:56 +0200
Subject: [PATCH] replace `xt::stddev` with custom made one to bypass bug with
rtensor
See https://gitlab.irstea.fr/HYCAR-Hydro/evalhyd/evalhyd-r/-/issues/1
for details on rtensor bug. Once this bug is fixed, `xt::stddev` should
be brought back.
---
src/determinist/evaluator.hpp | 23 +++++++++++++++++------
1 file changed, 17 insertions(+), 6 deletions(-)
diff --git a/src/determinist/evaluator.hpp b/src/determinist/evaluator.hpp
index da85e3f..ad130f9 100644
--- a/src/determinist/evaluator.hpp
+++ b/src/determinist/evaluator.hpp
@@ -6,6 +6,19 @@
namespace evalhyd
{
+
+ namespace detail {
+ // function to calculate standard deviation on last axis
+ //
+ // TODO: substitute with `xt::stddev` when fixed for `xt::rtensor`
+ // (see https://gitlab.irstea.fr/HYCAR-Hydro/evalhyd/evalhyd-r/-/issues/1)
+ template <class A1, class A2>
+ inline auto std_dev(A1&& arr, A2&& mean_arr)
+ {
+ return xt::sqrt(xt::mean(xt::square(arr - mean_arr), -1, xt::keep_dims));
+ }
+ }
+
namespace determinist
{
template <class A>
@@ -255,9 +268,8 @@ namespace evalhyd
void Evaluator<A>::calc_KGE()
{
// calculate error in spread of flow $alpha$
- auto alpha =
- xt::stddev(q_prd, {q_prd.dimension() - 1}, xt::keep_dims)
- / xt::stddev(q_obs, {q_obs.dimension() - 1}, xt::keep_dims);
+ auto alpha = detail::std_dev(q_prd, mean_prd)
+ / detail::std_dev(q_obs, mean_obs);
// compute KGE
KGE = 1 - xt::sqrt(
@@ -294,9 +306,8 @@ namespace evalhyd
void Evaluator<A>::calc_KGEPRIME()
{
// calculate error in spread of flow $gamma$
- auto gamma =
- (xt::stddev(q_prd, {q_prd.dimension() - 1}, xt::keep_dims) / mean_prd)
- / (xt::stddev(q_obs, {q_obs.dimension() - 1}, xt::keep_dims) / mean_obs);
+ auto gamma = (detail::std_dev(q_prd, mean_prd) / mean_prd)
+ / (detail::std_dev(q_obs, mean_obs) / mean_obs);
// compute KGE
KGEPRIME = 1 - xt::sqrt(
--
GitLab