#include <TauolappInterface.h>

Inheritance diagram for gen::TauolappInterface:

Public Member Functions
HepMC::GenEvent *	decay (HepMC::GenEvent *) override

void	disablePolarization () override

void	enablePolarization () override

void	init (const edm::EventSetup &) override

const std::vector< int > &	operatesOnParticles () override

void	SetLHE (lhef::LHEEvent *l) override

void	setRandomEngine (CLHEP::HepRandomEngine *v) override

void	statistics () override

	TauolappInterface (const edm::ParameterSet &)

	~TauolappInterface () override

Public Member Functions inherited from gen::TauolaInterfaceBase
virtual void	SetDecayRandomEngine (CLHEP::HepRandomEngine *decayRandomEngine)

	TauolaInterfaceBase ()

	TauolaInterfaceBase (const edm::ParameterSet &)

virtual	~TauolaInterfaceBase ()

Static Public Member Functions
static double	flat ()

Private Member Functions
void	BoostProdToLabLifeTimeInDecays (HepMC::GenParticle *p, TLorentzVector &lab, TLorentzVector &prod)

void	decodeMDTAU (int)

HepMC::GenParticle *	FirstTauInChain (HepMC::GenParticle *tau)

HepMC::GenParticle *	GetMother (HepMC::GenParticle *tau)

bool	isLastTauInChain (const HepMC::GenParticle *tau)

HepMC::GenEvent *	make_simple_tau_event (const TLorentzVector &l, int pdgid, int status)

double	MatchedLHESpinUp (HepMC::GenParticle *tau, std::vector< HepMC::GenParticle > &p, std::vector< double > &spinup, std::vector< int > &m_idx)

void	selectDecayByMDTAU ()

int	selectHadronic ()

int	selectLeptonic ()

void	update_particles (HepMC::GenParticle partHep, HepMC::GenEvent theEvent, HepMC::GenParticle *p, TVector3 &boost)

Private Attributes
double	dmMatch

bool	dolhe

bool	dolheBosonCorr

std::vector< int >	fHadronModes

bool	fIsInitialized

std::vector< int >	fLeptonModes

int	fMDTAU

std::vector< int >	fPDGs

edm::ESHandle< HepPDT::ParticleDataTable >	fPDGTable

bool	fPolarization

edm::ParameterSet *	fPSet

std::vector< double >	fScaledHadronBrRatios

std::vector< double >	fScaledLeptonBrRatios

bool	fSelectDecayByEvent

lhef::LHEEvent *	lhe

double	lifetime

int	ntries

Static Private Attributes
static CLHEP::HepRandomEngine *	fRandomEngine = nullptr

Additional Inherited Members
Protected Attributes inherited from gen::TauolaInterfaceBase
std::vector< int >	fPDGs

Detailed Description

Definition at line 27 of file TauolappInterface.h.

Constructor & Destructor Documentation

◆ TauolappInterface()

TauolappInterface::TauolappInterface ( const edm::ParameterSet & pset )

Definition at line 50 of file TauolappInterface.cc.

     : fPolarization(false),
       fPSet(nullptr),
       fIsInitialized(false),
       fMDTAU(-1),
       fSelectDecayByEvent(false),
       lhe(nullptr),
       dmMatch(0.5),
       dolhe(false),
       dolheBosonCorr(false),
       ntries(10) {
   if (fPSet != nullptr)
     throw cms::Exception("TauolappInterfaceError") << "Attempt to override Tauola an existing ParameterSet\n"
                                                    << std::endl;
   fPSet = new ParameterSet(pset);
 }

References Exception, fPSet, and muonDTDigis_cfi::pset.

◆ ~TauolappInterface()

TauolappInterface::~TauolappInterface ( )

override

Definition at line 67 of file TauolappInterface.cc.

                                       {
   if (fPSet != nullptr)
     delete fPSet;
 }

References fPSet.

Member Function Documentation

◆ BoostProdToLabLifeTimeInDecays()

void TauolappInterface::BoostProdToLabLifeTimeInDecays	(	HepMC::GenParticle *	p,
		TLorentzVector &	lab,
		TLorentzVector &	prod
	)

private

Definition at line 804 of file TauolappInterface.cc.

                                                                              {
   if (p->end_vertex() && p->production_vertex()) {
     HepMC::GenVertex* PGenVtx = p->production_vertex();
     HepMC::GenVertex* EGenVtx = p->end_vertex();
     double VxDec = PGenVtx->position().x() + lab.Px() / prod.Px() * (EGenVtx->position().x() - PGenVtx->position().x());
     double VyDec = PGenVtx->position().y() + lab.Py() / prod.Py() * (EGenVtx->position().y() - PGenVtx->position().y());
     double VzDec = PGenVtx->position().z() + lab.Pz() / prod.Pz() * (EGenVtx->position().z() - PGenVtx->position().z());
     double VtDec = PGenVtx->position().t() + lab.Pt() / prod.Pt() * (EGenVtx->position().t() - PGenVtx->position().t());
     EGenVtx->set_position(HepMC::FourVector(VxDec, VyDec, VzDec, VtDec));
     for (HepMC::GenVertex::particle_iterator dau = p->end_vertex()->particles_begin(HepMC::children);
          dau != p->end_vertex()->particles_end(HepMC::children);
          dau++) {
       BoostProdToLabLifeTimeInDecays((*dau), lab, prod);  //recursively modify everything in the decay chain
     }
   }
 }

References class-composition::children, gen::p, and dumpMFGeometry_cfg::prod.

Referenced by decay().

◆ decay()

HepMC::GenEvent * TauolappInterface::decay ( HepMC::GenEvent * evt )

overridevirtual

Reimplemented from gen::TauolaInterfaceBase.

