From 4c838decd2ad43d6d8758436b0af6f7fae118679 Mon Sep 17 00:00:00 2001
From: Maximilian Reininghaus <maximilian.reininghaus@tu-dortmund.de>
Date: Sun, 5 May 2019 15:39:32 -0300
Subject: [PATCH] use WeightedSum in GetInteractionLength
---
Environment/NuclearComposition.h | 6 +++---
Processes/Pythia/Interaction.cc | 26 ++++++--------------------
2 files changed, 9 insertions(+), 23 deletions(-)
diff --git a/Environment/NuclearComposition.h b/Environment/NuclearComposition.h
index e7f0d6d8..a3ebfe8f 100644
--- a/Environment/NuclearComposition.h
+++ b/Environment/NuclearComposition.h
@@ -16,6 +16,7 @@
#include <corsika/units/PhysicalUnits.h>
#include <cassert>
+#include <functional>
#include <numeric>
#include <random>
#include <stdexcept>
@@ -82,7 +83,6 @@ namespace corsika::environment {
auto WeightedSum(TFunction func) const {
using ResultQuantity = decltype(func(*fComponents.cbegin()));
- auto const sum = [](auto x, auto y) { return x + y; };
auto const prod = [&](auto const compID, auto const fraction) {
return func(compID) * fraction;
};
@@ -91,12 +91,12 @@ namespace corsika::environment {
return std::inner_product(
fComponents.cbegin(), fComponents.cend(), fNumberFractions.cbegin(),
ResultQuantity::zero(), // .zero() is defined for quantity types only
- sum, prod);
+ std::plus<ResultQuantity>(), prod);
} else {
return std::inner_product(
fComponents.cbegin(), fComponents.cend(), fNumberFractions.cbegin(),
ResultQuantity(0), // in other cases we have to use a bare 0
- sum, prod);
+ std::plus<ResultQuantity>(), prod);
}
}
diff --git a/Processes/Pythia/Interaction.cc b/Processes/Pythia/Interaction.cc
index ae810435..f4272b16 100644
--- a/Processes/Pythia/Interaction.cc
+++ b/Processes/Pythia/Interaction.cc
@@ -215,26 +215,12 @@ namespace corsika::process::pythia {
const auto mediumComposition =
currentNode->GetModelProperties().GetNuclearComposition();
// determine average interaction length
- // weighted sum
- int i = -1;
- si::CrossSectionType weightedProdCrossSection = 0_mb;
- // get weights of components from environment/medium
- const auto& w = mediumComposition.GetFractions();
- // loop over components in medium
- for (auto const targetId : mediumComposition.GetComponents()) {
- i++;
- cout << "Interaction: get interaction length for target: " << targetId << endl;
-
- auto const [productionCrossSection, elaCrossSection] =
- GetCrossSection(corsikaBeamId, targetId, ECoM);
- [[maybe_unused]] const auto& dummy_elaCrossSection = elaCrossSection;
-
- cout << "Interaction: IntLength: pythia return (mb): "
- << productionCrossSection / 1_mb << endl
- << "Interaction: IntLength: weight : " << w[i] << endl;
-
- weightedProdCrossSection += w[i] * productionCrossSection;
- }
+
+ auto const weightedProdCrossSection =
+ mediumComposition.WeightedSum([=](auto vTargetID) {
+ return std::get<0>(this->GetCrossSection(corsikaBeamId, vTargetID, ECoM));
+ });
+
cout << "Interaction: IntLength: weighted CrossSection (mb): "
<< weightedProdCrossSection / 1_mb << endl
<< "Interaction: IntLength: average mass number: "
--
GitLab