From d3e58cb790cb793a4c32140fe5e9d21907f9edcb Mon Sep 17 00:00:00 2001 From: ralfulrich <ralf.ulrich@kit.edu> Date: Sat, 20 Feb 2021 01:22:07 +0100 Subject: [PATCH] fixed, improvements and switch off debugging output --- corsika/corsika.hpp | 2 + corsika/detail/framework/core/Cascade.inl | 121 ++++++++++-------- .../detail/modules/LongitudinalProfile.inl | 10 +- corsika/detail/modules/StackInspector.inl | 10 +- corsika/detail/modules/pythia8/Decay.inl | 33 ++--- .../process/ContinuousProcessIndex.hpp | 2 + do-clang-format.py | 17 +++ examples/vertical_EAS.cpp | 4 +- 8 files changed, 114 insertions(+), 85 deletions(-) diff --git a/corsika/corsika.hpp b/corsika/corsika.hpp index 36fa90b63..6d9a01aba 100644 --- a/corsika/corsika.hpp +++ b/corsika/corsika.hpp @@ -1,3 +1,5 @@ + + /* * (c) Copyright 2020 CORSIKA Project, corsika-project@lists.kit.edu * diff --git a/corsika/detail/framework/core/Cascade.inl b/corsika/detail/framework/core/Cascade.inl index 3c73d38c8..acfd29c69 100644 --- a/corsika/detail/framework/core/Cascade.inl +++ b/corsika/detail/framework/core/Cascade.inl @@ -121,11 +121,17 @@ namespace corsika { ContinuousProcessStepLength const continuousMaxStep = sequence_.getMaxStepLength(vParticle, step); LengthType const continuous_max_dist = continuousMaxStep; - ContinuousProcessIndex const limitingId = continuousMaxStep; + ContinuousProcessIndex limitingId; // take minimum of geometry, interaction, decay for next step - auto const min_distance = - std::min({distance_interact, distance_decay, continuous_max_dist, geomMaxLength}); + LengthType const min_discrete = std::min(distance_interact, distance_decay); + LengthType const min_non_continuous = std::min(min_discrete, geomMaxLength); + LengthType const min_distance = std::min(min_non_continuous, continuous_max_dist); + + if (continuous_max_dist < min_non_continuous) { + limitingId = + continuousMaxStep; // the current step IS limited by known continuous process + } CORSIKA_LOG_DEBUG( "transport particle by : {} m " @@ -157,51 +163,43 @@ namespace corsika { } return; } + if (continuous_max_dist < min_non_continuous) { + return; // there is nothing further further + } CORSIKA_LOG_DEBUG("sth. happening before geometric limit ? {}", ((min_distance < geomMaxLength) ? "yes" : "no")); - if (min_distance < geomMaxLength) { // interaction to happen within geometric limit + if (geomMaxLength < min_discrete) { + // geometric / tracking limit - // check whether decay or interaction limits this step the - // outcome of decay or interaction MAY be a) new particles in - // secondaries, b) the projectile particle deleted (or - // changed) + if (nextVol != currentLogicalNode) { + // boundary crossing, step is limited by volume boundary - TStackView secondaries(vParticle); + CORSIKA_LOG_DEBUG("volume boundary crossing to {}", fmt::ptr(nextVol)); - if (min_distance < continuous_max_dist) { + if (nextVol == environment_.getUniverse().get()) { + CORSIKA_LOG_DEBUG( + "particle left physics world, is now in unknown space -> delete"); + vParticle.erase(); + } + vParticle.setNode(nextVol); /* - Create SecondaryView object on Stack. The data container - remains untouched and identical, and 'projectil' is identical - to 'vParticle' above this line. However, - projectil.AddSecondaries populate the SecondaryView, which can - then be used afterwards for further processing. Thus: it is - important to use projectle/view (and not vParticle) for Interaction, - and Decay! - */ + doBoundary may delete the particle (or not) - [[maybe_unused]] auto projectile = secondaries.getProjectile(); - - if (distance_interact < distance_decay) { - interaction(secondaries); - } else { - decay(secondaries); - // make sure particle actually did decay if it should have done so - if (secondaries.getSize() == 1 && - projectile.getPID() == secondaries.getNextParticle().getPID()) - throw std::runtime_error(fmt::format("Particle {} decays into itself!", - get_name(projectile.getPID()))); - } + caveat: any changes to vParticle, or even the production + of new secondaries is currently not passed to ParticleCut, + thus, particles outside the desired phase space may be produced. - sequence_.doSecondaries(secondaries); - vParticle.erase(); - } else { // step-length limitation within volume + \todo: this must be fixed. + */ - CORSIKA_LOG_DEBUG("step-length limitation"); - // no further physics happens here. just proceed to next step. + sequence_.doBoundaryCrossing(vParticle, *currentLogicalNode, *nextVol); + return; } + CORSIKA_LOG_DEBUG("step limit reached. nothing further happens."); + [[maybe_unused]] auto const assertion = [&] { auto const* numericalNodeAfterStep = environment_.getUniverse()->getContainingNode(vParticle.getPosition()); @@ -214,31 +212,42 @@ namespace corsika { assert(assertion()); // numerical and logical nodes should match, since // we did not cross any volume boundary - } else { // boundary crossing, step is limited by volume boundary - - if (nextVol != currentLogicalNode) { - - CORSIKA_LOG_DEBUG("volume boundary crossing to {}", fmt::ptr(nextVol)); - - if (nextVol == environment_.getUniverse().get()) { - CORSIKA_LOG_DEBUG( - "particle left physics world, is now in unknown space -> delete"); - vParticle.erase(); - } - vParticle.setNode(nextVol); - /* - doBoundary may delete the particle (or not) + // step length limit + return; + } - caveat: any changes to vParticle, or even the production - of new secondaries is currently not passed to ParticleCut, - thus, particles outside the desired phase space may be produced. + // interaction or decay to happen in this step + // the outcome of decay or interaction MAY be a) new particles in + // secondaries, b) the projectile particle deleted (or + // changed) - \todo: this must be fixed. - */ + TStackView secondaries(vParticle); - sequence_.doBoundaryCrossing(vParticle, *currentLogicalNode, *nextVol); - } + /* + Create SecondaryView object on Stack. The data container + remains untouched and identical, and 'projectil' is identical + to 'vParticle' above this line. However, + projectil.AddSecondaries populate the SecondaryView, which can + then be used afterwards for further processing. Thus: it is + important to use projectle/view (and not vParticle) for Interaction, + and Decay! + */ + + [[maybe_unused]] auto projectile = secondaries.getProjectile(); + + if (distance_interact < distance_decay) { + interaction(secondaries); + } else { + decay(secondaries); + // make sure particle actually did decay if it should have done so + if (secondaries.getSize() == 1 && + projectile.getPID() == secondaries.getNextParticle().getPID()) + throw std::runtime_error(fmt::format("Particle {} decays into itself!", + get_name(projectile.getPID()))); } + + sequence_.doSecondaries(secondaries); + vParticle.erase(); } template <typename TTracking, typename TProcessList, typename TStack, diff --git a/corsika/detail/modules/LongitudinalProfile.inl b/corsika/detail/modules/LongitudinalProfile.inl index 0a2ec2e77..48deafe6e 100644 --- a/corsika/detail/modules/LongitudinalProfile.inl +++ b/corsika/detail/modules/LongitudinalProfile.inl @@ -35,11 +35,11 @@ namespace corsika { GrammageType const grammageStart = shower_axis_.getProjectedX(vTrack.getPosition(0)); GrammageType const grammageEnd = shower_axis_.getProjectedX(vTrack.getPosition(1)); - CORSIKA_LOG_INFO("pos1={} m, pos2={}, X1={} g/cm2, X2={} g/cm2", - vTrack.getPosition(0).getCoordinates() / 1_m, - vTrack.getPosition(1).getCoordinates() / 1_m, - grammageStart / 1_g * square(1_cm), - grammageEnd / 1_g * square(1_cm)); + CORSIKA_LOG_TRACE("pos1={} m, pos2={}, X1={} g/cm2, X2={} g/cm2", + vTrack.getPosition(0).getCoordinates() / 1_m, + vTrack.getPosition(1).getCoordinates() / 1_m, + grammageStart / 1_g * square(1_cm), + grammageEnd / 1_g * square(1_cm)); // Note: particle may go also "upward", thus, grammageEnd<grammageStart const int binStart = std::ceil(grammageStart / dX_); diff --git a/corsika/detail/modules/StackInspector.inl b/corsika/detail/modules/StackInspector.inl index 10b0f1582..caae6c59c 100644 --- a/corsika/detail/modules/StackInspector.inl +++ b/corsika/detail/modules/StackInspector.inl @@ -24,8 +24,8 @@ namespace corsika { template <typename TStack> - inline StackInspector<TStack>::StackInspector(const int vNStep, const bool vReportStack, - const HEPEnergyType vE0) + inline StackInspector<TStack>::StackInspector(int const vNStep, bool const vReportStack, + HEPEnergyType const vE0) : StackProcess<StackInspector<TStack>>(vNStep) , ReportStack_(vReportStack) , E0_(vE0) @@ -35,12 +35,12 @@ namespace corsika { inline StackInspector<TStack>::~StackInspector() {} template <typename TStack> - inline void StackInspector<TStack>::doStack(const TStack& vS) { + inline void StackInspector<TStack>::doStack(TStack const& vS) { [[maybe_unused]] int i = 0; HEPEnergyType Etot = 0_GeV; - for (const auto& iterP : vS) { + for (auto const& iterP : vS) { HEPEnergyType E = iterP.getEnergy(); Etot += E; if (ReportStack_) { @@ -60,7 +60,7 @@ namespace corsika { } auto const now = std::chrono::system_clock::now(); - const std::chrono::duration<double> elapsed_seconds = now - StartTime_; + std::chrono::duration<double> const elapsed_seconds = now - StartTime_; std::time_t const now_time = std::chrono::system_clock::to_time_t(now); auto const dE = E0_ - Etot; if (dE < dE_threshold_) return; diff --git a/corsika/detail/modules/pythia8/Decay.inl b/corsika/detail/modules/pythia8/Decay.inl index d41c14902..d58ac21d8 100644 --- a/corsika/detail/modules/pythia8/Decay.inl +++ b/corsika/detail/modules/pythia8/Decay.inl @@ -82,7 +82,7 @@ namespace corsika::pythia8 { inline void Decay::setHandleDecay(Code const vParticleCode) { handleAllDecays_ = false; - CORSIKA_LOG_INFO("Pythia::Decay: set to handle decay of {} ", vParticleCode); + CORSIKA_LOG_DEBUG("Pythia::Decay: set to handle decay of {} ", vParticleCode); if (Decay::canHandleDecay(vParticleCode)) handledDecays_.insert(vParticleCode); else @@ -106,12 +106,12 @@ namespace corsika::pythia8 { } inline void Decay::setUnstable(Code const pCode) { - CORSIKA_LOG_INFO("Pythia::Decay: setting {} unstable..", pCode); + CORSIKA_LOG_DEBUG("Pythia::Decay: setting {} unstable..", pCode); Pythia8::Pythia::particleData.mayDecay(static_cast<int>(get_PDG(pCode)), true); } inline void Decay::setStable(Code const pCode) { - CORSIKA_LOG_INFO("Pythia::Decay: setting {} stable..", pCode); + CORSIKA_LOG_DEBUG("Pythia::Decay: setting {} stable..", pCode); Pythia8::Pythia::particleData.mayDecay(static_cast<int>(get_PDG(pCode)), false); } @@ -122,8 +122,8 @@ namespace corsika::pythia8 { inline bool Decay::canDecay(Code const pCode) { bool const ans = Pythia8::Pythia::particleData.canDecay(static_cast<int>(get_PDG(pCode))); - CORSIKA_LOG_INFO("Pythia::Decay: checking if particle: {} can decay in PYTHIA? {} ", - pCode, ans); + CORSIKA_LOG_DEBUG("Pythia::Decay: checking if particle: {} can decay in PYTHIA? {} ", + pCode, ans); return ans; } @@ -153,13 +153,13 @@ namespace corsika::pythia8 { TimeType const t0 = get_lifetime(pid); auto const lifetime = gamma * t0; - CORSIKA_LOG_INFO("Pythia::Decay: code: {}", particle.getPID()); - CORSIKA_LOG_INFO("Pythia::Decay: MinStep: t0: {}", t0); - CORSIKA_LOG_INFO("Pythia::Decay: MinStep: energy: {} GeV", E / 1_GeV); - CORSIKA_LOG_INFO("Pythia::Decay: momentum: {} GeV", - particle.getMomentum().getComponents() / 1_GeV); - CORSIKA_LOG_INFO("Pythia::Decay: MinStep: gamma: {}", gamma); - CORSIKA_LOG_INFO("Pythia::Decay: MinStep: tau: {} ", lifetime); + CORSIKA_LOG_TRACE("Pythia::Decay: code: {}", particle.getPID()); + CORSIKA_LOG_TRACE("Pythia::Decay: MinStep: t0: {}", t0); + CORSIKA_LOG_TRACE("Pythia::Decay: MinStep: energy: {} GeV", E / 1_GeV); + CORSIKA_LOG_TRACE("Pythia::Decay: momentum: {} GeV", + particle.getMomentum().getComponents() / 1_GeV); + CORSIKA_LOG_TRACE("Pythia::Decay: MinStep: gamma: {}", gamma); + CORSIKA_LOG_TRACE("Pythia::Decay: MinStep: tau: {} ", lifetime); return lifetime; } else @@ -209,7 +209,7 @@ namespace corsika::pythia8 { if (!Pythia8::Pythia::next()) throw std::runtime_error("Pythia::Decay: decay failed!"); else - CORSIKA_LOG_INFO("Pythia::Decay: particles after decay: {} ", event.size()); + CORSIKA_LOG_DEBUG("Pythia::Decay: particles after decay: {} ", event.size()); if (print_listing_) { // list final state @@ -227,9 +227,10 @@ namespace corsika::pythia8 { FourVector const fourMomRest{Erest, pRest}; auto const fourMomLab = boost.fromCoM(fourMomRest); - CORSIKA_LOG_INFO("particle: id={} momentum={} energy={} ", pyId, - fourMomLab.getSpaceLikeComponents().getComponents(labCS) / 1_GeV, - fourMomLab.getTimeLikeComponent()); + CORSIKA_LOG_TRACE( + "particle: id={} momentum={} energy={} ", pyId, + fourMomLab.getSpaceLikeComponents().getComponents(labCS) / 1_GeV, + fourMomLab.getTimeLikeComponent()); view.addSecondary(std::make_tuple(pyId, fourMomLab.getTimeLikeComponent(), fourMomLab.getSpaceLikeComponents(), decayPoint, diff --git a/corsika/framework/process/ContinuousProcessIndex.hpp b/corsika/framework/process/ContinuousProcessIndex.hpp index bc2617d2b..9d6e84d0f 100644 --- a/corsika/framework/process/ContinuousProcessIndex.hpp +++ b/corsika/framework/process/ContinuousProcessIndex.hpp @@ -18,6 +18,8 @@ namespace corsika { class ContinuousProcessIndex { public: + ContinuousProcessIndex() + : id_(-1) {} // default ContinuousProcessIndex(int const id) : id_(id) {} void setIndex(int const id) { id_ = id; } diff --git a/do-clang-format.py b/do-clang-format.py index 3bff0e619..d6ea59753 100755 --- a/do-clang-format.py +++ b/do-clang-format.py @@ -10,6 +10,14 @@ import os import sys import re +do_progress = False +try: + from progress.bar import ChargingBar + do_progress = True +except ImportError as e: + do_progress = False + # no progress bar + parser = argparse.ArgumentParser(description=__doc__) parser.add_argument('--apply', action="store_true", help="Apply clang-format to files which need changes.") @@ -85,17 +93,26 @@ version = subp.check_output(cmd.split() + ["--version"]).decode("utf-8") print (version) print ("Note: the clang-format version has an impact on the result. Make sure you are consistent with current CI. Consider \'--docker\' option.") + +bar = None +if do_progress: + bar = ChargingBar('Processing', max=len(filelist)) + if args.apply: for filename in filelist: + if bar: bar.next() subp.check_call(cmd.split() + ["-i", filename]) + if bar: bar.finish() else: # only print files which need formatting files_need_formatting = 0 for filename in filelist: + if bar: bar.next() a = open(filename, "rb").read() b = subp.check_output(cmd.split() + [filename]) if a != b: files_need_formatting += 1 print(filename) + if bar: bar.finish() sys.exit(1 if files_need_formatting > 0 else 0) diff --git a/examples/vertical_EAS.cpp b/examples/vertical_EAS.cpp index 59893527a..7738e957d 100644 --- a/examples/vertical_EAS.cpp +++ b/examples/vertical_EAS.cpp @@ -6,8 +6,6 @@ * the license. */ -#define DEBUG 1 - /* clang-format off */ // InteractionCounter used boost/histogram, which // fails if boost/type_traits have been included before. Thus, we have @@ -96,7 +94,7 @@ using MyExtraEnv = MediumPropertyModel<UniformMagneticField<T>>; int main(int argc, char** argv) { corsika_logger->set_pattern("[%n:%^%-8l%$] %s:%#: %v"); - logging::set_level(logging::level::trace); + logging::set_level(logging::level::info); CORSIKA_LOG_INFO("vertical_EAS"); -- GitLab