IAP GITLAB

Skip to content
Snippets Groups Projects
Commit fc513cd2 authored by ralfulrich's avatar ralfulrich
Browse files

fixe conex test

parent 0f934064
No related branches found
No related tags found
No related merge requests found
......@@ -152,10 +152,16 @@ namespace corsika {
"CONEX and Corsika8 dX grammage binning are not the same!");
}
int const bin = int((bend + bstart) / 2);
size_t const bin = size_t((bend + bstart) / 2);
CORSIKA_LOGGER_TRACE(TOutput::getLogger(),
"add binned energy loss {} {} bin={} dE={} GeV ", bstart, bend,
bin, dE / 1_GeV);
if (bin >= profile_.size()) {
CORSIKA_LOGGER_WARN(TOutput::getLogger(),
"Grammage bin {} outside of profile {}. skipping.", bin,
profile_.size());
return;
}
profile_[bin][static_cast<int>(dEdX_output::ProfileIndex::Total)] += dE;
}
......
......@@ -113,7 +113,13 @@ namespace corsika {
"CONEX and Corsika8 dX grammage binning are not the same!");
}
int const bin = int((bend + bstart) / 2);
size_t const bin = size_t((bend + bstart) / 2);
if (bin >= profile_.size()) {
CORSIKA_LOGGER_WARN(TOutput::getLogger(),
"Grammage bin {} outside of profile {}. skipping.", bin,
profile_.size());
return;
}
if (pid == Code::Photon) {
profile_.at(bin)[static_cast<int>(number_profile::ProfileIndex::Photon)] += weight;
......
......@@ -17,10 +17,7 @@
#include <corsika/framework/geometry/Vector.hpp>
#include <corsika/media/ShowerAxis.hpp>
#include <corsika/modules/writers/WriterOff.hpp>
#include <corsika/modules/writers/SubWriter.hpp>
#include <corsika/modules/writers/EnergyLossWriterParquet.hpp>
#include <corsika/modules/writers/LongitudinalProfileWriterParquet.hpp>
#include <corsika/modules/writers/EnergyLossWriter.hpp>
#include <corsika/modules/writers/LongitudinalWriter.hpp>
......@@ -52,7 +49,7 @@ namespace corsika {
* @tparam TOutputE -- Output writer for dEdX data.
* @tparam TOutputN -- Output writer for particle number profile data.
*/
template <typename TOutputE = WriterOff, typename TOutputN = WriterOff2>
template <typename TOutputE, typename TOutputN>
class CONEXhybrid : public CascadeEquationsProcess<CONEXhybrid<TOutputE, TOutputN>>,
public SecondariesProcess<CONEXhybrid<TOutputE, TOutputN>>,
public SubWriter<TOutputE>,
......
......@@ -128,6 +128,11 @@ namespace corsika {
*/
YAML::Node getConfig() const;
number_profile::ProfileData const& getProfile(
number_profile::ProfileIndex index) const {
return profile_.at(static_cast<int>(index));
}
private:
ShowerAxis const& showerAxis_; ///< conversion between geometry and grammage
GrammageType dX_; ///< binning of profile.
......
......@@ -21,38 +21,29 @@ namespace corsika {
class WriterOff : public BaseOutput {
public:
WriterOff() {}
virtual ~WriterOff() {}
void startOfLibrary(boost::filesystem::path const&) final override {}
void endOfShower(unsigned int const) final override {}
void endOfLibrary() final override {}
/**
* The purpose of this catch-all constructor is to discard all parameters.
* @tparam TArgs
*/
template <typename... TArgs>
void write(TArgs&&...) {}
WriterOff(TArgs&&...) {}
virtual YAML::Node getConfig() const override { return YAML::Node(); }
}; // class WriterOff
class WriterOff2 : public BaseOutput {
public:
WriterOff2() {}
virtual ~WriterOff2() {}
virtual ~WriterOff() {}
void startOfLibrary(boost::filesystem::path const&) final override {}
void startOfLibrary(boost::filesystem::path const&) override {}
void endOfShower(unsigned int const) final override {}
void endOfShower(unsigned int const) override {}
void endOfLibrary() final override {}
void endOfLibrary() override {}
/**
* The purpose of this catch-all method is to discard all data.
* @tparam TArgs
*/
template <typename... TArgs>
void write(TArgs&&...) {}
virtual YAML::Node getConfig() const override { return YAML::Node(); }
virtual YAML::Node getConfig() const { return YAML::Node(); }
}; // class WriterOff
......
......@@ -102,10 +102,16 @@ TEST_CASE("CONEX") {
corsika::sibyll::Interaction sibyll;
[[maybe_unused]] corsika::sibyll::NuclearInteractionModel sibyllNuc(sibyll, env);
WriterOff w1;
WriterOff2 w2;
CONEXhybrid<WriterOff, WriterOff2> conex(center, showerAxis, t, injectionHeight, E0,
get_PDG(Code::Proton), w1, w2);
EnergyLossWriter<WriterOff> w1(showerAxis);
LongitudinalWriter<WriterOff> w2(showerAxis);
CONEXhybrid<decltype(w1), decltype(w2)> conex(center, showerAxis, t, injectionHeight,
E0, get_PDG(Code::Proton), w1, w2);
// initialize writers
w1.startOfLibrary("test");
w1.startOfShower(0);
w2.startOfLibrary("test");
w2.startOfShower(0);
// init conex
conex.initCascadeEquations();
HEPEnergyType const Eem{1_PeV};
......@@ -123,45 +129,14 @@ TEST_CASE("CONEX") {
emPosition.getCoordinates(conex.getObserverCS()),
emPosition.getCoordinates(rootCS));
conex.addParticle(Code::Proton, Eem, 0_eV, emPosition, momentum.normalized(), 0_s);
conex.addParticle(Code::Proton, Eem, Proton::mass, emPosition, momentum.normalized(),
0_s);
// supperimpose a photon
auto const momentumPhoton = showerAxis.getDirection() * 1_TeV;
conex.addParticle(Code::Photon, 1_TeV, 0_eV, emPosition, momentumPhoton.normalized(),
0_s);
DummyStack stack;
conex.doCascadeEquations(stack);
}
#include <algorithm>
#include <iterator>
#include <string>
#include <fstream>
TEST_CASE("ConexOutput", "[output validation]") {
logging::set_level(logging::level::info);
auto file = GENERATE(as<std::string>{}, "conex_fit", "conex_output");
SECTION(std::string("check saved data, ") + file + ".txt") {
// compare to reference data
std::ifstream file1(file + ".txt");
std::ifstream file1ref(refDataDir + "/" + file + "_REF.txt");
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();
}
CHECK(w1.getEnergyLost() / 1_TeV == Approx(1.0).epsilon(0.1));
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment