// pbarpEnv class definition file. -*- C++ -*- // Copyright 2012 Bertram Kopf #include <getopt.h> #include <fstream> #include "pbarpUtils/pbarpEnv.hh" #include "PwaUtils/IsobarDecay.hh" #include "PwaUtils/IsobarDecayList.hh" #include "pbarpUtils/pbarpReaction.hh" #include "pbarpUtils/pbarpEventList.hh" #include "qft++/relativistic-quantum-mechanics/Utils.hh" #include "ErrLogger/ErrLogger.hh" #include "Particle/Particle.hh" #include "Particle/ParticleTable.hh" #include "Particle/PdtParser.hh" pbarpEnv* pbarpEnv::_instance=0; pbarpEnv* pbarpEnv::instance() { if (0==_instance) _instance = new pbarpEnv(); return _instance; } pbarpEnv::pbarpEnv() : _alreadySetUp(false) ,_lmax(0) ,_noFinalStateParticles(0) ,_decList(new IsobarDecayList()) ,_prodDecList(new IsobarDecayList()) { } pbarpEnv::~pbarpEnv(){ } void pbarpEnv::setup(pbarpParser& thePbarpParser){ if(_alreadySetUp){ Alert << " pbarpEnv already set up" << endmsg; exit(1); } _alreadySetUp=true; // common options (move to base class later) _outputFileNameSuffix = thePbarpParser.outputFileNameSuffix(); // pdtTable PdtParser pdtParser; std::string theSourcePath=getenv("CMAKE_SOURCE_DIR"); std::string pdtFileRelPath=thePbarpParser.pdgTableFile(); std::string pdtFile(theSourcePath+pdtFileRelPath); _particleTable = new ParticleTable; if (!pdtParser.parse(pdtFile, *_particleTable)) { Alert << "can not parse particle table " << pdtFile << endmsg; exit(1); } //Lmax _lmax=thePbarpParser.getLMax(); //final state particles const std::vector<std::string> finalStateParticleStr=thePbarpParser.finalStateParticles(); std::vector<std::string>::const_iterator itStr; for ( itStr = finalStateParticleStr.begin(); itStr != finalStateParticleStr.end(); ++itStr){ Particle* currentParticle = _particleTable->particle(*itStr); _finalStateParticles.push_back(currentParticle); } _noFinalStateParticles= (int) _finalStateParticles.size(); //decays std::vector<std::string> decaySystem= thePbarpParser.decaySystem(); for ( itStr = decaySystem.begin(); itStr != decaySystem.end(); ++itStr){ std::stringstream stringStr; stringStr << (*itStr); std::string motherStr; stringStr >> motherStr; Particle* motherParticle = _particleTable->particle(motherStr); if( 0==motherParticle){ Alert << "mother particle\t" << motherStr << "\tdoes not exist in pdtTable" << endmsg; exit(1); } std::string daughter1Str; stringStr >> daughter1Str; Particle* daughter1Particle = _particleTable->particle(daughter1Str); if( 0==daughter1Particle){ Alert << "first daughter particle\t" << daughter1Str << "\tdoes not exist in pdtTable" << endmsg; exit(1); } std::string daughter2Str; stringStr >> daughter2Str; Particle* daughter2Particle = _particleTable->particle(daughter2Str); if( 0==daughter2Particle){ Alert << "second daughter particle\t" << daughter2Str << "\tdoes not exist in pdtTable" << endmsg; exit(1); } boost::shared_ptr<IsobarDecay> tmpDec(new IsobarDecay(motherParticle, daughter1Particle, daughter2Particle)); _decList->addDecay(tmpDec); } //add dynamics std::vector<boost::shared_ptr<IsobarDecay> > isoDecList= _decList->getList(); std::vector<std::string> decDynVec = thePbarpParser.decayDynamics(); for ( itStr = decDynVec.begin(); itStr != decDynVec.end(); ++itStr){ std::stringstream stringStr; stringStr << (*itStr); std::string particleStr; stringStr >> particleStr; std::string dynStr; stringStr >> dynStr; std::vector<boost::shared_ptr<IsobarDecay> >::iterator itDyn; for (itDyn=isoDecList.begin(); itDyn!=isoDecList.end(); ++itDyn){ std::string theDecName=(*itDyn)->name(); std::string toFind=particleStr+"To"; size_t found; found=theDecName.find(toFind); if (found!=string::npos){ (*itDyn)->enableDynamics(dynStr); } } } //produced particle pairs std::vector<std::string> productionSystem = thePbarpParser.productionSystem(); for ( itStr = productionSystem.begin(); itStr != productionSystem.end(); ++itStr){ std::stringstream stringStr; stringStr << (*itStr); std::string firstParticleStr; stringStr >> firstParticleStr; std::cout << "first particle:\t" << firstParticleStr << std::endl; std::string secondParticleStr; stringStr >> secondParticleStr; std::cout << "second particle:\t" << secondParticleStr << std::endl; Particle* firstParticle = _particleTable->particle(firstParticleStr); if( 0==firstParticle){ Alert << "particle\t" << firstParticleStr << "\tdoes not exist in pdtTable" << endmsg; exit(1); } Particle* secondParticle = _particleTable->particle(secondParticleStr); if( 0==secondParticle){ Alert << "particle\t" << secondParticleStr << "\tdoes not exist in pdtTable" << endmsg; exit(1); } _producedParticlePairs.push_back(make_pair(firstParticle, secondParticle)); } _pbarpReaction=boost::shared_ptr<pbarpReaction>(new pbarpReaction(_producedParticlePairs, _lmax)); //fill prodDecayList std::vector< boost::shared_ptr<IsobarDecay> > prodDecs= _pbarpReaction->productionDecays(); std::vector< boost::shared_ptr<IsobarDecay> >::iterator itDec; for (itDec=prodDecs.begin(); itDec!=prodDecs.end(); ++itDec){ _prodDecList->addDecay(*itDec); } //set suffixes std::vector<std::string> suffixVec = thePbarpParser.replaceSuffixNames(); std::map<std::string, std::string> decSuffixNames; for ( itStr = suffixVec.begin(); itStr != suffixVec.end(); ++itStr){ std::stringstream stringStr; stringStr << (*itStr); std::string classStr; stringStr >> classStr; std::string suffixStr; stringStr >> suffixStr; decSuffixNames[classStr]=suffixStr; } //set suffixes for decay classes std::map<std::string, std::string>::iterator itMapStrStr; for (itMapStrStr=decSuffixNames.begin(); itMapStrStr!=decSuffixNames.end(); ++itMapStrStr){ _decList->replaceSuffix(itMapStrStr->first, itMapStrStr->second); _prodDecList->replaceSuffix(itMapStrStr->first, itMapStrStr->second); boost::shared_ptr<IsobarDecay> theDec=_decList->decay(itMapStrStr->first); } //replace mass key std::vector<std::string> replMassKeyVec = thePbarpParser.replaceMassKey(); std::map<std::string, std::string> decRepMassKeyNames; for ( itStr = replMassKeyVec.begin(); itStr != replMassKeyVec.end(); ++itStr){ std::stringstream stringStr; stringStr << (*itStr); std::string oldStr; stringStr >> oldStr; std::string newStr; stringStr >> newStr; decRepMassKeyNames[oldStr]=newStr; } for (itMapStrStr=decRepMassKeyNames.begin(); itMapStrStr!=decRepMassKeyNames.end(); ++itMapStrStr){ _decList->replaceMassKey(itMapStrStr->first, itMapStrStr->second); } std::vector<std::string> theHistMassNames=thePbarpParser.histMassNames(); //fill vector histMassSystems for ( itStr = theHistMassNames.begin(); itStr != theHistMassNames.end(); ++itStr){ std::stringstream stringStr; stringStr << (*itStr); std::string tmpName; std::vector<std::string> currentStringVec; while(stringStr >> tmpName){ currentStringVec.push_back(tmpName); } _histMassSystems.push_back(currentStringVec); } std::vector<std::string> theHistAngleNames=thePbarpParser.histAngleNames(); //fill vector histMassSystems for ( itStr = theHistAngleNames.begin(); itStr != theHistAngleNames.end(); ++itStr){ std::stringstream stringStr; stringStr << (*itStr); std::string tmpName; std::vector<std::string> currentStringDecVec; std::vector<std::string> currentStringMotherVec; bool isDecParticle=true; while(stringStr >> tmpName){ if(tmpName=="from") { isDecParticle=false; continue; } if(isDecParticle) currentStringDecVec.push_back(tmpName); else currentStringMotherVec.push_back(tmpName); } boost::shared_ptr<angleHistData> currentAngleHistData(new angleHistData(currentStringMotherVec, currentStringDecVec)); _angleHistDataVec.push_back(currentAngleHistData); } }