Definition at line 225 of file TauolappInterface.cc.

                                                           {
   if (!fIsInitialized)
     return evt;
   Tauolapp::Tauola::setRandomGenerator(
       gen::TauolappInterface::flat);  // rest tauola++ random number incase other modules use tauola++
   int NPartBefore = evt->particles_size();
   int NVtxBefore = evt->vertices_size();
  
   // what do we do if Hep::GenEvent size is larger than 10K ???
   // Tauola (& Photos, BTW) can only handle up to 10K via HEPEVT,
   // and in case of CMS, it's only up to 4K !!!
   // override decay mode if needs be
  
   if (fSelectDecayByEvent) {
     selectDecayByMDTAU();
   }
   if (dolhe && lhe != nullptr) {
     std::vector<HepMC::GenParticle> particles;
     std::vector<int> m_idx;
     std::vector<double> spinup = lhe->getHEPEUP()->SPINUP;
     std::vector<int> pdg = lhe->getHEPEUP()->IDUP;
     for (unsigned int i = 0; i < spinup.size(); i++) {
       particles.push_back(HepMC::GenParticle(HepMC::FourVector(lhe->getHEPEUP()->PUP.at(i)[0],
                                                                lhe->getHEPEUP()->PUP.at(i)[1],
                                                                lhe->getHEPEUP()->PUP.at(i)[2],
                                                                lhe->getHEPEUP()->PUP.at(i)[3]),
                                              lhe->getHEPEUP()->IDUP.at(i)));
       int status = lhe->getHEPEUP()->ISTUP.at(i);
       particles.at(particles.size() - 1).set_generated_mass(lhe->getHEPEUP()->PUP.at(i)[4]);
       particles.at(particles.size() - 1).set_status(status > 0 ? (status == 2 ? 3 : status) : 3);
       m_idx.push_back(lhe->getHEPEUP()->MOTHUP.at(i).first - 1);  // correct for fortran index offset
     }
     // match to taus in hepmc and identify mother of taus
     bool hastaus(false);
     std::vector<HepMC::GenParticle*> match;
     for (HepMC::GenEvent::particle_const_iterator iter = evt->particles_begin(); iter != evt->particles_end(); iter++) {
       if (abs((*iter)->pdg_id()) == 15) {
         hastaus = true;
         int mother_pid(0);
         // check imediate parent to avoid parent tau ie tau->taugamma
         for (HepMC::GenVertex::particle_iterator mother = (*iter)->production_vertex()->particles_begin(HepMC::parents);
              mother != (*iter)->production_vertex()->particles_end(HepMC::parents);
              mother++) {
           mother_pid = (*mother)->pdg_id();
           if (mother_pid != (*iter)->pdg_id()) {
             // match against lhe record
             if (abs(mother_pid) == 24 ||  // W
                 abs(mother_pid) == 37 ||  // H+/-
                 abs(mother_pid) == 23 ||  // Z
                 abs(mother_pid) == 22 ||  // gamma
                 abs(mother_pid) == 25 ||  // H0 SM
                 abs(mother_pid) == 35 ||  // H0
                 abs(mother_pid) == 36     // A0
             ) {
               bool isfound = false;
               for (unsigned int k = 0; k < match.size(); k++) {
                 if ((*mother) == match.at(k))
                   isfound = true;
               }
               if (!isfound)
                 match.push_back(*mother);
             }
           }
         }
       }
     }
     if (hastaus) {
       // if is single gauge boson decay and match helicities
       if (match.size() == 1 && dolheBosonCorr) {
         for (int i = 0; i < ntries; i++) {
           // re-decay taus then check if helicities match
           auto* t_event = new Tauolapp::TauolaHepMCEvent(evt);
           t_event->undecayTaus();
           t_event->decayTaus();
           bool ismatch = true;
           for (unsigned int j = 0; j < spinup.size(); j++) {
             if (abs(pdg.at(j)) == 15) {
               double diffhelminus = (-1.0 * (double)Tauolapp::Tauola::getHelMinus() -
                                      spinup.at(j));  // -1.0 to correct for tauola feature
               double diffhelplus = ((double)Tauolapp::Tauola::getHelPlus() - spinup.at(j));
               if (pdg.at(j) == 15 && diffhelminus > 0.5)
                 ismatch = false;
               if (pdg.at(j) == -15 && diffhelplus > 0.5)
                 ismatch = false;
             }
           }
           delete t_event;
           if (ismatch)
             break;
         }
       } else {
         // If the event does not contain a single gauge boson the code will be run with
         // remove all tau decays
         auto* t_event = new Tauolapp::TauolaHepMCEvent(evt);
         t_event->undecayTaus();
         delete t_event;
         // decay all taus manually based on the helicity
         for (HepMC::GenEvent::particle_const_iterator iter = evt->particles_begin(); iter != evt->particles_end();
              iter++) {
           if (abs((*iter)->pdg_id()) == 15 && isLastTauInChain(*iter)) {
             TLorentzVector ltau(
                 (*iter)->momentum().px(), (*iter)->momentum().py(), (*iter)->momentum().pz(), (*iter)->momentum().e());
             HepMC::GenParticle* m = GetMother(*iter);
             TLorentzVector mother(m->momentum().px(), m->momentum().py(), m->momentum().pz(), m->momentum().e());
             TVector3 boost = -1.0 * mother.BoostVector();  // boost into mother's CM frame
             TLorentzVector ltau_lab = ltau;
             ltau.Boost(boost);
             mother.Boost(boost);
             HepMC::GenEvent* tauevt = make_simple_tau_event(ltau, (*iter)->pdg_id(), (*iter)->status());
             HepMC::GenParticle* p = (*(tauevt->particles_begin()));
             Tauolapp::TauolaParticle* tp = new Tauolapp::TauolaHepMCParticle(p);
             double helicity = MatchedLHESpinUp(*iter, particles, spinup, m_idx);  // get helicity from lhe
             if ((*iter)->pdg_id() == 15)
               helicity *= -1.0;
             tp->undecay();
             // use |S_{tau}|=0.999999 to avoid issues with numerical roundoff
             Tauolapp::Tauola::decayOne(tp, true, 0, 0, ((double)helicity) * 0.999999);
             boost *= -1.0;  // boost back to lab frame
             mother.Boost(boost);
             update_particles((*iter), evt, p, boost);
             //correct tau liftetime for boost (change rest frame from mothers to taus)
             BoostProdToLabLifeTimeInDecays((*iter), ltau_lab, ltau);
             delete tauevt;
           }
         }
       }
     }
   } else {
     //construct tmp TAUOLA event
     auto* t_event = new Tauolapp::TauolaHepMCEvent(evt);
     //t_event->undecayTaus();
     t_event->decayTaus();
     delete t_event;
   }
  
   for (int iv = NVtxBefore + 1; iv <= evt->vertices_size(); iv++) {
     HepMC::GenVertex* GenVtx = evt->barcode_to_vertex(-iv);
     //
     // now find decay products with funky barcode, weed out and replace with clones of sensible barcode
     // we can NOT change the barcode while iterating, because iterators do depend on the barcoding
     // thus we have to take a 2-step procedure
     //
     std::vector<int> BCodes;
     BCodes.clear();
     for (HepMC::GenVertex::particle_iterator pitr = GenVtx->particles_begin(HepMC::children);
          pitr != GenVtx->particles_end(HepMC::children);
          ++pitr) {
       if ((*pitr)->barcode() > 10000) {
         BCodes.push_back((*pitr)->barcode());
       }
     }
     if (!BCodes.empty()) {
       for (size_t ibc = 0; ibc < BCodes.size(); ibc++) {
         HepMC::GenParticle* p1 = evt->barcode_to_particle(BCodes[ibc]);
         int nbc = p1->barcode() - 10000 + NPartBefore;
         p1->suggest_barcode(nbc);
       }
     }
   }
  
   for (HepMC::GenEvent::particle_const_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
     if ((*p)->end_vertex() && (*p)->status() == 1)
       (*p)->set_status(2);
     if ((*p)->end_vertex() && (*p)->end_vertex()->particles_out_size() == 0)
       edm::LogWarning("TauolappInterface::decay error: empty end vertex!");
   }
  
   return evt;
 }

