From b4b2ffb554f11f01f8ca45f25b917ee33ee6e379 Mon Sep 17 00:00:00 2001 From: Maximilian Reininghaus <maximilian.reininghaus@kit.edu> Date: Fri, 14 Sep 2018 19:40:45 +0200 Subject: [PATCH] reworked particle properties --- Framework/Particles/ParticleClassNames.xml | 20 ++-- Framework/Particles/ParticleProperties.h | 52 ++++++---- Framework/Particles/pdxml_reader.py | 105 +++++++++------------ Framework/Particles/testParticles.cc | 2 +- Framework/Units/PhysicalConstants.h | 5 +- Framework/Units/PhysicalUnits.h | 2 +- 6 files changed, 96 insertions(+), 90 deletions(-) diff --git a/Framework/Particles/ParticleClassNames.xml b/Framework/Particles/ParticleClassNames.xml index 637f6b695..2c909fb1a 100644 --- a/Framework/Particles/ParticleClassNames.xml +++ b/Framework/Particles/ParticleClassNames.xml @@ -3,20 +3,20 @@ <!-- For selected particles it is possible to specify the C++ class names for a specific unique PDG code --> - <particle pdgID="0" classname="unknown" /> <!-- VOID in Pythia8 --> + <particle pdgID="0" classname="unknown"/> <!-- VOID in Pythia8 --> - <particle pdgID="11" classname="Electron" /> - <particle pdgID="-11" classname="Positron" /> - <particle pdgID="22" classname="Gamma" /> + <particle pdgID="11" classname="Electron"/> + <particle pdgID="-11" classname="Positron"/> + <particle pdgID="22" classname="Gamma"/> - <particle id="130" classname="K0Long"> </particle> - <particle id="310" classname="K0Short"> </particle> + <particle pdgID="130" classname="K0Long"/> + <particle pdgID="310" classname="K0Short"/> - <particle id="2112" classname="Neutron"> </particle> - <particle id="-2112" classname="AntiNeutron"> </particle> + <particle pdgID="2112" classname="Neutron"/> + <particle pdgID="-2112" classname="AntiNeutron"/> - <particle id="2212" classname="Proton"> </particle> - <particle id="-2212" classname="AntiProton"> </particle> + <particle pdgID="2212" classname="Proton"/> + <particle pdgID="-2212" classname="AntiProton"/> </chapter> diff --git a/Framework/Particles/ParticleProperties.h b/Framework/Particles/ParticleProperties.h index 7df92b661..23ef23698 100644 --- a/Framework/Particles/ParticleProperties.h +++ b/Framework/Particles/ParticleProperties.h @@ -4,15 +4,16 @@ Interface to particle properties */ -#ifndef _include_Particle_h_ -#define _include_Particle_h_ +#ifndef _include_ParticleProperties_h_ +#define _include_ParticleProperties_h_ #include <array> #include <cstdint> #include <iostream> +#include <type_traits> +#include <corsika/units/PhysicalConstants.h> #include <corsika/units/PhysicalUnits.h> -#include <corsika/particles/GeneratedParticleProperties.inc> /** * @namespace particle @@ -23,31 +24,50 @@ */ namespace corsika::particles { + enum class Code : int16_t; + + using PDGCodeType = int16_t; + using CodeIntType = std::underlying_type<Code>::type; + + // forward declarations to be used in GeneratedParticleProperties + int16_t constexpr GetElectricChargeNumber(Code const); + corsika::units::ElectricChargeType constexpr GetElectricCharge(Code const); + corsika::units::MassType constexpr GetMass(Code const); + PDGCodeType constexpr GetPDG(Code const); + std::string const GetName(Code const); + +#include <corsika/particles/GeneratedParticleProperties.inc> - /** - * @function GetMass - * - * return mass of particle + /*! + * returns mass of particle */ - auto constexpr GetMass(Code const p) { return masses[static_cast<uint8_t const>(p)]; } + corsika::units::MassType constexpr GetMass(Code const p) { + return masses[static_cast<CodeIntType const>(p)]; + } - auto constexpr GetPDG(Code const p) { return pdg_codes[static_cast<uint8_t const>(p)]; } + PDGCodeType constexpr GetPDG(Code const p) { + return pdg_codes[static_cast<CodeIntType const>(p)]; + } - auto constexpr GetElectricChargeNumber(Code const p) { - return electric_charge[static_cast<uint8_t const>(p)] / 3; + /*! + * returns electric charge of particle / (e/3). + */ + int16_t constexpr GetElectricChargeNumber(Code const p) { + return electric_charges[static_cast<CodeIntType const>(p)]; } - auto constexpr GetElectricCharge(Code const p) { - return GetElectricChargeNumber(p) * (corsika::units::constants::e); + corsika::units::ElectricChargeType constexpr GetElectricCharge(Code const p) { + return GetElectricChargeNumber(p) * (corsika::units::constants::e / 3.); } - auto const GetName(Code const p) { return names[static_cast<uint8_t const>(p)]; } + std::string const GetName(Code const p) { + return names[static_cast<CodeIntType const>(p)]; + } namespace io { std::ostream& operator<<(std::ostream& stream, Code const p) { - stream << GetName(p); - return stream; + return stream << GetName(p); } } // namespace io diff --git a/Framework/Particles/pdxml_reader.py b/Framework/Particles/pdxml_reader.py index 46958ea05..4c271e8fb 100755 --- a/Framework/Particles/pdxml_reader.py +++ b/Framework/Particles/pdxml_reader.py @@ -51,7 +51,7 @@ def class_names(filename): for particle in root.iter("particle"): name = particle.attrib["classname"] - pdg_id = int(particle.attrib["id"]) + pdg_id = int(particle.attrib["pdgID"]) map[pdg_id] = name return map @@ -97,9 +97,10 @@ def c_identifier(name): # # This function produces names of type "DeltaPlusPlus" # -def c_identifier_cammel(name): +def c_identifier_camel(name): orig = name name = name[0].upper() + name[1:].lower() # all lower case + for c in "() /": # replace funny characters name = name.replace(c, "_") @@ -111,9 +112,9 @@ def c_identifier_cammel(name): # move "Bar" to end of name ibar = name.find('Bar') - if (ibar>0 and ibar<len(name)-3) : + if ibar > 0 and ibar < len(name)-3: name = name[:ibar] + name[ibar+3:] + str('Bar') - + # cleanup "_"s while True: tmp = name.replace("__", "_") @@ -127,18 +128,18 @@ def c_identifier_cammel(name): istart = 0 while True: i = name.find('_', istart) - if (i<1 or i>len(name)-1): + if i < 1 or i > len(name)-1: break istart = i - if (name[i-1].isdigit() and name[i+1].isdigit()): + if name[i-1].isdigit() and name[i+1].isdigit(): # there is a number on both sides break name = name[:i] + name[i+1:] # and last, for example: make NuE out of Nue - if (name[i-1].islower() and name[i].islower()) : - if (i<len(name)-1) : + if name[i-1].islower() and name[i].islower(): + if i < len(name)-1: name = name[:i] + name[i].upper() + name[i+1:] - else : + else: name = name[:i] + name[i].upper() # check if name is valid C++ identifier @@ -146,27 +147,25 @@ def c_identifier_cammel(name): if pattern.match(name): return name else: - raise Exception("could not generate C identifier for '{:s}' name='{:s}'".format(orig, name)) + raise Exception("could not generate C identifier for '{:s}': result '{:s}'".format(orig, name)) -# ######################################################## +########################################################## # # returns dict containing all data from pythia-xml input # def build_pythia_db(filename, classnames): - particle_db = OrderedDict() + particle_db = OrderedDict() for (pdg, name, mass, electric_charge, antiName) in parse(filename): c_id = "unknown" - if (pdg in classnames): + if pdg in classnames: c_id = classnames[pdg] else: - c_id = c_identifier_cammel(name) # the cammel case names + c_id = c_identifier_camel(name) # the camel case names - #~ print(name, c_id, sep='\t', file=sys.stderr) - #~ enums += "{:s} = {:d}, ".format(c_id, corsika_id) particle_db[c_id] = { "name" : name, "antiName" : antiName, @@ -186,17 +185,20 @@ def build_pythia_db(filename, classnames): # def gen_internal_enum(pythia_db): - string = "enum class Code : uint8_t {\n" - - string += " FirstParticle = 1, // if you want to loop over particles, you want to start with \"1\" \n" # identifier for eventual loops... - last_ngc_id = 0 + string = ("enum class Code : int16_t {\n" + " FirstParticle = 1, // if you want to loop over particles, you want to start with \"1\" \n") # identifier for eventual loops... + + for k in filter(lambda k: "ngc_code" in pythia_db[k], pythia_db): last_ngc_id = pythia_db[k]['ngc_code'] string += " {key:s} = {code:d},\n".format(key = k, code = last_ngc_id) - string += " LastParticle = " + str(last_ngc_id+1) + ",\n" # identifier for eventual loops... + string += (" LastParticle = {:d},\n" # identifier for eventual loops... + "}};").format(last_ngc_id + 1) + + if last_ngc_id > 0x7fff: # does not fit into int16_t + raise Exception("Integer overflow in internal particle code definition prevented!") - string += "};" return string @@ -214,13 +216,13 @@ def gen_properties(pythia_db): string += "\n" # particle masses table - string += "static constexpr std::array<quantity<energy_d> const, size> masses = {\n" + string += "static constexpr std::array<corsika::units::MassType const, size> masses = {\n" for p in pythia_db.values(): - string += " {mass:f}_GeV, // {name:s}\n".format(mass = p['mass'], name = p['name']) + string += " {mass:f} * (1e9 * corsika::units::constants::eV / corsika::units::constants::cSquared), // {name:s}\n".format(mass = p['mass'], name = p['name']) string += "};\n\n" # PDG code table - string += "static constexpr std::array<PDGCode const, size> pdg_codes = {\n" + string += "static constexpr std::array<PDGCodeType const, size> pdg_codes = {\n" for p in pythia_db.values(): string += " {pdg:d}, // {name:s}\n".format(pdg = p['pdg'], name = p['name']) string += "};\n" @@ -232,7 +234,7 @@ def gen_properties(pythia_db): string += "};\n" # electric charges table - string += "static constexpr std::array<int16_t, size> electric_charge = {\n" + string += "static constexpr std::array<int16_t, size> electric_charges = {\n" for p in pythia_db.values(): string += " {charge:d},\n".format(charge = p['electric_charge']) string += "};\n" @@ -269,23 +271,23 @@ def gen_classes(pythia_db): string += "/** @class " + cname + "\n\n" string += " * Particle properties are taken from the PYTHIA8 ParticleData.xml file:<br>\n" string += " * - pdg=" + str(pythia_db[cname]['pdg']) +"\n" - string += " * - mass=" + str(pythia_db[cname]['mass']) + " GeV \n" + string += " * - mass=" + str(pythia_db[cname]['mass']) + " GeV/c² \n" string += " * - charge= " + str(pythia_db[cname]['electric_charge']/3) + " \n" string += " * - name=" + str(cname) + "\n" string += " * - anti=" + str(antiP) + "\n" string += "*/\n\n" - string += "class " + cname + "{\n" + string += "class " + cname + " {\n" string += " public:\n" - string += " static Code GetCode() { return Type; }\n" - string += " static quantity<energy_d> GetMass() { return masses[TypeIndex]; }\n" - string += " static quantity<electric_charge_d> GetCharge() { return corsika::units::constants::e*electric_charge[TypeIndex]/3; }\n" - string += " static int GetChargeNumber() { return electric_charge[TypeIndex]/3; }\n" - string += " static std::string GetName() { return names[TypeIndex]; }\n" - string += " static Code GetAntiParticle() { return AntiType; }\n" - string += " static const Code Type = Code::" + cname + ";\n" - string += " static const Code AntiType = Code::" + antiP + ";\n" + string += " static constexpr Code GetCode() { return Type; }\n" + string += " static constexpr corsika::units::MassType GetMass() { return corsika::particles::GetMass(Type); }\n" + string += " static constexpr corsika::units::ElectricChargeType GetCharge() { return corsika::particles::GetElectricCharge(Type); }\n" + string += " static constexpr int16_t GetChargeNumber() { return corsika::particles::GetElectricChargeNumber(Type); }\n" + string += " static std::string const GetName() { return corsika::particles::GetName(Type); }\n" + string += " static constexpr Code GetAntiParticle() { return AntiType; }\n" + string += " static constexpr Code Type = Code::" + cname + ";\n" + string += " static constexpr Code AntiType = Code::" + antiP + ";\n" string += " private:\n" - string += " static const uint8_t TypeIndex = static_cast<uint8_t const>(Type);\n" + string += " static constexpr CodeIntType TypeIndex = static_cast<CodeIntType const>(Type);\n" string += "};\n" return string @@ -296,22 +298,8 @@ def gen_classes(pythia_db): # def inc_start(): - string = "" - string += "#ifndef _include_GeneratedParticleDataTable_h_\n" - string += "#define _include_GeneratedParticleDataTable_h_\n\n" - string += "#include <corsika/units/PhysicalUnits.h>\n" - string += "#include <corsika/units/PhysicalConstants.h>\n" - string += "#include <array>\n" - string += "#include <cstdint>\n" -# string += "#include <iostream>\n\n" - string += "namespace corsika { \n\n" -# string += "using namespace literals; \n" - string += "namespace particles { \n\n" - string += "using corsika::units::energy_d;\n" - string += "using corsika::units::electric_charge_d;\n" - string += "using corsika::units::quantity;\n" - string += "using corsika::units::operator\"\"_GeV;\n" - string += "typedef int16_t PDGCode;\n\n" + string = ('// generated by pdxml_reader.py\n' + '// MANUAL EDITS ON OWN RISK\n') return string @@ -321,9 +309,6 @@ def inc_start(): def inc_end(): string = "" - string += "\n}\n\n" - string += "\n}\n\n" - string += "#endif\n" return string @@ -334,14 +319,14 @@ def inc_end(): if __name__ == "__main__": - if (len(sys.argv)!=3) : - print ("pdxml_reader.py Pythia8.xml ClassNames.xml") - sys.exit(0) + if len(sys.argv) != 3: + print("usage: {:s} <Pythia8.xml> <ClassNames.xml>".format(sys.argv[0])) + sys.exit(1) names = class_names(sys.argv[2]) pythia_db = build_pythia_db(sys.argv[1], names) - print ("\n pdxml_reader.py: Automatically produce particle-properties from PYTHIA8 xml file\n") + print("\n pdxml_reader.py: Automatically produce particle-properties from PYTHIA8 xml file\n") counter = itertools.count(0) diff --git a/Framework/Particles/testParticles.cc b/Framework/Particles/testParticles.cc index 476264171..a04f36515 100644 --- a/Framework/Particles/testParticles.cc +++ b/Framework/Particles/testParticles.cc @@ -14,7 +14,7 @@ TEST_CASE("Particles", "[Particles]") { SECTION("Types") { REQUIRE(Electron::GetCode() == Code::Electron); } SECTION("Data") { - REQUIRE(Electron::GetMass() / 0.511_MeV == Approx(1)); + REQUIRE(Electron::GetMass() / (511_keV / constants::cSquared) == Approx(1)); REQUIRE(Electron::GetMass() / GetMass(Code::Electron) == Approx(1)); REQUIRE(Electron::GetCharge() / constants::e == Approx(-1)); REQUIRE(Positron::GetCharge() / constants::e == Approx(+1)); diff --git a/Framework/Units/PhysicalConstants.h b/Framework/Units/PhysicalConstants.h index a666f15ef..5c1b3a91c 100644 --- a/Framework/Units/PhysicalConstants.h +++ b/Framework/Units/PhysicalConstants.h @@ -39,16 +39,17 @@ namespace corsika::units::constants { // Avogadro constant constexpr quantity<dimensions<0, 0, 0, 0, 0, -1> > N_sub_A{Rep(6.02214199e+23L) / mole}; // electronvolt - constexpr quantity<energy_d> eV{Rep(1.60217733e-19L) * joule}; + constexpr quantity<energy_d> eV{Rep(1.6021766208e-19L) * joule}; // elementary charge - constexpr quantity<electric_charge_d> e{Rep(1.602176462e-19L) * coulomb}; + constexpr quantity<electric_charge_d> e{Rep(1.6021766208e-19L) * coulomb}; // Planck constant constexpr quantity<dimensions<2, 1, -1> > h{Rep(6.62606876e-34L) * joule * second}; // speed of light in a vacuum constexpr quantity<speed_d> c{Rep(299792458L) * meter / second}; + constexpr auto cSquared = c * c; // unified atomic mass unit constexpr quantity<mass_d> u{Rep(1.6605402e-27L) * kilogram}; diff --git a/Framework/Units/PhysicalUnits.h b/Framework/Units/PhysicalUnits.h index 9aeefb151..c1e6bd204 100644 --- a/Framework/Units/PhysicalUnits.h +++ b/Framework/Units/PhysicalUnits.h @@ -22,7 +22,6 @@ namespace phys { } // namespace phys namespace corsika::units { - using namespace phys::units; using namespace phys::units::literals; // namespace literals = phys::units::literals; @@ -34,6 +33,7 @@ namespace corsika::units { using ElectricChargeType = phys::units::quantity<phys::units::electric_charge_d, double>; using EnergyType = phys::units::quantity<phys::units::energy_d, double>; + using MassType = phys::units::quantity<phys::units::mass_d, double>; } // end namespace corsika::units -- GitLab