File indexing completed on 2024-04-06 12:15:38
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHDecayToTkpTknSymChargeBuilder.h"
0012
0013
0014
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
0030
0031 using namespace std;
0032
0033
0034
0035
0036
0037
0038
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
0069
0070 void BPHDecayToTkpTknSymChargeBuilder::fillRecList() {
0071
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
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 }