References funct::abs(), BoostProdToLabLifeTimeInDecays(), class-composition::children, dolhe, dolheBosonCorr, fIsInitialized, flat(), fSelectDecayByEvent, GenParticle::GenParticle, nanoDQM_cfi::GenVtx, lhef::LHEEvent::getHEPEUP(), GetMother(), mps_fire::i, lhef::HEPEUP::IDUP, isLastTauInChain(), lhef::HEPEUP::ISTUP, gpuVertexFinder::iv, dqmiolumiharvest::j, gen::k, lhe, visualization-live-secondInstance_cfg::m, make_simple_tau_event(), match(), MatchedLHESpinUp(), lhef::HEPEUP::MOTHUP, ntries, gen::p, p1, parents, ecalTrigSettings_cff::particles, lhef::HEPEUP::PUP, selectDecayByMDTAU(), lhef::HEPEUP::SPINUP, mps_update::status, cmsswSequenceInfo::tp, and update_particles().

◆ decodeMDTAU()

void TauolappInterface::decodeMDTAU ( int mdtau )

private

Definition at line 397 of file TauolappInterface.cc.

                                              {
   // Note-1:
   // I have to hack the common block directly because set<...>DecayMode(...)
   // only changes it in the Tauola++ instance but does NOT passes it over
   // to the Fortran core - this it does only one, via initialize() stuff...
   //
   // So I'll do both ways of settings, just for consistency...
   // but I probably need to communicate it to the Tauola(++) team...
   //
  
   // Note-2:
   // originally, the 1xx settings are meant for tau's from hard event,
   // and the 2xx settings are for any tau in the event record;
   //
   // later one, we'll have to take this into account...
   // but first I'll have to sort out what happens in the 1xx case
   // to tau's coming outside of hard event (if any in the record)
   //
  
   if (mdtau == 101 || mdtau == 201) {
     // override with electron mode for both tau's
     //
     Tauolapp::jaki_.jak1 = 1;
     Tauolapp::jaki_.jak2 = 1;
     Tauolapp::Tauola::setSameParticleDecayMode(1);
     Tauolapp::Tauola::setOppositeParticleDecayMode(1);
     return;
   }
  
   if (mdtau == 102 || mdtau == 202) {
     // override with muon mode for both tau's
     //
     Tauolapp::jaki_.jak1 = 2;
     Tauolapp::jaki_.jak2 = 2;
     Tauolapp::Tauola::setSameParticleDecayMode(2);
     Tauolapp::Tauola::setOppositeParticleDecayMode(2);
     return;
   }
  
   if (mdtau == 111 || mdtau == 211) {
     // override with electron mode for 1st tau
     // and any mode for 2nd tau
     //
     Tauolapp::jaki_.jak1 = 1;
     Tauolapp::jaki_.jak2 = 0;
     Tauolapp::Tauola::setSameParticleDecayMode(1);
     Tauolapp::Tauola::setOppositeParticleDecayMode(0);
     return;
   }
  
   if (mdtau == 112 || mdtau == 212) {
     // override with muon mode for the 1st tau
     // and any mode for the 2nd tau
     //
     Tauolapp::jaki_.jak1 = 2;
     Tauolapp::jaki_.jak2 = 0;
     Tauolapp::Tauola::setSameParticleDecayMode(2);
     Tauolapp::Tauola::setOppositeParticleDecayMode(0);
     return;
   }
  
   if (mdtau == 121 || mdtau == 221) {
     // override with any mode for the 1st tau
     // and electron mode for the 2nd tau
     //
     Tauolapp::jaki_.jak1 = 0;
     Tauolapp::jaki_.jak2 = 1;
     Tauolapp::Tauola::setSameParticleDecayMode(0);
     Tauolapp::Tauola::setOppositeParticleDecayMode(1);
     return;
   }
  
   if (mdtau == 122 || mdtau == 222) {
     // override with any mode for the 1st tau
     // and muon mode for the 2nd tau
     //
     Tauolapp::jaki_.jak1 = 0;
     Tauolapp::jaki_.jak2 = 2;
     Tauolapp::Tauola::setSameParticleDecayMode(0);
     Tauolapp::Tauola::setOppositeParticleDecayMode(2);
     return;
   }
  
   if (mdtau == 140 || mdtau == 240) {
     // override with pi+/- nutau mode for both tau's
     //
     Tauolapp::jaki_.jak1 = 3;
     Tauolapp::jaki_.jak2 = 3;
     Tauolapp::Tauola::setSameParticleDecayMode(3);
     Tauolapp::Tauola::setOppositeParticleDecayMode(3);
     return;
   }
  
   if (mdtau == 141 || mdtau == 241) {
     // override with pi+/- nutau mode for the 1st tau
     // and any mode for the 2nd tau
     //
     Tauolapp::jaki_.jak1 = 3;
     Tauolapp::jaki_.jak2 = 0;
     Tauolapp::Tauola::setSameParticleDecayMode(3);
     Tauolapp::Tauola::setOppositeParticleDecayMode(0);
     return;
   }
  
   if (mdtau == 142 || mdtau == 242) {
     // override with any mode for the 1st tau
     // and pi+/- nutau mode for 2nd tau
     //
     Tauolapp::jaki_.jak1 = 0;
     Tauolapp::jaki_.jak2 = 3;
     Tauolapp::Tauola::setSameParticleDecayMode(0);
     Tauolapp::Tauola::setOppositeParticleDecayMode(3);
     return;
   }
  
   // OK, we come here for semi-inclusive modes
   //
  
   // First of all, leptons and hadron modes sums
   //
   // re-scale branching ratios, just in case...
   //
   double sumBra = 0;
  
   // the number of decay modes is hardcoded at 22 because that's what it is right now in Tauola
   // in the future, perhaps an asscess method would be useful - communicate to Tauola team...
   //
  
   for (int i = 0; i < 22; i++) {
     sumBra += Tauolapp::taubra_.gamprt[i];
   }
   if (sumBra == 0.)
     return;  // perhaps need to throw ?
   for (int i = 0; i < 22; i++) {
     double newBra = Tauolapp::taubra_.gamprt[i] / sumBra;
     Tauolapp::Tauola::setTauBr(i + 1, newBra);
   }
   sumBra = 1.0;
  
   double sumLeptonBra = Tauolapp::taubra_.gamprt[0] + Tauolapp::taubra_.gamprt[1];
   double sumHadronBra = sumBra - sumLeptonBra;
  
   for (int i = 0; i < 2; i++) {
     fLeptonModes.push_back(i + 1);
     fScaledLeptonBrRatios.push_back((Tauolapp::taubra_.gamprt[i] / sumLeptonBra));
   }
   for (int i = 2; i < 22; i++) {
     fHadronModes.push_back(i + 1);
     fScaledHadronBrRatios.push_back((Tauolapp::taubra_.gamprt[i] / sumHadronBra));
   }
  
   fSelectDecayByEvent = true;
   return;
 }

