SpecificDecay/src/BPHDecayToTkpTknSymChargeBuilder.cc

0001 /*
0002  *  See header file for a description of this class.
0003  *
0004  *  \author Paolo Ronchese INFN Padova
0005  *
0006  */
0007
0008 //-----------------------
0009 // This Class' Header --
0010 //-----------------------
0011 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHDecayToTkpTknSymChargeBuilder.h"
0012
0013 //-------------------------------
0014 // Collaborating Class Headers --
0015 //-------------------------------
0016 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHParticlePtSelect.h"
0017 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHParticleEtaSelect.h"
0018 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHMassSelect.h"
0019 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHMassSymSelect.h"
0020 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHChi2Select.h"
0021 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHParticleMasses.h"
0022 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoBuilder.h"
0023 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHPlusMinusCandidate.h"
0024 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHTrackReference.h"
0025 #include "DataFormats/Candidate/interface/Candidate.h"
0026 #include "DataFormats/Math/interface/LorentzVector.h"
0027
0028 //---------------
0029 // C++ Headers --
0030 //---------------
0031 using namespace std;
0032
0033 //-------------------
0034 // Initializations --
0035 //-------------------
0036
0037 //----------------
0038 // Constructors --
0039 //----------------
0040 BPHDecayToTkpTknSymChargeBuilder::BPHDecayToTkpTknSymChargeBuilder(
0041     const BPHEventSetupWrapper& es,
0042     const string& daug1Name,
0043     double daug1Mass,
0044     double daug1Sigma,
0045     const string& daug2Name,
0046     double daug2Mass,
0047     double daug2Sigma,
0048     const BPHRecoBuilder::BPHGenericCollection* posCollection,
0049     const BPHRecoBuilder::BPHGenericCollection* negCollection,
0050     double expectedMass)
0051     : BPHDecayGenericBuilderBase(es),
0052       d1Name(daug1Name),
0053       d1Mass(daug1Mass),
0054       d1Sigma(daug1Sigma),
0055       d2Name(daug2Name),
0056       d2Mass(daug2Mass),
0057       d2Sigma(daug2Sigma),
0058       eMass(expectedMass),
0059       pCollection(posCollection),
0060       nCollection(negCollection),
0061       pt1Min(-1.0),
0062       pt2Min(-1.0),
0063       eta1Max(10.0),
0064       eta2Max(10.0),
0065       dzMax(1.0) {}
0066
0067 //--------------
0068 // Operations --
0069 //--------------
0070 void BPHDecayToTkpTknSymChargeBuilder::fillRecList() {
0071   // extract basic informations from input collections
0072
0073   vector<Particle*> pList;
0074   vector<Particle*> nList;
0075
0076   addParticle(pCollection, +1, pList);
0077   addParticle(nCollection, -1, nList);
0078   int iPos;
0079   int iNeg;
0080   int nPos = pList.size();
0081   int nNeg = nList.size();
0082   double massMin = getMassMin();
0083   double massMax = getMassMax();
0084   double mSqMin = massMin * massMin * 0.9;
0085   double mSqMax = massMax * massMax * 1.1;
0086   if (mSqMin < 0.0)
0087     mSqMin = 0.0;
0088
0089   for (iPos = 0; iPos < nPos; ++iPos) {
0090     Particle* pc = pList[iPos];
0091     const reco::Track* pt = pc->track;
0092     double px = pc->px;
0093     double py = pc->py;
0094     double pz = pc->pz;
0095     double p1 = pc->e1;
0096     double p2 = pc->e2;
0097     for (iNeg = 0; iNeg < nNeg; ++iNeg) {
0098       Particle* nc = nList[iNeg];
0099       const reco::Track* nt = nc->track;
0100       if (fabs(nt->dz() - pt->dz()) > dzMax)
0101         continue;
0102       double nx = nc->px;
0103       double ny = nc->py;
0104       double nz = nc->pz;
0105       double n1 = nc->e1;
0106       double n2 = nc->e2;
0107       const float tx = px + nx;
0108       const float ty = py + ny;
0109       const float tz = pz + nz;
0110       const float ta = ((p1 > 0.0) && (n2 > 0.0) ? p1 + n2 : -1.0);
0111       const float tb = ((p2 > 0.0) && (n1 > 0.0) ? p2 + n1 : -1.0);
0112       float ma = (ta > 0 ? (ta * ta) - ((tx * tx) + (ty * ty) + (tz * tz)) : -1.0);
0113       float mb = (tb > 0 ? (tb * tb) - ((tx * tx) + (ty * ty) + (tz * tz)) : -1.0);
0114       if (((ma < mSqMin) || (ma > mSqMax)) && ((mb < mSqMin) || (mb > mSqMax)))
0115         continue;
0116       BPHPlusMinusCandidatePtr rc(nullptr);
0117       float rcMass = -1.0;
0118       if (ma > 0.0) {
0119         rc = BPHPlusMinusCandidateWrap::create(evSetup);
0120         rc->add(d1Name, pc->cand, d1Mass, d1Sigma);
0121         rc->add(d2Name, nc->cand, d2Mass, d2Sigma);
0122         rcMass = rc->composite().mass();
0123       }
0124       BPHPlusMinusCandidatePtr rb(nullptr);
0125       float rbMass = -1.0;
0126       if (mb > 0.0) {
0127         rb = BPHPlusMinusCandidateWrap::create(evSetup);
0128         rb->add(d1Name, nc->cand, d1Mass, d1Sigma);
0129         rb->add(d2Name, pc->cand, d2Mass, d2Sigma);
0130         rbMass = rb->composite().mass();
0131       }
0132       BPHPlusMinusCandidatePtr* rp(nullptr);
0133       double mass = -1.0;
0134       if ((rc.get() != nullptr) && ((rb.get() == nullptr) || (fabs(rcMass - eMass) < fabs(rbMass - eMass)))) {
0135         mass = rcMass;
0136         rp = &rc;
0137       } else {
0138         mass = rbMass;
0139         rp = &rb;
0140       }
0141       BPHPlusMinusCandidate* rr = rp->get();
0142       if (mass < massMin)
0143         continue;
0144       if (mass > massMax)
0145         continue;
0146       if (!chi2Sel->accept(*rr))
0147         continue;
0148       recList.push_back(*rp);
0149     }
0150   }
0151
0152   for (iPos = 0; iPos < nPos; ++iPos)
0153     delete pList[iPos];
0154   for (iNeg = 0; iNeg < nNeg; ++iNeg)
0155     delete nList[iNeg];
0156
0157   return;
0158 }
0159
0160 /// set cuts
0161 void BPHDecayToTkpTknSymChargeBuilder::setTrk1PtMin(double pt) {
0162   outdated = true;
0163   pt1Min = pt;
0164   return;
0165 }
0166
0167 void BPHDecayToTkpTknSymChargeBuilder::setTrk2PtMin(double pt) {
0168   outdated = true;
0169   pt2Min = pt;
0170   return;
0171 }
0172
0173 void BPHDecayToTkpTknSymChargeBuilder::setTrk1EtaMax(double eta) {
0174   outdated = true;
0175   eta1Max = eta;
0176   return;
0177 }
0178
0179 void BPHDecayToTkpTknSymChargeBuilder::setTrk2EtaMax(double eta) {
0180   outdated = true;
0181   eta2Max = eta;
0182   return;
0183 }
0184
0185 void BPHDecayToTkpTknSymChargeBuilder::setDzMax(double dz) {
0186   outdated = true;
0187   dzMax = dz;
0188   return;
0189 }
0190
0191 void BPHDecayToTkpTknSymChargeBuilder::addParticle(const BPHRecoBuilder::BPHGenericCollection* collection,
0192                                                    int charge,
0193                                                    vector<Particle*>& list) {
0194   int i;
0195   int n = collection->size();
0196   list.reserve(n);
0197   for (i = 0; i < n; ++i) {
0198     const reco::Candidate& cand = collection->get(i);
0199     int q = cand.charge();
0200     if ((charge > 0) && (q <= 0))
0201       continue;
0202     if ((charge < 0) && (q >= 0))
0203       continue;
0204     const reco::Candidate::LorentzVector p4 = cand.p4();
0205     const reco::Track* tk = BPHTrackReference::getTrack(cand, "cfhp");
0206     if (tk == nullptr)
0207       continue;
0208     double px = p4.px();
0209     double py = p4.py();
0210     double pz = p4.pz();
0211     double p2 = (px * px) + (py * py) + (pz * pz);
0212     double e1 = ((p4.pt() >= pt1Min) && (p4.eta() <= eta1Max) ? sqrt(p2 + (d1Mass * d1Mass)) : -1.0);
0213     double e2 = ((p4.pt() >= pt2Min) && (p4.eta() <= eta2Max) ? sqrt(p2 + (d2Mass * d2Mass)) : -1.0);
0214     if ((e1 > 0.0) || (e2 > 0.0))
0215       list.push_back(new Particle(&cand, tk, px, py, pz, e1, e2));
0216   }
0217   return;
0218 }