diff --git a/Documentation/Examples/geometry_example.cc b/Documentation/Examples/geometry_example.cc index 4c38d1dbfce534385482f480804ec1cd4ed62ae5..da23d81e02baae5381358ff13b12f93d071c5dfe 100644 --- a/Documentation/Examples/geometry_example.cc +++ b/Documentation/Examples/geometry_example.cc @@ -22,8 +22,8 @@ using namespace corsika::units::si; int main() { - std::cout << "geometry_example" << std::endl; - + std::cout << "geometry_example" << std::endl; + // define the root coordinate system geometry::CoordinateSystem& root = geometry::RootCoordinateSystem::GetInstance().GetRootCoordinateSystem(); diff --git a/Documentation/Examples/helix_example.cc b/Documentation/Examples/helix_example.cc index d200cc99678942c427bc3e1a769430d9810ac833..9f3eba67e610c87f63f193b3db3aad5978eda4fc 100644 --- a/Documentation/Examples/helix_example.cc +++ b/Documentation/Examples/helix_example.cc @@ -22,7 +22,7 @@ using namespace corsika::units::si; int main() { std::cout << "helix_example" << std::endl; - + geometry::CoordinateSystem& root = geometry::RootCoordinateSystem::GetInstance().GetRootCoordinateSystem(); diff --git a/Documentation/Examples/particle_list_example.cc b/Documentation/Examples/particle_list_example.cc index 5aefdd41401a4b6ecd67e737f1984166534c3bf3..342c08e7ff6e7cceec63f2aa3583f5fb46199297 100644 --- a/Documentation/Examples/particle_list_example.cc +++ b/Documentation/Examples/particle_list_example.cc @@ -27,7 +27,7 @@ using namespace std; int main() { std::cout << "particle_list_example" << std::endl; - + cout << "------------------------------------------" << "particles in CORSIKA" << "------------------------------------------" << endl; diff --git a/Documentation/Examples/stack_example.cc b/Documentation/Examples/stack_example.cc index dd7de9274ab51f2d0b2d6af55deea35ded47be42..cc638dbe8e622f27db71210a109784a0f4de7eff 100644 --- a/Documentation/Examples/stack_example.cc +++ b/Documentation/Examples/stack_example.cc @@ -51,7 +51,7 @@ void read(corsika::stack::super_stupid::SuperStupidStack& s) { int main() { std::cout << "stack_example" << std::endl; - + corsika::stack::super_stupid::SuperStupidStack s; fill(s); read(s); diff --git a/Documentation/Examples/staticsequence_example.cc b/Documentation/Examples/staticsequence_example.cc index 5959ca5f8fb0c5bcb2dba8a371995cdcb7b4116b..ca75fe4b22172488675fddca7601748acf8e3214 100644 --- a/Documentation/Examples/staticsequence_example.cc +++ b/Documentation/Examples/staticsequence_example.cc @@ -102,7 +102,7 @@ void modular() { int main() { std::cout << "staticsequence_example" << std::endl; - + modular(); return 0; } diff --git a/Documentation/Examples/stopping_power.cc b/Documentation/Examples/stopping_power.cc index 6da3663dca95b830f20bd9ded7beb02b4be0b89f..100594f6a0060a5b58ea5561586a0ee13af4d678 100644 --- a/Documentation/Examples/stopping_power.cc +++ b/Documentation/Examples/stopping_power.cc @@ -35,7 +35,7 @@ using namespace corsika::units::si; int main() { std::cout << "stopping_power" << std::endl; - + feenableexcept(FE_INVALID); // setup environment, geometry diff --git a/Framework/Cascade/testCascade.cc b/Framework/Cascade/testCascade.cc index 9db2c73ccba340df465f7ac75c518d163fef6cd2..57d565d55f68948ad4cf9f0e4a4db4a757da687a 100644 --- a/Framework/Cascade/testCascade.cc +++ b/Framework/Cascade/testCascade.cc @@ -77,11 +77,13 @@ public: view.AddSecondary( std::tuple<particles::Code, units::si::HEPEnergyType, corsika::stack::MomentumVector, geometry::Point, units::si::TimeType>{ - projectile.GetPID(), E / 2, projectile.GetMomentum(), projectile.GetPosition(), projectile.GetTime()}); + projectile.GetPID(), E / 2, projectile.GetMomentum(), + projectile.GetPosition(), projectile.GetTime()}); view.AddSecondary( std::tuple<particles::Code, units::si::HEPEnergyType, corsika::stack::MomentumVector, geometry::Point, units::si::TimeType>{ - projectile.GetPID(), E / 2, projectile.GetMomentum(), projectile.GetPosition(), projectile.GetTime()}); + projectile.GetPID(), E / 2, projectile.GetMomentum(), + projectile.GetPosition(), projectile.GetTime()}); return EProcessReturn::eInteracted; } diff --git a/Framework/StackInterface/SecondaryView.h b/Framework/StackInterface/SecondaryView.h index 4a33b404c9b036b16f198e5c7443c5e8d90f5c2c..27f8dc6f316a818f5d22a7f97bbb855c29e48fb7 100644 --- a/Framework/StackInterface/SecondaryView.h +++ b/Framework/StackInterface/SecondaryView.h @@ -43,13 +43,13 @@ namespace corsika::stack { *Further information about implementation (for developers):* All data is stored in the original stack privided at construction time. The secondary particle (view) indices are stored in an - extra std::vector of SecondaryView class 'fIndices' referring to + extra std::vector of SecondaryView class 'indices_' referring to the original stack slot indices. The index of the primary projectle particle is also explicitly stored in - 'fProjectileIndex'. StackIterator indices + 'projectile_index_'. StackIterator indices 'i = StackIterator::GetIndex()' are referring to those numbers, - where 'i==0' refers to the 'fProjectileIndex', and - 'StackIterator::GetIndex()>0' to 'fIndices[i-1]', see function + where 'i==0' refers to the 'projectile_index_', and + 'StackIterator::GetIndex()>0' to 'indices_[i-1]', see function GetIndexFromIterator. */ @@ -62,8 +62,8 @@ namespace corsika::stack { /** * Helper type for inside this class */ - using InnerStackType = Stack<StackDataType&, ParticleInterface>; - using InnerStackType::getDeleted; + using InnerStackTypeRef = Stack<StackDataType&, ParticleInterface>; + using InnerStackTypeRef::getDeleted; /** * @name We need this "special" types with non-reference StackData for @@ -76,6 +76,9 @@ namespace corsika::stack { using StackIteratorValue = StackIteratorInterface<typename std::remove_reference<StackDataType>::type, ParticleInterface, InnerStackTypeValue>; + using ConstStackIteratorValue = + ConstStackIteratorInterface<typename std::remove_reference<StackDataType>::type, + ParticleInterface, InnerStackTypeValue>; /// @} using StackIterator = @@ -106,9 +109,9 @@ namespace corsika::stack { SecondaryView(Args... args) = delete; private: - InnerStackTypeValue& innerstack_; - unsigned int fProjectileIndex; - std::vector<unsigned int> fIndices; + InnerStackTypeValue& inner_stack_; + unsigned int projectile_index_; + std::vector<unsigned int> indices_; public: /** @@ -117,14 +120,29 @@ namespace corsika::stack { **/ SecondaryView(StackIteratorValue& vI) : Stack<StackDataType&, ParticleInterface>(vI.GetStackData()) - , innerstack_(vI.GetStack()) - , fProjectileIndex(vI.GetIndex()) {} + , inner_stack_(vI.GetStack()) + , projectile_index_(vI.GetIndex()) {} + /** + * This returns the projectile/parent in the original Stack, where this + * SecondaryView is derived from. This projectile should not be + * used to modify the Stack! + */ + + ConstStackIteratorValue parent() const { + return ConstStackIteratorValue(inner_stack_, projectile_index_); + } + + /** + * This return a projectile of this SecondaryView, which can be + * used to modify the SecondaryView + */ StackIterator GetProjectile() { // NOTE: 0 is special marker here for PROJECTILE, see GetIndexFromIterator return StackIterator(*this, 0); } + public: template <typename... Args> StackIterator AddSecondary(const Args... v) { StackIterator proj = GetProjectile(); @@ -134,13 +152,13 @@ namespace corsika::stack { template <typename... Args> StackIterator AddSecondary(StackIterator& proj, const Args... v) { // make space on stack - InnerStackType::GetStackData().IncrementSize(); - innerstack_.deleted_.push_back(false); + InnerStackTypeRef::GetStackData().IncrementSize(); + inner_stack_.deleted_.push_back(false); // get current number of secondaries on stack const unsigned int idSec = getSize(); // determine index on (inner) stack where new particle will be located - const unsigned int index = InnerStackType::GetStackData().GetSize() - 1; - fIndices.push_back(index); + const unsigned int index = InnerStackTypeRef::GetStackData().GetSize() - 1; + indices_.push_back(index); // NOTE: "+1" is since "0" is special marker here for PROJECTILE, see // GetIndexFromIterator return StackIterator(*this, idSec + 1, proj, v...); @@ -149,7 +167,7 @@ namespace corsika::stack { /** * overwrite Stack::GetSize to return actual number of secondaries */ - unsigned int getSize() const { return fIndices.size(); } + unsigned int getSize() const { return indices_.size(); } unsigned int getEntries() const { return getSize() - getDeleted(); } bool IsEmpty() const { return getEntries() == 0; } @@ -215,7 +233,7 @@ namespace corsika::stack { * SecondaryView::DeleteLast * * The particle is deleted on the underlying (internal) stack. The - * local references in SecondaryView in fIndices must be fixed, + * local references in SecondaryView in indices_ must be fixed, * too. The approach is to a) check if the particle 'p' is at the * very end of the internal stack, b) if not: move it there by * copying the last particle to the current particle location, c) @@ -229,8 +247,8 @@ namespace corsika::stack { if (isDeleted(p.GetIndex() - 1)) { /*error*/ throw std::runtime_error("Stack, cannot delete entry since already deleted"); } - innerstack_.Delete(GetIndexFromIterator(p.GetIndex())); - InnerStackType::nDeleted_++; // also count in SecondaryView + inner_stack_.Delete(GetIndexFromIterator(p.GetIndex())); + InnerStackTypeRef::nDeleted_++; // also count in SecondaryView } /** @@ -272,9 +290,9 @@ namespace corsika::stack { bool purgeLastIfDeleted() { if (!isDeleted(getSize() - 1)) return false; // the last particle is not marked for deletion. Do nothing. - innerstack_.purge(GetIndexFromIterator(getSize())); - InnerStackType::nDeleted_--; - fIndices.pop_back(); + inner_stack_.purge(GetIndexFromIterator(getSize())); + InnerStackTypeRef::nDeleted_--; + indices_.pop_back(); return true; } @@ -291,31 +309,31 @@ namespace corsika::stack { unsigned int size = getSize(); while (iStack < size) { if (isDeleted(iStack)) { - innerstack_.purge(iStack); - fIndices.erase(fIndices.begin() + iStack); + inner_stack_.purge(iStack); + indices_.erase(indices_.begin() + iStack); } size = getSize(); iStack++; } - InnerStackType::nDeleted_ = 0; + InnerStackTypeRef::nDeleted_ = 0; } protected: // forward to inner stack // this also checks the allowed bounds of 'i' bool isDeleted(unsigned int i) const { - if (i >= fIndices.size()) return false; - return innerstack_.isDeleted(GetIndexFromIterator(i + 1)); + if (i >= indices_.size()) return false; + return inner_stack_.isDeleted(GetIndexFromIterator(i + 1)); } /** - * We only want to 'see' secondaries indexed in fIndices. In this + * We only want to 'see' secondaries indexed in indices_. In this * function the conversion form iterator-index to stack-index is * performed. */ unsigned int GetIndexFromIterator(const unsigned int vI) const { - if (vI == 0) return fProjectileIndex; - return fIndices[vI - 1]; + if (vI == 0) return projectile_index_; + return indices_[vI - 1]; } }; diff --git a/Framework/StackInterface/StackIteratorInterface.h b/Framework/StackInterface/StackIteratorInterface.h index 5b1d25757256fbcdafeb616865894d5a41b00750..e02ab055c97819a1a9be60a8c88e4cd23e636eea 100644 --- a/Framework/StackInterface/StackIteratorInterface.h +++ b/Framework/StackInterface/StackIteratorInterface.h @@ -23,6 +23,10 @@ namespace corsika::stack { template <typename TStackData, template <typename> typename TParticleInterface> class SecondaryView; // forward decl + template <typename TStackData, template <typename> typename TParticleInterface, + typename StackType> + class ConstStackIteratorInterface; // forward decl + /** @class StackIteratorInterface @@ -81,11 +85,15 @@ namespace corsika::stack { friend class SecondaryView<TStackData, TParticleInterface>; // access for SecondaryView - template <typename T> - friend class corsika::history::HistorySecondaryView; + template <typename T> + friend class corsika::history::HistorySecondaryView; - private: + friend class ConstStackIteratorInterface<TStackData, TParticleInterface, StackType>; + + protected: unsigned int index_ = 0; + + private: StackType* data_ = 0; // info: Particles and StackIterators become invalid when parent // Stack is copied or deleted! @@ -164,11 +172,21 @@ namespace corsika::stack { GetStack().isDeleted(*this)); // this also check the allowed bounds of index_ return tmp; } - StackIteratorInterface operator+(int delta) { + StackIteratorInterface operator+(int delta) const { return StackIteratorInterface(*data_, index_ + delta); } - bool operator==(const StackIteratorInterface& rhs) { return index_ == rhs.index_; } - bool operator!=(const StackIteratorInterface& rhs) { return index_ != rhs.index_; } + bool operator==(const StackIteratorInterface& rhs) const { + return index_ == rhs.index_; + } + bool operator!=(const StackIteratorInterface& rhs) const { + return index_ != rhs.index_; + } + bool operator==( + const ConstStackIteratorInterface<TStackData, TParticleInterface, StackType>& rhs) + const; // implement below + bool operator!=( + const ConstStackIteratorInterface<TStackData, TParticleInterface, StackType>& rhs) + const; // implement below /** * Convert iterator to value type, where value type is the user-provided particle @@ -234,8 +252,12 @@ namespace corsika::stack { friend class SecondaryView<TStackData, TParticleInterface>; // access for SecondaryView - private: + friend class StackIteratorInterface<TStackData, TParticleInterface, StackType>; + + protected: unsigned int index_ = 0; + + private: const StackType* data_ = 0; // info: Particles and StackIterators become invalid when // parent Stack is copied or deleted! @@ -281,14 +303,21 @@ namespace corsika::stack { GetStack().isDeleted(*this)); // this also check the allowed bounds of index_ return tmp; } - ConstStackIteratorInterface operator+(const int delta) { - for (int i = 0; i < delta; ++i) operator++(); - return ConstStackIteratorInterface(*data_, index_); + ConstStackIteratorInterface operator+(const int delta) const { + return ConstStackIteratorInterface(*data_, index_ + delta); + } + bool operator==(const ConstStackIteratorInterface& rhs) const { + return index_ == rhs.index_; + } + bool operator!=(const ConstStackIteratorInterface& rhs) const { + return index_ != rhs.index_; } - bool operator==(const ConstStackIteratorInterface& rhs) { + bool operator==(const StackIteratorInterface<TStackData, TParticleInterface, + StackType>& rhs) const { return index_ == rhs.index_; } - bool operator!=(const ConstStackIteratorInterface& rhs) { + bool operator!=(const StackIteratorInterface<TStackData, TParticleInterface, + StackType>& rhs) const { return index_ != rhs.index_; } @@ -312,4 +341,20 @@ namespace corsika::stack { ///@} }; // end class ConstStackIterator + template <typename TStackData, template <typename> typename TParticleInterface, + typename StackType> + bool StackIteratorInterface<TStackData, TParticleInterface, StackType>::operator==( + const ConstStackIteratorInterface<TStackData, TParticleInterface, StackType>& rhs) + const { + return index_ == rhs.index_; + } + + template <typename TStackData, template <typename> typename TParticleInterface, + typename StackType> + bool StackIteratorInterface<TStackData, TParticleInterface, StackType>::operator!=( + const ConstStackIteratorInterface<TStackData, TParticleInterface, StackType>& rhs) + const { + return index_ != rhs.index_; + } + } // namespace corsika::stack diff --git a/Framework/StackInterface/testSecondaryView.cc b/Framework/StackInterface/testSecondaryView.cc index fb27467dd260630538680bbcc610ce8bce673cb0..e4c82d15bddf3eb0d9292559f9bb4ed7d5acb7be 100644 --- a/Framework/StackInterface/testSecondaryView.cc +++ b/Framework/StackInterface/testSecondaryView.cc @@ -125,6 +125,7 @@ TEST_CASE("SecondaryStack", "[stack]") { { auto proj = view.GetProjectile(); CHECK(proj.GetData() == particle.GetData()); + CHECK(particle == view.parent()); } } diff --git a/Processes/CONEXSourceCut/CONEXSourceCut.cc b/Processes/CONEXSourceCut/CONEXSourceCut.cc index 7a0da45a726bbf40e69aab6264baa288968d47ec..8039d8ed867053756f631a950815aaf50d9c380f 100644 --- a/Processes/CONEXSourceCut/CONEXSourceCut.cc +++ b/Processes/CONEXSourceCut/CONEXSourceCut.cc @@ -31,11 +31,11 @@ corsika::process::EProcessReturn CONEXSourceCut::DoSecondaries( [=](auto const& p) { return pid == p.first; }); if (it != egs_em_codes_.cend()) { // EM particle - + auto const egs_pid = it->second; addParticle(egs_pid, p.GetEnergy(), p.GetMass(), p.GetPosition(), - p.GetMomentum().normalized(), p.GetTime()); + p.GetMomentum().normalized(), p.GetTime()); p.Delete(); } diff --git a/Setup/SetupStack.h b/Setup/SetupStack.h index 88e9f77114bed928e4174c02d880297063616495..b984b70e2d7244716a970c15a05b36c639574671 100644 --- a/Setup/SetupStack.h +++ b/Setup/SetupStack.h @@ -21,7 +21,7 @@ namespace corsika::setup { // ------------------------------------------ // add geometry node tracking data to stack: - + // the GeometryNode stack needs to know the type of geometry-nodes from the // environment: template <typename TStackIter> @@ -40,26 +40,24 @@ namespace corsika::setup { corsika::stack::node::GeometryData<setup::SetupEnvironment>, StackWithGeometryInterface>; - // ------------------------------------------ // Add [optional] history data to stack, too: // combine dummy stack with geometry information for tracking template <typename TStackIter> using StackWithHistoryInterface = corsika::stack::CombinedParticleInterface< - StackWithGeometry::PIType, history::HistoryEventDataInterface, - TStackIter>; + StackWithGeometry::PIType, history::HistoryEventDataInterface, TStackIter>; + + using StackWithHistory = + corsika::stack::CombinedStack<typename StackWithGeometry::StackImpl, + history::HistoryEventData, + StackWithHistoryInterface>; - using StackWithHistory = corsika::stack::CombinedStack< - typename StackWithGeometry::StackImpl, - history::HistoryEventData, StackWithHistoryInterface>; - } // namespace detail - // --------------------------------------- // this is the FINAL stack we use in C8: - //using Stack = detail::StackWithGeometry; + // using Stack = detail::StackWithGeometry; using Stack = detail::StackWithHistory; /*