An error occurred while loading the file. Please try again.
-
fbourgin authoredec42e5d8
// Copyright (c) 2023, INRAE.
// Distributed under the terms of the GPL-3 Licence.
// The full licence is in the file LICENCE, distributed with this software.
#include <fstream>
#include <vector>
#include <tuple>
#include <array>
#include <string>
#include <unordered_map>
#include <gtest/gtest.h>
#include <xtl/xoptional.hpp>
#include <xtensor/xtensor.hpp>
#include <xtensor/xarray.hpp>
#include <xtensor/xview.hpp>
#include <xtensor/xmath.hpp>
#include <xtensor/xsort.hpp>
#include <xtensor/xmanipulation.hpp>
#include <xtensor/xcsv.hpp>
#include "evalhyd/evalp.hpp"
#ifndef EVALHYD_DATA_DIR
#error "need to define data directory"
#endif
using namespace xt::placeholders; // required for `_` to work
std::vector<std::string> all_metrics_p = {
"BS", "BSS", "BS_CRD", "BS_LBD", "REL_DIAG", "CRPS_FROM_BS",
"CRPS_FROM_ECDF",
"QS", "CRPS_FROM_QS",
"CONT_TBL", "POD", "POFD", "FAR", "CSI", "ROCSS",
"RANK_HIST", "DS", "AS",
"CR", "AW", "AWN", "WS",
"ES"
};
std::tuple<xt::xtensor<double, 1>, xt::xtensor<double, 2>> load_data_p()
{
// read in data
std::ifstream ifs;
ifs.open(EVALHYD_DATA_DIR "/data/q_obs.csv");
xt::xtensor<double, 1> observed = xt::squeeze(xt::load_csv<double>(ifs));
ifs.close();
ifs.open(EVALHYD_DATA_DIR "/data/q_prd.csv");
xt::xtensor<double, 2> predicted = xt::load_csv<double>(ifs);
ifs.close();
return std::make_tuple(observed, predicted);
}
std::unordered_map<std::string, xt::xarray<double>> load_expected_p()
{
// read in expected results
std::ifstream ifs;
std::unordered_map<std::string, xt::xarray<double>> expected;
for (const auto& metric : all_metrics_p)
{
ifs.open(EVALHYD_DATA_DIR "/expected/evalp/" + metric + ".csv");
expected[metric] = xt::view(
xt::squeeze(xt::load_csv<double>(ifs)),
xt::newaxis(), xt::newaxis(), xt::newaxis(),
xt::newaxis(), xt::all()
);
7172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140
ifs.close();
}
return expected;
}
TEST(ProbabilistTests, TestBrier)
{
// read in data
xt::xtensor<double, 1> observed;
xt::xtensor<double, 2> predicted;
std::tie(observed, predicted) = load_data_p();
// read in expected results
auto expected = load_expected_p();
// compute scores
xt::xtensor<double, 2> thresholds = {{690, 534, 445, NAN}};
std::vector<std::string> metrics = {"BS", "BSS", "BS_CRD", "BS_LBD", "REL_DIAG", "CRPS_FROM_BS"};
std::vector<xt::xarray<double>> results =
evalhyd::evalp(
// shape: (sites [1], time [t])
xt::eval(xt::view(observed, xt::newaxis(), xt::all())),
// shape: (sites [1], lead times [1], members [m], time [t])
xt::eval(xt::view(predicted, xt::newaxis(), xt::newaxis(), xt::all(), xt::all())),
metrics,
thresholds,
"high"
);
// check results
for (std::size_t m = 0; m < metrics.size(); m++)
{
if ( metrics[m] == "REL_DIAG" )
{
// /!\ stacked-up thresholds in CSV file because 7D metric,
// so need to resize array
expected[metrics[m]].resize(
{std::size_t {1}, std::size_t {1}, std::size_t {1},
std::size_t {1}, thresholds.shape(1),
predicted.shape(0) + 1, std::size_t {3}}
);
}
EXPECT_TRUE(xt::all(xt::isclose(
results[m], expected[metrics[m]], 1e-05, 1e-08, true
))) << "Failure for (" << metrics[m] << ")";
}
}
TEST(ProbabilistTests, TestCDF)
{
// read in data
xt::xtensor<double, 1> observed;
xt::xtensor<double, 2> predicted;
std::tie(observed, predicted) = load_data_p();
// read in expected results
auto expected = load_expected_p();
// compute scores
std::vector<std::string> metrics = {"CRPS_FROM_ECDF"};
std::vector<xt::xarray<double>> results =
evalhyd::evalp(
// shape: (sites [1], time [t])
xt::eval(xt::view(observed, xt::newaxis(), xt::all())),
// shape: (sites [1], lead times [1], members [m], time [t])
xt::eval(xt::view(predicted, xt::newaxis(), xt::newaxis(), xt::all(), xt::all())),
metrics
141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210
);
// check results
for (std::size_t m = 0; m < metrics.size(); m++)
{
EXPECT_TRUE(xt::all(xt::isclose(
results[m], expected[metrics[m]], 1e-05, 1e-08, true
))) << "Failure for (" << metrics[m] << ")";
}
}
TEST(ProbabilistTests, TestQuantiles)
{
// read in data
xt::xtensor<double, 1> observed;
xt::xtensor<double, 2> predicted;
std::tie(observed, predicted) = load_data_p();
// read in expected results
auto expected = load_expected_p();
// compute scores
std::vector<std::string> metrics = {"QS", "CRPS_FROM_QS"};
std::vector<xt::xarray<double>> results =
evalhyd::evalp(
// shape: (sites [1], time [t])
xt::eval(xt::view(observed, xt::newaxis(), xt::all())),
// shape: (sites [1], lead times [1], members [m], time [t])
xt::eval(xt::view(predicted, xt::newaxis(), xt::newaxis(), xt::all(), xt::all())),
metrics
);
// check results
for (std::size_t m = 0; m < metrics.size(); m++)
{
EXPECT_TRUE(xt::all(xt::isclose(
results[m], expected[metrics[m]], 1e-05, 1e-08, true
))) << "Failure for (" << metrics[m] << ")";
}
}
TEST(ProbabilistTests, TestContingency)
{
// read in data
xt::xtensor<double, 1> observed;
xt::xtensor<double, 2> predicted;
std::tie(observed, predicted) = load_data_p();
// read in expected results
auto expected = load_expected_p();
// compute scores
xt::xtensor<double, 2> thresholds = {{690, 534, 445, NAN}};
std::vector<std::string> metrics = {"CONT_TBL", "POD", "POFD", "FAR", "CSI", "ROCSS"};
std::vector<xt::xarray<double>> results =
evalhyd::evalp(
// shape: (sites [1], time [t])
xt::eval(xt::view(observed, xt::newaxis(), xt::all())),
// shape: (sites [1], lead times [1], members [m], time [t])
xt::eval(xt::view(predicted, xt::newaxis(), xt::newaxis(), xt::all(), xt::all())),
metrics,
thresholds,
"low"
);
// check results
for (std::size_t m = 0; m < metrics.size(); m++)
{
211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280
if (metrics[m] == "CONT_TBL")
{
// /!\ stacked-up thresholds and cells in CSV file because 7D metric,
// so need to resize array accordingly
expected[metrics[m]].resize(
{std::size_t {1}, std::size_t {1}, std::size_t {1}, std::size_t {1},
predicted.shape(0) + 1, thresholds.shape(1), std::size_t {4}}
);
}
EXPECT_TRUE(xt::all(xt::isclose(
results[m], expected[metrics[m]], 1e-05, 1e-08, true
))) << "Failure for (" << metrics[m] << ")";
}
}
TEST(ProbabilistTests, TestRanks)
{
// read in data
xt::xtensor<double, 1> observed;
xt::xtensor<double, 2> predicted;
std::tie(observed, predicted) = load_data_p();
// read in expected results
auto expected = load_expected_p();
std::vector<std::string> metrics = {"RANK_HIST", "DS", "AS"};
// compute scores
std::vector<xt::xarray<double>> results =
evalhyd::evalp(
// shape: (sites [1], time [t])
xt::eval(xt::view(observed, xt::newaxis(), xt::all())),
// shape: (sites [1], lead times [1], members [m], time [t])
xt::eval(xt::view(predicted, xt::newaxis(), xt::newaxis(), xt::all(), xt::all())),
{"RANK_HIST", "DS", "AS"},
xt::xtensor<double, 2>({}),
"high", // events
{}, // c_lvl
{}, // q_lvl
xt::xtensor<bool, 4>({}), // t_msk
xt::xtensor<std::array<char, 32>, 2>({}), // m_cdt
xtl::missing<const std::unordered_map<std::string, int>>(), // bootstrap
{}, // dts
7 // seed
);
// check results
for (std::size_t m = 0; m < metrics.size(); m++)
{
EXPECT_TRUE(xt::all(xt::isclose(
results[m], expected[metrics[m]], 1e-05, 1e-08, true
))) << "Failure for (" << metrics[m] << ")";
}
}
TEST(ProbabilistTests, TestIntervals)
{
// read in data
xt::xtensor<double, 1> observed;
xt::xtensor<double, 2> predicted;
std::tie(observed, predicted) = load_data_p();
// read in expected results
auto expected = load_expected_p();
// compute scores
std::vector<std::string> metrics = {"CR", "AW", "AWN", "WS"};
std::vector<xt::xarray<double>> results =
evalhyd::evalp(