diff --git a/include/SimCore/KaonPhysics.h b/include/SimCore/KaonPhysics.h index dddaf85..1ced0b7 100644 --- a/include/SimCore/KaonPhysics.h +++ b/include/SimCore/KaonPhysics.h @@ -9,6 +9,8 @@ #include #include #include +#include +#include #include #include @@ -38,22 +40,25 @@ class KaonPhysics : public G4VPhysicsConstructor { * K^+ -> \mu^+ + \nu_\mu * K^+ -> \pi^+ + \pi^0 * K^+ -> \pi^+ + \pi^- + \pi^+ - * K^+ -> \pi^+ + \pi^0 + \pi^0 * K^+ -> \pi^0 + e^+ + \nu_e * K^+ -> \pi^0 + \mu^+ + \nu_\mu + * K^+ -> \pi^+ + \pi^0 + \pi^0 * * And vice versa for K^-. * The indices here correspond to the position of the branching ratio for * that process in the corresponding parameter as well as the position in * the decay table. + * + * @note: The order in the the decay table is sorted by the branching ratios + * of the default physics settings! */ enum ChargedKaonDecayChannel { mu_nu = 0, pi_pi0 = 1, pi_pi_pi = 2, - pi_pi0_pi0 = 3, - pi0_e_nu = 4, - pi0_mu_nu = 5 + pi0_e_nu = 3, + pi0_mu_nu = 4, + pi_pi0_pi0 = 5, }; /** * @@ -62,26 +67,29 @@ class KaonPhysics : public G4VPhysicsConstructor { * * The processes are * - * K^0_L -> \pi^0 + \pi^0 + \pi^0 - * K^0_L -> \pi^0 + \pi^+ + \pi^- * K^0_L -> \pi^- + e^+ + \nu_e * K^0_L -> \pi^+ + e^- + \nu_e + * K^0_L -> \pi^0 + \pi^0 + \pi^0 * K^0_L -> \pi^- + \mu^+ + \nu_\mu * K^0_L -> \pi^+ + \mu^- + \nu_\mu + * K^0_L -> \pi^0 + \pi^+ + \pi^- * * and * * K^0_S -> \pi^+ + \pi^- * K^0_S -> \pi^0 + \pi^0 * + * @note: The order in the the decay table is sorted by the branching ratios + * of the default physics settings! + * **/ enum KaonZeroLongDecayChannel { - pi0_pi0_pi0 = 0, - pi0_pip_pim = 1, - pip_e_nu = 2, - pim_e_nu = 3, - pim_mu_nu = 4, - pip_mu_nu = 5, + pim_e_nu = 0, + pip_e_nu = 1, + pi0_pi0_pi0 = 2, + pim_mu_nu = 3, + pip_mu_nu = 4, + pi0_pip_pim = 5, }; enum KaonZeroShortDecayChannel { pip_pim = 0, @@ -102,6 +110,10 @@ class KaonPhysics : public G4VPhysicsConstructor { std::vector k0l_branching_ratios; std::vector k0s_branching_ratios; + // If > 0, dump details about what was changed + // If > 1, dump details about the initial branching ratios + int verbosity; + public: KaonPhysics(const G4String& name, const framework::config::Parameters& parameters); @@ -127,6 +139,8 @@ class KaonPhysics : public G4VPhysicsConstructor { void ConstructParticle() override; + void DumpDecayDetails(const G4ParticleDefinition* kaon) const; + /** * Construct processes * diff --git a/python/kaon_physics.py b/python/kaon_physics.py index 3ead464..a35aa3b 100644 --- a/python/kaon_physics.py +++ b/python/kaon_physics.py @@ -43,31 +43,37 @@ class KaonPhysics(): Multiplicative factor to scale the lifetime of kaons by. Default is to scale by 1 (no scaling) + verbosity : int + + If > 0: Dump details about the decay after the update + If > 1: Also dump details about the decay before the update to show the + difference + """ def __init__(self): self.kplus_branching_ratios = [ 0.6355, # K^+ -> mu^+ + nu_mu 0.2066, # K^+ -> pi^+ + pi^0 0.0559, # K^+ -> pi^+ + pi^- + pi^+ - 0.01761, # K^+ -> pi^+ + pi^0 + pi^0 0.0507, # K^+ -> pi^0 + e^+ + nu_e 0.0335, # K^+ -> pi^0 + mu^+ + nu_mu + 0.01761, # K^+ -> pi^+ + pi^0 + pi^0 ] self.kminus_branching_ratios = [ 0.6355, # K^- -> mu^- + anti_nu_mu 0.2066, # K^- -> pi^- + pi^0 0.0559, # K^- -> pi^- + pi^+ + pi^- - 0.01761, # K-+ -> pi^- + pi^0 + pi^0 0.0507, # K-+ -> pi^0 + e^- + anti_nu_e 0.0335, # K-+ -> pi^0 + mu^- + anti_nu_mu + 0.01761, # K-+ -> pi^- + pi^0 + pi^0 ] self.k0l_branching_ratios = [ - 0.1952, # K^0_L -> pi^0 + pi^0 + pi^0 - 0.1254, # K^0_L -> pi^0 + pi^+ + pi^- 0.2027, # K^0_L -> pi^- + e^+ + nu_e 0.2027, # K^0_L -> pi^+ + e^- + anti_nu_e + 0.1952, # K^0_L -> pi^0 + pi^0 + pi^0 0.1352, # K^0_L -> pi^- + mu^+ + nu_mu 0.1352, # K^0_L -> pi^+ + mu^- + anti_nu_mu + 0.1254, # K^0_L -> pi^0 + pi^+ + pi^- ] self.k0s_branching_ratios = [ 0.6920, # K^0_S -> pi^+ + pi^- @@ -78,6 +84,8 @@ def __init__(self): self.k0l_lifetime_factor = 1. self.k0s_lifetime_factor = 1. + self.verbosity=0 + def upKaons(): """Returns a configuration of the kaon physics corresponding @@ -93,22 +101,23 @@ def upKaons(): kaon_physics = KaonPhysics() kaon_physics.kplus_branching_ratios = [ 0.8831, # K^+ -> mu^+ + nu_mu - 0., # K^+ -> pi^+ + pi^0 - 0., # K^+ -> pi^+ + pi^- + pi^+ - 0., # K^+ -> pi^+ + pi^0 + pi^0 + 0., # K^+ -> pi^+ + pi^0 + 0., # K^+ -> pi^+ + pi^- + pi^+ 0.0704, # K^+ -> pi^0 + e^+ + nu_e 0.0465, # K^+ -> pi^0 + mu^+ + nu_mu + 0., # K^+ -> pi^+ + pi^0 + pi^0 ] kaon_physics.kminus_branching_ratios = [ 0.8831, # K^- -> mu^- + anti_nu_mu - 0., # K^- -> pi^- + pi^0 - 0., # K^- -> pi^- + pi^+ + pi^- - 0., # K-+ -> pi^- + pi^0 + pi^0 - 0.0704, # K-+ -> pi^0 + e^- + anti_nu_e - 0.0464, # K-+ -> pi^0 + mu^- + anti_nu_mu + 0., # K^- -> pi^- + pi^0 + 0., # K^- -> pi^- + pi^+ + pi^- + 0.0704, # K^- -> pi^0 + e^- + anti_nu_e + 0.0464, # K^- -> pi^0 + mu^- + anti_nu_mu + 0., # K^- -> pi^- + pi^0 + pi^0 ] kaon_physics.kplus_lifetime_factor = 1/50. kaon_physics.kminus_lifetime_factor = 1/50. + kaon_physics.verbosity = 2 return kaon_physics def __setattr__(self, key, value): diff --git a/src/SimCore/KaonPhysics.cxx b/src/SimCore/KaonPhysics.cxx index 0cbc18b..84546b8 100644 --- a/src/SimCore/KaonPhysics.cxx +++ b/src/SimCore/KaonPhysics.cxx @@ -18,16 +18,22 @@ KaonPhysics::KaonPhysics(const G4String& name, parameters.getParameter("kminus_lifetime_factor"); k0l_lifetime_factor = parameters.getParameter("k0l_lifetime_factor"); k0s_lifetime_factor = parameters.getParameter("k0s_lifetime_factor"); + verbosity = parameters.getParameter("verbosity"); } void KaonPhysics::setDecayProperties( G4ParticleDefinition* kaon, const std::vector& branching_ratios, double lifetime_factor) const { - kaon->SetPDGLifeTime(kaon->GetPDGLifeTime() * lifetime_factor); auto table{kaon->GetDecayTable()}; if (!table) { EXCEPTION_RAISE("KaonPhysics", "Unable to get the decay table from " + kaon->GetParticleName()); } + if (verbosity > 1) { + std::cout << "Decay details after setting branching ratios and lifetimes" + << std::endl; + DumpDecayDetails(kaon); + } + kaon->SetPDGLifeTime(kaon->GetPDGLifeTime() * lifetime_factor); if (kaon == G4KaonZeroLong::Definition()) { (*table)[KaonZeroLongDecayChannel::pi0_pi0_pi0]->SetBR( branching_ratios[KaonZeroLongDecayChannel::pi0_pi0_pi0]); @@ -60,6 +66,11 @@ void KaonPhysics::setDecayProperties( (*table)[ChargedKaonDecayChannel::pi0_mu_nu]->SetBR( branching_ratios[ChargedKaonDecayChannel::pi0_mu_nu]); } + if (verbosity > 0) { + std::cout << "Decay details after setting branching ratios and lifetimes" + << std::endl; + DumpDecayDetails(kaon); + } } void KaonPhysics::ConstructParticle() { auto kaonPlus{G4KaonPlus::Definition()}; @@ -79,4 +90,27 @@ void KaonPhysics::ConstructParticle() { setDecayProperties(kaonShort, k0s_branching_ratios, k0s_lifetime_factor); } +void KaonPhysics::DumpDecayDetails(const G4ParticleDefinition* kaon) const { + std::cout << "Decay table details for " << kaon->GetParticleName() + << std::endl + << std::scientific << std::setprecision(15); + std::cout << "PDG Lifetime " << kaon->GetPDGLifeTime() << std::endl; + const auto table{kaon->GetDecayTable()}; + const int entries{table->entries()}; + for (auto i{0}; i < entries; ++i) { + const auto channel{(*table)[i]}; + std::cout << "Channel " << i << " Kinematics type " + << channel->GetKinematicsName() << " with BR " << channel->GetBR() + << std::endl; + std::cout << kaon->GetParticleName() << " -> "; + const auto daughters{channel->GetNumberOfDaughters()}; + for (auto j{0}; j < daughters - 1; ++j) { + std::cout << channel->GetDaughter(j)->GetParticleName() << " + "; + } + // Special formatting for last one :) + std::cout << channel->GetDaughter(daughters - 1)->GetParticleName() + << std::endl; + } +} + } // namespace simcore