diff --git a/corsika/detail/media/NuclearComposition.inl b/corsika/detail/media/NuclearComposition.inl index 85582510e5cb0e115c7e97133ee215de1193b2b2..37d8da57d6cd0a21696c8d0c2748e0b916fd163c 100644 --- a/corsika/detail/media/NuclearComposition.inl +++ b/corsika/detail/media/NuclearComposition.inl @@ -97,10 +97,14 @@ namespace corsika { return components_[iChannel]; } - // Note: when this class ever modifies its internal data, the hash + // Note: when this class ever modifies its internal data, the hash // must be updated, too! inline size_t NuclearComposition::getHash() const { return hash_; } + inline bool NuclearComposition::operator==(NuclearComposition const& v) const { + return v.hash_ == hash_; + } + inline void NuclearComposition::updateHash() { std::vector<std::size_t> hashes; for (float ifrac : this->getFractions()) hashes.push_back(std::hash<float>{}(ifrac)); diff --git a/corsika/detail/modules/qgsjetII/Interaction.inl b/corsika/detail/modules/qgsjetII/Interaction.inl index c41b436cc64d60360984f66d59d0cc74ad21a305..828b9fcfeba8ae193ac622249a74ef317eaacc13 100644 --- a/corsika/detail/modules/qgsjetII/Interaction.inl +++ b/corsika/detail/modules/qgsjetII/Interaction.inl @@ -77,7 +77,7 @@ namespace corsika::qgsjetII { std::ostringstream txt; txt << "QgsjetII projectile outside range. Aprojectile=" << iProjectile; throw std::runtime_error(txt.str().c_str()); - } + } } CORSIKA_LOG_DEBUG( @@ -237,7 +237,7 @@ namespace corsika::qgsjetII { int targetMassNumber = 1; // proton if (is_nucleus(targetCode)) { // nucleus - targetMassNumber = get_nucleus_A(targetCode); + targetMassNumber = get_nucleus_A(targetCode); if (targetMassNumber > maxMassNumber_) throw std::runtime_error("QgsjetII target mass outside range."); } else { diff --git a/corsika/media/NuclearComposition.hpp b/corsika/media/NuclearComposition.hpp index 5341554cd3f132d7c6415112a6eca72058de025b..9e435adf60d2b6c1325ac65df56005e875142d37 100644 --- a/corsika/media/NuclearComposition.hpp +++ b/corsika/media/NuclearComposition.hpp @@ -53,9 +53,9 @@ namespace corsika { **/ size_t getSize() const; - /// Returns a const reference to the fraction + //! Returns a const reference to the fraction std::vector<float> const& getFractions() const; - /// Returns a const reference to the fraction + //! Returns a const reference to the fraction std::vector<Code> const& getComponents() const; double const getAverageMassNumber() const; @@ -63,10 +63,13 @@ namespace corsika { Code sampleTarget(std::vector<CrossSectionType> const& sigma, TRNG& randomStream) const; - // Note: when this class ever modifies its internal data, the hash + // Note: when this class ever modifies its internal data, the hash // must be updated, too! size_t getHash() const; + //! based on hash value + bool operator==(NuclearComposition const& v) const; + private: void updateHash(); diff --git a/tests/media/testEnvironment.cpp b/tests/media/testEnvironment.cpp index ca16d9f589d0e84932e004b4df440335be90fd9f..27f99a61f6179902e0e45606ace1d1f44939722d 100644 --- a/tests/media/testEnvironment.cpp +++ b/tests/media/testEnvironment.cpp @@ -46,6 +46,9 @@ TEST_CASE("HomogeneousMedium") { NuclearComposition const protonComposition(std::vector<Code>{Code::Proton}, std::vector<float>{1.f}); HomogeneousMedium<IMediumModel> const medium(19.2_g / cube(1_cm), protonComposition); + + CHECK_THROWS(NuclearComposition({Code::Proton}, {1.1})); + CHECK_THROWS(NuclearComposition({Code::Proton}, {0.99})); } TEST_CASE("FlatExponential") { @@ -216,6 +219,8 @@ TEST_CASE("InhomogeneousMedium") { inhMedium.getIntegratedGrammage(trajectory, l)); CHECK(rho.getArclengthFromGrammage(trajectory, 20_g / (1_cm * 1_cm)) == inhMedium.getArclengthFromGrammage(trajectory, 20_g / (1_cm * 1_cm))); + CHECK(inhMedium.getNuclearComposition() == composition); + CHECK(inhMedium.getMassDensity({gCS, {0_m, 0_m, 0_m}}) == 1_kg/static_pow<3>(1_m)); } } @@ -234,6 +239,8 @@ TEST_CASE("LayeredSphericalAtmosphereBuilder") { builder.addLinearLayer(2_km, 20_km); builder.addExponentialLayer(540.1778_g / (1_cm * 1_cm), 772170.16_cm, 30_km); + CHECK_THROWS(builder.addLinearLayer(0.5_km, 5_km)); + CHECK(builder.getSize() == 3); auto const builtEnv = builder.assemble(); @@ -300,7 +307,7 @@ TEST_CASE("LayeredSphericalAtmosphereBuilder w/ magnetic field") { ->getModelProperties() .getMagneticField(pTest2) .getComponents(gCS)); -} + } TEST_CASE("media", "LayeredSphericalAtmosphereBuilder USStd") { // setup environment, geometry diff --git a/tests/media/testShowerAxis.cpp b/tests/media/testShowerAxis.cpp index 8ba454c342eb3ba1f972c4c6129ab0876fb1bcf9..bd21f6805768604fd4dbbfc3d59ec3b48eb0631c 100644 --- a/tests/media/testShowerAxis.cpp +++ b/tests/media/testShowerAxis.cpp @@ -62,7 +62,7 @@ TEST_CASE("Homogeneous Density") { Point const injectionPos = showerCore + Vector<dimensionless_d>{cs, {0, 0, 1}} * t; ShowerAxis const showerAxis{injectionPos, (showerCore - injectionPos), *env, - false, // -> do not throw exceptions + true, // -> do not throw exceptions 20}; // -> number of bins CHECK(showerAxis.getSteplength() == 500_m); @@ -80,4 +80,7 @@ TEST_CASE("Homogeneous Density") { const Vector<dimensionless_d> dir{cs, {0, 0, -1}}; CHECK(showerAxis.getDirection().getComponents(cs) == dir.getComponents(cs)); CHECK(showerAxis.getStart().getCoordinates() == injectionPos.getCoordinates()); + + CHECK_THROWS(showerAxis.getX(-1_m)); + CHECK_THROWS(showerAxis.getX((injectionPos-showerCore).getNorm()+1_m)); } diff --git a/tests/modules/testQGSJetII.cpp b/tests/modules/testQGSJetII.cpp index 09a28979c07e7e0a1875e3090a2ed5edec15ee61..477dedc6ca15a123b26f5524cbeb888fb2700cd7 100644 --- a/tests/modules/testQGSJetII.cpp +++ b/tests/modules/testQGSJetII.cpp @@ -79,6 +79,8 @@ TEST_CASE("QgsjetII", "[processes]") { SECTION("QgsjetII -> Corsika") { CHECK(Code::PiPlus == corsika::qgsjetII::convertFromQgsjetII( corsika::qgsjetII::QgsjetIICode::PiPlus)); + CHECK_THROWS( + corsika::qgsjetII::convertFromQgsjetII(corsika::qgsjetII::QgsjetIICode::Unknown)); } SECTION("Corsika -> QgsjetII") { @@ -169,12 +171,12 @@ TEST_CASE("QgsjetIIInterface", "[processes]") { CHECK((secMomSum - projectileMomentum).getNorm() / projectileMomentum.getNorm() == Approx(0).margin(1e-2)); } - + SECTION("InteractionInterface Nuclei") { auto [stackPtr, secViewPtr] = setup::testing::setup_stack( - Code::Nucleus, 20, 10, 10100_GeV, (setup::Environment::BaseNodeType* const)nodePtr, - *csPtr); + Code::Nucleus, 20, 10, 10100_GeV, + (setup::Environment::BaseNodeType* const)nodePtr, *csPtr); setup::StackView& view = *(secViewPtr.get()); auto particle = stackPtr->first(); auto projectile = secViewPtr->getProjectile(); @@ -188,37 +190,37 @@ TEST_CASE("QgsjetIIInterface", "[processes]") { } SECTION("Heavy nuclei") { - + auto [stackPtr, secViewPtr] = setup::testing::setup_stack( - Code::Nucleus, 1000, 1000, 1100_GeV, (setup::Environment::BaseNodeType* const)nodePtr, - *csPtr); + Code::Nucleus, 1000, 1000, 1100_GeV, + (setup::Environment::BaseNodeType* const)nodePtr, *csPtr); setup::StackView& view = *(secViewPtr.get()); auto particle = stackPtr->first(); auto projectile = secViewPtr->getProjectile(); auto const projectileMomentum = projectile.getMomentum(); - + corsika::qgsjetII::Interaction model; - - CHECK_THROWS(model.getCrossSection(Code::Nucleus, Code::Nucleus, 100_GeV, 10., 1000.)); - CHECK_THROWS(model.getCrossSection(Code::Nucleus, Code::Nucleus, 100_GeV, 1000., 10.)); + + CHECK_THROWS( + model.getCrossSection(Code::Nucleus, Code::Nucleus, 100_GeV, 10., 1000.)); + CHECK_THROWS( + model.getCrossSection(Code::Nucleus, Code::Nucleus, 100_GeV, 1000., 10.)); CHECK_THROWS(model.doInteraction(view)); CHECK_THROWS(model.getInteractionLength(particle)); } SECTION("Allowed Particles") { - + auto [stackPtr, secViewPtr] = setup::testing::setup_stack( Code::Electron, 0, 0, 1100_GeV, (setup::Environment::BaseNodeType* const)nodePtr, *csPtr); - setup::StackView& view = *(secViewPtr.get()); auto particle = stackPtr->first(); auto projectile = secViewPtr->getProjectile(); auto const projectileMomentum = projectile.getMomentum(); - + corsika::qgsjetII::Interaction model; - - GrammageType const length = model.getInteractionLength(particle); - CHECK(length / (1_g / square(1_cm)) == std::numeric_limits<double>::infinity()); - } + GrammageType const length = model.getInteractionLength(particle); + CHECK(length / (1_g / square(1_cm)) == std::numeric_limits<double>::infinity()); + } }