diff --git a/corsika/detail/modules/epos/InteractionModel.inl b/corsika/detail/modules/epos/InteractionModel.inl index 35ed15a884c6e608056fcb5b562ddfd1d0c5e50e..f9e8925b7b2da6d94744a6d2afcad1c4fefdfa6a 100644 --- a/corsika/detail/modules/epos/InteractionModel.inl +++ b/corsika/detail/modules/epos/InteractionModel.inl @@ -430,11 +430,6 @@ namespace corsika::epos { MomentumVector P_final(originalCS, {0.0_GeV, 0.0_GeV, 0.0_GeV}); HEPEnergyType E_final = 0_GeV; - // position and time of interaction, not used in QgsjetII - auto const& projectile = view.getProjectile(); - Point const& pOrig = projectile.getPosition(); - TimeType const tOrig = projectile.getTime(); - // secondaries EposStack es; CORSIKA_LOGGER_DEBUG(logger_, "number of entries on Epos stack: {}", es.getSize()); @@ -449,12 +444,14 @@ namespace corsika::epos { EposCode const eposId = psec.getPID(); Code const pid = epos::convertFromEpos(eposId); + HEPEnergyType const mass = get_mass(pid); + HEPEnergyType const Ekin = sqrt(p3output.getSquaredNorm() + mass * mass) - mass; CORSIKA_LOGGER_TRACE(logger_, " id= {}" " p= {}", pid, p3output.getComponents() / 1_GeV); - auto pnew = view.addSecondary(std::make_tuple(pid, p3output, pOrig, tOrig)); + auto pnew = view.addSecondary(std::make_tuple(pid, Ekin, p3output.normalized())); P_final += pnew.getMomentum(); E_final += pnew.getEnergy(); } diff --git a/corsika/detail/modules/pythia8/Decay.inl b/corsika/detail/modules/pythia8/Decay.inl index 81306df4509492448ba4425b1fd8e874f307ad41..3da0ebdbf3869072b97b40198cb53a10f0ec16b8 100644 --- a/corsika/detail/modules/pythia8/Decay.inl +++ b/corsika/detail/modules/pythia8/Decay.inl @@ -170,9 +170,6 @@ namespace corsika::pythia8 { auto projectile = view.getProjectile(); - auto const& decayPoint = projectile.getPosition(); - auto const t0 = projectile.getTime(); - auto const& labMomentum = projectile.getMomentum(); [[maybe_unused]] CoordinateSystemPtr const& labCS = labMomentum.getCoordinateSystem(); @@ -230,14 +227,17 @@ namespace corsika::pythia8 { {event[i].px() * 1_GeV, event[i].py() * 1_GeV, event[i].pz() * 1_GeV}); FourVector const fourMomRest{Erest, pRest}; auto const fourMomLab = boost.fromCoM(fourMomRest); + auto const p3 = fourMomLab.getSpaceLikeComponents(); + + HEPEnergyType const mass = get_mass(pyId); + HEPEnergyType const Ekin = sqrt(p3.getSquaredNorm() + mass * mass) - mass; CORSIKA_LOG_TRACE( "particle: id={} momentum={} energy={} ", pyId, fourMomLab.getSpaceLikeComponents().getComponents(labCS) / 1_GeV, fourMomLab.getTimeLikeComponent()); - view.addSecondary( - std::make_tuple(pyId, fourMomLab.getSpaceLikeComponents(), decayPoint, t0)); + view.addSecondary(std::make_tuple(pyId, Ekin, p3.normalized())); } // set particle stable diff --git a/corsika/detail/modules/pythia8/Interaction.inl b/corsika/detail/modules/pythia8/Interaction.inl index fe55d587412537f86dff4ca4f0ca0f9dd5ef99ab..2144380e7480f3325375a45052c0822880c9e036 100644 --- a/corsika/detail/modules/pythia8/Interaction.inl +++ b/corsika/detail/modules/pythia8/Interaction.inl @@ -247,8 +247,12 @@ namespace corsika::pythia8 { MomentumVector const pyPlab(labCS, {p8p.px() * 1_GeV, p8p.py() * 1_GeV, p8p.pz() * 1_GeV}); + HEPEnergyType const mass = get_mass(pyId); + HEPEnergyType const Ekin = sqrt(pyPlab.getSquaredNorm() + mass * mass) - mass; + // add to corsika stack - auto pnew = projectile.addSecondary(std::make_tuple(pyId, pyPlab, pOrig, tOrig)); + auto pnew = + projectile.addSecondary(std::make_tuple(pyId, Ekin, pyPlab.normalized())); Plab_final += pnew.getMomentum(); Elab_final += pnew.getEnergy(); diff --git a/corsika/detail/modules/qgsjetII/InteractionModel.inl b/corsika/detail/modules/qgsjetII/InteractionModel.inl index 8203faa49d1c9c6c3fe16de453f01dcde138e3bb..81b62894b26508a00d83ae8a366629c3c7cc08aa 100644 --- a/corsika/detail/modules/qgsjetII/InteractionModel.inl +++ b/corsika/detail/modules/qgsjetII/InteractionModel.inl @@ -178,11 +178,6 @@ namespace corsika::qgsjetII { // internal QGSJetII system COMBoost boostInternal({Elab, pLab}, get_mass(targetId)); - // position and time of interaction, not used in QgsjetII - auto const projectile = view.getProjectile(); - Point const pOrig = projectile.getPosition(); - TimeType const tOrig = projectile.getTime(); - // fragments QGSJetIIFragmentsStack qfs; for (auto& fragm : qfs) { @@ -208,11 +203,16 @@ namespace corsika::qgsjetII { auto p3output = P4output.getSpaceLikeComponents(); p3output.rebase(originalCS); // transform back into standard lab frame + HEPEnergyType const mass = get_mass(idFragm); + HEPEnergyType const Ekin = sqrt(p3output.getSquaredNorm() + mass * mass) - mass; + CORSIKA_LOG_DEBUG( "secondary fragment> id= {}" " p={}", idFragm, p3output.getComponents()); - auto pnew = view.addSecondary(std::make_tuple(idFragm, p3output, pOrig, tOrig)); + + auto pnew = + view.addSecondary(std::make_tuple(idFragm, Ekin, p3output.normalized())); Plab_final += pnew.getMomentum(); Elab_final += pnew.getEnergy(); @@ -252,15 +252,19 @@ namespace corsika::qgsjetII { auto p3output = P4output.getSpaceLikeComponents(); p3output.rebase(originalCS); // transform back into standard lab frame + Code const idFragm = get_nucleus_code(A, Z); + HEPEnergyType const mass = get_mass(idFragm); + HEPEnergyType const Ekin = sqrt(p3output.getSquaredNorm() + mass * mass) - mass; + CORSIKA_LOG_DEBUG( "secondary fragment> id={}" " p={}" " A={}" " Z={}", - get_nucleus_code(A, Z), p3output.getComponents(), A, Z); + idFragm, p3output.getComponents(), A, Z); - auto pnew = view.addSecondary( - std::make_tuple(get_nucleus_code(A, Z), p3output, pOrig, tOrig)); + auto pnew = + view.addSecondary(std::make_tuple(idFragm, Ekin, p3output.normalized())); Plab_final += pnew.getMomentum(); Elab_final += pnew.getEnergy(); } @@ -278,11 +282,12 @@ namespace corsika::qgsjetII { auto p3output = P4output.getSpaceLikeComponents(); p3output.rebase(originalCS); // transform back into standard lab frame - CORSIKA_LOG_DEBUG("secondary> id= {}, p= {}", - corsika::qgsjetII::convertFromQgsjetII(psec.getPID()), - p3output.getComponents()); - auto pnew = view.addSecondary(std::make_tuple( - corsika::qgsjetII::convertFromQgsjetII(psec.getPID()), p3output, pOrig, tOrig)); + Code const pid = corsika::qgsjetII::convertFromQgsjetII(psec.getPID()); + HEPEnergyType const mass = get_mass(pid); + HEPEnergyType const Ekin = sqrt(p3output.getSquaredNorm() + mass * mass) - mass; + + CORSIKA_LOG_DEBUG("secondary> id= {}, p= {}", pid, p3output.getComponents()); + auto pnew = view.addSecondary(std::make_tuple(pid, Ekin, p3output.normalized())); Plab_final += pnew.getMomentum(); Elab_final += pnew.getEnergy(); } diff --git a/corsika/detail/modules/sibyll/Decay.inl b/corsika/detail/modules/sibyll/Decay.inl index ab0114182b41f5db321f19ec6d91458cc9b5000d..b804feb534aedea966eddd1c5f378c7c4ba0329d 100644 --- a/corsika/detail/modules/sibyll/Decay.inl +++ b/corsika/detail/modules/sibyll/Decay.inl @@ -173,9 +173,7 @@ namespace corsika::sibyll { throw std::runtime_error("STOP! Sibyll not configured to execute this decay!"); count_++; - // remember position - Point const& decayPoint = projectile.getPosition(); - TimeType const t0 = projectile.getTime(); + // switch on decay for this particle setUnstable(pCode); printDecayConfig(pCode); @@ -222,8 +220,11 @@ namespace corsika::sibyll { CORSIKA_LOG_TRACE("Sibyll::Decay: i={} id={} p={} GeV", i, pid, components / 1_GeV); - projectile.addSecondary( - std::make_tuple(pid, MomentumVector(rootCS, components), decayPoint, t0)); + auto const p3 = MomentumVector(rootCS, components); + HEPEnergyType const mass = get_mass(pid); + HEPEnergyType const Ekin = sqrt(p3.getSquaredNorm() + mass * mass) - mass; + + projectile.addSecondary(std::make_tuple(pid, Ekin, p3.normalized())); } } diff --git a/corsika/detail/modules/sibyll/InteractionModel.inl b/corsika/detail/modules/sibyll/InteractionModel.inl index 99f078673c2d80df3c6251bd6d22eed25a077c4f..dd00e7606aabf439c677bf2dbaf74b34e1ddab82 100644 --- a/corsika/detail/modules/sibyll/InteractionModel.inl +++ b/corsika/detail/modules/sibyll/InteractionModel.inl @@ -138,11 +138,6 @@ namespace corsika::sibyll { // add particles from sibyll to stack - // position and time of interaction, not used in Sibyll - auto const& projectile = secondaries.parent(); - Point const& pOrig = projectile.getPosition(); - TimeType const tOrig = projectile.getTime(); // no time in sibyll - // link to sibyll stack SibStack ss; @@ -163,9 +158,13 @@ namespace corsika::sibyll { auto const P4lab = boost.fromCoM(FourVector{eCoM, pCoM}); auto const p3lab = P4lab.getSpaceLikeComponents(); + Code const pid = corsika::sibyll::convertFromSibyll(psib.getPID()); + HEPEnergyType const mass = get_mass(pid); + HEPEnergyType const Ekin = sqrt(p3lab.getSquaredNorm() + mass * mass) - mass; + // add to corsika stack - auto pnew = secondaries.addSecondary(std::make_tuple( - corsika::sibyll::convertFromSibyll(psib.getPID()), p3lab, pOrig, tOrig)); + auto pnew = + secondaries.addSecondary(std::make_tuple(pid, Ekin, p3lab.normalized())); Plab_final += pnew.getMomentum(); Elab_final += pnew.getEnergy(); diff --git a/corsika/detail/modules/sibyll/NuclearInteractionModel.inl b/corsika/detail/modules/sibyll/NuclearInteractionModel.inl index 4a564c7e96bcc64e059dfec58b43eb3f5b8ed7ac..06de518c44dd1880b170f9434566aebced499a4f 100644 --- a/corsika/detail/modules/sibyll/NuclearInteractionModel.inl +++ b/corsika/detail/modules/sibyll/NuclearInteractionModel.inl @@ -292,14 +292,6 @@ namespace corsika::sibyll { } // (LCOV_EXCL_STOP) - // position and time of interaction, not used in NUCLIB - auto const& projectile = view.parent(); - // position and time of interaction, not used in NUCLI - Point const& pOrig = projectile.getPosition(); - TimeType const delay = projectile.getTime(); - - CORSIKA_LOG_DEBUG("Interaction: position of interaction: {}, {} ns", - pOrig.getCoordinates(), delay / 1_ns); CORSIKA_LOG_DEBUG("number of fragments: {}", nFragments); CORSIKA_LOG_DEBUG("adding nuclear fragments to particle stack.."); // put nuclear fragments on corsika stack @@ -324,8 +316,11 @@ namespace corsika::sibyll { // CORSIKA 7 way // spectators inherit momentum from original projectile auto const p3lab = p3NucleonLab * nuclA; + + HEPEnergyType const Ekin = sqrt(p3lab.getSquaredNorm() + mass * mass) - mass; + CORSIKA_LOG_DEBUG("fragment momentum {}", p3lab.getComponents() / 1_GeV); - view.addSecondary(std::make_tuple(specCode, p3lab, pOrig, delay)); + view.addSecondary(std::make_tuple(specCode, Ekin, p3lab.normalized())); } // add elastic nucleons to corsika stack @@ -340,7 +335,11 @@ namespace corsika::sibyll { // elastic nucleons inherit momentum from original projectile // neglecting momentum transfer in interaction auto const p3lab = p3NucleonLab; - view.addSecondary(std::make_tuple(elaNucCode, p3lab, pOrig, delay)); + + HEPEnergyType const mass = get_mass(elaNucCode); + HEPEnergyType const Ekin = sqrt(p3lab.getSquaredNorm() + mass * mass) - mass; + + view.addSecondary(std::make_tuple(elaNucCode, Ekin, p3lab.normalized())); } // add inelastic interactions @@ -348,12 +347,18 @@ namespace corsika::sibyll { for (int j = 0; j < nInelNucleons; ++j) { // TODO: sample neutron or proton auto const pCode = Code::Proton; + HEPEnergyType const mass = get_mass(pCode); + HEPEnergyType const Ekin = sqrt(p3NucleonLab.getSquaredNorm() + mass * mass) - mass; + // temporarily add to stack, will be removed after interaction in DoInteraction CORSIKA_LOG_DEBUG("inelastic interaction no. {}", j); typename TSecondaryView::inner_stack_value_type nucleonStack; - auto inelasticNucleon = - nucleonStack.addParticle(std::make_tuple(pCode, p3NucleonLab, pOrig, delay)); + Point const pDummy(boost.getOriginalCS(), {0_m, 0_m, 0_m}); + TimeType const tDummy = 0_ns; + auto inelasticNucleon = nucleonStack.addParticle( + std::make_tuple(pCode, Ekin, p3NucleonLab.normalized(), pDummy, tDummy)); inelasticNucleon.setNode(view.getProjectile().getNode()); + // create inelastic interaction for each nucleon CORSIKA_LOG_TRACE("calling HadronicInteraction..."); // create new StackView for each of the nucleons @@ -362,8 +367,13 @@ namespace corsika::sibyll { hadronicInteraction_.doInteraction(nucleon_secondaries, pCode, targetId, nucleonP4, targetP4); for (const auto& pSec : nucleon_secondaries) { - view.addSecondary(std::make_tuple(pSec.getPID(), pSec.getMomentum(), - pSec.getPosition(), pSec.getTime())); + + auto const p3lab = pSec.getMomentum(); + Code const pid = pSec.getPID(); + HEPEnergyType const mass = get_mass(pid); + HEPEnergyType const Ekin = sqrt(p3lab.getSquaredNorm() + mass * mass) - mass; + + view.addSecondary(std::make_tuple(pid, Ekin, p3lab.normalized())); } } } diff --git a/corsika/detail/modules/urqmd/UrQMD.inl b/corsika/detail/modules/urqmd/UrQMD.inl index 1bf0f0d59a42428df231717f2df869f72c30bd45..9e070f1c62df8860f9bbb4d89df26f9c87728e09 100644 --- a/corsika/detail/modules/urqmd/UrQMD.inl +++ b/corsika/detail/modules/urqmd/UrQMD.inl @@ -259,10 +259,6 @@ namespace corsika::urqmd { iflb_; // flag for retrying interaction in case of empty event, 0 means retry ::urqmd::urqmd_(iflb); - auto projectile = view.getProjectile(); - auto const& projectilePosition = projectile.getPosition(); - auto const projectileTime = projectile.getTime(); - // now retrieve secondaries from UrQMD COMBoost const boost(projectileP4, targetP4); auto const& originalCS = boost.getOriginalCS(); @@ -283,8 +279,15 @@ namespace corsika::urqmd { momentum.rebase(originalCS); // transform back into standard lab frame CORSIKA_LOG_DEBUG(" {} {} {} ", i, code, momentum.getComponents()); - view.addSecondary( - std::make_tuple(code, momentum, projectilePosition, projectileTime)); + HEPEnergyType const mass = get_mass(code); + HEPEnergyType Ekin = sqrt(momentum.getSquaredNorm() + mass * mass) - mass; + if (Ekin <= 0_GeV) { + CORSIKA_LOG_WARN("Negative kinetic energy {} {}. Skipping.", code, Ekin); + view.addSecondary( + std::make_tuple(code, 0_eV, DirectionVector{originalCS, {0, 0, 0}})); + } else { + view.addSecondary(std::make_tuple(code, Ekin, momentum.normalized())); + } } CORSIKA_LOG_DEBUG("UrQMD generated {} secondaries!", ::urqmd::sys_.npart); } diff --git a/corsika/detail/stack/VectorStack.inl b/corsika/detail/stack/VectorStack.inl index 57bc094a438dc6f7383baa39215434099777690b..fc62767e1779aabbc4d0caf8f8b4ae6201aa302d 100644 --- a/corsika/detail/stack/VectorStack.inl +++ b/corsika/detail/stack/VectorStack.inl @@ -22,41 +22,6 @@ namespace corsika { - template <typename StackIteratorInterface> - inline void ParticleInterface<StackIteratorInterface>::setParticleData( - particle_data_momentum_type const& v) { - this->setPID(std::get<0>(v)); - MomentumVector const p = std::get<1>(v); - auto const P2 = p.getSquaredNorm(); - HEPMassType const M = getMass(); - this->setKineticEnergy(sqrt(P2 + square(M)) - M); - if (P2 == static_pow<2>(0_eV)) { - this->setDirection(DirectionVector(p.getCoordinateSystem(), {0, 0, 0})); - } else { - this->setDirection(p / sqrt(P2)); - } - this->setPosition(std::get<2>(v)); - this->setTime(std::get<3>(v)); - } - - template <typename StackIteratorInterface> - inline void ParticleInterface<StackIteratorInterface>::setParticleData( - ParticleInterface<StackIteratorInterface> const& parent, - particle_data_momentum_type const& v) { - this->setPID(std::get<0>(v)); - MomentumVector const p = std::get<1>(v); - auto const P2 = p.getSquaredNorm(); - HEPMassType const M = getMass(); - this->setKineticEnergy(sqrt(P2 + square(M)) - M); - if (P2 == static_pow<2>(0_eV)) { - this->setDirection(DirectionVector(p.getCoordinateSystem(), {0, 0, 0})); - } else { - this->setDirection(p / sqrt(P2)); - } - this->setPosition(std::get<2>(v)); - this->setTime(std::get<3>(v) + parent.getTime()); // parent time is added - } - template <typename StackIteratorInterface> inline void ParticleInterface<StackIteratorInterface>::setParticleData( particle_data_type const& v) { @@ -70,12 +35,12 @@ namespace corsika { template <typename StackIteratorInterface> inline void ParticleInterface<StackIteratorInterface>::setParticleData( ParticleInterface<StackIteratorInterface> const& parent, - particle_data_type const& v) { + secondary_data_type const& v) { this->setPID(std::get<0>(v)); this->setKineticEnergy(std::get<1>(v)); this->setDirection(std::get<2>(v)); - this->setPosition(std::get<3>(v)); - this->setTime(std::get<4>(v) + parent.getTime()); // parent time is added + this->setPosition(parent.getPosition()); // position + this->setTime(parent.getTime()); // parent time } template <typename StackIteratorInterface> diff --git a/corsika/stack/VectorStack.hpp b/corsika/stack/VectorStack.hpp index 22873790c22c1e478a0c6883008ed95e54a0562b..029f75490c56acbd4e7c7a25bb6e3aadb4674015 100644 --- a/corsika/stack/VectorStack.hpp +++ b/corsika/stack/VectorStack.hpp @@ -33,50 +33,41 @@ namespace corsika { typedef ParticleBase<TStackIterator> super_type; public: + /** + * particle data information content. + * + * PID, Ekin, direction, position, time. + */ typedef std::tuple<Code, HEPEnergyType, DirectionVector, Point, TimeType> particle_data_type; - typedef std::tuple<Code, MomentumVector, Point, TimeType> particle_data_momentum_type; - - std::string asString() const; - /** - * Set data of new particle. - * - * @param v tuple containing: PID, Momentum Vector, Position, Time + * secondary particle data information content. * - * MomentumVector is only used to determine the DirectionVector, the normalization - * is lost. + * PID, Ekin, direction. */ - void setParticleData(particle_data_type const& v); + typedef std::tuple<Code, HEPEnergyType, DirectionVector> secondary_data_type; + + std::string asString() const; /** * Set data of new particle. * - * @param parent parent particle - * @param v tuple containing: PID, Momentum Vector, Position, Time + * @param v tuple containing of type particle_data_type * - * MomentumVector is only used to determine the DirectionVector, the normalization + * MomentumVector is only used to determine the DirectionVector, the normalization * is lost. */ - void setParticleData(ParticleInterface<TStackIterator> const& parent, - particle_data_type const& v); - - /** - * Set data of new particle. - * - * @param v tuple containing: PID, kinetic Energy, Direction Vector, Position, Time - */ - void setParticleData(particle_data_momentum_type const& v); + void setParticleData(particle_data_type const& v); /** * Set data of new particle. * * @param parent parent particle - * @param v tuple containing: PID, kinetic Energy, Direction Vector, Position, Time + * @param v tuple containing of type secondary_data_type. */ void setParticleData(ParticleInterface<TStackIterator> const& parent, - particle_data_momentum_type const& v); + secondary_data_type const& v); ///! Set particle corsika::Code void setPID(Code const id) { diff --git a/tests/common/SetupTestStack.hpp b/tests/common/SetupTestStack.hpp index a6141e582262b77e4062eb020bfdc7e3a19adece..476c99a54bcd2a5487b9d4cdee591e88607a3409 100644 --- a/tests/common/SetupTestStack.hpp +++ b/tests/common/SetupTestStack.hpp @@ -27,9 +27,6 @@ namespace corsika::setup::testing { * * standard stack setup for unit tests. * - * - * - * * \return a tuple with element 0 being a Stack object filled with * one particle, and element 1 the StackView on it. */ @@ -43,9 +40,11 @@ namespace corsika::setup::testing { Point const origin(cs, {0_m, 0_m, 0_m}); MomentumVector const pLab(cs, {vMomentum, 0_GeV, 0_GeV}); + HEPMassType const mass = get_mass(vProjectileType); + HEPEnergyType const Ekin = sqrt(vMomentum * vMomentum + mass * mass) - mass; - auto particle = - stack->addParticle(std::make_tuple(vProjectileType, pLab, origin, 0_ns)); + auto particle = stack->addParticle( + std::make_tuple(vProjectileType, Ekin, pLab.normalized(), origin, 0_ns)); particle.setNode(vNodePtr); return std::make_tuple(std::move(stack), std::make_unique<setup::StackView>(particle)); diff --git a/tests/framework/testCascade.cpp b/tests/framework/testCascade.cpp index 18732cb1f5faaa646fd9333a89faa1e9e2861973..624dc041575d1304c4e7bdcb012810e7ca3cc478 100644 --- a/tests/framework/testCascade.cpp +++ b/tests/framework/testCascade.cpp @@ -98,10 +98,10 @@ public: ++calls_; auto vP = view.getProjectile(); const HEPEnergyType Ekin = vP.getKineticEnergy(); - vP.addSecondary(std::make_tuple(vP.getPID(), Ekin / 2, vP.getMomentum().normalized(), - vP.getPosition(), vP.getTime())); - vP.addSecondary(std::make_tuple(vP.getPID(), Ekin / 2, vP.getMomentum().normalized(), - vP.getPosition(), vP.getTime())); + vP.addSecondary( + std::make_tuple(vP.getPID(), Ekin / 2, vP.getMomentum().normalized())); + vP.addSecondary( + std::make_tuple(vP.getPID(), Ekin / 2, vP.getMomentum().normalized())); } int getCalls() const { return calls_; } @@ -162,11 +162,10 @@ TEST_CASE("Cascade", "[Cascade]") { auto sequence = make_sequence(nullModel, stackInspect, split, cut); TestCascadeStack stack; stack.clear(); - stack.addParticle(std::make_tuple( - Code::Electron, - MomentumVector(rootCS, {0_GeV, 0_GeV, - -sqrt(E0 * E0 - static_pow<2>(get_mass(Code::Electron)))}), - Point(rootCS, {0_m, 0_m, 10_km}), 0_ns)); + stack.addParticle(std::make_tuple(Code::Electron, + E0 - get_mass(Code::Electron), // Ekin + DirectionVector(rootCS, {0, 0, -1}), + Point(rootCS, {0_m, 0_m, 10_km}), 0_ns)); DummyTracking tracking; DummyOutputManager output; diff --git a/tests/modules/testParticleCut.cpp b/tests/modules/testParticleCut.cpp index b293d5eb04a98c2dca8f375f094c22dbcb711ec2..fc546ebf332a80ca3417c3873b38fc9891a067c2 100644 --- a/tests/modules/testParticleCut.cpp +++ b/tests/modules/testParticleCut.cpp @@ -67,8 +67,8 @@ TEST_CASE("ParticleCut", "process,continuous,secondary") { // add secondaries, all with energies above the threshold // only cut is by species for (auto proType : particleList) - projectile.addSecondary(std::make_tuple( - proType, Eabove, DirectionVector(rootCS, {1, 0, 0}), point0, 0_ns)); + projectile.addSecondary( + std::make_tuple(proType, Eabove, DirectionVector(rootCS, {1, 0, 0}))); CHECK(view.getEntries() == 11); CHECK(stack.getEntries() == 12); @@ -94,8 +94,8 @@ TEST_CASE("ParticleCut", "process,continuous,secondary") { // add secondaries, all with energies above the threshold // only cut is by species for (auto proType : particleList) { - projectile.addSecondary(std::make_tuple( - proType, Eabove, DirectionVector(rootCS, {1, 0, 0}), point0, 0_ns)); + projectile.addSecondary( + std::make_tuple(proType, Eabove, DirectionVector(rootCS, {1, 0, 0}))); } cut.doSecondaries(view); @@ -118,16 +118,14 @@ TEST_CASE("ParticleCut", "process,continuous,secondary") { // add secondaries, all with energies below the threshold // only cut is by species for (auto proType : particleList) - projectile.addSecondary(std::make_tuple( - proType, Ebelow, DirectionVector(rootCS, {1, 0, 0}), point0, 0_ns)); + projectile.addSecondary( + std::make_tuple(proType, Ebelow, DirectionVector(rootCS, {1, 0, 0}))); unsigned short A = 18; unsigned short Z = 8; projectile.addSecondary(std::make_tuple(get_nucleus_code(A, Z), Eabove * A, - DirectionVector(rootCS, {1, 0, 0}), point0, - 0_ns)); + DirectionVector(rootCS, {1, 0, 0}))); projectile.addSecondary(std::make_tuple(get_nucleus_code(A, Z), Ebelow * A, - DirectionVector(rootCS, {1, 0, 0}), point0, - 0_ns)); + DirectionVector(rootCS, {1, 0, 0}))); cut.doSecondaries(view); @@ -148,21 +146,19 @@ TEST_CASE("ParticleCut", "process,continuous,secondary") { // only this way the secondary view is populated auto projectile = view.getProjectile(); // add secondaries - projectile.addSecondary(std::make_tuple( - Code::Photon, 3_MeV, DirectionVector(rootCS, {1, 0, 0}), point0, 0_ns)); - projectile.addSecondary(std::make_tuple( - Code::Electron, 3_MeV, DirectionVector(rootCS, {1, 0, 0}), point0, 0_ns)); - projectile.addSecondary(std::make_tuple( - Code::PiPlus, 4_GeV, DirectionVector(rootCS, {1, 0, 0}), point0, 0_ns)); + projectile.addSecondary( + std::make_tuple(Code::Photon, 3_MeV, DirectionVector(rootCS, {1, 0, 0}))); + projectile.addSecondary( + std::make_tuple(Code::Electron, 3_MeV, DirectionVector(rootCS, {1, 0, 0}))); + projectile.addSecondary( + std::make_tuple(Code::PiPlus, 4_GeV, DirectionVector(rootCS, {1, 0, 0}))); unsigned short A = 18; unsigned short Z = 8; projectile.addSecondary(std::make_tuple(get_nucleus_code(A, Z), 4_GeV * A, - DirectionVector(rootCS, {1, 0, 0}), point0, - 0_ns)); + DirectionVector(rootCS, {1, 0, 0}))); projectile.addSecondary(std::make_tuple(get_nucleus_code(A, Z), 6_GeV * A, - DirectionVector(rootCS, {1, 0, 0}), point0, - 0_ns)); + DirectionVector(rootCS, {1, 0, 0}))); cut.doSecondaries(view); @@ -185,7 +181,7 @@ TEST_CASE("ParticleCut", "process,continuous,secondary") { // add primary particle to stack auto particle = stack.addParticle(std::make_tuple( - Code::Proton, Eabove, DirectionVector(rootCS, {1, 0, 0}), point0, 1_ns)); + Code::Proton, Eabove, DirectionVector(rootCS, {1, 0, 0}), point0, too_late)); // view on secondary particles setup::StackView view(particle); // ref. to primary particle through the secondary view. @@ -194,8 +190,8 @@ TEST_CASE("ParticleCut", "process,continuous,secondary") { // add secondaries, all with energies above the threshold // only cut is by time for (auto proType : particleList) { - projectile.addSecondary(std::make_tuple( - proType, Eabove, DirectionVector(rootCS, {1, 0, 0}), point0, too_late)); + projectile.addSecondary( + std::make_tuple(proType, Eabove, DirectionVector(rootCS, {1, 0, 0}))); } cut.doSecondaries(view); diff --git a/tests/modules/testStackInspector.cpp b/tests/modules/testStackInspector.cpp index ec2157e3e79cff3a1680692b482bf6f5eac159f3..f76eab9a93f7177ee20cfc4bf6b531e153aba0e1 100644 --- a/tests/modules/testStackInspector.cpp +++ b/tests/modules/testStackInspector.cpp @@ -35,11 +35,11 @@ TEST_CASE("StackInspector", "modules") { TestCascadeStack stack; stack.clear(); HEPEnergyType E0 = 100_GeV; - stack.addParticle(std::make_tuple(Code::Electron, - MomentumVector(rootCS, {0_GeV, 0_GeV, -1_GeV}), + stack.addParticle(std::make_tuple(Code::Electron, 1_GeV, + DirectionVector(rootCS, {0, 0, -1}), Point(rootCS, {0_m, 0_m, 10_km}), 0_ns)); - stack.addParticle(std::make_tuple(get_nucleus_code(16, 8), - MomentumVector(rootCS, {0_GeV, 0_GeV, -1_GeV}), + stack.addParticle(std::make_tuple(get_nucleus_code(16, 8), 1_GeV, + DirectionVector(rootCS, {0, 0, -1}), Point(rootCS, {0_m, 0_m, 10_km}), 0_ns)); SECTION("interface") { diff --git a/tests/modules/testUrQMD.cpp b/tests/modules/testUrQMD.cpp index f75d59b75b44a41d25376d013c9d9e0d813efa32..10f5954330887a57a09c7b4095311c99849feb31 100644 --- a/tests/modules/testUrQMD.cpp +++ b/tests/modules/testUrQMD.cpp @@ -121,8 +121,6 @@ TEST_CASE("UrQMD") { auto [stackPtr, secViewPtr] = setup::testing::setup_stack( Code::PiPlus, 40_GeV, (setup::Environment::BaseNodeType* const)nodePtr, *csPtr); - CHECK(stackPtr->getEntries() == 1); - CHECK(secViewPtr->getEntries() == 0); // must be assigned to variable, cannot be used as rvalue?! auto projectile = secViewPtr->getProjectile(); diff --git a/tests/stack/testVectorStack.cpp b/tests/stack/testVectorStack.cpp index b44e3d55896e1c5d259642210f9b4dfa8de19a68..2adb1d49caeb2f9500141394af7387f83ae51133 100644 --- a/tests/stack/testVectorStack.cpp +++ b/tests/stack/testVectorStack.cpp @@ -40,12 +40,6 @@ TEST_CASE("VectorStack", "stack") { CHECK(pout.getTime() == 100_s); CHECK(pout.getChargeNumber() == -1); - // particle with no momentum, has no direction - auto const p0 = s.addParticle( - std::make_tuple(Code::Proton, MomentumVector(dummyCS, {0_GeV, 0_GeV, 0_GeV}), - Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s)); - CHECK(p0.getDirection().getNorm() == 0); - s.clear(); CHECK(s.getEntries() == 0); CHECK(s.getSize() == 0); @@ -57,7 +51,7 @@ TEST_CASE("VectorStack", "stack") { for (int i = 0; i < 99; ++i) s.addParticle( - std::make_tuple(Code::Electron, MomentumVector(dummyCS, {1_GeV, 1_GeV, 1_GeV}), + std::make_tuple(Code::Electron, 1_GeV, DirectionVector(dummyCS, {1, 0, 0}), Point(dummyCS, {1 * meter, 1 * meter, 1 * meter}), 100_s)); CHECK(s.getSize() == 99);