Commit df97a059 authored by Thibault Hallouin's avatar Thibault Hallouin
Browse files

document template parameters

1 merge request!3release v0.1.0
Pipeline #42926 passed with stage
in 2 minutes and 6 seconds
Showing with 122 additions and 98 deletions
+122 -98
...@@ -33,9 +33,9 @@ namespace evalhyd ...@@ -33,9 +33,9 @@ namespace evalhyd
// \return // \return
// Mean observed streamflow. // Mean observed streamflow.
// shape: (subsets, samples, series, 1) // shape: (subsets, samples, series, 1)
template <class D2> template <class XD2>
inline xt::xtensor<double, 4> calc_mean_obs( inline xt::xtensor<double, 4> calc_mean_obs(
const D2& q_obs, const XD2& q_obs,
const xt::xtensor<bool, 3>& t_msk, const xt::xtensor<bool, 3>& t_msk,
const std::vector<xt::xkeep_slice<int>>& b_exp, const std::vector<xt::xkeep_slice<int>>& b_exp,
std::size_t n_srs, std::size_t n_srs,
...@@ -85,9 +85,9 @@ namespace evalhyd ...@@ -85,9 +85,9 @@ namespace evalhyd
// \return // \return
// Mean predicted streamflow. // Mean predicted streamflow.
// shape: (subsets, samples, series, 1) // shape: (subsets, samples, series, 1)
template <class D2> template <class XD2>
inline xt::xtensor<double, 4> calc_mean_prd( inline xt::xtensor<double, 4> calc_mean_prd(
const D2& q_prd, const XD2& q_prd,
const xt::xtensor<bool, 3>& t_msk, const xt::xtensor<bool, 3>& t_msk,
const std::vector<xt::xkeep_slice<int>>& b_exp, const std::vector<xt::xkeep_slice<int>>& b_exp,
std::size_t n_srs, std::size_t n_srs,
...@@ -128,10 +128,10 @@ namespace evalhyd ...@@ -128,10 +128,10 @@ namespace evalhyd
// \return // \return
// Quadratic errors between observations and predictions. // Quadratic errors between observations and predictions.
// shape: (series, time) // shape: (series, time)
template <class D2> template <class XD2>
inline xt::xtensor<double, 2> calc_quad_err( inline xt::xtensor<double, 2> calc_quad_err(
const D2& q_obs, const XD2& q_obs,
const D2& q_prd const XD2& q_prd
) )
{ {
return xt::square(q_obs - q_prd); return xt::square(q_obs - q_prd);
...@@ -159,9 +159,9 @@ namespace evalhyd ...@@ -159,9 +159,9 @@ namespace evalhyd
// \return // \return
// Quadratic errors between observations and mean observation. // Quadratic errors between observations and mean observation.
// shape: (subsets, samples, series, time) // shape: (subsets, samples, series, time)
template <class D2> template <class XD2>
inline xt::xtensor<double, 4> calc_quad_obs( inline xt::xtensor<double, 4> calc_quad_obs(
const D2& q_obs, const XD2& q_obs,
const xt::xtensor<double, 4>& mean_obs, const xt::xtensor<double, 4>& mean_obs,
const xt::xtensor<bool, 3>& t_msk, const xt::xtensor<bool, 3>& t_msk,
std::size_t n_srs, std::size_t n_srs,
...@@ -212,9 +212,9 @@ namespace evalhyd ...@@ -212,9 +212,9 @@ namespace evalhyd
// \return // \return
// Quadratic errors between predictions and mean prediction. // Quadratic errors between predictions and mean prediction.
// shape: (subsets, samples, series, time) // shape: (subsets, samples, series, time)
template <class D2> template <class XD2>
inline xt::xtensor<double, 4> calc_quad_prd( inline xt::xtensor<double, 4> calc_quad_prd(
const D2& q_prd, const XD2& q_prd,
const xt::xtensor<double, 4>& mean_prd, const xt::xtensor<double, 4>& mean_prd,
const xt::xtensor<bool, 3>& t_msk, const xt::xtensor<bool, 3>& t_msk,
std::size_t n_srs, std::size_t n_srs,
...@@ -278,10 +278,10 @@ namespace evalhyd ...@@ -278,10 +278,10 @@ namespace evalhyd
// \return // \return
// Pearson correlation coefficients. // Pearson correlation coefficients.
// shape: (subsets, samples, series) // shape: (subsets, samples, series)
template <class D2> template <class XD2>
inline xt::xtensor<double, 3> calc_r_pearson( inline xt::xtensor<double, 3> calc_r_pearson(
const D2& q_obs, const XD2& q_obs,
const D2& q_prd, const XD2& q_prd,
const xt::xtensor<double, 4>& mean_obs, const xt::xtensor<double, 4>& mean_obs,
const xt::xtensor<double, 4>& mean_prd, const xt::xtensor<double, 4>& mean_prd,
const xt::xtensor<double, 4>& quad_obs, const xt::xtensor<double, 4>& quad_obs,
...@@ -358,10 +358,10 @@ namespace evalhyd ...@@ -358,10 +358,10 @@ namespace evalhyd
// \return // \return
// Alphas, ratios of standard deviations. // Alphas, ratios of standard deviations.
// shape: (subsets, samples, series) // shape: (subsets, samples, series)
template <class D2> template <class XD2>
inline xt::xtensor<double, 3> calc_alpha( inline xt::xtensor<double, 3> calc_alpha(
const D2& q_obs, const XD2& q_obs,
const D2& q_prd, const XD2& q_prd,
const xt::xtensor<double, 4>& mean_obs, const xt::xtensor<double, 4>& mean_obs,
const xt::xtensor<double, 4>& mean_prd, const xt::xtensor<double, 4>& mean_prd,
const xt::xtensor<bool, 3>& t_msk, const xt::xtensor<bool, 3>& t_msk,
...@@ -419,10 +419,10 @@ namespace evalhyd ...@@ -419,10 +419,10 @@ namespace evalhyd
// \return // \return
// Biases. // Biases.
// shape: (subsets, samples, series) // shape: (subsets, samples, series)
template <class D2> template <class XD2>
inline xt::xtensor<double, 3> calc_bias( inline xt::xtensor<double, 3> calc_bias(
const D2& q_obs, const XD2& q_obs,
const D2& q_prd, const XD2& q_prd,
const xt::xtensor<bool, 3>& t_msk, const xt::xtensor<bool, 3>& t_msk,
const std::vector<xt::xkeep_slice<int>>& b_exp, const std::vector<xt::xkeep_slice<int>>& b_exp,
std::size_t n_srs, std::size_t n_srs,
......
...@@ -15,13 +15,13 @@ namespace evalhyd ...@@ -15,13 +15,13 @@ namespace evalhyd
{ {
namespace determinist namespace determinist
{ {
template <class D2, class B2> template <class XD2, class XB2>
class Evaluator class Evaluator
{ {
private: private:
// members for input data // members for input data
const D2& q_obs; const XD2& q_obs;
const D2& q_prd; const XD2& q_prd;
xt::xtensor<bool, 3> t_msk; xt::xtensor<bool, 3> t_msk;
const std::vector<xt::xkeep_slice<int>>& b_exp; const std::vector<xt::xkeep_slice<int>>& b_exp;
...@@ -52,7 +52,7 @@ namespace evalhyd ...@@ -52,7 +52,7 @@ namespace evalhyd
{ {
if (!mean_obs.has_value()) if (!mean_obs.has_value())
{ {
mean_obs = elements::calc_mean_obs<D2>( mean_obs = elements::calc_mean_obs<XD2>(
q_obs, t_msk, b_exp, n_srs, n_msk, n_exp q_obs, t_msk, b_exp, n_srs, n_msk, n_exp
); );
} }
...@@ -63,7 +63,7 @@ namespace evalhyd ...@@ -63,7 +63,7 @@ namespace evalhyd
{ {
if (!mean_prd.has_value()) if (!mean_prd.has_value())
{ {
mean_prd = elements::calc_mean_prd<D2>( mean_prd = elements::calc_mean_prd<XD2>(
q_prd, t_msk, b_exp, n_srs, n_msk, n_exp q_prd, t_msk, b_exp, n_srs, n_msk, n_exp
); );
} }
...@@ -74,7 +74,7 @@ namespace evalhyd ...@@ -74,7 +74,7 @@ namespace evalhyd
{ {
if (!quad_err.has_value()) if (!quad_err.has_value())
{ {
quad_err = elements::calc_quad_err<D2>( quad_err = elements::calc_quad_err<XD2>(
q_obs, q_prd q_obs, q_prd
); );
} }
...@@ -85,7 +85,7 @@ namespace evalhyd ...@@ -85,7 +85,7 @@ namespace evalhyd
{ {
if (!quad_obs.has_value()) if (!quad_obs.has_value())
{ {
quad_obs = elements::calc_quad_obs<D2>( quad_obs = elements::calc_quad_obs<XD2>(
q_obs, get_mean_obs(), t_msk, q_obs, get_mean_obs(), t_msk,
n_srs, n_tim, n_msk, n_exp n_srs, n_tim, n_msk, n_exp
); );
...@@ -97,7 +97,7 @@ namespace evalhyd ...@@ -97,7 +97,7 @@ namespace evalhyd
{ {
if (!quad_prd.has_value()) if (!quad_prd.has_value())
{ {
quad_prd = elements::calc_quad_prd<D2>( quad_prd = elements::calc_quad_prd<XD2>(
q_prd, get_mean_prd(), t_msk, q_prd, get_mean_prd(), t_msk,
n_srs, n_tim, n_msk, n_exp n_srs, n_tim, n_msk, n_exp
); );
...@@ -109,7 +109,7 @@ namespace evalhyd ...@@ -109,7 +109,7 @@ namespace evalhyd
{ {
if (!r_pearson.has_value()) if (!r_pearson.has_value())
{ {
r_pearson = elements::calc_r_pearson<D2>( r_pearson = elements::calc_r_pearson<XD2>(
q_obs, q_prd, get_mean_obs(), get_mean_prd(), q_obs, q_prd, get_mean_obs(), get_mean_prd(),
get_quad_obs(), get_quad_prd(), t_msk, b_exp, get_quad_obs(), get_quad_prd(), t_msk, b_exp,
n_srs, n_msk, n_exp n_srs, n_msk, n_exp
...@@ -122,7 +122,7 @@ namespace evalhyd ...@@ -122,7 +122,7 @@ namespace evalhyd
{ {
if (!alpha.has_value()) if (!alpha.has_value())
{ {
alpha = elements::calc_alpha<D2>( alpha = elements::calc_alpha<XD2>(
q_obs, q_prd, get_mean_obs(), get_mean_prd(), q_obs, q_prd, get_mean_obs(), get_mean_prd(),
t_msk, b_exp, n_srs, n_msk, n_exp t_msk, b_exp, n_srs, n_msk, n_exp
); );
...@@ -134,7 +134,7 @@ namespace evalhyd ...@@ -134,7 +134,7 @@ namespace evalhyd
{ {
if (!bias.has_value()) if (!bias.has_value())
{ {
bias = elements::calc_bias<D2>( bias = elements::calc_bias<XD2>(
q_obs, q_prd, t_msk, b_exp, n_srs, n_msk, n_exp q_obs, q_prd, t_msk, b_exp, n_srs, n_msk, n_exp
); );
} }
...@@ -143,9 +143,9 @@ namespace evalhyd ...@@ -143,9 +143,9 @@ namespace evalhyd
public: public:
// constructor method // constructor method
Evaluator(const D2& obs, Evaluator(const XD2& obs,
const D2& prd, const XD2& prd,
const B2& msk, const XB2& msk,
const std::vector<xt::xkeep_slice<int>>& exp) : const std::vector<xt::xkeep_slice<int>>& exp) :
q_obs{obs}, q_prd{prd}, b_exp{exp} q_obs{obs}, q_prd{prd}, b_exp{exp}
{ {
......
...@@ -23,10 +23,10 @@ namespace evalhyd ...@@ -23,10 +23,10 @@ namespace evalhyd
// \return // \return
// Event observed outcome. // Event observed outcome.
// shape: (thresholds, time) // shape: (thresholds, time)
template<class V1D2> template<class XV1D2>
inline xt::xtensor<double, 2> calc_o_k( inline xt::xtensor<double, 2> calc_o_k(
const V1D2& q_obs, const XV1D2& q_obs,
const V1D2& q_thr const XV1D2& q_thr
) )
{ {
// determine observed realisation of threshold(s) exceedance // determine observed realisation of threshold(s) exceedance
...@@ -100,10 +100,10 @@ namespace evalhyd ...@@ -100,10 +100,10 @@ namespace evalhyd
// \return // \return
// Event probability forecast. // Event probability forecast.
// shape: (thresholds, time) // shape: (thresholds, time)
template<class V2D4, class V1D2> template<class XV2D4, class XV1D2>
inline xt::xtensor<double, 2> calc_y_k( inline xt::xtensor<double, 2> calc_y_k(
const V2D4& q_prd, const XV2D4& q_prd,
const V1D2& q_thr, const XV1D2& q_thr,
std::size_t n_mbr std::size_t n_mbr
) )
{ {
...@@ -171,10 +171,10 @@ namespace evalhyd ...@@ -171,10 +171,10 @@ namespace evalhyd
// \return // \return
// Brier score for each subset and for each threshold. // Brier score for each subset and for each threshold.
// shape: (subsets, samples, thresholds) // shape: (subsets, samples, thresholds)
template <class V1D2> template <class XV1D2>
inline xt::xtensor<double, 3> calc_BS( inline xt::xtensor<double, 3> calc_BS(
const xt::xtensor<double, 2>& bs, const xt::xtensor<double, 2>& bs,
const V1D2& q_thr, const XV1D2& q_thr,
const xt::xtensor<bool, 2>& t_msk, const xt::xtensor<bool, 2>& t_msk,
const std::vector<xt::xkeep_slice<int>>& b_exp, const std::vector<xt::xkeep_slice<int>>& b_exp,
std::size_t n_thr, std::size_t n_thr,
...@@ -253,9 +253,9 @@ namespace evalhyd ...@@ -253,9 +253,9 @@ namespace evalhyd
// Brier score components (reliability, resolution, uncertainty) // Brier score components (reliability, resolution, uncertainty)
// for each subset and for each threshold. // for each subset and for each threshold.
// shape: (subsets, samples, thresholds, components) // shape: (subsets, samples, thresholds, components)
template <class V1D2> template <class XV1D2>
inline xt::xtensor<double, 4> calc_BS_CRD( inline xt::xtensor<double, 4> calc_BS_CRD(
const V1D2& q_thr, const XV1D2& q_thr,
const xt::xtensor<double, 2>& o_k, const xt::xtensor<double, 2>& o_k,
const xt::xtensor<double, 2>& y_k, const xt::xtensor<double, 2>& y_k,
const xt::xtensor<double, 3>& bar_o, const xt::xtensor<double, 3>& bar_o,
...@@ -402,9 +402,9 @@ namespace evalhyd ...@@ -402,9 +402,9 @@ namespace evalhyd
// Brier score components (type 2 bias, discrimination, sharpness) // Brier score components (type 2 bias, discrimination, sharpness)
// for each subset and for each threshold. // for each subset and for each threshold.
// shape: (subsets, samples, thresholds, components) // shape: (subsets, samples, thresholds, components)
template <class V1D2> template <class XV1D2>
inline xt::xtensor<double, 4> calc_BS_LBD( inline xt::xtensor<double, 4> calc_BS_LBD(
const V1D2& q_thr, const XV1D2& q_thr,
const xt::xtensor<double, 2>& o_k, const xt::xtensor<double, 2>& o_k,
const xt::xtensor<double, 2>& y_k, const xt::xtensor<double, 2>& y_k,
const xt::xtensor<bool, 2>& t_msk, const xt::xtensor<bool, 2>& t_msk,
...@@ -560,10 +560,10 @@ namespace evalhyd ...@@ -560,10 +560,10 @@ namespace evalhyd
// \return // \return
// Brier skill score for each subset and for each threshold. // Brier skill score for each subset and for each threshold.
// shape: (subsets, samples, thresholds) // shape: (subsets, samples, thresholds)
template <class V1D2> template <class XV1D2>
inline xt::xtensor<double, 3> calc_BSS( inline xt::xtensor<double, 3> calc_BSS(
const xt::xtensor<double, 2>& bs, const xt::xtensor<double, 2>& bs,
const V1D2& q_thr, const XV1D2& q_thr,
const xt::xtensor<double, 2>& o_k, const xt::xtensor<double, 2>& o_k,
const xt::xtensor<double, 3>& bar_o, const xt::xtensor<double, 3>& bar_o,
const xt::xtensor<bool, 2>& t_msk, const xt::xtensor<bool, 2>& t_msk,
......
...@@ -16,13 +16,13 @@ namespace evalhyd ...@@ -16,13 +16,13 @@ namespace evalhyd
{ {
namespace probabilist namespace probabilist
{ {
template <class D2, class D4, class B4> template <class XD2, class XD4, class XB4>
class Evaluator class Evaluator
{ {
private: private:
using view1d_xtensor2d_double_type = decltype( using view1d_xtensor2d_double_type = decltype(
xt::view( xt::view(
std::declval<const D2&>(), std::declval<const XD2&>(),
std::declval<std::size_t>(), std::declval<std::size_t>(),
xt::all() xt::all()
) )
...@@ -30,7 +30,7 @@ namespace evalhyd ...@@ -30,7 +30,7 @@ namespace evalhyd
using view2d_xtensor4d_double_type = decltype( using view2d_xtensor4d_double_type = decltype(
xt::view( xt::view(
std::declval<const D4&>(), std::declval<const XD4&>(),
std::declval<std::size_t>(), std::declval<std::size_t>(),
std::declval<std::size_t>(), std::declval<std::size_t>(),
xt::all(), xt::all(),
...@@ -40,7 +40,7 @@ namespace evalhyd ...@@ -40,7 +40,7 @@ namespace evalhyd
using view2d_xtensor4d_bool_type = decltype( using view2d_xtensor4d_bool_type = decltype(
xt::view( xt::view(
std::declval<const B4&>(), std::declval<const XB4&>(),
std::declval<std::size_t>(), std::declval<std::size_t>(),
std::declval<std::size_t>(), std::declval<std::size_t>(),
xt::all(), xt::all(),
......
...@@ -20,9 +20,9 @@ namespace evalhyd ...@@ -20,9 +20,9 @@ namespace evalhyd
// \return // \return
// Streamflow forecast quantiles. // Streamflow forecast quantiles.
// shape: (quantiles, time) // shape: (quantiles, time)
template <class V2D4> template <class XV2D4>
inline xt::xtensor<double, 2> calc_q_qnt( inline xt::xtensor<double, 2> calc_q_qnt(
const V2D4& q_prd const XV2D4& q_prd
) )
{ {
return xt::sort(q_prd, 0); return xt::sort(q_prd, 0);
...@@ -42,9 +42,9 @@ namespace evalhyd ...@@ -42,9 +42,9 @@ namespace evalhyd
// \return // \return
// Quantile scores for each time step. // Quantile scores for each time step.
// shape: (quantiles, time) // shape: (quantiles, time)
template <class V1D2> template <class XV1D2>
inline xt::xtensor<double, 2> calc_qs( inline xt::xtensor<double, 2> calc_qs(
const V1D2 &q_obs, const XV1D2 &q_obs,
const xt::xtensor<double, 2>& q_qnt, const xt::xtensor<double, 2>& q_qnt,
std::size_t n_mbr std::size_t n_mbr
) )
......
...@@ -21,16 +21,26 @@ namespace evalhyd ...@@ -21,16 +21,26 @@ namespace evalhyd
/// ///
/// \rst /// \rst
/// ///
/// :Template Parameters:
///
/// XD2: Any 2-dimensional container class storing numeric elements
/// (e.g. ``xt::xtensor<double, 2>``, ``xt::pytensor<double, 2>``,
/// ``xt::rtensor<double, 2>``, etc.).
///
/// XB2: Any 2-dimensional container class storing boolean elements
/// (e.g. ``xt::xtensor<bool, 4>``, ``xt::pytensor<bool, 4>``,
/// ``xt::rtensor<bool, 4>``, etc.).
///
/// :Parameters: /// :Parameters:
/// ///
/// q_obs: ``xt::xtensor<double, 2>`` /// q_obs: ``XD2``
/// Streamflow observations. Time steps with missing observations /// Streamflow observations. Time steps with missing observations
/// must be assigned `NAN` values. Those time steps will be ignored /// must be assigned `NAN` values. Those time steps will be ignored
/// both in the observations and the predictions before the /// both in the observations and the predictions before the
/// *metrics* are computed. /// *metrics* are computed.
/// shape: (1, time) /// shape: (1, time)
/// ///
/// q_prd: ``xt::xtensor<double, 2>`` /// q_prd: ``XD2``
/// Streamflow predictions. Time steps with missing predictions /// Streamflow predictions. Time steps with missing predictions
/// must be assigned `NAN` values. Those time steps will be ignored /// must be assigned `NAN` values. Those time steps will be ignored
/// both in the observations and the predictions before the /// both in the observations and the predictions before the
...@@ -62,7 +72,7 @@ namespace evalhyd ...@@ -62,7 +72,7 @@ namespace evalhyd
/// `Pushpalatha et al. (2012) /// `Pushpalatha et al. (2012)
/// <https://doi.org/10.1016/j.jhydrol.2011.11.055>`_. /// <https://doi.org/10.1016/j.jhydrol.2011.11.055>`_.
/// ///
/// t_msk: ``xt::xtensor<bool, 2>``, optional /// t_msk: ``XB2``, optional
/// Mask used to temporally subset of the whole streamflow time series /// Mask used to temporally subset of the whole streamflow time series
/// (where True/False is used for the time steps to include/discard in /// (where True/False is used for the time steps to include/discard in
/// the subset). If provided, masks must feature the same number of /// the subset). If provided, masks must feature the same number of
...@@ -144,15 +154,15 @@ namespace evalhyd ...@@ -144,15 +154,15 @@ namespace evalhyd
/// evalhyd::evald(obs, prd, {"NSE"}, "none", 1, -9, msk); /// evalhyd::evald(obs, prd, {"NSE"}, "none", 1, -9, msk);
/// ///
/// \endrst /// \endrst
template <class D2, class B2> template <class XD2, class XB2>
std::vector<xt::xarray<double>> evald( std::vector<xt::xarray<double>> evald(
const xt::xexpression<D2>& q_obs, const xt::xexpression<XD2>& q_obs,
const xt::xexpression<D2>& q_prd, const xt::xexpression<XD2>& q_prd,
const std::vector<std::string>& metrics, const std::vector<std::string>& metrics,
const std::string& transform = "none", const std::string& transform = "none",
double exponent = 1, double exponent = 1,
double epsilon = -9, double epsilon = -9,
const xt::xexpression<B2>& t_msk = B2({}), const xt::xexpression<XB2>& t_msk = XB2({}),
const xt::xtensor<std::array<char, 32>, 1>& m_cdt = {}, const xt::xtensor<std::array<char, 32>, 1>& m_cdt = {},
const std::unordered_map<std::string, int>& bootstrap = const std::unordered_map<std::string, int>& bootstrap =
{{"n_samples", -9}, {"len_sample", -9}, {"summary", 0}}, {{"n_samples", -9}, {"len_sample", -9}, {"summary", 0}},
...@@ -160,13 +170,13 @@ namespace evalhyd ...@@ -160,13 +170,13 @@ namespace evalhyd
) )
{ {
// check ranks of tensors // check ranks of tensors
if (xt::get_rank<D2>::value != 2) if (xt::get_rank<XD2>::value != 2)
{ {
throw std::runtime_error( throw std::runtime_error(
"observations and/or predictions are not two-dimensional" "observations and/or predictions are not two-dimensional"
); );
} }
if (xt::get_rank<B2>::value != 2) if (xt::get_rank<XB2>::value != 2)
{ {
throw std::runtime_error( throw std::runtime_error(
"temporal masks are not two-dimensional" "temporal masks are not two-dimensional"
...@@ -174,10 +184,10 @@ namespace evalhyd ...@@ -174,10 +184,10 @@ namespace evalhyd
} }
// retrieve real types of the expressions // retrieve real types of the expressions
const D2& q_obs_ = q_obs.derived_cast(); const XD2& q_obs_ = q_obs.derived_cast();
const D2& q_prd_ = q_prd.derived_cast(); const XD2& q_prd_ = q_prd.derived_cast();
const B2& t_msk_ = t_msk.derived_cast(); const XB2& t_msk_ = t_msk.derived_cast();
// check that the metrics to be computed are valid // check that the metrics to be computed are valid
utils::check_metrics( utils::check_metrics(
...@@ -243,7 +253,7 @@ namespace evalhyd ...@@ -243,7 +253,7 @@ namespace evalhyd
// else if m_cdt provided, use them to generate t_msk // else if m_cdt provided, use them to generate t_msk
else if (m_cdt.size() > 0) else if (m_cdt.size() > 0)
{ {
B2 c_msk = xt::zeros<bool>({n_msk, n_tim}); XB2 c_msk = xt::zeros<bool>({n_msk, n_tim});
for (std::size_t m = 0; m < n_msk; m++) for (std::size_t m = 0; m < n_msk; m++)
{ {
...@@ -258,14 +268,14 @@ namespace evalhyd ...@@ -258,14 +268,14 @@ namespace evalhyd
// if neither t_msk nor m_cdt provided, generate dummy mask // if neither t_msk nor m_cdt provided, generate dummy mask
else else
{ {
return B2({xt::ones<bool>({std::size_t{1}, n_tim})}); return XB2({xt::ones<bool>({std::size_t{1}, n_tim})});
} }
}; };
auto msk = gen_msk(); auto msk = gen_msk();
// apply streamflow transformation if requested // apply streamflow transformation if requested
auto q_transform = [&](const D2& q) auto q_transform = [&](const XD2& q)
{ {
if ( transform == "none" || (transform == "pow" && exponent == 1)) if ( transform == "none" || (transform == "pow" && exponent == 1))
{ {
...@@ -273,7 +283,7 @@ namespace evalhyd ...@@ -273,7 +283,7 @@ namespace evalhyd
} }
else if ( transform == "sqrt" ) else if ( transform == "sqrt" )
{ {
return D2(xt::sqrt(q)); return XD2(xt::sqrt(q));
} }
else if ( transform == "inv" ) else if ( transform == "inv" )
{ {
...@@ -283,7 +293,7 @@ namespace evalhyd ...@@ -283,7 +293,7 @@ namespace evalhyd
epsilon = xt::mean(q_obs_)() * 0.01; epsilon = xt::mean(q_obs_)() * 0.01;
} }
return D2(1. / (q + epsilon)); return XD2(1. / (q + epsilon));
} }
else if ( transform == "log" ) else if ( transform == "log" )
{ {
...@@ -293,7 +303,7 @@ namespace evalhyd ...@@ -293,7 +303,7 @@ namespace evalhyd
epsilon = xt::mean(q_obs_)() * 0.01; epsilon = xt::mean(q_obs_)() * 0.01;
} }
return D2(xt::log(q + epsilon)); return XD2(xt::log(q + epsilon));
} }
else if ( transform == "pow" ) else if ( transform == "pow" )
{ {
...@@ -305,11 +315,11 @@ namespace evalhyd ...@@ -305,11 +315,11 @@ namespace evalhyd
epsilon = xt::mean(q_obs_)() * 0.01; epsilon = xt::mean(q_obs_)() * 0.01;
} }
return D2(xt::pow(q + epsilon, exponent)); return XD2(xt::pow(q + epsilon, exponent));
} }
else else
{ {
return D2(xt::pow(q, exponent)); return XD2(xt::pow(q, exponent));
} }
} }
else else
...@@ -320,8 +330,8 @@ namespace evalhyd ...@@ -320,8 +330,8 @@ namespace evalhyd
} }
}; };
const D2 obs = q_transform(q_obs_); const XD2 obs = q_transform(q_obs_);
const D2 prd = q_transform(q_prd_); const XD2 prd = q_transform(q_prd_);
// generate bootstrap experiment if requested // generate bootstrap experiment if requested
std::vector<xt::xkeep_slice<int>> exp; std::vector<xt::xkeep_slice<int>> exp;
...@@ -349,7 +359,7 @@ namespace evalhyd ...@@ -349,7 +359,7 @@ namespace evalhyd
} }
// instantiate determinist evaluator // instantiate determinist evaluator
determinist::Evaluator<D2, B2> evaluator(obs, prd, msk, exp); determinist::Evaluator<XD2, XB2> evaluator(obs, prd, msk, exp);
// retrieve or compute requested metrics // retrieve or compute requested metrics
std::vector<xt::xarray<double>> r; std::vector<xt::xarray<double>> r;
......
...@@ -21,16 +21,30 @@ namespace evalhyd ...@@ -21,16 +21,30 @@ namespace evalhyd
/// ///
/// \rst /// \rst
/// ///
/// :Template Parameters:
///
/// XD2: Any 2-dimensional container class storing numeric elements
/// (e.g. ``xt::xtensor<double, 2>``, ``xt::pytensor<double, 2>``,
/// ``xt::rtensor<double, 2>``, etc.).
///
/// XD4: Any 4-dimensional container class storing numeric elements
/// (e.g. ``xt::xtensor<double, 4>``, ``xt::pytensor<double, 4>``,
/// ``xt::rtensor<double, 4>``, etc.).
///
/// XB4: Any 4-dimensional container class storing boolean elements
/// (e.g. ``xt::xtensor<bool, 4>``, ``xt::pytensor<bool, 4>``,
/// ``xt::rtensor<bool, 4>``, etc.).
///
/// :Parameters: /// :Parameters:
/// ///
/// q_obs: ``xt::xtensor<double, 2>`` /// q_obs: ``XD2``
/// Streamflow observations. Time steps with missing observations /// Streamflow observations. Time steps with missing observations
/// must be assigned `NAN` values. Those time steps will be ignored /// must be assigned `NAN` values. Those time steps will be ignored
/// both in the observations and the predictions before the /// both in the observations and the predictions before the
/// *metrics* are computed. /// *metrics* are computed.
/// shape: (sites, time) /// shape: (sites, time)
/// ///
/// q_prd: ``xt::xtensor<double, 4>`` /// q_prd: ``XD4``
/// Streamflow predictions. Time steps with missing predictions /// Streamflow predictions. Time steps with missing predictions
/// must be assigned `NAN` values. Those time steps will be ignored /// must be assigned `NAN` values. Those time steps will be ignored
/// both in the observations and the predictions before the /// both in the observations and the predictions before the
...@@ -40,11 +54,11 @@ namespace evalhyd ...@@ -40,11 +54,11 @@ namespace evalhyd
/// metrics: ``std::vector<std::string>`` /// metrics: ``std::vector<std::string>``
/// The sequence of evaluation metrics to be computed. /// The sequence of evaluation metrics to be computed.
/// ///
/// q_thr: ``xt::xtensor<double, 2>``, optional /// q_thr: ``XD2``, optional
/// Streamflow exceedance threshold(s). /// Streamflow exceedance threshold(s).
/// shape: (sites, thresholds) /// shape: (sites, thresholds)
/// ///
/// t_msk: ``xt::xtensor<bool, 4>``, optional /// t_msk: ``XB4``, optional
/// Mask(s) used to generate temporal subsets of the whole streamflow /// Mask(s) used to generate temporal subsets of the whole streamflow
/// time series (where True/False is used for the time steps to /// time series (where True/False is used for the time steps to
/// include/discard in a given subset). If not provided and neither /// include/discard in a given subset). If not provided and neither
...@@ -122,13 +136,13 @@ namespace evalhyd ...@@ -122,13 +136,13 @@ namespace evalhyd
/// evalhyd::evalp(obs, prd, {"CRPS"}); /// evalhyd::evalp(obs, prd, {"CRPS"});
/// ///
/// \endrst /// \endrst
template <class D2, class D4, class B4> template <class XD2, class XD4, class XB4>
std::vector<xt::xarray<double>> evalp( std::vector<xt::xarray<double>> evalp(
const xt::xexpression<D2>& q_obs, const xt::xexpression<XD2>& q_obs,
const xt::xexpression<D4>& q_prd, const xt::xexpression<XD4>& q_prd,
const std::vector<std::string>& metrics, const std::vector<std::string>& metrics,
const xt::xexpression<D2>& q_thr = D2({}), const xt::xexpression<XD2>& q_thr = XD2({}),
const xt::xexpression<B4>& t_msk = B4({}), const xt::xexpression<XB4>& t_msk = XB4({}),
const xt::xtensor<std::array<char, 32>, 2>& m_cdt = {}, const xt::xtensor<std::array<char, 32>, 2>& m_cdt = {},
const std::unordered_map<std::string, int>& bootstrap = const std::unordered_map<std::string, int>& bootstrap =
{{"n_samples", -9}, {"len_sample", -9}, {"summary", 0}}, {{"n_samples", -9}, {"len_sample", -9}, {"summary", 0}},
...@@ -136,19 +150,19 @@ namespace evalhyd ...@@ -136,19 +150,19 @@ namespace evalhyd
) )
{ {
// check ranks of tensors // check ranks of tensors
if (xt::get_rank<D2>::value != 2) if (xt::get_rank<XD2>::value != 2)
{ {
throw std::runtime_error( throw std::runtime_error(
"observations and/or thresholds are not two-dimensional" "observations and/or thresholds are not two-dimensional"
); );
} }
if (xt::get_rank<D4>::value != 4) if (xt::get_rank<XD4>::value != 4)
{ {
throw std::runtime_error( throw std::runtime_error(
"predictions are not four-dimensional" "predictions are not four-dimensional"
); );
} }
if (xt::get_rank<B4>::value != 4) if (xt::get_rank<XB4>::value != 4)
{ {
throw std::runtime_error( throw std::runtime_error(
"temporal masks are not four-dimensional" "temporal masks are not four-dimensional"
...@@ -156,11 +170,11 @@ namespace evalhyd ...@@ -156,11 +170,11 @@ namespace evalhyd
} }
// retrieve real types of the expressions // retrieve real types of the expressions
const D2& q_obs_ = q_obs.derived_cast(); const XD2& q_obs_ = q_obs.derived_cast();
const D4& q_prd_ = q_prd.derived_cast(); const XD4& q_prd_ = q_prd.derived_cast();
const D2& q_thr_ = q_thr.derived_cast(); const XD2& q_thr_ = q_thr.derived_cast();
const B4& t_msk_ = t_msk.derived_cast(); const XB4& t_msk_ = t_msk.derived_cast();
// check that the metrics to be computed are valid // check that the metrics to be computed are valid
utils::check_metrics( utils::check_metrics(
...@@ -268,7 +282,7 @@ namespace evalhyd ...@@ -268,7 +282,7 @@ namespace evalhyd
// generate masks from conditions if provided // generate masks from conditions if provided
auto gen_msk = [&]() auto gen_msk = [&]()
{ {
B4 c_msk = xt::zeros<bool>({n_sit, n_ltm, n_msk, n_tim}); XB4 c_msk = xt::zeros<bool>({n_sit, n_ltm, n_msk, n_tim});
if (m_cdt.size() > 0) if (m_cdt.size() > 0)
{ {
for (std::size_t s = 0; s < n_sit; s++) for (std::size_t s = 0; s < n_sit; s++)
...@@ -290,7 +304,7 @@ namespace evalhyd ...@@ -290,7 +304,7 @@ namespace evalhyd
return c_msk; return c_msk;
}; };
const B4 c_msk = gen_msk(); const XB4 c_msk = gen_msk();
// generate bootstrap experiment if requested // generate bootstrap experiment if requested
std::vector<xt::xkeep_slice<int>> b_exp; std::vector<xt::xkeep_slice<int>> b_exp;
...@@ -341,7 +355,7 @@ namespace evalhyd ...@@ -341,7 +355,7 @@ namespace evalhyd
xt::view(c_msk, s, l, xt::all(), xt::all()) : xt::view(c_msk, s, l, xt::all(), xt::all()) :
xt::view(t_msk_, s, l, xt::all(), xt::all())); xt::view(t_msk_, s, l, xt::all(), xt::all()));
probabilist::Evaluator<D2, D4, B4> evaluator( probabilist::Evaluator<XD2, XD4, XB4> evaluator(
q_obs_v, q_prd_v, q_thr_v, t_msk_v, b_exp q_obs_v, q_prd_v, q_thr_v, t_msk_v, b_exp
); );
......
Supports Markdown
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