diff --git a/Framework/Utilities/SaveBoostHistogram.hpp b/Framework/Utilities/SaveBoostHistogram.hpp
index c6ffbe064df6b0d707a37177bb1d811d40534ffa..508c982aa4b0e06444eb0d95186d3e5138381a9d 100644
--- a/Framework/Utilities/SaveBoostHistogram.hpp
+++ b/Framework/Utilities/SaveBoostHistogram.hpp
@@ -15,10 +15,14 @@
#include <numeric>
#include <utility>
#include <vector>
+#include <string>
#pragma once
namespace corsika::utl {
+
+ enum class SaveMode { overwrite, append };
+
/**
* This functions saves a boost::histogram into a numpy file. Only rather basic axis
* types are supported: regular, variable, integer, category<int>. Only "ordinary" bin
@@ -29,9 +33,12 @@ namespace corsika::utl {
*/
template <class Axes, class Storage>
inline void save_hist(boost::histogram::histogram<Axes, Storage> const& h,
- std::string const& filename) {
+ std::string const& filename, SaveMode mode = SaveMode::append) {
unsigned const rank = h.rank();
+ // append vs. overwrite
+ const std::string mode_str = (mode == SaveMode::append ? "a" : "w");
+
std::vector<size_t> axes_dims;
axes_dims.reserve(rank);
@@ -61,7 +68,7 @@ namespace corsika::utl {
ax_edges.push_back(ax.bin(ax.size() - 1).upper());
cnpy::npz_save(filename, std::string{"binedges_"} + std::to_string(i),
- ax_edges.data(), {ax_edges.size()}, "a");
+ ax_edges.data(), {ax_edges.size()}, mode_str);
} else {
axis_types.push_back('d');
std::vector<int64_t> bins; // we assume that discrete axes have integer bins
@@ -70,13 +77,14 @@ namespace corsika::utl {
for (int j = 0; j < ax.size(); ++j) { bins.push_back(ax.bin(j).lower()); }
cnpy::npz_save(filename, std::string{"bins_"} + std::to_string(i), bins.data(),
- {bins.size()}, "a");
+ {bins.size()}, mode_str);
}
}
- cnpy::npz_save(filename, std::string{"axistypes"}, axis_types.data(), {rank}, "a");
- cnpy::npz_save(filename, std::string{"overflow"}, overflow.get(), {rank}, "a");
- cnpy::npz_save(filename, std::string{"underflow"}, underflow.get(), {rank}, "a");
+ cnpy::npz_save(filename, std::string{"axistypes"}, axis_types.data(), {rank},
+ mode_str);
+ cnpy::npz_save(filename, std::string{"overflow"}, overflow.get(), {rank}, mode_str);
+ cnpy::npz_save(filename, std::string{"underflow"}, underflow.get(), {rank}, mode_str);
auto const prod_axis_size = std::accumulate(axes_dims.cbegin(), axes_dims.cend(),
unsigned{1}, std::multiplies<>());
@@ -100,7 +108,7 @@ namespace corsika::utl {
temp[p] = *x;
}
- cnpy::npz_save(filename, "data", temp.get(), axes_dims, "a");
+ cnpy::npz_save(filename, "data", temp.get(), axes_dims, mode_str);
// In Python this array can directly be assigned to a histogram view if that
// histogram has its axes correspondingly: hist.view(flow=True)[:] = file['data']
}
diff --git a/Processes/InteractionCounter/CMakeLists.txt b/Processes/InteractionCounter/CMakeLists.txt
index e7f374dbf345b888659c7b035b03ec1f24cf7494..a8c0cf5dda53fd966c7bce56098cb1d479f3ea1f 100644
--- a/Processes/InteractionCounter/CMakeLists.txt
+++ b/Processes/InteractionCounter/CMakeLists.txt
@@ -50,6 +50,11 @@ CORSIKA_ADD_TEST(testInteractionCounter SOURCES
${MODEL_HEADERS}
)
+target_compile_definitions (
+ testInteractionCounter
+ PRIVATE
+ REFDATADIR="${CMAKE_CURRENT_SOURCE_DIR}"
+ )
target_link_libraries (
testInteractionCounter
ProcessInteractionCounter
diff --git a/Processes/InteractionCounter/InteractionHistogram.cc b/Processes/InteractionCounter/InteractionHistogram.cc
index 8ae626ffe57d1b0dd1cc016667b409a26b697f9a..45d04469b4c8fa7f9dd7ed1f22b41817ab4f3a2a 100644
--- a/Processes/InteractionCounter/InteractionHistogram.cc
+++ b/Processes/InteractionCounter/InteractionHistogram.cc
@@ -41,12 +41,14 @@ void InteractionHistogram::fill(particles::Code projectile_id,
}
}
-void InteractionHistogram::saveLab(std::string const& filename) const {
- corsika::utl::save_hist(inthist_lab_, filename);
+void InteractionHistogram::saveLab(std::string const& filename,
+ utl::SaveMode mode) const {
+ corsika::utl::save_hist(inthist_lab_, filename, mode);
}
-void InteractionHistogram::saveCMS(std::string const& filename) const {
- corsika::utl::save_hist(inthist_cms_, filename);
+void InteractionHistogram::saveCMS(std::string const& filename,
+ utl::SaveMode mode) const {
+ corsika::utl::save_hist(inthist_cms_, filename, mode);
}
InteractionHistogram& InteractionHistogram::operator+=(
diff --git a/Processes/InteractionCounter/InteractionHistogram.hpp b/Processes/InteractionCounter/InteractionHistogram.hpp
index 96d243f1223505498d211d675fc8cba6ec620fb6..d374eb8d1778379134222b8de0beeb0dcfc181fb 100644
--- a/Processes/InteractionCounter/InteractionHistogram.hpp
+++ b/Processes/InteractionCounter/InteractionHistogram.hpp
@@ -54,9 +54,11 @@ namespace corsika::process::interaction_counter {
hist_type const& labHist() const { return inthist_lab_; }
- void saveLab(std::string const& filename) const;
+ void saveLab(std::string const& filename,
+ utl::SaveMode mode = utl::SaveMode::append) const;
- void saveCMS(std::string const& filename) const;
+ void saveCMS(std::string const& filename,
+ utl::SaveMode mode = utl::SaveMode::append) const;
InteractionHistogram& operator+=(InteractionHistogram const& other);
diff --git a/Processes/InteractionCounter/testInteractionCounter.cc b/Processes/InteractionCounter/testInteractionCounter.cc
index 95b8925cc9ca4fce1160b188568b59cc8f559e09..ad88405b5ae0717f45f89e2c1d37f5b67a20da44 100644
--- a/Processes/InteractionCounter/testInteractionCounter.cc
+++ b/Processes/InteractionCounter/testInteractionCounter.cc
@@ -27,6 +27,8 @@ using namespace corsika::process::interaction_counter;
using namespace corsika::units;
using namespace corsika::units::si;
+const std::string refDataDir = std::string(REFDATADIR);
+
struct DummyProcess {
template <typename TParticle>
GrammageType GetInteractionLength([[maybe_unused]] TParticle const& particle) {
@@ -39,7 +41,7 @@ struct DummyProcess {
}
};
-TEST_CASE("InteractionCounter") {
+TEST_CASE("InteractionCounter", "[process]") {
logging::SetLevel(logging::level::debug);
@@ -69,15 +71,14 @@ TEST_CASE("InteractionCounter") {
auto const& h2 = countedProcess.GetHistogram().CMSHist();
REQUIRE(h2.at(h2.axis(0).index(1'000'070'140), h2.axis(1).index(1.6e12)) == 1);
- // REQUIRE(h2.at(1'000'070'140, 92) == 1); // bin 1.584 .. 1.995 TeV √s
REQUIRE(std::accumulate(h2.cbegin(), h2.cend(), 0) == 1);
- countedProcess.GetHistogram().saveLab("testInteractionCounter_file1.npz");
- countedProcess.GetHistogram().saveCMS("testInteractionCounter_file2.npz");
+ countedProcess.GetHistogram().saveLab("testInteractionCounter_file1.npz",
+ utl::SaveMode::overwrite);
+ countedProcess.GetHistogram().saveCMS("testInteractionCounter_file2.npz",
+ utl::SaveMode::overwrite);
}
- SECTION("check saving") {}
-
SECTION("DoInteraction Lambda") {
auto constexpr code = particles::Code::Lambda0;
auto [stackPtr, secViewPtr] =
@@ -97,3 +98,38 @@ TEST_CASE("InteractionCounter") {
REQUIRE(std::accumulate(h2.cbegin(), h2.cend(), 0) == 1);
}
}
+
+#include <algorithm>
+#include <iterator>
+#include <string>
+#include <fstream>
+
+TEST_CASE("InteractionCounterOutput", "[output validation]") {
+
+ auto file = GENERATE(as<std::string>{}, "testInteractionCounter_file1",
+ "testInteractionCounter_file2");
+
+ SECTION(std::string("check saved data, ") + file + ".npz") {
+
+ std::cout << file + ".npz vs " << refDataDir + "/" + file + "_REF.npz" << std::endl;
+
+ // compare to binary reference data
+ std::ifstream file1(file + ".npz");
+ std::ifstream file1ref(refDataDir + "/" + file + "_REF.npz");
+
+ std::istreambuf_iterator<char> begin1(file1);
+ std::istreambuf_iterator<char> begin1ref(file1ref);
+
+ std::istreambuf_iterator<char> end;
+
+ while (begin1 != end && begin1ref != end) {
+ CHECK(*begin1 == *begin1ref);
+ ++begin1;
+ ++begin1ref;
+ }
+ CHECK(begin1 == end);
+ CHECK(begin1ref == end);
+ file1.close();
+ file1ref.close();
+ }
+}
diff --git a/Processes/InteractionCounter/testInteractionCounter_file1_REF.npz b/Processes/InteractionCounter/testInteractionCounter_file1_REF.npz
new file mode 100644
index 0000000000000000000000000000000000000000..13e2d5d49e16b6a2e7eb6b95dc1a26bccab18ae8
Binary files /dev/null and b/Processes/InteractionCounter/testInteractionCounter_file1_REF.npz differ
diff --git a/Processes/InteractionCounter/testInteractionCounter_file2_REF.npz b/Processes/InteractionCounter/testInteractionCounter_file2_REF.npz
new file mode 100644
index 0000000000000000000000000000000000000000..8509b27e4d690a386c2710a24ceda6a1143bf757
Binary files /dev/null and b/Processes/InteractionCounter/testInteractionCounter_file2_REF.npz differ