File indexing completed on 2021-02-14 12:49:20
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoBuilder.h"
0012
0013
0014
0015
0016 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoSelect.h"
0017 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHMomentumSelect.h"
0018 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHVertexSelect.h"
0019 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHFitSelect.h"
0020 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoCandidate.h"
0021 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHTrackReference.h"
0022 #include "DataFormats/RecoCandidate/interface/RecoCandidate.h"
0023 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0024
0025
0026
0027
0028 using namespace std;
0029
0030
0031
0032
0033
0034
0035
0036
0037 BPHRecoBuilder::BPHRecoBuilder(const edm::EventSetup& es) : evSetup(&es), minPDiff(-1.0) {
0038 msList.reserve(5);
0039 vsList.reserve(5);
0040 }
0041
0042
0043
0044
0045 BPHRecoBuilder::~BPHRecoBuilder() {
0046 int n = sourceList.size();
0047 while (n--) {
0048 delete sourceList[n]->collection;
0049 delete sourceList[n];
0050 }
0051 int m = srCompList.size();
0052 while (m--)
0053 delete srCompList[m];
0054 while (!compCollectList.empty()) {
0055 const vector<BPHRecoConstCandPtr>* cCollection = *compCollectList.begin();
0056 delete cCollection;
0057 compCollectList.erase(cCollection);
0058 }
0059 }
0060
0061
0062
0063
0064 BPHRecoBuilder::BPHGenericCollection* BPHRecoBuilder::createCollection(const vector<const reco::Candidate*>& candList,
0065 const string& list) {
0066 return new BPHSpecificCollection<vector<const reco::Candidate*>>(candList, list);
0067 }
0068
0069 void BPHRecoBuilder::add(const string& name, const BPHGenericCollection* collection, double mass, double msig) {
0070 BPHRecoSource* rs;
0071 if (sourceId.find(name) != sourceId.end()) {
0072 edm::LogPrint("TooManyParticles") << "BPHRecoBuilder::add: "
0073 << "Decay product already inserted with name " << name << " , skipped";
0074 return;
0075 }
0076 rs = new BPHRecoSource;
0077 rs->name = &sourceId.insert(make_pair(name, sourceList.size())).first->first;
0078 rs->collection = collection;
0079 rs->selector.reserve(5);
0080 rs->mass = mass;
0081 rs->msig = msig;
0082 sourceList.push_back(rs);
0083 return;
0084 }
0085
0086 void BPHRecoBuilder::add(const string& name, const vector<BPHRecoConstCandPtr>& collection) {
0087 BPHCompSource* cs;
0088 if (srCompId.find(name) != srCompId.end()) {
0089 edm::LogPrint("TooManyParticles") << "BPHRecoBuilder::add: "
0090 << "Decay product already inserted with name " << name << " , skipped";
0091 return;
0092 }
0093 cs = new BPHCompSource;
0094 cs->name = &srCompId.insert(make_pair(name, srCompList.size())).first->first;
0095 cs->collection = &collection;
0096 srCompList.push_back(cs);
0097 return;
0098 }
0099
0100 void BPHRecoBuilder::filter(const string& name, const BPHRecoSelect& sel) const {
0101 map<string, int>::const_iterator iter = sourceId.find(name);
0102 if (iter == sourceId.end())
0103 return;
0104 BPHRecoSource* rs = sourceList[iter->second];
0105 rs->selector.push_back(&sel);
0106 return;
0107 }
0108
0109 void BPHRecoBuilder::filter(const string& name, const BPHMomentumSelect& sel) const {
0110 map<string, int>::const_iterator iter = srCompId.find(name);
0111 if (iter == sourceId.end())
0112 return;
0113 BPHCompSource* cs = srCompList[iter->second];
0114 cs->momSelector.push_back(&sel);
0115 return;
0116 }
0117
0118 void BPHRecoBuilder::filter(const string& name, const BPHVertexSelect& sel) const {
0119 map<string, int>::const_iterator iter = srCompId.find(name);
0120 if (iter == sourceId.end())
0121 return;
0122 BPHCompSource* cs = srCompList[iter->second];
0123 cs->vtxSelector.push_back(&sel);
0124 return;
0125 }
0126
0127 void BPHRecoBuilder::filter(const string& name, const BPHFitSelect& sel) const {
0128 map<string, int>::const_iterator iter = srCompId.find(name);
0129 if (iter == sourceId.end())
0130 return;
0131 BPHCompSource* cs = srCompList[iter->second];
0132 cs->fitSelector.push_back(&sel);
0133 return;
0134 }
0135
0136 void BPHRecoBuilder::filter(const BPHMomentumSelect& sel) {
0137 msList.push_back(&sel);
0138 return;
0139 }
0140
0141 void BPHRecoBuilder::filter(const BPHVertexSelect& sel) {
0142 vsList.push_back(&sel);
0143 return;
0144 }
0145
0146 void BPHRecoBuilder::filter(const BPHFitSelect& sel) {
0147 fsList.push_back(&sel);
0148 return;
0149 }
0150
0151 bool BPHRecoBuilder::accept(const BPHRecoCandidate& cand) const {
0152 int i;
0153 int n;
0154 n = msList.size();
0155 for (i = 0; i < n; ++i) {
0156 if (!msList[i]->accept(cand))
0157 return false;
0158 }
0159 n = vsList.size();
0160 for (i = 0; i < n; ++i) {
0161 if (!vsList[i]->accept(cand))
0162 return false;
0163 }
0164 n = fsList.size();
0165 for (i = 0; i < n; ++i) {
0166 if (!fsList[i]->accept(cand))
0167 return false;
0168 }
0169 return true;
0170 }
0171
0172 void BPHRecoBuilder::setMinPDiffererence(double pMin) {
0173 minPDiff = pMin;
0174 return;
0175 }
0176
0177 vector<BPHRecoBuilder::ComponentSet> BPHRecoBuilder::build() const {
0178 daugMap.clear();
0179 compMap.clear();
0180 vector<ComponentSet> candList;
0181 ComponentSet compSet;
0182 build(candList, compSet, sourceList.begin(), sourceList.end(), srCompList.begin(), srCompList.end());
0183 return candList;
0184 }
0185
0186 const edm::EventSetup* BPHRecoBuilder::eventSetup() const { return evSetup; }
0187
0188 bool BPHRecoBuilder::sameTrack(const reco::Candidate* lCand, const reco::Candidate* rCand, double minPDifference) {
0189 const reco::Track* lrcTrack = BPHTrackReference::getFromRC(*lCand);
0190 const reco::Track* rrcTrack = BPHTrackReference::getFromRC(*rCand);
0191 const reco::Track* lpfTrack = BPHTrackReference::getFromPF(*lCand);
0192 const reco::Track* rpfTrack = BPHTrackReference::getFromPF(*rCand);
0193 if ((lrcTrack != nullptr) && ((lrcTrack == rrcTrack) || (lrcTrack == rpfTrack)))
0194 return true;
0195 if ((lpfTrack != nullptr) && ((lpfTrack == rrcTrack) || (lpfTrack == rpfTrack)))
0196 return true;
0197 reco::Candidate::Vector pDiff = (lCand->momentum() - rCand->momentum());
0198 reco::Candidate::Vector pMean = (lCand->momentum() + rCand->momentum());
0199 double pDMod = pDiff.mag2();
0200 double pMMod = pMean.mag2();
0201 if (((pDMod / pMMod) < minPDifference) && (lCand->charge() == rCand->charge()))
0202 return true;
0203 return false;
0204 }
0205
0206 void BPHRecoBuilder::build(vector<ComponentSet>& compList,
0207 ComponentSet& compSet,
0208 vector<BPHRecoSource*>::const_iterator r_iter,
0209 vector<BPHRecoSource*>::const_iterator r_iend,
0210 vector<BPHCompSource*>::const_iterator c_iter,
0211 vector<BPHCompSource*>::const_iterator c_iend) const {
0212 if (r_iter == r_iend) {
0213 if (c_iter == c_iend) {
0214 compSet.compMap = compMap;
0215 compList.push_back(compSet);
0216 return;
0217 }
0218 BPHCompSource* source = *c_iter++;
0219 const vector<BPHRecoConstCandPtr>* collection = source->collection;
0220 vector<const BPHMomentumSelect*> momSelector = source->momSelector;
0221 vector<const BPHVertexSelect*> vtxSelector = source->vtxSelector;
0222 vector<const BPHFitSelect*> fitSelector = source->fitSelector;
0223 int i;
0224 int j;
0225 int n = collection->size();
0226 int m;
0227 bool skip;
0228 for (i = 0; i < n; ++i) {
0229 skip = false;
0230 BPHRecoConstCandPtr cand = collection->at(i);
0231 if (contained(compSet, cand))
0232 continue;
0233 m = momSelector.size();
0234 for (j = 0; j < m; ++j) {
0235 if (!momSelector[j]->accept(*cand)) {
0236 skip = true;
0237 break;
0238 }
0239 }
0240 if (skip)
0241 continue;
0242 m = vtxSelector.size();
0243 for (j = 0; j < m; ++j) {
0244 if (!vtxSelector[j]->accept(*cand)) {
0245 skip = true;
0246 break;
0247 }
0248 }
0249 if (skip)
0250 continue;
0251 m = fitSelector.size();
0252 for (j = 0; j < m; ++j) {
0253 if (!fitSelector[j]->accept(*cand)) {
0254 skip = true;
0255 break;
0256 }
0257 }
0258 if (skip)
0259 continue;
0260 compMap[*source->name] = cand;
0261 build(compList, compSet, r_iter, r_iend, c_iter, c_iend);
0262 compMap.erase(*source->name);
0263 }
0264 return;
0265 }
0266 BPHRecoSource* source = *r_iter++;
0267 const BPHGenericCollection* collection = source->collection;
0268 vector<const BPHRecoSelect*>& selector = source->selector;
0269 int i;
0270 int j;
0271 int n = collection->size();
0272 int m = selector.size();
0273 bool skip;
0274 for (i = 0; i < n; ++i) {
0275 const reco::Candidate& cand = collection->get(i);
0276 if (contained(compSet, &cand))
0277 continue;
0278 skip = false;
0279 for (j = 0; j < m; ++j) {
0280 if (!selector[j]->accept(cand, this)) {
0281 skip = true;
0282 break;
0283 }
0284 }
0285 if (skip)
0286 continue;
0287 BPHDecayMomentum::Component comp;
0288 comp.cand = &cand;
0289 comp.mass = source->mass;
0290 comp.msig = source->msig;
0291 comp.searchList = collection->searchList();
0292 compSet.daugMap[*source->name] = comp;
0293 daugMap[*source->name] = &cand;
0294 build(compList, compSet, r_iter, r_iend, c_iter, c_iend);
0295 daugMap.erase(*source->name);
0296 compSet.daugMap.erase(*source->name);
0297 }
0298 return;
0299 }
0300
0301 bool BPHRecoBuilder::contained(ComponentSet& compSet, const reco::Candidate* cand) const {
0302 map<string, BPHDecayMomentum::Component>& dMap = compSet.daugMap;
0303 map<string, BPHDecayMomentum::Component>::const_iterator d_iter;
0304 map<string, BPHDecayMomentum::Component>::const_iterator d_iend = dMap.end();
0305 for (d_iter = dMap.begin(); d_iter != d_iend; ++d_iter) {
0306 const reco::Candidate* cChk = d_iter->second.cand;
0307 if (cand == cChk)
0308 return true;
0309 if (sameTrack(cand, cChk))
0310 return true;
0311 }
0312 return false;
0313 }
0314
0315 bool BPHRecoBuilder::contained(ComponentSet& compSet, BPHRecoConstCandPtr cCand) const {
0316 map<string, BPHRecoConstCandPtr>::const_iterator c_iter;
0317 map<string, BPHRecoConstCandPtr>::const_iterator c_iend = compMap.end();
0318 const vector<const reco::Candidate*>& dCand = cCand->daughFull();
0319 int j;
0320 int m = dCand.size();
0321 int k;
0322 int l;
0323 for (j = 0; j < m; ++j) {
0324 const reco::Candidate* cand = cCand->originalReco(dCand[j]);
0325 map<string, BPHDecayMomentum::Component>& dMap = compSet.daugMap;
0326 map<string, BPHDecayMomentum::Component>::const_iterator d_iter;
0327 map<string, BPHDecayMomentum::Component>::const_iterator d_iend = dMap.end();
0328 for (d_iter = dMap.begin(); d_iter != d_iend; ++d_iter) {
0329 const reco::Candidate* cChk = d_iter->second.cand;
0330 if (cand == cChk)
0331 return true;
0332 if (sameTrack(cand, cChk))
0333 return true;
0334 }
0335
0336 for (c_iter = compMap.begin(); c_iter != c_iend; ++c_iter) {
0337 const pair<string, BPHRecoConstCandPtr>& entry = *c_iter;
0338 BPHRecoConstCandPtr cCChk = entry.second;
0339 const vector<const reco::Candidate*>& dCChk = cCChk->daughFull();
0340 l = dCChk.size();
0341 for (k = 0; k < l; ++k) {
0342 const reco::Candidate* cChk = cCChk->originalReco(dCChk[k]);
0343 if (cand == cChk)
0344 return true;
0345 if (sameTrack(cand, cChk))
0346 return true;
0347 }
0348 }
0349 }
0350
0351 return false;
0352 }
0353
0354 bool BPHRecoBuilder::sameTrack(const reco::Candidate* lCand, const reco::Candidate* rCand) const {
0355 return sameTrack(lCand, rCand, minPDiff);
0356 }