From c6eaf499d7c05f07075c998ff32286604ccfbcab Mon Sep 17 00:00:00 2001
From: Alan Coleman <alanc@udel.edu>
Date: Wed, 6 Dec 2023 22:25:40 +0100
Subject: [PATCH] auto-calculate the writer binning
---
applications/c8_air_shower.cpp | 6 ++++--
examples/cascade_examples/em_shower.cpp | 6 ++++--
examples/cascade_examples/hybrid_MC.cpp | 6 ++++--
examples/cascade_examples/mars.cpp | 6 ++++--
examples/cascade_examples/radio_em_shower.cpp | 6 ++++--
examples/cascade_examples/water.cpp | 6 ++++--
examples/physics_examples/stopping_power.cpp | 1 -
7 files changed, 24 insertions(+), 13 deletions(-)
diff --git a/applications/c8_air_shower.cpp b/applications/c8_air_shower.cpp
index cf6b72f2b..0837d31c4 100644
--- a/applications/c8_air_shower.cpp
+++ b/applications/c8_air_shower.cpp
@@ -352,6 +352,8 @@ int main(int argc, char** argv) {
// we make the axis much longer than the inj-core distance since the
// profile will go beyond the core, depending on zenith angle
ShowerAxis const showerAxis{injectionPos, (showerCore - injectionPos) * 1.2, env};
+ auto const dX = 10_g / square(1_cm); //Binning of the writers along the shower axis
+ uint const nAxisBins = showerAxis.getMaximumX() / dX + 1; //Get maximum number of bins
/* === END: CONSTRUCT GEOMETRY === */
double const emthinfrac = app["--emthin"]->as<double>();
@@ -370,7 +372,7 @@ int main(int argc, char** argv) {
OutputManager output(app["--filename"]->as<std::string>(), seed, args.str(), outputDir);
// register energy losses as output
- EnergyLossWriter dEdX{showerAxis, 10_g / square(1_cm), 200};
+ EnergyLossWriter dEdX{showerAxis, dX, nAxisBins};
output.add("energyloss", dEdX);
DynamicInteractionProcess<StackType> heModel;
@@ -432,7 +434,7 @@ int main(int argc, char** argv) {
auto emContinuous =
make_select(EMHadronSwitch(), emContinuousBethe, emContinuousProposal);
- LongitudinalWriter profile{showerAxis, 200, 10_g / square(1_cm)};
+ LongitudinalWriter profile{showerAxis, nAxisBins, dX};
output.add("profile", profile);
LongitudinalProfile<SubWriter<decltype(profile)>> longprof{profile};
diff --git a/examples/cascade_examples/em_shower.cpp b/examples/cascade_examples/em_shower.cpp
index a685d70f1..614572bc9 100644
--- a/examples/cascade_examples/em_shower.cpp
+++ b/examples/cascade_examples/em_shower.cpp
@@ -133,6 +133,8 @@ int main(int argc, char** argv) {
ShowerAxis const showerAxis{injectionPos, (showerCore - injectionPos) * 1.02, env,
false, 1000};
+ auto const dX = 10_g / square(1_cm); //Binning of the writers along the shower axis
+ uint const nAxisBins = showerAxis.getMaximumX() / dX + 1; //Get maximum number of bins
CORSIKA_LOG_INFO("Primary particle: {}", beamCode);
CORSIKA_LOG_INFO("Zenith angle: {} (rad)", theta);
@@ -143,7 +145,7 @@ int main(int argc, char** argv) {
(showerCore - injectionPos).getNorm() * 1.02);
// setup processes, decays and interactions
- EnergyLossWriter energyloss{showerAxis, 10_g / square(1_cm), 200};
+ EnergyLossWriter energyloss{showerAxis, dX, nAxisBins};
ParticleCut<SubWriter<decltype(energyloss)>> cut(5_MeV, 5_MeV, 100_GeV, 100_GeV, true,
energyloss);
@@ -165,7 +167,7 @@ int main(int argc, char** argv) {
TrackWriter tracks;
output.add("tracks", tracks);
- LongitudinalWriter profile{showerAxis, 10_g / square(1_cm)};
+ LongitudinalWriter profile{showerAxis, nAxisBins, dX};
output.add("profile", profile);
LongitudinalProfile<SubWriter<decltype(profile)>> longprof{profile};
diff --git a/examples/cascade_examples/hybrid_MC.cpp b/examples/cascade_examples/hybrid_MC.cpp
index b3241ffa9..6d6841f56 100644
--- a/examples/cascade_examples/hybrid_MC.cpp
+++ b/examples/cascade_examples/hybrid_MC.cpp
@@ -206,6 +206,8 @@ int main(int argc, char** argv) {
ShowerAxis const showerAxis{injectionPos, (showerCore - injectionPos) * 1.02, env,
false, 1000};
+ auto const dX = 10_g / square(1_cm); //Binning of the writers along the shower axis
+ uint const nAxisBins = showerAxis.getMaximumX() / dX + 1; //Get maximum number of bins
CORSIKA_LOG_INFO("Primary particle: {}", beamCode);
CORSIKA_LOG_INFO("Zenith angle: {} (rad)", theta);
@@ -220,7 +222,7 @@ int main(int argc, char** argv) {
OutputManager output("hybrid_MC_outputs");
// register energy losses as output
- EnergyLossWriter dEdX{showerAxis, 10_g / square(1_cm), 200};
+ EnergyLossWriter dEdX{showerAxis, dX, nAxisBins};
output.add("energyloss", dEdX);
// create a track writer and register it with the output manager
@@ -230,7 +232,7 @@ int main(int argc, char** argv) {
ParticleCut<SubWriter<decltype(dEdX)>> cut(3_GeV, false, dEdX);
BetheBlochPDG<SubWriter<decltype(dEdX)>> eLoss(dEdX);
- LongitudinalWriter profile{showerAxis, 10_g / square(1_cm)};
+ LongitudinalWriter profile{showerAxis, nAxisBins, dX};
output.add("profile", profile);
LongitudinalProfile<SubWriter<decltype(profile)>> longprof{profile};
diff --git a/examples/cascade_examples/mars.cpp b/examples/cascade_examples/mars.cpp
index f382ef119..e4c6945a8 100644
--- a/examples/cascade_examples/mars.cpp
+++ b/examples/cascade_examples/mars.cpp
@@ -303,8 +303,10 @@ int main(int argc, char** argv) {
OutputManager output(app["--filename"]->as<std::string>());
ShowerAxis const showerAxis{injectionPos, (showerCore - injectionPos) * 1.2, env};
+ auto const dX = 10_g / square(1_cm); //Binning of the writers along the shower axis
+ uint const nAxisBins = showerAxis.getMaximumX() / dX + 1; //Get maximum number of bins
- EnergyLossWriter dEdX{showerAxis};
+ EnergyLossWriter dEdX{showerAxis, dX, nAxisBins};
output.add("energyloss", dEdX);
HEPEnergyType const emcut = 1_GeV;
@@ -336,7 +338,7 @@ int main(int argc, char** argv) {
auto emContinuous =
make_select(EMHadronSwitch(), emContinuousBethe, emContinuousProposal);
- LongitudinalWriter longprof{showerAxis};
+ LongitudinalWriter longprof{showerAxis, nAxisBins, dX};
output.add("profile", longprof);
LongitudinalProfile<SubWriter<decltype(longprof)>> profile{longprof};
diff --git a/examples/cascade_examples/radio_em_shower.cpp b/examples/cascade_examples/radio_em_shower.cpp
index 9631ce895..199946217 100644
--- a/examples/cascade_examples/radio_em_shower.cpp
+++ b/examples/cascade_examples/radio_em_shower.cpp
@@ -157,6 +157,8 @@ int main(int argc, char** argv) {
ShowerAxis const showerAxis{injectionPos, (showerCore - injectionPos) * 1.02, env,
false, 1000};
+ auto const dX = 10_g / square(1_cm); //Binning of the writers along the shower axis
+ uint const nAxisBins = showerAxis.getMaximumX() / dX + 1; //Get maximum number of bins
// setup the radio antennas
TimeType const groundHitTime{(showerCore - injectionPos).getNorm() / constants::c};
@@ -227,7 +229,7 @@ int main(int argc, char** argv) {
// }
// setup processes, decays and interactions
- EnergyLossWriter energyloss{showerAxis, 10_g / square(1_cm), 200};
+ EnergyLossWriter energyloss{showerAxis, dX, nAxisBins};
ParticleCut<SubWriter<decltype(energyloss)>> cut(5_MeV, 5_MeV, 100_GeV, 100_GeV, true,
energyloss);
@@ -249,7 +251,7 @@ int main(int argc, char** argv) {
TrackWriter tracks;
output.add("tracks", tracks);
- LongitudinalWriter profile{showerAxis, 10_g / square(1_cm)};
+ LongitudinalWriter profile{showerAxis, nAxisBins, dX};
output.add("profile", profile);
LongitudinalProfile<SubWriter<decltype(profile)>> longprof{profile};
diff --git a/examples/cascade_examples/water.cpp b/examples/cascade_examples/water.cpp
index c0f995f63..5e398d5b5 100644
--- a/examples/cascade_examples/water.cpp
+++ b/examples/cascade_examples/water.cpp
@@ -199,11 +199,13 @@ int main(int argc, char** argv) {
// * longitutional profile
ShowerAxis const showerAxis{injectionPos, 1.2 * injectorLength * downVec, env};
- LongitudinalWriter longiWriter{showerAxis, 5500, 1_g / square(1_cm)};
+ auto const dX = 1_g / square(1_cm); //Binning of the writers along the shower axis
+ uint const nAxisBins = showerAxis.getMaximumX() / dX + 1; //Get maximum number of bins
+ LongitudinalWriter longiWriter{showerAxis, nAxisBins, dX};
LongitudinalProfile<SubWriter<decltype(longiWriter)>> longprof{longiWriter};
// * energy loss profile
- EnergyLossWriter dEdX{showerAxis, 1_g / square(1_cm), 5500};
+ EnergyLossWriter dEdX{showerAxis, dX, nAxisBins};
// * physical process list
// particle production threshold
diff --git a/examples/physics_examples/stopping_power.cpp b/examples/physics_examples/stopping_power.cpp
index cbf1733ff..28a9ff7f0 100644
--- a/examples/physics_examples/stopping_power.cpp
+++ b/examples/physics_examples/stopping_power.cpp
@@ -9,7 +9,6 @@
#include <corsika/media/Environment.hpp>
#include <corsika/media/HomogeneousMedium.hpp>
#include <corsika/media/IMediumModel.hpp>
-#include <corsika/media/ShowerAxis.hpp>
#include <corsika/framework/geometry/Sphere.hpp>
#include <corsika/modules/BetheBlochPDG.hpp>
--
GitLab