References fHadronModes, fLeptonModes, fScaledHadronBrRatios, fScaledLeptonBrRatios, fSelectDecayByEvent, gamprt, mps_fire::i, mdtau, and taubra_.

Referenced by init().

◆ disablePolarization()

void gen::TauolappInterface::disablePolarization ( )

inlineoverridevirtual

Reimplemented from gen::TauolaInterfaceBase.

Definition at line 37 of file TauolappInterface.h.

                                         {
       fPolarization = false;
       return;
     }

References fPolarization.

◆ enablePolarization()

void gen::TauolappInterface::enablePolarization ( )

inlineoverridevirtual

Reimplemented from gen::TauolaInterfaceBase.

Definition at line 33 of file TauolappInterface.h.

                                        {
       fPolarization = true;
       return;
     }

References fPolarization.

◆ FirstTauInChain()

HepMC::GenParticle * TauolappInterface::FirstTauInChain ( HepMC::GenParticle * tau )

private

Definition at line 777 of file TauolappInterface.cc.

                                                                           {
   if (tau->production_vertex()) {
     HepMC::GenVertex::particle_iterator mother;
     for (mother = tau->production_vertex()->particles_begin(HepMC::parents);
          mother != tau->production_vertex()->particles_end(HepMC::parents);
          mother++) {
       if ((*mother)->pdg_id() == tau->pdg_id())
         return FirstTauInChain(*mother);  // recursive call to get mother with different pdgid
     }
   }
   return tau;
 }

References parents, and metsig::tau.

Referenced by MatchedLHESpinUp().

◆ flat()

double TauolappInterface::flat ( void )

static

Definition at line 215 of file TauolappInterface.cc.

                                {
   if (!fRandomEngine) {
     throw cms::Exception("LogicError")
         << "TauolaInterface::flat: Attempt to generate random number when engine pointer is null\n"
         << "This might mean that the code was modified to generate a random number outside the\n"
         << "event and beginLuminosityBlock methods, which is not allowed.\n";
   }
   return fRandomEngine->flat();
 }

References Exception, and fRandomEngine.

Referenced by decay(), init(), gen::ranmar_(), and selectLeptonic().

◆ GetMother()

HepMC::GenParticle * TauolappInterface::GetMother ( HepMC::GenParticle * tau )

private

Definition at line 790 of file TauolappInterface.cc.

                                                                     {
   if (tau->production_vertex()) {
     HepMC::GenVertex::particle_iterator mother;
     for (mother = tau->production_vertex()->particles_begin(HepMC::parents);
          mother != tau->production_vertex()->particles_end(HepMC::parents);
          mother++) {
       if ((*mother)->pdg_id() == tau->pdg_id())
         return GetMother(*mother);  // recursive call to get mother with different pdgid
       return (*mother);
     }
   }
   return tau;
 }

References parents, and metsig::tau.

Referenced by decay(), and MatchedLHESpinUp().

◆ init()

void TauolappInterface::init ( const edm::EventSetup & es )

overridevirtual

Reimplemented from gen::TauolaInterfaceBase.

Definition at line 72 of file TauolappInterface.cc.

                                                     {
   if (fIsInitialized)
     return;  // do init only once
   if (fPSet == nullptr)
     throw cms::Exception("TauolappInterfaceError") << "Attempt to initialize Tauola with an empty ParameterSet\n"
                                                    << std::endl;
  
   fIsInitialized = true;
  
   es.getData(fPDGTable);
  
   Tauolapp::Tauola::setDecayingParticle(15);
  
   // LHE Information
   dmMatch = fPSet->getUntrackedParameter<double>("dmMatch", 0.5);
   dolhe = fPSet->getUntrackedParameter<bool>("dolhe", false);
   dolheBosonCorr = fPSet->getUntrackedParameter<bool>("dolheBosonCorr", true);
   ntries = fPSet->getUntrackedParameter<int>("ntries", 10);
  
   // polarization switch
   // fPolarization = fPSet->getParameter<bool>("UseTauolaPolarization") ? 1 : 0 ;
   fPolarization = fPSet->getParameter<bool>("UseTauolaPolarization");
  
   // read tau decay mode switches
   //
   ParameterSet cards = fPSet->getParameter<ParameterSet>("InputCards");
  
   fMDTAU = cards.getParameter<int>("mdtau");
  
   if (fMDTAU == 0 || fMDTAU == 1) {
     Tauolapp::Tauola::setSameParticleDecayMode(cards.getParameter<int>("pjak1"));
     Tauolapp::Tauola::setOppositeParticleDecayMode(cards.getParameter<int>("pjak2"));
   }
  
   Tauolapp::Tauola::spin_correlation.setAll(fPolarization);
  
   // some more options, copied over from an example
   // Default values
   //Tauola::setEtaK0sPi(0,0,0); // switches to decay eta K0_S and pi0 1/0 on/off.
  
   const HepPDT::ParticleData* PData =
       fPDGTable->particle(HepPDT::ParticleID(abs(Tauolapp::Tauola::getDecayingParticle())));
   double lifetime = PData->lifetime().value();
   Tauolapp::Tauola::setTauLifetime(lifetime);
  
   fPDGs.push_back(Tauolapp::Tauola::getDecayingParticle());
  
   Tauolapp::Tauola::setRandomGenerator(gen::TauolappInterface::flat);
  
   if (fPSet->exists("parameterSets")) {
     std::vector<std::string> par = fPSet->getParameter<std::vector<std::string> >("parameterSets");
     for (unsigned int ip = 0; ip < par.size(); ++ip) {
       std::string curSet = par[ip];
       if (curSet == "setNewCurrents")
         Tauolapp::Tauola::setNewCurrents(fPSet->getParameter<int>(curSet));
     }
   }
  
   Tauolapp::Tauola::initialize();
  
   Tauolapp::Tauola::spin_correlation.setAll(
       fPolarization);  // Tauola switches this on during Tauola::initialise(); so we add this here to keep it on/off
  
   if (fPSet->exists("parameterSets")) {
     std::vector<std::string> par = fPSet->getParameter<std::vector<std::string> >("parameterSets");
     for (unsigned int ip = 0; ip < par.size(); ++ip) {
       std::string curSet = par[ip];
       if (curSet == "spinCorrelationSetAll")
         Tauolapp::Tauola::spin_correlation.setAll(fPSet->getParameter<bool>(curSet));
       if (curSet == "spinCorrelationGAMMA")
         Tauolapp::Tauola::spin_correlation.GAMMA = fPSet->getParameter<bool>(curSet);
       if (curSet == "spinCorrelationZ0")
         Tauolapp::Tauola::spin_correlation.Z0 = fPSet->getParameter<bool>(curSet);
       if (curSet == "spinCorrelationHIGGS")
         Tauolapp::Tauola::spin_correlation.HIGGS = fPSet->getParameter<bool>(curSet);
       if (curSet == "spinCorrelationHIGGSH")
         Tauolapp::Tauola::spin_correlation.HIGGS_H = fPSet->getParameter<bool>(curSet);
       if (curSet == "spinCorrelationHIGGSA")
         Tauolapp::Tauola::spin_correlation.HIGGS_A = fPSet->getParameter<bool>(curSet);
       if (curSet == "spinCorrelationHIGGSPLUS")
         Tauolapp::Tauola::spin_correlation.HIGGS_PLUS = fPSet->getParameter<bool>(curSet);
       if (curSet == "spinCorrelationHIGGSMINUS")
         Tauolapp::Tauola::spin_correlation.HIGGS_MINUS = fPSet->getParameter<bool>(curSet);
       if (curSet == "spinCorrelationWPLUS")
         Tauolapp::Tauola::spin_correlation.W_PLUS = fPSet->getParameter<bool>(curSet);
       if (curSet == "spinCorrelationWMINUS")
         Tauolapp::Tauola::spin_correlation.W_MINUS = fPSet->getParameter<bool>(curSet);
  
       if (curSet == "setHiggsScalarPseudoscalarPDG")
         Tauolapp::Tauola::setHiggsScalarPseudoscalarPDG(fPSet->getParameter<int>(curSet));
       if (curSet == "setHiggsScalarPseudoscalarMixingAngle")
         Tauolapp::Tauola::setHiggsScalarPseudoscalarMixingAngle(fPSet->getParameter<double>(curSet));
  
       if (curSet == "setRadiation")
         Tauolapp::Tauola::setRadiation(fPSet->getParameter<bool>(curSet));
       if (curSet == "setRadiationCutOff")
         Tauolapp::Tauola::setRadiationCutOff(fPSet->getParameter<double>(curSet));
  
       if (curSet == "setEtaK0sPi") {
         std::vector<int> vpar = fPSet->getParameter<std::vector<int> >(curSet);
         if (vpar.size() == 3)
           Tauolapp::Tauola::setEtaK0sPi(vpar[0], vpar[1], vpar[2]);
         else {
           std::cout << "WARNING invalid size for setEtaK0sPi: " << vpar.size() << " Require 3 elements " << std::endl;
         }
       }
  
       if (curSet == "setTaukle") {
         std::vector<double> vpar = fPSet->getParameter<std::vector<double> >(curSet);
         if (vpar.size() == 4)
           Tauolapp::Tauola::setTaukle(vpar[0], vpar[1], vpar[2], vpar[3]);
         else {
           std::cout << "WARNING invalid size for setTaukle: " << vpar.size() << " Require 4 elements " << std::endl;
         }
       }
  
       if (curSet == "setTauBr") {
         edm::ParameterSet cfg = fPSet->getParameter<edm::ParameterSet>(curSet);
         std::vector<int> vJAK = cfg.getParameter<std::vector<int> >("JAK");
         std::vector<double> vBR = cfg.getParameter<std::vector<double> >("BR");
         if (vJAK.size() == vBR.size()) {
           for (unsigned int i = 0; i < vJAK.size(); i++)
             Tauolapp::Tauola::setTauBr(vJAK[i], vBR[i]);
         } else {
           std::cout << "WARNING invalid size for setTauBr - JAK: " << vJAK.size() << " BR: " << vBR.size() << std::endl;
         }
       }
     }
   }
  
   // override decay modes if needs be
   //
   // we have to do it AFTER init because otherwises branching ratios are NOT filled in
   //
   if (fMDTAU != 0 && fMDTAU != 1) {
     decodeMDTAU(fMDTAU);
   }
  
   Tauolapp::Log::LogWarning(false);
  
   return;
 }

