#ifndef __JAANET__JEVTWEIGHTTOOLKIT__ #define __JAANET__JEVTWEIGHTTOOLKIT__ #include #include #include "km3net-dataformat/offline/Head.hh" #include "km3net-dataformat/offline/Evt.hh" #include "JLang/JException.hh" #include "JLang/JSharedPointer.hh" #include "Jeep/JProperties.hh" #include "JAAnet/JHead.hh" #include "JAAnet/JHeadToolkit.hh" #include "JAAnet/JEvtWeight.hh" #include "JAAnet/JEvtWeightDAQ.hh" #include "JAAnet/JEvtWeightMupage.hh" #include "JAAnet/JEvtWeightCorsika.hh" #include "JAAnet/JEvtWeightGSeaGen.hh" #include "JAAnet/JEvtWeightKM3BUU.hh" #include "JAAnet/JEvtWeightMiscellaneous.hh" #include "JAAnet/JFlatFlux.hh" #include "JAAnet/JPowerLawFlux.hh" #include "JAAnet/JAtmosphericNeutrinoFlux.hh" #include "JAAnet/JEvtWeightFactorMultiParticle.hh" #include "JAAnet/JParticleTypes.hh" #include "JAAnet/JAAnetToolkit.hh" /** * \author bjung */ namespace JAANET {} namespace JPP { using namespace JAANET; } namespace JAANET { /** * Look-up table for event weighters. */ struct JEvtWeighter : public std::vector< JLANG::JSharedPointer > { /** * Constructor */ JEvtWeighter() { this->push_back(new JEvtWeightGSeaGen()); this->push_back(new JEvtWeightKM3BUU()); this->push_back(new JEvtWeightCorsika()); this->push_back(new JEvtWeightMupage()); this->push_back(new JEvtWeightDAQ()); this->push_back(new JEvtWeightMiscellaneous()); } /** * Get event weighter corresponding to given header. * * \param header header * \return event weighter */ const JEvtWeight& operator()(const JHead& header) const { for (const_iterator i = this->begin(); i != this->end(); ++i) { if ((*i)->check(header)) { return *(*i); } } THROW(JValueOutOfRange, "JEvtWeighter::operator(): No event weighter found for given header."); } }; /** * Auxiliary class for parsing a vector of neutrino PDG identifiers. */ struct JNeutrinoTypeCollection : std::vector { /** * Default constructor. */ JNeutrinoTypeCollection() {} /** * Stream input. * * \param in input stream * \param collection collection of neutrino PDG types * \return input stream */ friend inline std::istream& operator>>(std::istream& in, JNeutrinoTypeCollection& collection) { for (int type; in >> type; ) { if (abs(type) == TRACK_TYPE_NUE || abs(type) == TRACK_TYPE_NUMU || abs(type) == TRACK_TYPE_NUTAU) { collection.push_back(type); } else { THROW(JValueOutOfRange, "JNeutrinoTypeCollection::operator>>(): Invalid particle type: " << type); } } return in; } /** * Stream output. * * \param out output stream * \param collection collection of neutrino PDG types * \return output stream */ friend inline std::ostream& operator<<(std::ostream& out, const JNeutrinoTypeCollection& collection) { for (std::vector::const_iterator i = collection.cbegin(); i != collection.cend(); ++i) { out << ' ' << *i; } return out; } }; /** * Auxiliary class for parsing multiparticle fluxes. */ struct JFluxMapParser : public JEEP::JProperties, public JOscProbHelper { /** * Constructor. */ JFluxMapParser() { this->put("flat", flatFluxes); this->put("powerlaw", powerLawFluxes); this->put("atmospheric", atmosphericFluxes); } /** * Get multiparticle flux function. * * \return multiparticle flux function */ JFluxMultiParticle getMultiParticleFlux() const { using namespace std; using namespace JPP; JFluxMultiParticle multiFlux; for (map::const_iterator i = flatFluxes.cbegin(); i != flatFluxes.cend(); ++i) { multiFlux.insert(i->first, make_fluxFunction(i->second)); } for (map::const_iterator i = powerLawFluxes.cbegin(); i != powerLawFluxes.cend(); ++i) { multiFlux.insert(i->first, make_fluxFunction(i->second)); } if (is_valid()) { const JOscProb& oscProbInterface = *(this->get()); const JAtmosphericNeutrinoFlux atmFlux(oscProbInterface); const JEvtWeightFactorFunction atmFluxFunction = make_fluxFunction(atmFlux); for (vector::const_iterator i = atmosphericFluxes.cbegin(); i != atmosphericFluxes.cend(); ++i) { multiFlux.insert(*i, atmFluxFunction); } } return multiFlux; } /** * Conversion operator. * * \return multiparticle flux function */ operator JFluxMultiParticle() const { return getMultiParticleFlux(); } /** * Stream input. * * \param in input stream * \param fluxMapParser flux map parser * \return input stream */ friend inline std::istream& operator>>(std::istream& in, JFluxMapParser& fluxMapParser) { using namespace std; using namespace JPP; JStringStream is(in); if (getFileStatus(is.str().c_str())) { is.load(); } is >> static_cast(fluxMapParser); return in; } /** * Stream output. * * \param out output stream * \param fluxMapParser flux map parser * \return output stream */ friend inline std::ostream& operator<<(std::ostream& out, const JFluxMapParser& fluxMapParser) { return out << static_cast(fluxMapParser); } std::map flatFluxes; //!< Uniform flux functions std::map powerLawFluxes; //!< Power-law flux functions JNeutrinoTypeCollection atmosphericFluxes; //!< Atmospheric neutrino flux functions }; extern JEvtWeighter getEventWeighter; //!< Function object for mapping header to event weighter. } #endif