Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 12:49:24

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/BPHDecayToV0DiffMassBuilder.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/SpecificDecay/interface/BPHDecayToTkpTknSymChargeBuilder.h"
0019 #include "DataFormats/Candidate/interface/Candidate.h"
0020 
0021 //---------------
0022 // C++ Headers --
0023 //---------------
0024 #include <cmath>
0025 using namespace std;
0026 
0027 //-------------------
0028 // Initializations --
0029 //-------------------
0030 
0031 //----------------
0032 // Constructors --
0033 //----------------
0034 BPHDecayToV0DiffMassBuilder::BPHDecayToV0DiffMassBuilder(const edm::EventSetup& es,
0035                                                          const std::string& d1Name,
0036                                                          double d1Mass,
0037                                                          double d1Sigma,
0038                                                          const std::string& d2Name,
0039                                                          double d2Mass,
0040                                                          double d2Sigma,
0041                                                          const BPHRecoBuilder::BPHGenericCollection* posCollection,
0042                                                          const BPHRecoBuilder::BPHGenericCollection* negCollection,
0043                                                          double expectedMass)
0044     : BPHDecayToV0Builder(es, d1Name, d2Name, posCollection, negCollection),
0045       p1Mass(d1Mass),
0046       p2Mass(d2Mass),
0047       p1Sigma(d1Sigma),
0048       p2Sigma(d2Sigma),
0049       expMass(expectedMass) {}
0050 
0051 BPHDecayToV0DiffMassBuilder::BPHDecayToV0DiffMassBuilder(const edm::EventSetup& es,
0052                                                          const std::string& d1Name,
0053                                                          double d1Mass,
0054                                                          double d1Sigma,
0055                                                          const std::string& d2Name,
0056                                                          double d2Mass,
0057                                                          double d2Sigma,
0058                                                          const std::vector<reco::VertexCompositeCandidate>* v0Collection,
0059                                                          double expectedMass,
0060                                                          const std::string& searchList)
0061     : BPHDecayToV0Builder(es, d1Name, d2Name, v0Collection, searchList),
0062       p1Mass(d1Mass),
0063       p2Mass(d2Mass),
0064       p1Sigma(d1Sigma),
0065       p2Sigma(d2Sigma),
0066       expMass(expectedMass) {}
0067 
0068 BPHDecayToV0DiffMassBuilder::BPHDecayToV0DiffMassBuilder(
0069     const edm::EventSetup& es,
0070     const std::string& d1Name,
0071     double d1Mass,
0072     double d1Sigma,
0073     const std::string& d2Name,
0074     double d2Mass,
0075     double d2Sigma,
0076     const std::vector<reco::VertexCompositePtrCandidate>* vpCollection,
0077     double expectedMass,
0078     const std::string& searchList)
0079     : BPHDecayToV0Builder(es, d1Name, d2Name, vpCollection, searchList),
0080       p1Mass(d1Mass),
0081       p2Mass(d2Mass),
0082       p1Sigma(d1Sigma),
0083       p2Sigma(d2Sigma),
0084       expMass(expectedMass) {}
0085 
0086 //--------------
0087 // Destructor --
0088 //--------------
0089 BPHDecayToV0DiffMassBuilder::~BPHDecayToV0DiffMassBuilder() {}
0090 
0091 //--------------
0092 // Operations --
0093 //--------------
0094 void BPHDecayToV0DiffMassBuilder::buildFromBPHGenericCollection() {
0095   BPHDecayToTkpTknSymChargeBuilder b(
0096       *evSetup, p1Name, p1Mass, p1Sigma, p2Name, p2Mass, p2Sigma, p1Collection, p2Collection, expMass);
0097 
0098   b.setTrk1PtMin(ptMin);
0099   b.setTrk2PtMin(ptMin);
0100   b.setTrk1EtaMax(etaMax);
0101   b.setTrk2EtaMax(etaMax);
0102   b.setMassRange(getMassMin(), getMassMax());
0103   b.setProbMin(getProbMin());
0104 
0105   cList = b.build();
0106 
0107   return;
0108 }
0109 
0110 BPHPlusMinusCandidatePtr BPHDecayToV0DiffMassBuilder::buildCandidate(const reco::Candidate* c1,
0111                                                                      const reco::Candidate* c2,
0112                                                                      const void* v0,
0113                                                                      v0Type type) {
0114   BPHPlusMinusCandidatePtr candX = BPHPlusMinusCandidateWrap::create(evSetup);
0115   BPHPlusMinusCandidatePtr candY = BPHPlusMinusCandidateWrap::create(evSetup);
0116   BPHPlusMinusCandidate* cptrX = candX.get();
0117   BPHPlusMinusCandidate* cptrY = candY.get();
0118   cptrX->add(p1Name, c1, sList, p1Mass, p1Sigma);
0119   cptrX->add(p2Name, c2, sList, p2Mass, p2Sigma);
0120   cptrY->add(p1Name, c2, sList, p1Mass, p1Sigma);
0121   cptrY->add(p2Name, c1, sList, p2Mass, p2Sigma);
0122   double mv0 = 0.0;
0123   switch (type) {
0124     case VertexCompositeCandidate:
0125       mv0 = static_cast<const reco::VertexCompositeCandidate*>(v0)->mass();
0126       break;
0127     case VertexCompositePtrCandidate:
0128       mv0 = static_cast<const reco::VertexCompositePtrCandidate*>(v0)->mass();
0129       break;
0130     default:
0131       mv0 = expMass;
0132       break;
0133   }
0134   double m1 = 0.0;
0135   double m2 = 0.0;
0136   if (p1Mass > p2Mass) {
0137     m1 = c1->mass();
0138     m2 = c2->mass();
0139   } else {
0140     m1 = c2->mass();
0141     m2 = c1->mass();
0142   }
0143   // check daughter masses in V0 CompositeCandidate
0144   double mcut = (p1Mass + p2Mass) / 2;
0145   if ((m1 > mcut) && (m2 < mcut))
0146     return candX;
0147   if ((m1 < mcut) && (m2 > mcut))
0148     return candY;
0149   // choose combination having the best invariant mass
0150   return (fabs(mv0 - cptrX->mass()) < fabs(mv0 - cptrY->mass()) ? candX : candY);
0151 }