References funct::abs(), looper::cfg, gather_cfg::cout, decodeMDTAU(), dmMatch, dolhe, dolheBosonCorr, Exception, edm::ParameterSet::exists(), fIsInitialized, flat(), fMDTAU, fPDGs, fPDGTable, fPolarization, fPSet, edm::EventSetup::getData(), edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), mps_fire::i, initialize(), lifetime, ntries, LHEGenericFilter_cfi::ParticleID, and AlCaHLTBitMon_QueryRunRegistry::string.

◆ isLastTauInChain()

bool TauolappInterface::isLastTauInChain ( const HepMC::GenParticle * tau )

private

Definition at line 740 of file TauolappInterface.cc.

                                                                     {
   if (tau->end_vertex()) {
     HepMC::GenVertex::particle_iterator dau;
     for (dau = tau->end_vertex()->particles_begin(HepMC::children);
          dau != tau->end_vertex()->particles_end(HepMC::children);
          dau++) {
       int dau_pid = (*dau)->pdg_id();
       if (dau_pid == tau->pdg_id())
         return false;
     }
   }
   return true;
 }

References class-composition::children, and metsig::tau.

Referenced by decay().

◆ make_simple_tau_event()

HepMC::GenEvent * TauolappInterface::make_simple_tau_event	(	const TLorentzVector &	l,
		int	pdgid,
		int	status
	)

private

Definition at line 698 of file TauolappInterface.cc.

                                                                                                       {
   HepMC::GenEvent* event = new HepMC::GenEvent();
   // make tau's four vector
   HepMC::FourVector momentum_tau1(l.Px(), l.Py(), l.Pz(), l.E());
   // make particles
   HepMC::GenParticle* tau1 = new HepMC::GenParticle(momentum_tau1, pdgid, status);
   // make the vertex
   HepMC::GenVertex* vertex = new HepMC::GenVertex();
   vertex->add_particle_out(tau1);
   event->add_vertex(vertex);
   return event;
 }

References edmPickEvents::event, GenParticle::GenParticle, cmsLHEtoEOSManager::l, EgammaValidation_cff::pdgid, mps_update::status, jets_cff::tau1, and bphysicsOniaDQM_cfi::vertex.

Referenced by decay().

◆ MatchedLHESpinUp()

double TauolappInterface::MatchedLHESpinUp	(	HepMC::GenParticle *	tau,
		std::vector< HepMC::GenParticle > &	p,
		std::vector< double > &	spinup,
		std::vector< int > &	m_idx
	)

private

Definition at line 754 of file TauolappInterface.cc.

                                                                   {
   HepMC::GenParticle* Tau = FirstTauInChain(tau);
   HepMC::GenParticle* mother = GetMother(Tau);
   TLorentzVector t(tau->momentum().px(), tau->momentum().py(), tau->momentum().pz(), tau->momentum().e());
   TLorentzVector m(mother->momentum().px(), mother->momentum().py(), mother->momentum().pz(), mother->momentum().e());
   for (unsigned int i = 0; i < p.size(); i++) {
     if (tau->pdg_id() == p.at(i).pdg_id()) {
       if (mother->pdg_id() == p.at(m_idx.at(i)).pdg_id()) {
         TLorentzVector pm(p.at(m_idx.at(i)).momentum().px(),
                           p.at(m_idx.at(i)).momentum().py(),
                           p.at(m_idx.at(i)).momentum().pz(),
                           p.at(m_idx.at(i)).momentum().e());
         if (fabs(m.M() - pm.M()) < dmMatch)
           return spinup.at(i);
       }
     }
   }
   return 0;
 }

References dmMatch, FirstTauInChain(), GenParticle::GenParticle, GetMother(), mps_fire::i, visualization-live-secondInstance_cfg::m, gen::p, HiggsValidation_cfi::pdg_id, submitPVValidationJobs::t, and metsig::tau.

Referenced by decay().

◆ operatesOnParticles()

const std::vector<int>& gen::TauolappInterface::operatesOnParticles ( )

inlineoverridevirtual

Reimplemented from gen::TauolaInterfaceBase.

Definition at line 42 of file TauolappInterface.h.

42 { return fPDGs; }

References fPDGs.

◆ selectDecayByMDTAU()

void TauolappInterface::selectDecayByMDTAU ( )

private

