SpecificDecay/src/BPHX3872ToJPsiPiPiBuilder.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/BPHX3872ToJPsiPiPiBuilder.h"
0012
0013 //-------------------------------
0014 // Collaborating Class Headers --
0015 //-------------------------------
0016 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoBuilder.h"
0017 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHPlusMinusCandidate.h"
0018 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHTrackReference.h"
0019 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHMassSelect.h"
0020 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHChi2Select.h"
0021 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHMassFitSelect.h"
0022 #include "HeavyFlavorAnalysis/SpecificDecay/interface/BPHParticleMasses.h"
0023
0024 //---------------
0025 // C++ Headers --
0026 //---------------
0027 using namespace std;
0028
0029 //-------------------
0030 // Initializations --
0031 //-------------------
0032
0033 //----------------
0034 // Constructors --
0035 //----------------
0036 BPHX3872ToJPsiPiPiBuilder::BPHX3872ToJPsiPiPiBuilder(const edm::EventSetup& es,
0037                                                      const std::vector<BPHPlusMinusConstCandPtr>& jpsiCollection,
0038                                                      const BPHRecoBuilder::BPHGenericCollection* posCollection,
0039                                                      const BPHRecoBuilder::BPHGenericCollection* negCollection)
0040     : jPsiName("JPsi"),
0041       pionPosName("PionPos"),
0042       pionNegName("PionNeg"),
0043       evSetup(&es),
0044       jCollection(&jpsiCollection),
0045       pCollection(posCollection),
0046       nCollection(negCollection) {
0047   jpsiSel = new BPHMassSelect(2.80, 3.40);
0048   ptMin = 1.0;
0049   etaMax = 10.0;
0050   massSel = new BPHMassSelect(3.00, 4.50);
0051   chi2Sel = new BPHChi2Select(0.02);
0052   mFitSel = new BPHMassFitSelect(jPsiName, BPHParticleMasses::jPsiMass, BPHParticleMasses::jPsiMWidth, 3.50, 4.20);
0053   massConstr = true;
0054   minPDiff = 1.0e-4;
0055   updated = false;
0056 }
0057
0058 //--------------
0059 // Destructor --
0060 //--------------
0061 BPHX3872ToJPsiPiPiBuilder::~BPHX3872ToJPsiPiPiBuilder() {
0062   delete jpsiSel;
0063   delete massSel;
0064   delete chi2Sel;
0065   delete mFitSel;
0066 }
0067
0068 //--------------
0069 // Operations --
0070 //--------------
0071 vector<BPHRecoConstCandPtr> BPHX3872ToJPsiPiPiBuilder::build() {
0072   if (updated)
0073     return x3872List;
0074
0075   x3872List.clear();
0076
0077   // extract basic informations from input collections
0078
0079   class Particle {
0080   public:
0081     Particle(const reco::Candidate* c, const reco::Track* tk, double x, double y, double z, double p)
0082         : cand(c), track(tk), px(x), py(y), pz(z), ePion(p) {}
0083     const reco::Candidate* cand;
0084     const reco::Track* track;
0085     double px;
0086     double py;
0087     double pz;
0088     double ePion;
0089   };
0090
0091   vector<Particle> pList;
0092   vector<Particle> nList;
0093
0094   int nPos = pCollection->size();
0095   int nNeg = nCollection->size();
0096
0097   pList.reserve(nPos);
0098   nList.reserve(nNeg);
0099
0100   // filter input collections
0101
0102   int iPos;
0103   int iNeg;
0104
0105   for (iPos = 0; iPos < nPos; ++iPos) {
0106     const reco::Candidate& cand = pCollection->get(iPos);
0107     if (cand.charge() != +1)
0108       continue;
0109     const reco::Candidate::LorentzVector p4 = cand.p4();
0110     if (p4.pt() < ptMin)
0111       continue;
0112     if (p4.eta() > etaMax)
0113       continue;
0114     const reco::Track* tk = BPHTrackReference::getTrack(cand, "cfhp");
0115     if (tk == nullptr)
0116       continue;
0117     double px = p4.px();
0118     double py = p4.py();
0119     double pz = p4.pz();
0120     double p2 = (px * px) + (py * py) + (pz * pz);
0121     pList.push_back(
0122         Particle(&cand, tk, px, py, pz, sqrt(p2 + (BPHParticleMasses::pionMass * BPHParticleMasses::pionMass))));
0123   }
0124
0125   for (iNeg = 0; iNeg < nNeg; ++iNeg) {
0126     const reco::Candidate& cand = nCollection->get(iNeg);
0127     if (cand.charge() != -1)
0128       continue;
0129     const reco::Candidate::LorentzVector p4 = cand.p4();
0130     if (p4.pt() < ptMin)
0131       continue;
0132     if (p4.eta() > etaMax)
0133       continue;
0134     const reco::Track* tk = BPHTrackReference::getTrack(cand, "cfhp");
0135     if (tk == nullptr)
0136       continue;
0137     double px = p4.px();
0138     double py = p4.py();
0139     double pz = p4.pz();
0140     double p2 = (px * px) + (py * py) + (pz * pz);
0141     nList.push_back(
0142         Particle(&cand, tk, px, py, pz, sqrt(p2 + (BPHParticleMasses::pionMass * BPHParticleMasses::pionMass))));
0143   }
0144
0145   // filter basic candidates
0146
0147   nPos = pList.size();
0148   nNeg = nList.size();
0149
0150   double mMax = getMassMax() - (0.8 * BPHParticleMasses::jPsiMass);
0151
0152   struct PionPair {
0153     const reco::Candidate* posPion;
0154     const reco::Candidate* negPion;
0155   };
0156   vector<const PionPair*> pionPairs;
0157   pionPairs.reserve(nPos * nNeg);
0158
0159   for (iPos = 0; iPos < nPos; ++iPos) {
0160     Particle& pc = pList[iPos];
0161     const reco::Track* pt = pc.track;
0162     double px = pc.px;
0163     double py = pc.py;
0164     double pz = pc.pz;
0165     double pe = pc.ePion;
0166     for (iNeg = 0; iNeg < nNeg; ++iNeg) {
0167       Particle& nc = nList[iNeg];
0168       const reco::Track* nt = nc.track;
0169       if (fabs(nt->dz() - pt->dz()) > 1.0)
0170         continue;
0171       double nx = nc.px;
0172       double ny = nc.py;
0173       double nz = nc.pz;
0174       double ne = nc.ePion;
0175       const float tx = px + nx;
0176       const float ty = py + ny;
0177       const float tz = pz + nz;
0178       const float te = pe + ne;
0179       float mass = (te * te) - ((tx * tx) + (ty * ty) + (tz * tz));
0180       if (mass > mMax)
0181         continue;
0182       PionPair* pp = new PionPair;
0183       pp->posPion = pc.cand;
0184       pp->negPion = nc.cand;
0185       pionPairs.push_back(pp);
0186     }
0187   }
0188
0189   vector<BPHPlusMinusConstCandPtr> jPsi;
0190   int nJPsi = jCollection->size();
0191   int iJPsi;
0192   jPsi.reserve(nJPsi);
0193   for (iJPsi = 0; iJPsi < nJPsi; ++iJPsi) {
0194     const BPHPlusMinusConstCandPtr& jpCand = jCollection->at(iJPsi);
0195     if (jpsiSel->accept(*jpCand))
0196       jPsi.push_back(jpCand);
0197   }
0198   nJPsi = jPsi.size();
0199
0200   int nPair = pionPairs.size();
0201   int iPair;
0202   for (iPair = 0; iPair < nPair; ++iPair) {
0203     const PionPair* pp = pionPairs[iPair];
0204     for (iJPsi = 0; iJPsi < nJPsi; ++iJPsi) {
0205       BPHRecoCandidate* x3872 = new BPHRecoCandidate(evSetup);
0206       BPHRecoCandidatePtr xPtr(x3872);
0207       x3872->add(jPsiName, jPsi[iJPsi]);
0208       x3872->add(pionPosName, pp->posPion, BPHParticleMasses::pionMass, BPHParticleMasses::pionMSigma);
0209       x3872->add(pionNegName, pp->negPion, BPHParticleMasses::pionMass, BPHParticleMasses::pionMSigma);
0210       if (!massSel->accept(*x3872))
0211         continue;
0212       if ((chi2Sel != nullptr) && !chi2Sel->accept(*x3872))
0213         continue;
0214       if (!mFitSel->accept(*x3872))
0215         continue;
0216       x3872List.push_back(xPtr);
0217     }
0218     delete pp;
0219   }
0220
0221   updated = true;
0222   return x3872List;
0223 }
0224
0225 /// set cuts
0226 void BPHX3872ToJPsiPiPiBuilder::setJPsiMassMin(double m) {
0227   updated = false;
0228   jpsiSel->setMassMin(m);
0229   return;
0230 }
0231
0232 void BPHX3872ToJPsiPiPiBuilder::setJPsiMassMax(double m) {
0233   updated = false;
0234   jpsiSel->setMassMax(m);
0235   return;
0236 }
0237
0238 void BPHX3872ToJPsiPiPiBuilder::setPiPtMin(double pt) {
0239   updated = false;
0240   ptMin = pt;
0241   return;
0242 }
0243
0244 void BPHX3872ToJPsiPiPiBuilder::setPiEtaMax(double eta) {
0245   updated = false;
0246   etaMax = eta;
0247   return;
0248 }
0249
0250 void BPHX3872ToJPsiPiPiBuilder::setMassMin(double m) {
0251   updated = false;
0252   massSel->setMassMin(m);
0253   return;
0254 }
0255
0256 void BPHX3872ToJPsiPiPiBuilder::setMassMax(double m) {
0257   updated = false;
0258   massSel->setMassMax(m);
0259   return;
0260 }
0261
0262 void BPHX3872ToJPsiPiPiBuilder::setProbMin(double p) {
0263   updated = false;
0264   delete chi2Sel;
0265   chi2Sel = (p < 0.0 ? nullptr : new BPHChi2Select(p));
0266   return;
0267 }
0268
0269 void BPHX3872ToJPsiPiPiBuilder::setMassFitMin(double m) {
0270   updated = false;
0271   mFitSel->setMassMin(m);
0272   return;
0273 }
0274
0275 void BPHX3872ToJPsiPiPiBuilder::setMassFitMax(double m) {
0276   updated = false;
0277   mFitSel->setMassMax(m);
0278   return;
0279 }
0280
0281 void BPHX3872ToJPsiPiPiBuilder::setConstr(bool flag) {
0282   updated = false;
0283   massConstr = flag;
0284   return;
0285 }
0286
0287 /// get current cuts
0288 double BPHX3872ToJPsiPiPiBuilder::getJPsiMassMin() const { return jpsiSel->getMassMin(); }
0289
0290 double BPHX3872ToJPsiPiPiBuilder::getJPsiMassMax() const { return jpsiSel->getMassMax(); }
0291
0292 double BPHX3872ToJPsiPiPiBuilder::getPiPtMin() const { return ptMin; }
0293
0294 double BPHX3872ToJPsiPiPiBuilder::getPiEtaMax() const { return etaMax; }
0295
0296 double BPHX3872ToJPsiPiPiBuilder::getMassMin() const { return massSel->getMassMin(); }
0297
0298 double BPHX3872ToJPsiPiPiBuilder::getMassMax() const { return massSel->getMassMax(); }
0299
0300 double BPHX3872ToJPsiPiPiBuilder::getProbMin() const { return (chi2Sel == nullptr ? -1.0 : chi2Sel->getProbMin()); }
0301
0302 double BPHX3872ToJPsiPiPiBuilder::getMassFitMin() const { return mFitSel->getMassMin(); }
0303
0304 double BPHX3872ToJPsiPiPiBuilder::getMassFitMax() const { return mFitSel->getMassMax(); }
0305
0306 bool BPHX3872ToJPsiPiPiBuilder::getConstr() const { return massConstr; }