From dd4079cb1b7deb2ca82f3448be39cd67292ac797 Mon Sep 17 00:00:00 2001
From: Maximilian Reininghaus <maximilian.reininghaus@kit.edu>
Date: Wed, 23 Sep 2020 19:10:20 +0200
Subject: [PATCH] simple HistoryObservationPlane
---
Stack/History/Event.hpp | 3 +++
Stack/History/HistoryObservationPlane.cc | 21 +++++++++++++++++----
Stack/History/HistoryObservationPlane.hpp | 14 ++++++++++----
3 files changed, 30 insertions(+), 8 deletions(-)
diff --git a/Stack/History/Event.hpp b/Stack/History/Event.hpp
index e47eaf2ec..c89f8064e 100644
--- a/Stack/History/Event.hpp
+++ b/Stack/History/Event.hpp
@@ -42,6 +42,9 @@ namespace corsika::history {
template <typename TStackIterator>
TStackIterator projectile(TStackIterator begin) {
+ // MR: This is dangerous. You can pass any iterator though it must
+ // be stack.begin() to yield the correct projectile
+
return begin + projectileIndex_;
}
diff --git a/Stack/History/HistoryObservationPlane.cc b/Stack/History/HistoryObservationPlane.cc
index ebc8ab11d..f65ce8993 100644
--- a/Stack/History/HistoryObservationPlane.cc
+++ b/Stack/History/HistoryObservationPlane.cc
@@ -14,10 +14,12 @@ using namespace corsika::units::si;
using namespace corsika::history;
using namespace corsika;
-HistoryObservationPlane::HistoryObservationPlane(geometry::Plane const& obsPlane,
+HistoryObservationPlane::HistoryObservationPlane(setup::Stack const& stack,
+ geometry::Plane const& obsPlane,
bool deleteOnHit)
- : plane_(obsPlane)
- , deleteOnHit_(deleteOnHit) {}
+ : stack_{stack}
+ , plane_{obsPlane}
+ , deleteOnHit_{deleteOnHit} {}
corsika::process::EProcessReturn HistoryObservationPlane::DoContinuous(
setup::Stack::ParticleType const& particle, setup::Trajectory const& trajectory) {
@@ -55,7 +57,18 @@ LengthType HistoryObservationPlane::MaxStepLength(setup::Stack::ParticleType con
return (trajectory.GetR0() - pointOfIntersection).norm() * 1.0001;
}
-void fillHistoryHistogram(setup::Stack::ParticleType const& muon) {
+void HistoryObservationPlane::fillHistoryHistogram(
+ setup::Stack::ParticleType const& muon) {
double const muonEnergy = muon.GetEnergy() / 1_eV;
+ auto parent = stack_.begin() + muon.GetEvent()->projectileIndex();
+
+ auto* event = muon.GetEvent().get();
+
+ int intCounter = 0;
+ while (event) {
+ event = event->parentEvent().get();
+ intCounter++;
+ }
+ histogram_(intCounter);
}
diff --git a/Stack/History/HistoryObservationPlane.hpp b/Stack/History/HistoryObservationPlane.hpp
index 8c2d5700d..a603963db 100644
--- a/Stack/History/HistoryObservationPlane.hpp
+++ b/Stack/History/HistoryObservationPlane.hpp
@@ -21,7 +21,7 @@
namespace corsika::history {
namespace detail {
auto hist_factory() {
- auto h = boost::histogram::make_histogram(
+ /*auto h = boost::histogram::make_histogram(
boost::histogram::axis::regular<double, boost::histogram::axis::transform::log>{
130, 1e8, 1e21, "muon energy/eV"},
boost::histogram::axis::integer<int, boost::histogram::use_default,
@@ -30,7 +30,12 @@ namespace corsika::history {
boost::histogram::axis::regular<double, boost::histogram::axis::transform::log>{
130, 1e8, 1e21, "hadronic energy/eV"},
boost::histogram::axis::category<int, boost::histogram::use_default,
- boost::histogram::axis::option::growth_t>{});
+ boost::histogram::axis::option::growth_t>{});*/
+
+ auto h = boost::histogram::make_histogram(
+ boost::histogram::axis::integer<int, boost::histogram::use_default,
+ boost::histogram::axis::option::growth_t>{
+ 0, 10, "hadronic generation"});
return h;
}
} // namespace detail
@@ -38,9 +43,9 @@ namespace corsika::history {
class HistoryObservationPlane
: public corsika::process::ContinuousProcess<HistoryObservationPlane> {
public:
- HistoryObservationPlane(geometry::Plane const&, bool = true);
+ HistoryObservationPlane(setup::Stack const&, geometry::Plane const&, bool = true);
- void save(std::string const&);
+ //~ void save(std::string const&);
corsika::units::si::LengthType MaxStepLength(
corsika::setup::Stack::ParticleType const&,
@@ -53,6 +58,7 @@ namespace corsika::history {
private:
void fillHistoryHistogram(setup::Stack::ParticleType const&);
+ setup::Stack const& stack_;
geometry::Plane const plane_;
bool const deleteOnHit_;
--
GitLab