Back to home page

Project CMSSW displayed by LXR

 
 

    


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

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/RecoDecay/interface/BPHDecayVertex.h"
0012 
0013 //-------------------------------
0014 // Collaborating Class Headers --
0015 //-------------------------------
0016 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoCandidate.h"
0017 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoBuilder.h"
0018 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHTrackReference.h"
0019 #include "FWCore/Framework/interface/EventSetup.h"
0020 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0021 #include "DataFormats/TrackReco/interface/Track.h"
0022 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0023 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
0024 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
0025 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
0026 #include "RecoVertex/VertexPrimitives/interface/TransientVertex.h"
0027 #include "RecoVertex/KalmanVertexFit/interface/KalmanVertexFitter.h"
0028 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0029 
0030 //---------------
0031 // C++ Headers --
0032 //---------------
0033 using namespace std;
0034 
0035 //-------------------
0036 // Initializations --
0037 //-------------------
0038 
0039 //----------------
0040 // Constructors --
0041 //----------------
0042 BPHDecayVertex::BPHDecayVertex(const edm::EventSetup* es)
0043     : evSetup(es),
0044       oldTracks(true),
0045       oldVertex(true),
0046       validTks(false),
0047       savedFitter(nullptr),
0048       savedBS(nullptr),
0049       savedPP(nullptr),
0050       savedPE(nullptr) {}
0051 
0052 BPHDecayVertex::BPHDecayVertex(const BPHDecayVertex* ptr, const edm::EventSetup* es)
0053     : evSetup(es),
0054       oldTracks(true),
0055       oldVertex(true),
0056       validTks(false),
0057       savedFitter(nullptr),
0058       savedBS(nullptr),
0059       savedPP(nullptr),
0060       savedPE(nullptr) {
0061   map<const reco::Candidate*, const reco::Candidate*> iMap;
0062   const vector<const reco::Candidate*>& daug = daughters();
0063   const vector<Component>& list = ptr->componentList();
0064   int i;
0065   int n = daug.size();
0066   for (i = 0; i < n; ++i) {
0067     const reco::Candidate* cand = daug[i];
0068     iMap[originalReco(cand)] = cand;
0069   }
0070   for (i = 0; i < n; ++i) {
0071     const Component& c = list[i];
0072     searchMap[iMap[c.cand]] = c.searchList;
0073   }
0074   const vector<BPHRecoConstCandPtr>& dComp = daughComp();
0075   int j;
0076   int m = dComp.size();
0077   for (j = 0; j < m; ++j) {
0078     const map<const reco::Candidate*, string>& dMap = dComp[j]->searchMap;
0079     searchMap.insert(dMap.begin(), dMap.end());
0080   }
0081 }
0082 
0083 //--------------
0084 // Destructor --
0085 //--------------
0086 BPHDecayVertex::~BPHDecayVertex() {}
0087 
0088 //--------------
0089 // Operations --
0090 //--------------
0091 bool BPHDecayVertex::validTracks() const {
0092   if (oldTracks)
0093     tTracks();
0094   return validTks;
0095 }
0096 
0097 bool BPHDecayVertex::validVertex() const {
0098   vertex();
0099   return validTks && fittedVertex.isValid();
0100 }
0101 
0102 const reco::Vertex& BPHDecayVertex::vertex(VertexFitter<5>* fitter,
0103                                            const reco::BeamSpot* bs,
0104                                            const GlobalPoint* priorPos,
0105                                            const GlobalError* priorError) const {
0106   if ((fitter == nullptr) && (bs == nullptr) && (priorPos == nullptr) && (priorError == nullptr)) {
0107     fitter = savedFitter;
0108     bs = savedBS;
0109     priorPos = savedPP;
0110     priorError = savedPE;
0111   }
0112   if (oldVertex || (fitter != savedFitter) || (bs != savedBS) || (priorPos != savedPP) || (priorError != savedPE)) {
0113     if (fitter != nullptr) {
0114       fitVertex(fitter, bs, priorPos, priorError);
0115     } else {
0116       KalmanVertexFitter kvf(true);
0117       fitVertex(&kvf, bs, priorPos, priorError);
0118     }
0119   }
0120   return fittedVertex;
0121 }
0122 
0123 const vector<const reco::Track*>& BPHDecayVertex::tracks() const {
0124   if (oldTracks)
0125     tTracks();
0126   return rTracks;
0127 }
0128 
0129 const reco::Track* BPHDecayVertex::getTrack(const reco::Candidate* cand) const {
0130   if (oldTracks)
0131     tTracks();
0132   map<const reco::Candidate*, const reco::Track*>::const_iterator iter = tkMap.find(cand);
0133   map<const reco::Candidate*, const reco::Track*>::const_iterator iend = tkMap.end();
0134   return (iter != iend ? iter->second : nullptr);
0135 }
0136 
0137 const vector<reco::TransientTrack>& BPHDecayVertex::transientTracks() const {
0138   if (oldTracks)
0139     tTracks();
0140   return trTracks;
0141 }
0142 
0143 reco::TransientTrack* BPHDecayVertex::getTransientTrack(const reco::Candidate* cand) const {
0144   if (oldTracks)
0145     tTracks();
0146   map<const reco::Candidate*, reco::TransientTrack*>::const_iterator iter = ttMap.find(cand);
0147   map<const reco::Candidate*, reco::TransientTrack*>::const_iterator iend = ttMap.end();
0148   return (iter != iend ? iter->second : nullptr);
0149 }
0150 
0151 /// retrieve EventSetup
0152 const edm::EventSetup* BPHDecayVertex::getEventSetup() const { return evSetup; }
0153 
0154 const string& BPHDecayVertex::getTrackSearchList(const reco::Candidate* cand) const {
0155   static string dum = "";
0156   map<const reco::Candidate*, string>::const_iterator iter = searchMap.find(cand);
0157   if (iter != searchMap.end())
0158     return iter->second;
0159   return dum;
0160 }
0161 
0162 void BPHDecayVertex::addV(const string& name, const reco::Candidate* daug, const string& searchList, double mass) {
0163   addP(name, daug, mass);
0164   searchMap[daughters().back()] = searchList;
0165   return;
0166 }
0167 
0168 void BPHDecayVertex::addV(const string& name, const BPHRecoConstCandPtr& comp) {
0169   addP(name, comp);
0170   const map<const reco::Candidate*, string>& dMap = comp->searchMap;
0171   searchMap.insert(dMap.begin(), dMap.end());
0172   return;
0173 }
0174 
0175 void BPHDecayVertex::setNotUpdated() const {
0176   BPHDecayMomentum::setNotUpdated();
0177   oldTracks = oldVertex = true;
0178   validTks = false;
0179   return;
0180 }
0181 
0182 void BPHDecayVertex::tTracks() const {
0183   oldTracks = false;
0184   rTracks.clear();
0185   trTracks.clear();
0186   tkMap.clear();
0187   ttMap.clear();
0188   edm::ESHandle<TransientTrackBuilder> ttB;
0189   evSetup->get<TransientTrackRecord>().get("TransientTrackBuilder", ttB);
0190   const vector<const reco::Candidate*>& dL = daughFull();
0191   int n = dL.size();
0192   trTracks.reserve(n);
0193   validTks = true;
0194   while (n--) {
0195     const reco::Candidate* rp = dL[n];
0196     tkMap[rp] = nullptr;
0197     ttMap[rp] = nullptr;
0198     if (!rp->charge())
0199       continue;
0200     const reco::Track* tp;
0201     const char* searchList = "cfhp";
0202     map<const reco::Candidate*, string>::const_iterator iter = searchMap.find(rp);
0203     if (iter != searchMap.end())
0204       searchList = iter->second.c_str();
0205     tp = BPHTrackReference::getTrack(*originalReco(rp), searchList);
0206     if (tp == nullptr) {
0207       edm::LogPrint("DataNotFound") << "BPHDecayVertex::tTracks: "
0208                                     << "no track for reco::(PF)Candidate";
0209       validTks = false;
0210       continue;
0211     }
0212     rTracks.push_back(tp);
0213     trTracks.push_back(ttB->build(tp));
0214     reco::TransientTrack* ttp = &trTracks.back();
0215     tkMap[rp] = tp;
0216     ttMap[rp] = ttp;
0217   }
0218   return;
0219 }
0220 
0221 void BPHDecayVertex::fitVertex(VertexFitter<5>* fitter,
0222                                const reco::BeamSpot* bs,
0223                                const GlobalPoint* priorPos,
0224                                const GlobalError* priorError) const {
0225   oldVertex = false;
0226   savedFitter = fitter;
0227   savedBS = bs;
0228   savedPP = priorPos;
0229   savedPE = priorError;
0230   if (oldTracks)
0231     tTracks();
0232   if (trTracks.size() < 2)
0233     return;
0234   try {
0235     if (bs == nullptr) {
0236       if (priorPos == nullptr) {
0237         TransientVertex tv = fitter->vertex(trTracks);
0238         fittedVertex = tv;
0239       } else {
0240         if (priorError == nullptr) {
0241           TransientVertex tv = fitter->vertex(trTracks, *priorPos);
0242           fittedVertex = tv;
0243         } else {
0244           TransientVertex tv = fitter->vertex(trTracks, *priorPos, *priorError);
0245           fittedVertex = tv;
0246         }
0247       }
0248     } else {
0249       TransientVertex tv = fitter->vertex(trTracks, *bs);
0250       fittedVertex = tv;
0251     }
0252   } catch (std::exception const&) {
0253     reco::Vertex tv;
0254     fittedVertex = tv;
0255     edm::LogPrint("FitFailed") << "BPHDecayVertex::fitVertex: "
0256                                << "vertex fit failed";
0257   }
0258   return;
0259 }