diff --git a/Documentation/Examples/vertical_EAS.cc b/Documentation/Examples/vertical_EAS.cc index f7b1868962883510b5c1730609384053bf5d40c6..be02d308c7b42deea043eae2413cc89ad843802f 100644 --- a/Documentation/Examples/vertical_EAS.cc +++ b/Documentation/Examples/vertical_EAS.cc @@ -74,8 +74,6 @@ void registerRandomStreams(const int seed) { int main(int argc, char** argv) { - logging::SetLevel(logging::level::debug); - C8LOG_INFO("vertical_EAS"); if (argc < 4) { diff --git a/Framework/Cascade/Cascade.h b/Framework/Cascade/Cascade.h index de033d190b037c5e30b713be5d9d6b39eedfeef5..f5f9040b05e343d83f019b1853f21ca53749b03c 100644 --- a/Framework/Cascade/Cascade.h +++ b/Framework/Cascade/Cascade.h @@ -17,6 +17,8 @@ #include <corsika/stack/SecondaryView.h> #include <corsika/units/PhysicalUnits.h> +#include <corsika/logging/Logging.h> + #include <corsika/setup/SetupTrajectory.h> /* see Issue 161, we need to include SetupStack only because we need @@ -61,8 +63,9 @@ namespace corsika::cascade { template <typename TTracking, typename TProcessList, typename TStack, /* - TStackView is needed as template parameter because of issue 161 and the - inability of clang to understand "MakeView" so far. + TStackView is needed as explicit template parameter because + of issue 161 and the + inability of clang to understand "stack::MakeView" so far. */ typename TStackView = corsika::setup::StackView> class Cascade { @@ -96,9 +99,8 @@ namespace corsika::cascade { , fTracking(tr) , fProcessSequence(pl) , fStack(stack) - , energy_cut_(0 * corsika::units::si::electronvolt) { - - std::cout << c8_ascii_ << std::endl; + , count_(0) { + C8LOG_INFO(c8_ascii_); } corsika::units::si::HEPEnergyType GetEnergyCut() const { return energy_cut_; } @@ -124,13 +126,17 @@ namespace corsika::cascade { while (!fStack.IsEmpty()) { while (!fStack.IsEmpty()) { + C8LOG_TRACE(fmt::format("Stack: {}", fStack.as_string())); + for (const auto& tmp_p : fStack) + std::cout << " test " << tmp_p.GetPID() << std::endl; count_++; auto pNext = fStack.GetNextParticle(); - std::cout << "========= next: count=" << count_ << ", pid=" << pNext.GetPID() - << ", stack entries=" << fStack.getEntries() - << ", stack deleted=" << fStack.getDeleted() << std::endl; + C8LOG_DEBUG(fmt::format( + "============== next particle : count={}, pid={}, " + ", stack entries={}" + ", stack deleted={}", + count_, pNext.GetPID(), fStack.getEntries(), fStack.getDeleted())); Step(pNext); - std::cout << "========= stack ============" << std::endl; fProcessSequence.DoStack(fStack); } // do cascade equations, which can put new particles on Stack, @@ -145,7 +151,7 @@ namespace corsika::cascade { * want to call forceInteraction() for the primary interaction. */ void forceInteraction() { - std::cout << "forced interaction!" << std::endl; + C8LOG_DEBUG("forced interaction!"); auto vParticle = fStack.GetNextParticle(); TStackView secondaries(vParticle); interaction(vParticle, secondaries); @@ -180,8 +186,10 @@ namespace corsika::cascade { corsika::random::ExponentialDistribution expDist(1 / total_inv_lambda); GrammageType const next_interact = expDist(fRNG); - std::cout << "total_inv_lambda=" << total_inv_lambda - << ", next_interact=" << next_interact << std::endl; + C8LOG_DEBUG( + "total_lambda={} g/cm2, " + ", next_interact={} g/cm2", + double((1. / total_inv_lambda) / 1_g * 1_cm * 1_cm), double(next_interact / 1_g * 1_cm * 1_cm)); auto const* currentLogicalNode = vParticle.GetNode(); @@ -197,7 +205,7 @@ namespace corsika::cascade { // determine the maximum geometric step length from continuous processes LengthType const distance_max = fProcessSequence.MaxStepLength(vParticle, step); - std::cout << "distance_max=" << distance_max << std::endl; + C8LOG_DEBUG("distance_max={} m", distance_max/1_m); // determine combined total inverse decay time InverseTimeType const total_inv_lifetime = @@ -206,8 +214,10 @@ namespace corsika::cascade { // sample random exponential decay time corsika::random::ExponentialDistribution expDistDecay(1 / total_inv_lifetime); TimeType const next_decay = expDistDecay(fRNG); - std::cout << "total_inv_lifetime=" << total_inv_lifetime - << ", next_decay=" << next_decay << std::endl; + C8LOG_DEBUG( + "total_lifetime={} s" + ", next_decay={} s", + (1/total_inv_lifetime)/1_s, next_decay/1_s); // convert next_decay from time to length [m] LengthType const distance_decay = next_decay * vParticle.GetMomentum().norm() / @@ -217,7 +227,7 @@ namespace corsika::cascade { auto const min_distance = std::min({distance_interact, distance_decay, distance_max, geomMaxLength}); - std::cout << " move particle by : " << min_distance << std::endl; + C8LOG_DEBUG("transport particle by : {} m", min_distance/1_m); // here the particle is actually moved along the trajectory to new position: // std::visit(setup::ParticleUpdate<Particle>{vParticle}, step); @@ -231,15 +241,14 @@ namespace corsika::cascade { process::EProcessReturn status = fProcessSequence.DoContinuous(vParticle, step); if (status == process::EProcessReturn::eParticleAbsorbed) { - std::cout << "Cascade: delete absorbed particle " << vParticle.GetPID() << " " - << vParticle.GetEnergy() / 1_GeV << "GeV" << std::endl; - energy_cut_ += vParticle.GetEnergy(); + C8LOG_DEBUG("Cascade: delete absorbed particle PID={} E={} GeV", + vParticle.GetPID(), vParticle.GetEnergy() / 1_GeV); vParticle.Delete(); return; } - std::cout << "sth. happening before geometric limit ? " - << ((min_distance < geomMaxLength) ? "yes" : "no") << std::endl; + C8LOG_DEBUG("sth. happening before geometric limit ? {}", + ((min_distance < geomMaxLength) ? "yes" : "no")); if (min_distance < geomMaxLength) { // interaction to happen within geometric limit @@ -257,7 +266,7 @@ namespace corsika::cascade { 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 (and not vParticle) for Interaction, + important to use projectle/view (and not vParticle) for Interaction, and Decay! */ @@ -277,15 +286,10 @@ namespace corsika::cascade { } fProcessSequence.DoSecondaries(secondaries); - vParticle.Delete(); // todo: this should be reviewed. Where - // exactly are particles best deleted, and - // where they should NOT be - // deleted... maybe Delete function should - // be "protected" and not accessible to physics + vParticle.Delete(); } else { // step-length limitation within volume - - std::cout << "step-length limitation" << std::endl; + C8LOG_DEBUG("step-length limitation"); // no extra physics happens here. just proceed to next step. } @@ -299,10 +303,8 @@ namespace corsika::cascade { }; assert(assertion()); // numerical and logical nodes don't match - - } else { // boundary crossing, step is limited by volume boundary - - std::cout << "boundary crossing! next node = " << nextVol << std::endl; + } else { // boundary crossing, step is limited by volume boundary + // C8LOG_DEBUG("boundary crossing! next node = {}", int(nextVol)); vParticle.SetNode(nextVol); /* DoBoundary may delete the particle (or not) @@ -316,7 +318,7 @@ namespace corsika::cascade { } auto decay(Particle& particle, TStackView& view) { - std::cout << "decay" << std::endl; + C8LOG_DEBUG("decay"); units::si::InverseTimeType const actual_decay_time = fProcessSequence.GetTotalInverseLifetime(particle); @@ -329,7 +331,7 @@ namespace corsika::cascade { } auto interaction(Particle& particle, TStackView& view) { - std::cout << "collide" << std::endl; + C8LOG_DEBUG("collide"); units::si::InverseGrammageType const current_inv_length = fProcessSequence.GetTotalInverseInteractionLength(particle); diff --git a/Framework/StackInterface/CombinedStack.h b/Framework/StackInterface/CombinedStack.h index 2534ced4e89e3b5fc6fc231c57e0d0887d760ee0..733f56b9a43b08a1045bb51eea44e0ffc251eccf 100644 --- a/Framework/StackInterface/CombinedStack.h +++ b/Framework/StackInterface/CombinedStack.h @@ -11,6 +11,7 @@ #include <corsika/particles/ParticleProperties.h> #include <corsika/stack/Stack.h> #include <corsika/units/PhysicalUnits.h> +#include <corsika/logging/Logging.h> namespace corsika::stack { @@ -37,10 +38,6 @@ namespace corsika::stack { class CombinedParticleInterface : public ParticleInterfaceB<ParticleInterfaceA<StackIterator>> { - // template<template <typename> typename _PI> - // template <typename StackDataType, template <typename> typename ParticleInterface> - // template<typename T1, template <typename> typename T2> friend class Stack<T1, T2>; - using PI_C = CombinedParticleInterface<ParticleInterfaceA, ParticleInterfaceB, StackIterator>; using PI_A = ParticleInterfaceA<StackIterator>; @@ -91,6 +88,11 @@ namespace corsika::stack { PI_B::SetParticleData(static_cast<PI_B&>(p), vB); } ///@} + + std::string as_string() const { + return fmt::format("[[{}][{}]]",PI_A::as_string(), PI_B::as_string()); + } + }; /** @@ -113,7 +115,7 @@ namespace corsika::stack { unsigned int GetSize() const { return Stack1Impl::GetSize(); } unsigned int GetCapacity() const { return Stack1Impl::GetCapacity(); } - + /** * Function to copy particle at location i1 in stack to i2 */ diff --git a/Framework/StackInterface/SecondaryView.h b/Framework/StackInterface/SecondaryView.h index 732779a2473a94d2d93557ef479bb30df375a9ca..d69e6aabbddd241ab39136cdfa50032d077be10a 100644 --- a/Framework/StackInterface/SecondaryView.h +++ b/Framework/StackInterface/SecondaryView.h @@ -276,6 +276,7 @@ namespace corsika::stack { * */ void Delete(StackIterator p) { + C8LOG_TRACE("SecondaryView::Delete"); if (IsEmpty()) { /*error*/ throw std::runtime_error("Stack, cannot delete entry since size is zero"); } @@ -323,6 +324,7 @@ namespace corsika::stack { * the function will just return false and do nothing. */ bool purgeLastIfDeleted() { + C8LOG_TRACE("SecondaryView::purgeLastIfDeleted"); if (!isDeleted(getSize() - 1)) return false; // the last particle is not marked for deletion. Do nothing. inner_stack_.purge(GetIndexFromIterator(getSize())); @@ -367,6 +369,7 @@ namespace corsika::stack { * performed. */ unsigned int GetIndexFromIterator(const unsigned int vI) const { + // this is too much: C8LOG_TRACE("SecondaryView::GetIndexFromIterator({})={}", vI, (vI?indices_[vI-1]:projectile_index_)); if (vI == 0) return projectile_index_; return indices_[vI - 1]; } diff --git a/Framework/StackInterface/Stack.h b/Framework/StackInterface/Stack.h index d1c10ca97b3bed7c88d14605b983c928370585d3..2158c0fd6ebdea20adcc02b09aed980c55e4c70e 100644 --- a/Framework/StackInterface/Stack.h +++ b/Framework/StackInterface/Stack.h @@ -10,10 +10,12 @@ #include <corsika/stack/StackIteratorInterface.h> #include <corsika/utl/MetaProgramming.h> +#include <corsika/logging/Logging.h> #include <stdexcept> #include <type_traits> #include <vector> +#include <string> /** All classes around management of particles on a stack. @@ -85,9 +87,9 @@ namespace corsika::stack { * stacks, where the inner data container is always a reference * and cannot be initialized here. */ - template <typename... Args, typename _ = TStackData, + template <typename... TArgs, typename _ = TStackData, typename = utl::disable_if<std::is_reference<_>>> - Stack(Args... args) + Stack(TArgs... args) : data_(args...) , deleted_(std::vector<bool>(data_.GetSize(), false)) , nDeleted_(0) {} @@ -125,7 +127,6 @@ namespace corsika::stack { friend class StackIteratorInterface<StackDataValueType, MParticleInterface, Stack>; friend class ConstStackIteratorInterface<StackDataValueType, MParticleInterface, Stack>; - // friend class SecondaryView<StackDataValueType, MParticleInterface>; template <typename T1, //=TStackData, template <typename> typename M1, //=MParticleInterface, @@ -145,8 +146,8 @@ namespace corsika::stack { unsigned int getDeleted() const { return nDeleted_; } unsigned int getEntries() const { return getSize() - getDeleted(); } - template <typename... Args> - void Clear(Args... args) { + template <typename... TArgs> + void Clear(TArgs... args) { data_.Clear(args...); deleted_ = std::vector<bool>(data_.GetSize(), false); nDeleted_ = 0; @@ -195,7 +196,7 @@ namespace corsika::stack { for (; i < getSize(); ++i) { if (!deleted_[i]) break; } - return ConstStackIterator(*this, 0); + return ConstStackIterator(*this, i); } ConstStackIterator cend() const { return ConstStackIterator(*this, getSize()); } ConstStackIterator clast() const { @@ -215,8 +216,9 @@ namespace corsika::stack { /** * increase stack size, create new particle at end of stack */ - template <typename... Args> - StackIterator AddParticle(const Args... v) { + template <typename... TArgs> + StackIterator AddParticle(const TArgs... v) { + C8LOG_TRACE("Stack::AddParticle"); data_.IncrementSize(); deleted_.push_back(false); return StackIterator(*this, getSize() - 1, v...); @@ -230,8 +232,9 @@ namespace corsika::stack { * This should only get internally called from a * StackIterator::AddSecondary via ParticleBase */ - template <typename... Args> - StackIterator AddSecondary(StackIterator& parent, const Args... v) { + template <typename... TArgs> + StackIterator AddSecondary(StackIterator& parent, const TArgs... v) { + C8LOG_TRACE("Stack::AddSecondary"); data_.IncrementSize(); deleted_.push_back(false); return StackIterator(*this, getSize() - 1, parent, v...); @@ -239,27 +242,43 @@ namespace corsika::stack { public: void Swap(StackIterator a, StackIterator b) { + C8LOG_TRACE("Stack::Swap"); data_.Swap(a.GetIndex(), b.GetIndex()); std::swap(deleted_[a.GetIndex()], deleted_[b.GetIndex()]); } void Copy(StackIterator a, StackIterator b) { + C8LOG_TRACE("Stack::Copy"); data_.Copy(a.GetIndex(), b.GetIndex()); if (deleted_[b.GetIndex()] && !deleted_[a.GetIndex()]) nDeleted_--; if (!deleted_[b.GetIndex()] && deleted_[a.GetIndex()]) nDeleted_++; deleted_[b.GetIndex()] = deleted_[a.GetIndex()]; } void Copy(ConstStackIterator a, StackIterator b) { + C8LOG_TRACE("Stack::Copy"); data_.Copy(a.GetIndex(), b.GetIndex()); if (deleted_[b.GetIndex()] && !deleted_[a.GetIndex()]) nDeleted_--; if (!deleted_[b.GetIndex()] && deleted_[a.GetIndex()]) nDeleted_++; deleted_[b.GetIndex()] = deleted_[a.GetIndex()]; } + std::string as_string() const { + std::string str(fmt::format("size {}, entries {}, deleted {} \n", getSize(), getEntries(), getDeleted())); + // we make our own begin/end since we want ALL entries + std::string new_line = " "; + for (unsigned int iPart = 0; iPart!=getSize(); ++iPart) { + ConstStackIterator itPart(*this, iPart); + str += fmt::format("{}{}{}", new_line, itPart.as_string(), (deleted_[itPart.GetIndex()]?" [deleted]":"")); + new_line = "\n "; + } + return str; + } + /** * delete this particle */ public: void Delete(StackIterator p) { + C8LOG_TRACE("Stack::Delete"); if (IsEmpty()) { /*error*/ throw std::runtime_error("Stack, cannot delete entry since size is zero"); } @@ -294,6 +313,7 @@ namespace corsika::stack { bool purgeLastIfDeleted() { if (!deleted_.back()) return false; // the last particle is not marked for deletion. Do nothing. + C8LOG_TRACE("Stack::purgeLastIfDeleted: yes"); data_.DecrementSize(); nDeleted_--; deleted_.pop_back(); @@ -359,7 +379,10 @@ namespace corsika::stack { * TStackData data_. By default (and in almost all cases) this * should just be identiy. See class SecondaryView for an alternative implementation. */ - unsigned int GetIndexFromIterator(const unsigned int vI) const { return vI; } + unsigned int GetIndexFromIterator(const unsigned int vI) const { + // this is too much: C8LOG_TRACE("Stack::GetIndexFromIterator({})={}", vI, vI); + return vI; + } /** * @name Return reference to TStackData object data_ for data access diff --git a/Setup/SetupStack.h b/Setup/SetupStack.h index c58d8f0d3fa9bc4dd3519f8ecd6e10b76573f6fa..b03e70abb5a954dd9f0c31a70586a71f6c508190 100644 --- a/Setup/SetupStack.h +++ b/Setup/SetupStack.h @@ -86,7 +86,7 @@ namespace corsika::setup { typename corsika::setup::Stack::StackImpl, // CHECK with CLANG: corsika::setup::Stack::MPIType>; // corsika::setup::detail::StackWithGeometryInterface>; - corsika::setup::detail::StackWithHistoryInterface>; + corsika::setup::detail::StackWithHistoryInterface, StackViewProducer>; #elif defined(__GNUC__) || defined(__GNUG__) using TheStackView = corsika::stack::MakeView<corsika::setup::Stack, StackViewProducer>::type; diff --git a/Stack/DummyStack/DummyStack.h b/Stack/DummyStack/DummyStack.h index 8427475abd443584c74362716c956bdb9eb98261..a6cfb02be6cba78e1544e0ab9c3cf1e542ee6d6d 100644 --- a/Stack/DummyStack/DummyStack.h +++ b/Stack/DummyStack/DummyStack.h @@ -9,10 +9,12 @@ #pragma once #include <corsika/particles/ParticleProperties.h> +#include <corsika/logging/Logging.h> #include <corsika/stack/Stack.h> #include <corsika/units/PhysicalUnits.h> #include <tuple> +#include <string> namespace corsika::stack { @@ -46,6 +48,8 @@ namespace corsika::stack { void SetParticleData(const std::tuple<NoData>& /*v*/) {} void SetParticleData(ParticleInterface<StackIteratorInterface>& /*parent*/, const std::tuple<NoData>& /*v*/) {} + + std::string as_string() const { return "dummy-data"; } }; /** diff --git a/Stack/GeometryNodeStackExtension/GeometryNodeStackExtension.h b/Stack/GeometryNodeStackExtension/GeometryNodeStackExtension.h index d33997b8eecddca3d53be9fd4f1f7b1b95ad4835..4b87872f2ffc5f4114a3140c60369b41422a5b22 100644 --- a/Stack/GeometryNodeStackExtension/GeometryNodeStackExtension.h +++ b/Stack/GeometryNodeStackExtension/GeometryNodeStackExtension.h @@ -10,6 +10,7 @@ #include <corsika/logging/Logging.h> #include <corsika/stack/Stack.h> +#include <corsika/logging/Logging.h> #include <tuple> #include <utility> diff --git a/Stack/History/Event.hpp b/Stack/History/Event.hpp index c89f8064eb8e4a527384ba87c8aee9c3601546f5..327fdf670cec232b65e5c98e0a001c557a11d0ad 100644 --- a/Stack/History/Event.hpp +++ b/Stack/History/Event.hpp @@ -10,6 +10,7 @@ #include <corsika/particles/ParticleProperties.h> #include <corsika/history/SecondaryParticle.hpp> +#include <corsika/logging/Logging.h> #include <iostream> #include <memory> @@ -20,7 +21,7 @@ namespace corsika::history { class Event; using EventPtr = std::shared_ptr<history::Event>; - + class Event { size_t projectileIndex_ = 0; //!< index of projectile on stack @@ -28,20 +29,23 @@ namespace corsika::history { EventPtr parent_event_; std::optional<corsika::particles::Code> - targetCode_; // cannot be const, value set only after construction + targetCode_; public: Event() = default; void setParentEvent(EventPtr const& evt) { parent_event_ = evt; } - EventPtr parentEvent() { return parent_event_; } - + bool hasParentEvent() const { return bool(parent_event_); } + EventPtr& parentEvent() { return parent_event_; } + const EventPtr& parentEvent() const { return parent_event_; } + void setProjectileIndex(size_t i) { projectileIndex_ = i; } size_t projectileIndex() const { return projectileIndex_; } template <typename TStackIterator> TStackIterator projectile(TStackIterator begin) { + // todo: change this // MR: This is dangerous. You can pass any iterator though it must // be stack.begin() to yield the correct projectile @@ -58,6 +62,10 @@ namespace corsika::history { std::vector<SecondaryParticle>& secondaries() { return secondaries_; } void setTargetCode(const particles::Code t) { targetCode_ = t; } + + std::string as_string() const { + return fmt::format("hasParent={}, projIndex={}, Nsec={}", hasParentEvent(), projectileIndex_, secondaries_.size()); + } }; } // namespace corsika::history diff --git a/Stack/History/HistoryObservationPlane.cc b/Stack/History/HistoryObservationPlane.cc index d9374eee7e00ee8177cd3027b86a9fac7bc6b9cd..16a3ae3d593b5e7de7d11eb079b806542b7f2df9 100644 --- a/Stack/History/HistoryObservationPlane.cc +++ b/Stack/History/HistoryObservationPlane.cc @@ -8,6 +8,7 @@ #include <corsika/logging/Logging.h> #include <corsika/history/HistoryObservationPlane.hpp> +#include <corsika/logging/Logging.h> #include <boost/histogram/ostream.hpp> diff --git a/Stack/History/HistoryStackExtension.h b/Stack/History/HistoryStackExtension.h index 76f136bcbbfaac621126160288a1d76e74b0f6e9..5a4881a6063092da95bdf95f72fd42c429161a8b 100644 --- a/Stack/History/HistoryStackExtension.h +++ b/Stack/History/HistoryStackExtension.h @@ -43,7 +43,8 @@ namespace corsika::history { // custom data access function void SetEvent(const int i, EventPtr v) { historyData_[i].first = std::move(v); } - EventPtr GetEvent(const int i) const { return historyData_[i].first; } + const EventPtr& GetEvent(const int i) const { return historyData_[i].first; } + EventPtr& GetEvent(const int i) { return historyData_[i].first; } void SetParentEventIndex(const int i, ParentEventIndex v) { historyData_[i].second = std::move(v); @@ -108,6 +109,10 @@ namespace corsika::history { int GetParentEventIndex() const { return GetStackData().GetParentEventIndex(GetIndex()); } + + std::string as_string() const { + return fmt::format("i_parent={}, [evt: {}]", GetParentEventIndex(), (bool(GetEvent())?GetEvent()->as_string():"n/a")); + } }; template <typename T, typename TEvent> diff --git a/Stack/History/testHistoryView.cc b/Stack/History/testHistoryView.cc index 6461ed1df759807e3fa476fb4521c94d965f1ea3..7e090e0b5921efaa6303ec25c0b8899796db0b62 100644 --- a/Stack/History/testHistoryView.cc +++ b/Stack/History/testHistoryView.cc @@ -14,6 +14,8 @@ #include <corsika/stack/dummy/DummyStack.h> #include <corsika/stack/nuclear_extension/NuclearStackExtension.h> +#include <corsika/logging/Logging.h> + #include <catch2/catch.hpp> #include <iostream> diff --git a/Stack/NuclearStackExtension/NuclearStackExtension.h b/Stack/NuclearStackExtension/NuclearStackExtension.h index 61c663e7e4b1a9580932b86e7d94a47cbcbfa6e0..4409a215eed0c331e922b311235f37e10ded24d8 100644 --- a/Stack/NuclearStackExtension/NuclearStackExtension.h +++ b/Stack/NuclearStackExtension/NuclearStackExtension.h @@ -18,6 +18,8 @@ #include <corsika/geometry/Point.h> #include <corsika/geometry/Vector.h> +#include <corsika/logging/Logging.h> + #include <algorithm> #include <tuple> #include <vector> @@ -60,6 +62,7 @@ namespace corsika::stack { protected: using InnerParticleInterface<StackIteratorInterface>::GetStackData; using InnerParticleInterface<StackIteratorInterface>::GetIndex; + using InnerParticleInterface<StackIteratorInterface>::as_string; public: void SetParticleData( @@ -141,6 +144,13 @@ namespace corsika::stack { std::get<4>(v)}); } + std::string as_string() const { + return fmt::format( + "{}, nuc({})", InnerParticleInterface<StackIteratorInterface>::as_string(), + (isNucleus() ? fmt::format("A={}, Z={}", GetNuclearA(), GetNuclearZ()) + : "n/a")); + } + /** * @name individual setters * @{ @@ -184,6 +194,7 @@ namespace corsika::stack { protected: void SetNucleusRef(const int vR) { GetStackData().SetNucleusRef(GetIndex(), vR); } + bool isNucleus() const { return GetStackData().isNucleus(GetIndex()); } }; /** @@ -240,6 +251,8 @@ namespace corsika::stack { throw std::runtime_error(err.str()); } + bool isNucleus(const unsigned int i) const { return fNucleusRef[i] >= 0; } + /** * Function to copy particle at location i1 in stack to i2 */ diff --git a/Stack/SuperStupidStack/CMakeLists.txt b/Stack/SuperStupidStack/CMakeLists.txt index 4edc8001c9276908ea569df2c0c601e4eeaff0e1..55b88c69c8c232a96c111e441426f4f23b8aca6a 100644 --- a/Stack/SuperStupidStack/CMakeLists.txt +++ b/Stack/SuperStupidStack/CMakeLists.txt @@ -12,6 +12,7 @@ target_link_libraries ( CORSIKAunits CORSIKAparticles CORSIKAgeometry + CORSIKAlogging ) target_include_directories ( @@ -33,8 +34,6 @@ install ( CORSIKA_ADD_TEST(testSuperStupidStack) target_link_libraries ( testSuperStupidStack - CORSIKAgeometry - CORSIKAparticles - CORSIKAunits + SuperStupidStack CORSIKAtesting ) diff --git a/Stack/SuperStupidStack/SuperStupidStack.h b/Stack/SuperStupidStack/SuperStupidStack.h index 7c88628b6c95bd18f7654b2e5215e3f184d4f5c8..eae265cc72d5329cc8d18cad7b762609bee06e9e 100644 --- a/Stack/SuperStupidStack/SuperStupidStack.h +++ b/Stack/SuperStupidStack/SuperStupidStack.h @@ -16,6 +16,7 @@ #include <corsika/geometry/RootCoordinateSystem.h> // remove #include <corsika/geometry/Vector.h> +#include <string> #include <tuple> #include <vector> @@ -40,6 +41,12 @@ namespace corsika::stack { using corsika::stack::ParticleBase<StackIteratorInterface>::GetIndex; public: + std::string as_string() const { + using namespace corsika::units::si; + return fmt::format("particle: i={}, PID={}, E={}GeV", GetIndex(), + particles::GetName(GetPID()), GetEnergy() / 1_GeV); + } + void SetParticleData( const std::tuple<corsika::particles::Code, corsika::units::si::HEPEnergyType, MomentumVector, corsika::geometry::Point,