diff --git a/Framework/Particles/CMakeLists.txt b/Framework/Particles/CMakeLists.txt index 8ea8654fad658692dbbc1819dbc2ab957cb64de5..3301139d597d1bc513e8b226db4708ca4c20f6ed 100644 --- a/Framework/Particles/CMakeLists.txt +++ b/Framework/Particles/CMakeLists.txt @@ -1,12 +1,14 @@ add_custom_command ( OUTPUT ${PROJECT_BINARY_DIR}/Framework/Particles/GeneratedParticleProperties.inc - ${PROJECT_BINARY_DIR}/Framework/Particles/pythia_db.pkl + ${PROJECT_BINARY_DIR}/Framework/Particles/particle_db.pkl COMMAND ${PROJECT_SOURCE_DIR}/Framework/Particles/pdxml_reader.py ${PROJECT_SOURCE_DIR}/Framework/Particles/ParticleData.xml + ${PROJECT_SOURCE_DIR}/Framework/Particles/NuclearData.xml ${PROJECT_SOURCE_DIR}/Framework/Particles/ParticleClassNames.xml DEPENDS pdxml_reader.py ParticleData.xml + NuclearData.xml ParticleClassNames.xml WORKING_DIRECTORY ${PROJECT_BINARY_DIR}/Framework/Particles/ diff --git a/Framework/Particles/NuclearData.xml b/Framework/Particles/NuclearData.xml new file mode 100644 index 0000000000000000000000000000000000000000..d6d2a46dee4c0ac3beb53e8ec5ca5db94bf20680 --- /dev/null +++ b/Framework/Particles/NuclearData.xml @@ -0,0 +1,27 @@ +<chapter name="Nuclear Data"> + +<particle id="1000010010" name="hydrogen" A="1" Z="1" > +</particle> + +<particle id="1000010020" name="deuterium" A="2" Z="1" > +</particle> + +<particle id="1000010030" name="tritium" A="3" Z="1" > +</particle> + +<particle id="1000020040" name="helium" A="4" Z="2" > +</particle> + +<particle id="1000060120" name="carbon" A="12" Z="6" > +</particle> + +<particle id="1000070140" name="nitrogen" A="14" Z="7" > +</particle> + +<particle id="1000080160" name="oxygen" A="16" Z="8" > +</particle> + +<particle id="1000180400" name="argon" A="40" Z="18" > +</particle> + +</chapter> diff --git a/Framework/Particles/ParticleProperties.h b/Framework/Particles/ParticleProperties.h index e7d48547f87779c8da84ea2c318de55ee8b163a8..91bb74eefd485b3f8afbd67fc0eaeabdbf4f7b42 100644 --- a/Framework/Particles/ParticleProperties.h +++ b/Framework/Particles/ParticleProperties.h @@ -40,7 +40,7 @@ namespace corsika::particles { enum class Code : int16_t; - using PDGCodeType = int16_t; + using PDGCodeType = int32_t; using CodeIntType = std::underlying_type<Code>::type; // forward declarations to be used in GeneratedParticleProperties @@ -51,24 +51,28 @@ namespace corsika::particles { constexpr std::string const& GetName(Code const); corsika::units::si::TimeType constexpr GetLifetime(Code const); + bool constexpr IsNucleus(Code const); + int constexpr GetNucleusA(Code const); + int constexpr GetNucleusZ(Code const); + #include <corsika/particles/GeneratedParticleProperties.inc> /*! * returns mass of particle */ corsika::units::hep::MassType constexpr GetMass(Code const p) { - return masses[static_cast<CodeIntType const>(p)]; + return detail::masses[static_cast<CodeIntType const>(p)]; } PDGCodeType constexpr GetPDG(Code const p) { - return pdg_codes[static_cast<CodeIntType const>(p)]; + return detail::pdg_codes[static_cast<CodeIntType const>(p)]; } /*! * returns electric charge of particle / (e/3). */ int16_t constexpr GetElectricChargeNumber(Code const p) { - return electric_charges[static_cast<CodeIntType const>(p)]; + return detail::electric_charges[static_cast<CodeIntType const>(p)]; } corsika::units::si::ElectricChargeType constexpr GetElectricCharge(Code const p) { @@ -76,13 +80,26 @@ namespace corsika::particles { } constexpr std::string const& GetName(Code const p) { - return names[static_cast<CodeIntType const>(p)]; + return detail::names[static_cast<CodeIntType const>(p)]; } corsika::units::si::TimeType constexpr GetLifetime(Code const p) { - return lifetime[static_cast<CodeIntType const>(p)]; + return detail::lifetime[static_cast<CodeIntType const>(p)] * corsika::units::si::second; } + bool constexpr IsNucleus(Code const p) { + return detail::isNucleus[static_cast<CodeIntType const>(p)]; + } + + int constexpr GetNucleusA(Code const p) { + return detail::nucleusA[static_cast<CodeIntType const>(p)]; + } + + int constexpr GetNucleusZ(Code const p) { + return detail::nucleusZ[static_cast<CodeIntType const>(p)]; + } + + namespace io { std::ostream& operator<<(std::ostream& stream, Code const p); diff --git a/Framework/Particles/pdxml_reader.py b/Framework/Particles/pdxml_reader.py index aa77d0c5a4e954f713c4cac0aeb938bcf00d8bff..ca7701d40dbb467781a63ffb304adfe8bdb2a7e3 100755 --- a/Framework/Particles/pdxml_reader.py +++ b/Framework/Particles/pdxml_reader.py @@ -5,21 +5,23 @@ import xml.etree.ElementTree as ET from collections import OrderedDict import pickle +GeVfm = 0.19732696312541853 +c_speed_of_light = 29.9792458e10 # mm / s +# for nuclear masses +mneutron = 0.9395654133 # GeV +mproton = 0.9382720813 # GeV + ############################################################## # # reading xml input data, return line by line particle data # - -def parse(filename): +def parsePythia(filename): tree = ET.parse(filename) root = tree.getroot() - - GeVfm = 0.19732696312541853 - c_speed_of_light = 29.9792458e10 # mm / s for particle in root.iter("particle"): - name = particle.attrib["name"] + name = particle.attrib["name"] antiName = "Unknown" if ("antiName" in particle.attrib): antiName = particle.attrib["antiName"] @@ -47,12 +49,39 @@ def parse(filename): yield (-pdg_id, antiName, mass, -electric_charge, name, ctau/c_speed_of_light) - ############################################################## # -# returns dict with particle codes and class names +# reading xml input data, return line by line particle data # +def parseNuclei(filename): + tree = ET.parse(filename) + root = tree.getroot() + + for particle in root.iter("particle"): + name = particle.attrib["name"] + antiName = "Unknown" + if ("antiName" in particle.attrib): + antiName = particle.attrib["antiName"] + pdg_id = int(particle.attrib["id"]) + A = int(particle.attrib["A"]) + Z = int(particle.attrib["Z"]) + # mass in GeV + if ("mass" in particle.attrib): + mass = particle.attrib["mass"] + else: + mass = (A-Z)*mneutron + Z*mproton + + electric_charge = Z*3 # in units of e/3 + ctau = float('Inf') + + yield (pdg_id, name, mass, electric_charge, antiName, ctau/c_speed_of_light, A, Z) + + +############################################################## +# +# returns dict with particle codes and class names +# def class_names(filename): tree = ET.parse(filename) root = tree.getroot() @@ -164,14 +193,12 @@ def c_identifier_camel(name): ########################################################## # # returns dict containing all data from pythia-xml input -# +# +def read_pythia_db(filename, particle_db, classnames): -def build_pythia_db(filename, classnames): - particle_db = OrderedDict() + counter = itertools.count(len(particle_db)) - counter = itertools.count(0) - - for (pdg, name, mass, electric_charge, antiName, lifetime) in parse(filename): + for (pdg, name, mass, electric_charge, antiName, lifetime) in parsePythia(filename): c_id = "Unknown" if pdg in classnames: @@ -186,7 +213,41 @@ def build_pythia_db(filename, classnames): "mass" : mass, # in GeV "electric_charge" : electric_charge, # in e/3 "lifetime" : lifetime, - "ngc_code" : next(counter) + "ngc_code" : next(counter), + "isNucleus" : False, + } + + return particle_db + + + +########################################################## +# +# returns dict containing all data from pythia-xml input +# +def read_nuclei_db(filename, particle_db, classnames): + + counter = itertools.count(len(particle_db)) + + for (pdg, name, mass, electric_charge, antiName, lifetime, A, Z) in parseNuclei(filename): + + c_id = "Unknown" + if pdg in classnames: + c_id = classnames[pdg] + else: + c_id = c_identifier_camel(name) + + particle_db[c_id] = { + "name" : name, + "antiName" : antiName, + "pdg" : pdg, + "mass" : mass, # in GeV + "electric_charge" : electric_charge, # in e/3 + "lifetime" : lifetime, + "ngc_code" : next(counter), + "A" : A, + "Z" : Z, + "isNucleus" : True, } return particle_db @@ -198,14 +259,13 @@ def build_pythia_db(filename, classnames): # # return string with enum of all internal particle codes # - -def gen_internal_enum(pythia_db): +def gen_internal_enum(particle_db): 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'] + for k in filter(lambda k: "ngc_code" in particle_db[k], particle_db): + last_ngc_id = particle_db[k]['ngc_code'] string += " {key:s} = {code:d},\n".format(key = k, code = last_ngc_id) string += (" LastParticle = {:d},\n" # identifier for eventual loops... @@ -223,52 +283,83 @@ def gen_internal_enum(pythia_db): # # return string with all data arrays # - -def gen_properties(pythia_db): +def gen_properties(particle_db): # number of particles, size of tables - string = "static constexpr std::size_t size = {size:d};\n".format(size = len(pythia_db)) + string = "static constexpr std::size_t size = {size:d};\n".format(size = len(particle_db)) string += "\n" # particle masses table string += "static constexpr std::array<corsika::units::hep::MassType const, size> masses = {\n" - for p in pythia_db.values(): - string += " {mass:f} * (1e9 * corsika::units::si::constants::eV), // {name:s}\n".format(mass = p['mass'], name = p['name']) + for p in particle_db.values(): + string += " {mass:e} * (1e9 * corsika::units::si::constants::eV), // {name:s}\n".format(mass = p['mass'], name = p['name']) string += "};\n\n" # PDG code table string += "static constexpr std::array<PDGCodeType const, size> pdg_codes = {\n" - for p in pythia_db.values(): + for p in particle_db.values(): string += " {pdg:d}, // {name:s}\n".format(pdg = p['pdg'], name = p['name']) string += "};\n" # name string table string += "static const std::array<std::string const, size> names = {\n" - for p in pythia_db.values(): + for p in particle_db.values(): string += " \"{name:s}\",\n".format(name = p['name']) string += "};\n" # electric charges table string += "static constexpr std::array<int16_t, size> electric_charges = {\n" - for p in pythia_db.values(): + for p in particle_db.values(): string += " {charge:d},\n".format(charge = p['electric_charge']) string += "};\n" # anti-particle table -# string += "static constexpr std::array<size, size> anti_particle = {\n" -# for p in pythia_db.values(): -# string += " {anti:d},\n".format(charge = p['anti_particle']) -# string += "};\n" + # string += "static constexpr std::array<size, size> anti_particle = {\n" + # for p in particle_db.values(): + # string += " {anti:d},\n".format(charge = p['anti_particle']) + # string += "};\n" # lifetime - string += "static constexpr std::array<corsika::units::si::TimeType const, size> lifetime = {\n" - for p in pythia_db.values(): + #string += "static constexpr std::array<corsika::units::si::TimeType const, size> lifetime = {\n" + string += "static constexpr std::array<double const, size> lifetime = {\n" + for p in particle_db.values(): if p['lifetime'] == float("Inf") : - string += " std::numeric_limits<double>::infinity() * corsika::units::si::second, \n" + string += " std::numeric_limits<double>::infinity(), \n" # * corsika::units::si::second, \n" else : - string += " {tau:f} * corsika::units::si::second, \n".format(tau = p['lifetime']) + string += " {tau:e}, \n".format(tau = p['lifetime']) + #string += " {tau:e} * corsika::units::si::second, \n".format(tau = p['lifetime']) + string += "};\n" + + + ### nuclear data ### + + # is nucleus flag + string += "static constexpr std::array<bool, size> isNucleus = {\n" + for p in particle_db.values(): + value = 'false' + if p['isNucleus']: + value = 'true' + string += " {val},\n".format(val = value) + string += "};\n" + + # nucleus mass number A + string += "static constexpr std::array<int, size> nucleusA = {\n" + for p in particle_db.values(): + A = 0 + if p['isNucleus']: + A = p['A'] + string += " {val},\n".format(val = A) string += "};\n" + # nucleus charge number Z + string += "static constexpr std::array<int, size> nucleusZ = {\n" + for p in particle_db.values(): + Z = 0 + if p['isNucleus']: + Z = p['Z'] + string += " {val},\n".format(val = Z) + string += "};\n" + return string @@ -278,27 +369,29 @@ def gen_properties(pythia_db): # # return string with a list of classes for all particles # - -def gen_classes(pythia_db): +def gen_classes(particle_db): string = "// list of C++ classes to access particle properties" - for cname in pythia_db: + for cname in particle_db: antiP = 'Unknown' - for cname_anti in pythia_db: - if (pythia_db[cname_anti]['name'] == pythia_db[cname]['antiName']): + for cname_anti in particle_db: + if (particle_db[cname_anti]['name'] == particle_db[cname]['antiName']): antiP = cname_anti break string += "\n"; 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 += " * - charge= " + str(pythia_db[cname]['electric_charge']/3) + " \n" + string += " * - pdg=" + str(particle_db[cname]['pdg']) +"\n" + string += " * - mass=" + str(particle_db[cname]['mass']) + " GeV \n" + string += " * - charge= " + str(particle_db[cname]['electric_charge']/3) + " \n" string += " * - name=" + str(cname) + "\n" string += " * - anti=" + str(antiP) + "\n" + if (particle_db[cname]['isNucleus']): + string += " * - nuclear A=" + str(particle_db[cname]['A']) + "\n" + string += " * - nuclear Z=" + str(particle_db[cname]['Z']) + "\n" string += "*/\n\n" string += "class " + cname + " {\n" string += " public:\n" @@ -308,6 +401,9 @@ def gen_classes(pythia_db): 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 bool IsNucleus() { return corsika::particles::IsNucleus(Type); }\n" + string += " static constexpr int GetNucleusA() { return corsika::particles::GetNucleusA(Type); }\n" + string += " static constexpr int GetNucleusZ() { return corsika::particles::GetNucleusZ(Type); }\n" string += " static constexpr Code Type = Code::" + cname + ";\n" string += " static constexpr Code AntiType = Code::" + antiP + ";\n" string += " private:\n" @@ -320,17 +416,30 @@ def gen_classes(pythia_db): ############################################################### # # - def inc_start(): string = ('// generated by pdxml_reader.py\n' - '// MANUAL EDITS ON OWN RISK\n') + '// MANUAL EDITS ON OWN RISK. THEY WILL BE OVERWRITTEN. \n') return string ############################################################### # # +def detail_start(): + string = ('namespace detail {\n\n') + return string + + +############################################################### +# +# +def detail_end(): + string = "\n}//end namespace detail\n" + return string +############################################################### +# +# def inc_end(): string = "" return string @@ -338,36 +447,38 @@ def inc_end(): ################################################################### # -# Serialize pythia_db into file +# Serialize particle_db into file # +def serialize_particle_db(particle_db, file): + pickle.dump(particle_db, file) -def serialize_pythia_db(pythia_db, file): - pickle.dump(pythia_db, file) - + ################################################################### # # Main function # - if __name__ == "__main__": - - if len(sys.argv) != 3: - print("usage: {:s} <Pythia8.xml> <ClassNames.xml>".format(sys.argv[0]), file=sys.stderr) + + if len(sys.argv) != 4: + print("usage: {:s} <Pythia8.xml> <Nuclei.xml> <ClassNames.xml>".format(sys.argv[0]), file=sys.stderr) sys.exit(1) - print("\n pdxml_reader.py: Automatically produce particle-properties from PYTHIA8 xml file\n") - - names = class_names(sys.argv[2]) - pythia_db = build_pythia_db(sys.argv[1], names) - - print (str(pythia_db)) + print("\n pdxml_reader.py: Automatically produce particle-properties from input files\n") + + names = class_names(sys.argv[3]) + particle_db = OrderedDict() + read_pythia_db(sys.argv[1], particle_db, names) + read_nuclei_db(sys.argv[2], particle_db, names) with open("GeneratedParticleProperties.inc", "w") as f: - print(inc_start(), file=f) - print(gen_internal_enum(pythia_db), file=f) - print(gen_properties(pythia_db), file=f) - print(gen_classes(pythia_db), file=f) + print(inc_start(), file=f) + print(gen_internal_enum(particle_db), file=f) + print(detail_start(), file=f) + print(gen_properties(particle_db), file=f) + print(detail_end(), file=f) + print(gen_classes(particle_db), file=f) print(inc_end(), file=f) - with open("pythia_db.pkl", "wb") as f: - serialize_pythia_db(pythia_db, f) + with open("particle_db.pkl", "wb") as f: + serialize_particle_db(particle_db, f) + diff --git a/Framework/Particles/testParticles.cc b/Framework/Particles/testParticles.cc index 78f34d0beb262166e51839bbedf3b83521355d6a..bf307798b57d4998115a54321de965e6e53d8bb3 100644 --- a/Framework/Particles/testParticles.cc +++ b/Framework/Particles/testParticles.cc @@ -58,11 +58,28 @@ TEST_CASE("ParticleProperties", "[Particles]") { REQUIRE(GetLifetime(Code::Electron) == std::numeric_limits<double>::infinity() * corsika::units::si::second); REQUIRE(GetLifetime(Code::DPlus) < GetLifetime(Code::Gamma)); - // REQUIRE(GetLifetime(Code::RhoPlus)/corsika::units::si::second == - // (Approx(4.414566727909413e-24).epsilon(1e-3))); - // REQUIRE(GetLifetime(Code::SigmaMinusBar)/corsika::units::si::second == - // (Approx(8.018880848563575e-11).epsilon(1e-5))); - // REQUIRE(GetLifetime(Code::MuPlus)/corsika::units::si::second == - // (Approx(2.1970332555864364e-06).epsilon(1e-5))); + REQUIRE(GetLifetime(Code::RhoPlus)/corsika::units::si::second == + (Approx(4.414566727909413e-24).epsilon(1e-3))); + REQUIRE(GetLifetime(Code::SigmaMinusBar)/corsika::units::si::second == + (Approx(8.018880848563575e-11).epsilon(1e-5))); + REQUIRE(GetLifetime(Code::MuPlus)/corsika::units::si::second == + (Approx(2.1970332555864364e-06).epsilon(1e-5))); } + + SECTION("Nuclei") { + REQUIRE(IsNucleus(Code::Gamma) == false); + REQUIRE(IsNucleus(Code::Argon) == true); + REQUIRE(IsNucleus(Code::Proton) == false); + REQUIRE(IsNucleus(Code::Hydrogen) == true); + REQUIRE(Argon::IsNucleus() == true); + REQUIRE(EtaC::IsNucleus() == false); + + REQUIRE(GetNucleusA(Code::Hydrogen) == 1); + REQUIRE(GetNucleusA(Code::Tritium) == 3); + REQUIRE(Hydrogen::GetNucleusZ() == 1); + REQUIRE(Tritium::GetNucleusA() == 3); + + + } + } diff --git a/Processes/Sibyll/CMakeLists.txt b/Processes/Sibyll/CMakeLists.txt index 31d204714dba4aa72301e97ab3f5402e1abdeb73..880ed76d5ed5c2da19b6228b39680984309f1147 100644 --- a/Processes/Sibyll/CMakeLists.txt +++ b/Processes/Sibyll/CMakeLists.txt @@ -1,11 +1,11 @@ add_custom_command ( OUTPUT ${PROJECT_BINARY_DIR}/Processes/Sibyll/Generated.inc COMMAND ${PROJECT_SOURCE_DIR}/Processes/Sibyll/code_generator.py - ${PROJECT_BINARY_DIR}/Framework/Particles/pythia_db.pkl + ${PROJECT_BINARY_DIR}/Framework/Particles/particle_db.pkl ${PROJECT_SOURCE_DIR}/Processes/Sibyll/sibyll_codes.dat DEPENDS code_generator.py sibyll_codes.dat - ${PROJECT_BINARY_DIR}/Framework/Particles/pythia_db.pkl + ${PROJECT_BINARY_DIR}/Framework/Particles/particle_db.pkl WORKING_DIRECTORY ${PROJECT_BINARY_DIR}/Processes/Sibyll/ COMMENT "Generate conversion tables for particle codes SIBYLL <-> CORSIKA" diff --git a/Processes/Sibyll/code_generator.py b/Processes/Sibyll/code_generator.py index 388668e1de0d5b04c8b5fa920308d20ab9982d39..e1cfcbfaf84e4d04aff4ed35b324daacf1fc66c1 100755 --- a/Processes/Sibyll/code_generator.py +++ b/Processes/Sibyll/code_generator.py @@ -13,16 +13,16 @@ import pickle, sys, itertools -# loads the pickled pythia_db (which is an OrderedDict) -def load_pythiadb(filename): +# loads the pickled particle_db (which is an OrderedDict) +def load_particledb(filename): with open(filename, "rb") as f: - pythia_db = pickle.load(f) - return pythia_db + particle_db = pickle.load(f) + return particle_db # -def read_sibyll_codes(filename, pythia_db): +def read_sibyll_codes(filename, particle_db): with open(filename) as f: for line in f: line = line.strip() @@ -30,19 +30,19 @@ def read_sibyll_codes(filename, pythia_db): continue identifier, sib_code, canInteractFlag, xsType = line.split() try: - pythia_db[identifier]["sibyll_code"] = int(sib_code) - pythia_db[identifier]["sibyll_canInteract"] = int(canInteractFlag) - pythia_db[identifier]["sibyll_xsType"] = int(xsType) + particle_db[identifier]["sibyll_code"] = int(sib_code) + particle_db[identifier]["sibyll_canInteract"] = int(canInteractFlag) + particle_db[identifier]["sibyll_xsType"] = int(xsType) except KeyError as e: - raise Exception("Identifier '{:s}' not found in pythia_db".format(identifier)) + raise Exception("Identifier '{:s}' not found in particle_db".format(identifier)) # generates the enum to access sibyll particles by readable names -def generate_sibyll_enum(pythia_db): +def generate_sibyll_enum(particle_db): output = "enum class SibyllCode : int8_t {\n" - for identifier, pData in pythia_db.items(): + for identifier, pData in particle_db.items(): if pData.get('sibyll_code') != None: output += " {:s} = {:d},\n".format(identifier, pData['sibyll_code']) output += "};\n" @@ -51,9 +51,9 @@ def generate_sibyll_enum(pythia_db): # generates the look-up table to convert corsika codes to sibyll codes -def generate_corsika2sibyll(pythia_db): - string = "std::array<SibyllCodeIntType, {:d}> constexpr corsika2sibyll = {{\n".format(len(pythia_db)) - for identifier, pData in pythia_db.items(): +def generate_corsika2sibyll(particle_db): + string = "std::array<SibyllCodeIntType, {:d}> constexpr corsika2sibyll = {{\n".format(len(particle_db)) + for identifier, pData in particle_db.items(): sibCode = pData.get("sibyll_code", 0) string += " {:d}, // {:s}\n".format(sibCode, identifier if sibCode else identifier + " (not implemented in SIBYLL)") string += "};\n" @@ -62,9 +62,9 @@ def generate_corsika2sibyll(pythia_db): # generates the look-up table to convert corsika codes to sibyll codes -def generate_corsika2sibyll_xsType(pythia_db): - string = "std::array<int, {:d}> constexpr corsika2sibyllXStype = {{\n".format(len(pythia_db)) - for identifier, pData in pythia_db.items(): +def generate_corsika2sibyll_xsType(particle_db): + string = "std::array<int, {:d}> constexpr corsika2sibyllXStype = {{\n".format(len(particle_db)) + for identifier, pData in particle_db.items(): sibCodeXS = pData.get("sibyll_xsType", -1) string += " {:d}, // {:s}\n".format(sibCodeXS, identifier if sibCodeXS else identifier + " (not implemented in SIBYLL)") string += "};\n" @@ -72,18 +72,18 @@ def generate_corsika2sibyll_xsType(pythia_db): # generates the look-up table to convert sibyll codes to corsika codes -def generate_sibyll2corsika(pythia_db) : +def generate_sibyll2corsika(particle_db) : string = "" minID = 0 - for identifier, pData in pythia_db.items() : + for identifier, pData in particle_db.items() : if 'sibyll_code' in pData: minID = min(minID, pData['sibyll_code']) string += "SibyllCodeIntType constexpr minSibyll = {:d};\n\n".format(minID) pDict = {} - for identifier, pData in pythia_db.items() : + for identifier, pData in particle_db.items() : if 'sibyll_code' in pData: sib_code = pData['sibyll_code'] - minID pDict[sib_code] = identifier @@ -104,13 +104,13 @@ def generate_sibyll2corsika(pythia_db) : # generates the bitset for the flag whether Sibyll knows the particle -def generate_known_particle(pythia_db): - num_particles = len(pythia_db) +def generate_known_particle(particle_db): + num_particles = len(particle_db) num_bytes = num_particles // 32 + 1 string = "Bitset2::bitset2<{:d}> constexpr isKnown{{ std::array<uint32_t, {:d}>{{{{\n".format(num_particles, num_bytes) numeric = 0 - for identifier, pData in reversed(pythia_db.items()): + for identifier, pData in reversed(particle_db.items()): handledBySibyll = int("sibyll_code" in pData) & 0x1 numeric = (numeric << 1) | handledBySibyll @@ -125,14 +125,14 @@ def generate_known_particle(pythia_db): # generates the bitset for the flag whether Sibyll can use particle as projectile -def generate_interacting_particle(pythia_db): - num_particles = len(pythia_db) +def generate_interacting_particle(particle_db): + num_particles = len(particle_db) num_bytes = num_particles // 32 + 1 string = "Bitset2::bitset2<{:d}> constexpr canInteract{{ std::array<uint32_t, {:d}>{{{{\n".format(num_particles, num_bytes) #string = "std::array<bool, {:d}> constexpr corsika2sibyll = {{\n".format(num_particles) numeric = 0 - for identifier, pData in reversed(pythia_db.items()): + for identifier, pData in reversed(particle_db.items()): can = 0 if 'sibyll_canInteract' in pData: if pData['sibyll_canInteract'] > 0: @@ -151,19 +151,19 @@ def generate_interacting_particle(pythia_db): if __name__ == "__main__": if len(sys.argv) != 3: - print("usage: {:s} <pythia_db.pkl> <sibyll_codes.dat>".format(sys.argv[0]), file=sys.stderr) + print("usage: {:s} <particle_db.pkl> <sibyll_codes.dat>".format(sys.argv[0]), file=sys.stderr) sys.exit(1) print("code_generator.py for SIBYLL") - pythia_db = load_pythiadb(sys.argv[1]) - read_sibyll_codes(sys.argv[2], pythia_db) + particle_db = load_particledb(sys.argv[1]) + read_sibyll_codes(sys.argv[2], particle_db) with open("Generated.inc", "w") as f: print("// this file is automatically generated\n// edit at your own risk!\n", file=f) - print(generate_sibyll_enum(pythia_db), file=f) - print(generate_corsika2sibyll(pythia_db), file=f) - print(generate_known_particle(pythia_db), file=f) - print(generate_sibyll2corsika(pythia_db), file=f) - print(generate_interacting_particle(pythia_db), file=f) - print(generate_corsika2sibyll_xsType(pythia_db), file=f) + print(generate_sibyll_enum(particle_db), file=f) + print(generate_corsika2sibyll(particle_db), file=f) + print(generate_known_particle(particle_db), file=f) + print(generate_sibyll2corsika(particle_db), file=f) + print(generate_interacting_particle(particle_db), file=f) + print(generate_corsika2sibyll_xsType(particle_db), file=f)