Definition at line 552 of file TauolappInterface.cc.

                                            {
   if (fMDTAU == 100 || fMDTAU == 200) {
     int mode = selectLeptonic();
     Tauolapp::jaki_.jak1 = mode;
     Tauolapp::Tauola::setSameParticleDecayMode(mode);
     mode = selectLeptonic();
     Tauolapp::jaki_.jak2 = mode;
     Tauolapp::Tauola::setOppositeParticleDecayMode(mode);
     return;
   }
  
   int modeL = selectLeptonic();
   int modeH = selectHadronic();
  
   if (fMDTAU == 110 || fMDTAU == 210) {
     Tauolapp::jaki_.jak1 = modeL;
     Tauolapp::jaki_.jak2 = 0;
     Tauolapp::Tauola::setSameParticleDecayMode(modeL);
     Tauolapp::Tauola::setOppositeParticleDecayMode(0);
     return;
   }
  
   if (fMDTAU == 120 || fMDTAU == 22) {
     Tauolapp::jaki_.jak1 = 0;
     Tauolapp::jaki_.jak2 = modeL;
     Tauolapp::Tauola::setSameParticleDecayMode(0);
     Tauolapp::Tauola::setOppositeParticleDecayMode(modeL);
     return;
   }
  
   if (fMDTAU == 114 || fMDTAU == 214) {
     Tauolapp::jaki_.jak1 = modeL;
     Tauolapp::jaki_.jak2 = modeH;
     Tauolapp::Tauola::setSameParticleDecayMode(modeL);
     Tauolapp::Tauola::setOppositeParticleDecayMode(modeH);
     return;
   }
  
   if (fMDTAU == 124 || fMDTAU == 224) {
     Tauolapp::jaki_.jak1 = modeH;
     Tauolapp::jaki_.jak2 = modeL;
     Tauolapp::Tauola::setSameParticleDecayMode(modeH);
     Tauolapp::Tauola::setOppositeParticleDecayMode(modeL);
     return;
   }
  
   if (fMDTAU == 115 || fMDTAU == 215) {
     Tauolapp::jaki_.jak1 = 1;
     Tauolapp::jaki_.jak2 = modeH;
     Tauolapp::Tauola::setSameParticleDecayMode(1);
     Tauolapp::Tauola::setOppositeParticleDecayMode(modeH);
     return;
   }
  
   if (fMDTAU == 125 || fMDTAU == 225) {
     Tauolapp::jaki_.jak1 = modeH;
     Tauolapp::jaki_.jak2 = 1;
     Tauolapp::Tauola::setSameParticleDecayMode(modeH);
     Tauolapp::Tauola::setOppositeParticleDecayMode(1);
     return;
   }
  
   if (fMDTAU == 116 || fMDTAU == 216) {
     Tauolapp::jaki_.jak1 = 2;
     Tauolapp::jaki_.jak2 = modeH;
     Tauolapp::Tauola::setSameParticleDecayMode(2);
     Tauolapp::Tauola::setOppositeParticleDecayMode(modeH);
     return;
   }
  
   if (fMDTAU == 126 || fMDTAU == 226) {
     Tauolapp::jaki_.jak1 = modeH;
     Tauolapp::jaki_.jak2 = 2;
     Tauolapp::Tauola::setSameParticleDecayMode(modeH);
     Tauolapp::Tauola::setOppositeParticleDecayMode(2);
     return;
   }
  
   if (fMDTAU == 130 || fMDTAU == 230) {
     Tauolapp::jaki_.jak1 = modeH;
     Tauolapp::jaki_.jak2 = selectHadronic();
     Tauolapp::Tauola::setSameParticleDecayMode(modeH);
     Tauolapp::Tauola::setOppositeParticleDecayMode(Tauolapp::jaki_.jak2);
     return;
   }
  
   if (fMDTAU == 131 || fMDTAU == 231) {
     Tauolapp::jaki_.jak1 = modeH;
     Tauolapp::jaki_.jak2 = 0;
     Tauolapp::Tauola::setSameParticleDecayMode(modeH);
     Tauolapp::Tauola::setOppositeParticleDecayMode(0);
     return;
   }
  
   if (fMDTAU == 132 || fMDTAU == 232) {
     Tauolapp::jaki_.jak1 = 0;
     Tauolapp::jaki_.jak2 = modeH;
     Tauolapp::Tauola::setSameParticleDecayMode(0);
     Tauolapp::Tauola::setOppositeParticleDecayMode(modeH);
     return;
   }
  
   // unlikely that we get here on unknown mdtau
   // - there's a protection earlier
   // but if we do, just set defaults
   // probably need to spit a warning...
   //
   Tauolapp::Tauola::setSameParticleDecayMode(0);
   Tauolapp::Tauola::setOppositeParticleDecayMode(0);
  
   return;
 }

References fMDTAU, jak2, ALCARECOPromptCalibProdSiPixelAli0T_cff::mode, selectHadronic(), and selectLeptonic().

Referenced by decay().

◆ selectHadronic()

int TauolappInterface::selectHadronic ( )

private

Definition at line 677 of file TauolappInterface.cc.

                                       {
   float prob = 0.;
   int len = 1;
   ranmar_(&prob, &len);
  
   double sumBra = fScaledHadronBrRatios[0];
   if (prob > 0. && prob <= sumBra) {
     return fHadronModes[0];
   } else {
     int NN = fScaledHadronBrRatios.size();
     for (int i = 1; i < NN; i++) {
       if (prob > sumBra && prob <= (sumBra + fScaledHadronBrRatios[i])) {
         return fHadronModes[i];
       }
       sumBra += fScaledHadronBrRatios[i];
     }
   }
  
   return 0;
 }

References fHadronModes, fScaledHadronBrRatios, mps_fire::i, TtFullHadEvtBuilder_cfi::prob, and gen::ranmar_().

Referenced by selectDecayByMDTAU().

◆ selectLeptonic()

int TauolappInterface::selectLeptonic ( )

private

Definition at line 665 of file TauolappInterface.cc.

                                       {
   float prob = flat();
  
   if (prob > 0. && prob <= fScaledLeptonBrRatios[0]) {
     return 1;
   } else if (prob > fScaledLeptonBrRatios[1] && prob <= 1.) {
     return 2;
   }
  
   return 0;
 }

References flat(), fScaledLeptonBrRatios, and TtFullHadEvtBuilder_cfi::prob.

Referenced by selectDecayByMDTAU().

◆ SetLHE()

void gen::TauolappInterface::SetLHE ( lhef::LHEEvent * l )

inlineoverridevirtual

Reimplemented from gen::TauolaInterfaceBase.

Definition at line 45 of file TauolappInterface.h.

45 { lhe = l; }

References cmsLHEtoEOSManager::l, and lhe.

◆ setRandomEngine()

void gen::TauolappInterface::setRandomEngine ( CLHEP::HepRandomEngine * v )

inlineoverridevirtual

Implements gen::TauolaInterfaceBase.

Definition at line 46 of file TauolappInterface.h.

46 { fRandomEngine = v; }

References fRandomEngine, and gen::v.

◆ statistics()

void TauolappInterface::statistics ( )

overridevirtual

Reimplemented from gen::TauolaInterfaceBase.

Definition at line 395 of file TauolappInterface.cc.

395 { return; }

◆ update_particles()

void TauolappInterface::update_particles	(	HepMC::GenParticle *	partHep,
		HepMC::GenEvent *	theEvent,
		HepMC::GenParticle *	p,
		TVector3 &	boost
	)

private

