IAP GITLAB

Skip to content
Snippets Groups Projects
Commit 3e14b07b authored by Ralf M Ulrich's avatar Ralf M Ulrich
Browse files

Fixed the unknown id problem in qgsjetII

parent ee9c472e
No related branches found
No related tags found
1 merge request!186Resolve "particle list example"
Pipeline #1255 passed
......@@ -58,8 +58,7 @@ namespace corsika::process::qgsjetII {
#include <corsika/process/qgsjetII/Generated.inc>
QgsjetIICode constexpr ConvertToQgsjetII(corsika::particles::Code pCode) {
return static_cast<QgsjetIICode>(
corsika2qgsjetII[static_cast<corsika::particles::CodeIntType>(pCode)]);
return corsika2qgsjetII[static_cast<corsika::particles::CodeIntType>(pCode)];
}
corsika::particles::Code constexpr ConvertFromQgsjetII(QgsjetIICode pCode) {
......
......@@ -13,15 +13,26 @@ import pickle, sys, itertools
# loads the pickled particle_db (which is an OrderedDict)
# definition of particle_db dict is: "name", "antiName", "pdg", "mass", "electric_charge", "lifetime", "ngc_code", "isNucleus", "isHadron"
def load_particledb(filename):
'''
loads the pickled particle_db (which is an OrderedDict)
definition of particle_db dict is: "name", "antiName", "pdg", "mass", "electric_charge", "lifetime", "ngc_code", "isNucleus", "isHadron"
'''
with open(filename, "rb") as f:
particle_db = pickle.load(f)
return particle_db
def set_default_qgsjetII_definition(particle_db):
'''
Also particles not explicitly known by QGSJetII may in fact interact via mapping
to cross section types (xsType) and hadron type (hadronType)
This is achieved here.
The function return nothing, but modified the input particle_db by adding the
fields 'xsType' and 'hadronType'
'''
for identifier, pData in particle_db.items():
# the cross-section types
xsType = "CannotInteract"
......@@ -43,7 +54,7 @@ def set_default_qgsjetII_definition(particle_db):
hadronType = "PiPlusType"
else:
hadronType = "PiMinusType"
elif ((pdg>=300 and pdg<400) or pdg in [130, 10313, 10323]): # kanos
elif ((pdg>=300 and pdg<400) or pdg in [130, 10313, 10323]): # kaons
xsType = "Kaons"
if (charge>0):
hadronType = "KaonPlusType"
......@@ -69,13 +80,15 @@ def set_default_qgsjetII_definition(particle_db):
hadronType = "AntiProtonType"
# all othe not-captured cased are hopefully irrelevant
pData['qgsjetII_xsType'] = xsType
pData['qgsjetII_hadronType'] = hadronType
# read the qgsjet-codes data file
def read_qgsjetII_codes(filename, particle_db):
'''
reads the qgsjet-codes data file. For particles known to QGSJetII the 'qgsjetII_code' is set in the particle_db, as
well as the 'xsType' is updated in case it is different from its default value set above.
'''
with open(filename) as f:
for line in f:
line = line.strip()
......@@ -91,28 +104,36 @@ def read_qgsjetII_codes(filename, particle_db):
raise Exception("Identifier '{:s}' not found in particle_db".format(identifier))
# generates the enum to access qgsjetII particles by readable names
def generate_qgsjetII_enum(particle_db):
'''
generates the enum to access qgsjetII particles by readable names
'''
output = "enum class QgsjetIICode : int8_t {\n"
for identifier, pData in particle_db.items():
if pData.get('qgsjetII_code') != None:
if 'qgsjetII_code' in pData:
output += " {:s} = {:d},\n".format(identifier, pData['qgsjetII_code'])
output += "};\n"
return output
# generates the look-up table to convert corsika codes to qgsjetII codes
def generate_corsika2qgsjetII(particle_db):
string = "std::array<QgsjetIICodeIntType, {:d}> constexpr corsika2qgsjetII = {{\n".format(len(particle_db))
'''
generates the look-up table to convert corsika codes to qgsjetII codes
'''
string = "std::array<QgsjetIICode, {:d}> constexpr corsika2qgsjetII = {{\n".format(len(particle_db))
for identifier, pData in particle_db.items():
modelCode = pData.get("qgsjetII_code", 0)
string += " {:d}, // {:s}\n".format(modelCode, identifier if modelCode else identifier + " (not implemented in QGSJETII)")
if 'qgsjetII_code' in pData:
string += " QgsjetIICode::{:s}, \n".format(identifier)
else:
string += " QgsjetIICode::Unknown, // {:s}\n".format(identifier + ' not implemented in QGSJetII')
string += "};\n"
return string
# generates the look-up table to convert corsika codes to qgsjetII codes
def generate_corsika2qgsjetII_xsType(particle_db):
'''
generates the look-up table to convert corsika codes to qgsjetII codes
'''
string = "std::array<QgsjetIIXSClass, {:d}> constexpr corsika2qgsjetIIXStype = {{\n".format(len(particle_db))
for identifier, pData in particle_db.items():
modelCodeXS = pData.get("qgsjetII_xsType", "CannotInteract")
......@@ -121,8 +142,10 @@ def generate_corsika2qgsjetII_xsType(particle_db):
return string
# generates the look-up table to convert corsika codes to qgsjetII codes
def generate_corsika2qgsjetII_hadronType(particle_db):
'''
generates the look-up table to convert corsika codes to qgsjetII codes
'''
string = "std::array<QgsjetIIHadronType, {:d}> constexpr corsika2qgsjetIIHadronType = {{\n".format(len(particle_db))
for identifier, pData in particle_db.items():
modelCode = pData.get("qgsjetII_hadronType", "UndefinedType")
......@@ -131,8 +154,10 @@ def generate_corsika2qgsjetII_hadronType(particle_db):
return string
# generates the look-up table to convert qgsjetII codes to corsika codes
def generate_qgsjetII2corsika(particle_db) :
'''
generates the look-up table to convert qgsjetII codes to corsika codes
'''
minID = 0
for identifier, pData in particle_db.items() :
if 'qgsjetII_code' in pData:
......
# input file for particle conversion to/from QGSJet
# the format of this file is: "corsika-identifier" "qgsjet-id" "xs-class"
# The 'Unknown' particle is needed to mark all particles qgsjetII does
# not know
# IMPORTANT Note: the code "20" MAY NOT BE USED by qgsjetII. Change to
# another positive integer if a conflict arises. Since we use std::array
# to store PID data keep the number as low as reasonable.
Unknown 20 CannotInteract
# Note, we list here only the particles, which are produced by
# QGSJetII as secondaries. There is additional mapping of corsika
# particles on QGSJetII "xs-class"es and "hadron-type"s, which are
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment