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