Definition at line 711 of file TauolappInterface.cc.

                                                           {
   partHep->set_status(p->status());
   if (p->end_vertex()) {
     if (!partHep->end_vertex()) {
       HepMC::GenVertex* vtx = new HepMC::GenVertex(p->end_vertex()->position());
       theEvent->add_vertex(vtx);
       vtx->add_particle_in(partHep);
     }
     if (p->end_vertex()->particles_out_size() != 0) {
       for (HepMC::GenVertex::particles_out_const_iterator d = p->end_vertex()->particles_out_const_begin();
            d != p->end_vertex()->particles_out_const_end();
            d++) {
         // Create daughter and add to event
         TLorentzVector l((*d)->momentum().px(), (*d)->momentum().py(), (*d)->momentum().pz(), (*d)->momentum().e());
         l.Boost(boost);
         HepMC::FourVector momentum(l.Px(), l.Py(), l.Pz(), l.E());
         HepMC::GenParticle* daughter = new HepMC::GenParticle(momentum, (*d)->pdg_id(), (*d)->status());
         daughter->suggest_barcode(theEvent->particles_size() + 1);
         partHep->end_vertex()->add_particle_out(daughter);
         if ((*d)->end_vertex())
           update_particles(daughter, theEvent, (*d), boost);
       }
     }
   }
 }

References ztail::d, GenParticle::GenParticle, cmsLHEtoEOSManager::l, gen::p, and extraflags_cff::vtx.

Referenced by decay().

Member Data Documentation

◆ dmMatch

double gen::TauolappInterface::dmMatch

private

Definition at line 86 of file TauolappInterface.h.

Referenced by init(), and MatchedLHESpinUp().

◆ dolhe

bool gen::TauolappInterface::dolhe

private

Definition at line 87 of file TauolappInterface.h.

Referenced by decay(), and init().

◆ dolheBosonCorr

bool gen::TauolappInterface::dolheBosonCorr

private

Definition at line 88 of file TauolappInterface.h.

Referenced by decay(), and init().

◆ fHadronModes

std::vector<int> gen::TauolappInterface::fHadronModes

private

Definition at line 81 of file TauolappInterface.h.

Referenced by decodeMDTAU(), and selectHadronic().

◆ fIsInitialized

bool gen::TauolappInterface::fIsInitialized

private

Definition at line 76 of file TauolappInterface.h.

Referenced by decay(), and init().

◆ fLeptonModes

std::vector<int> gen::TauolappInterface::fLeptonModes

private

Definition at line 80 of file TauolappInterface.h.

Referenced by decodeMDTAU().

◆ fMDTAU

int gen::TauolappInterface::fMDTAU

private

Definition at line 78 of file TauolappInterface.h.

Referenced by init(), and selectDecayByMDTAU().

◆ fPDGs

std::vector<int> gen::TauolappInterface::fPDGs

private

Definition at line 72 of file TauolappInterface.h.

Referenced by init(), and operatesOnParticles().

◆ fPDGTable

edm::ESHandle<HepPDT::ParticleDataTable> gen::TauolappInterface::fPDGTable

private

Definition at line 74 of file TauolappInterface.h.

Referenced by init().

◆ fPolarization

bool gen::TauolappInterface::fPolarization

private

Definition at line 73 of file TauolappInterface.h.

Referenced by disablePolarization(), enablePolarization(), and init().

◆ fPSet

edm::ParameterSet* gen::TauolappInterface::fPSet

private

Definition at line 75 of file TauolappInterface.h.

Referenced by init(), TauolappInterface(), and ~TauolappInterface().

◆ fRandomEngine

CLHEP::HepRandomEngine * TauolappInterface::fRandomEngine = nullptr

staticprivate

Definition at line 71 of file TauolappInterface.h.

Referenced by flat(), and setRandomEngine().

◆ fScaledHadronBrRatios

std::vector<double> gen::TauolappInterface::fScaledHadronBrRatios

private

Definition at line 83 of file TauolappInterface.h.

Referenced by decodeMDTAU(), and selectHadronic().

◆ fScaledLeptonBrRatios

std::vector<double> gen::TauolappInterface::fScaledLeptonBrRatios

private

Definition at line 82 of file TauolappInterface.h.

Referenced by decodeMDTAU(), and selectLeptonic().

◆ fSelectDecayByEvent

bool gen::TauolappInterface::fSelectDecayByEvent

private

Definition at line 79 of file TauolappInterface.h.

Referenced by decay(), and decodeMDTAU().

◆ lhe

lhef::LHEEvent* gen::TauolappInterface::lhe

private

Definition at line 84 of file TauolappInterface.h.

Referenced by decay(), and SetLHE().

◆ lifetime

double gen::TauolappInterface::lifetime

private

Definition at line 90 of file TauolappInterface.h.

Referenced by init().

◆ ntries

int gen::TauolappInterface::ntries

private

Definition at line 89 of file TauolappInterface.h.

Referenced by decay(), and init().

Public Member Functions

Static Public Member Functions

Private Member Functions

Private Attributes

Static Private Attributes

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

◆ TauolappInterface()

◆ ~TauolappInterface()

Member Function Documentation

◆ BoostProdToLabLifeTimeInDecays()

◆ decay()

◆ decodeMDTAU()

◆ disablePolarization()

◆ enablePolarization()

◆ FirstTauInChain()

◆ flat()

◆ GetMother()

◆ init()

◆ isLastTauInChain()

◆ make_simple_tau_event()

◆ MatchedLHESpinUp()

◆ operatesOnParticles()

◆ selectDecayByMDTAU()

◆ selectHadronic()

◆ selectLeptonic()

◆ SetLHE()

◆ setRandomEngine()

◆ statistics()

◆ update_particles()

Member Data Documentation

◆ dmMatch

◆ dolhe

◆ dolheBosonCorr

◆ fHadronModes

◆ fIsInitialized

◆ fLeptonModes

◆ fMDTAU

◆ fPDGs

◆ fPDGTable

◆ fPolarization

◆ fPSet

◆ fRandomEngine

◆ fScaledHadronBrRatios

◆ fScaledLeptonBrRatios

◆ fSelectDecayByEvent

◆ lhe

◆ lifetime

◆ ntries