EventGenerator/src/TauDecay_GenParticle.cc

0001 #include "Validation/EventGenerator/interface/TauDecay_GenParticle.h"
0002 #include "Validation/EventGenerator/interface/PdtPdgMini.h"
0003
0004 #include <iomanip>
0005 #include <cstdlib>
0006
0007 TauDecay_GenParticle::TauDecay_GenParticle() : TauDecay() {}
0008
0009 TauDecay_GenParticle::~TauDecay_GenParticle() {}
0010
0011 bool TauDecay_GenParticle::AnalyzeTau(
0012     const reco::GenParticle *Tau, unsigned int &MODE_ID, unsigned int &TauBitMask, bool dores, bool dopi0) {
0013   Reset();
0014   MotherIdx.clear();
0015   TauDecayProducts.clear();
0016   if (abs(Tau->pdgId()) == PdtPdgMini::tau_minus) {  // check that it is a tau
0017     unsigned int Tauidx = TauDecayProducts.size();
0018     TauDecayProducts.push_back(Tau);
0019     MotherIdx.push_back(Tauidx);
0020     for (unsigned int i = 0; i < Tau->numberOfDaughters(); i++) {
0021       const reco::Candidate *dau = Tau->daughter(i);
0022       Analyze(static_cast<const reco::GenParticle *>(dau), Tauidx, dores, dopi0);
0023     }
0024     ClassifyDecayMode(MODE_ID, TauBitMask);
0025     return true;
0026   }
0027   return false;
0028 }
0029
0030 void TauDecay_GenParticle::Analyze(const reco::GenParticle *Particle, unsigned int midx, bool dores, bool dopi0) {
0031   unsigned int pdgid = abs(Particle->pdgId());
0032   if (isTauFinalStateParticle(pdgid)) {
0033     if (!isTauParticleCounter(pdgid))
0034       std::cout << "TauDecay_GenParticle::Analyze WARNING: Unknow Final State Particle in Tau Decay... " << std::endl;
0035
0036     if (!AddRadInfo(static_cast<const reco::GenParticle *>(Particle), midx)) {
0037       TauDecayProducts.push_back(Particle);
0038       MotherIdx.push_back(midx);
0039     }
0040     // store information on pi0 decay products even though a pi0 is a finsal state particle (for 3PiPi0 studies)
0041     if (pdgid == PdtPdgMini::pi0 && dopi0) {
0042       midx = MotherIdx.size() - 1;
0043       for (unsigned int i = 0; i < Particle->numberOfDaughters(); i++) {
0044         const reco::Candidate *dau = Particle->daughter(i);
0045         AddPi0Info(static_cast<const reco::GenParticle *>(dau), midx);
0046       }
0047     }
0048     return;
0049   }
0050   if (Particle->status() == 1 || (isTauResonanceCounter(pdgid) && dores)) {
0051     TauDecayProducts.push_back(Particle);
0052     MotherIdx.push_back(midx);
0053     midx = MotherIdx.size() - 1;
0054   }
0055   for (unsigned int i = 0; i < Particle->numberOfDaughters(); i++) {
0056     const reco::Candidate *dau = Particle->daughter(i);
0057     Analyze(static_cast<const reco::GenParticle *>(dau), midx, dores, dopi0);
0058   }
0059 }
0060
0061 void TauDecay_GenParticle::AddPi0Info(const reco::GenParticle *Particle, unsigned int midx) {
0062   if (Particle->status() == 1) {
0063     TauDecayProducts.push_back(Particle);
0064     MotherIdx.push_back(midx);
0065     return;
0066   }
0067   for (unsigned int i = 0; i < Particle->numberOfDaughters(); i++) {
0068     const reco::Candidate *dau = Particle->daughter(i);
0069     AddPi0Info(static_cast<const reco::GenParticle *>(dau), midx);
0070   }
0071 }
0072
0073 bool TauDecay_GenParticle::AddRadInfo(const reco::GenParticle *Particle, unsigned int midx) {
0074   bool same(false), photon(false);
0075   int n = 0;
0076   for (unsigned int i = 0; i < Particle->numberOfDaughters(); i++) {
0077     const reco::Candidate *dau = Particle->daughter(i);
0078     if (Particle->pdgId() == dau->pdgId()) {
0079       same = true;
0080     }
0081     if (dau->pdgId() == PdtPdgMini::gamma)
0082       photon = true;
0083     n++;
0084   }
0085   if (same && photon && n == 2) {
0086     for (unsigned int i = 0; i < Particle->numberOfDaughters(); i++) {
0087       const reco::Candidate *dau = Particle->daughter(i);
0088       if (Particle->pdgId() == dau->pdgId()) {
0089         if (!AddRadInfo(static_cast<const reco::GenParticle *>(dau), midx)) {
0090           TauDecayProducts.push_back(static_cast<const reco::GenParticle *>(dau));
0091           MotherIdx.push_back(midx);
0092         }
0093       }
0094       if (dau->pdgId() == PdtPdgMini::gamma) {
0095         TauDecayProducts.push_back(static_cast<const reco::GenParticle *>(dau));
0096         MotherIdx.push_back(midx);
0097       }
0098     }
0099     return true;
0100   }
0101   return false;
0102 }