File indexing completed on 2021-02-14 12:49:20
0001 #include "HeavyFlavorAnalysis/RecoDecay/test/stubs/TestBPHRecoDecay.h"
0002
0003 #include "FWCore/Framework/interface/MakerMacros.h"
0004 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0006
0007 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoBuilder.h"
0008 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoSelect.h"
0009 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHRecoCandidate.h"
0010 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHPlusMinusCandidate.h"
0011 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHMomentumSelect.h"
0012 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHVertexSelect.h"
0013 #include "HeavyFlavorAnalysis/RecoDecay/interface/BPHTrackReference.h"
0014
0015 #include "DataFormats/PatCandidates/interface/Muon.h"
0016 #include "DataFormats/TrackReco/interface/Track.h"
0017
0018 #include "DataFormats/PatCandidates/interface/GenericParticle.h"
0019 #include "DataFormats/PatCandidates/interface/CompositeCandidate.h"
0020
0021 #include "RecoVertex/KinematicFit/interface/TwoTrackMassKinematicConstraint.h"
0022 #include <TH1.h>
0023 #include <TFile.h>
0024 #include <TMath.h>
0025
0026 #include <set>
0027 #include <string>
0028 #include <iostream>
0029 #include <fstream>
0030
0031 using namespace std;
0032
0033 #define SET_LABEL(NAME, PSET) (NAME = PSET.getParameter<string>(#NAME))
0034
0035
0036
0037
0038 TestBPHRecoDecay::TestBPHRecoDecay(const edm::ParameterSet& ps) {
0039 usePM = (!SET_LABEL(patMuonLabel, ps).empty());
0040 useCC = (!SET_LABEL(ccCandsLabel, ps).empty());
0041 usePF = (!SET_LABEL(pfCandsLabel, ps).empty());
0042 usePC = (!SET_LABEL(pcCandsLabel, ps).empty());
0043 useGP = (!SET_LABEL(gpCandsLabel, ps).empty());
0044
0045 if (usePM)
0046 consume<pat::MuonCollection>(patMuonToken, patMuonLabel);
0047 if (useCC)
0048 consume<vector<pat::CompositeCandidate>>(ccCandsToken, ccCandsLabel);
0049 if (usePF)
0050 consume<vector<reco::PFCandidate>>(pfCandsToken, pfCandsLabel);
0051 if (usePC)
0052 consume<vector<BPHTrackReference::candidate>>(pcCandsToken, pcCandsLabel);
0053 if (useGP)
0054 consume<vector<pat::GenericParticle>>(gpCandsToken, gpCandsLabel);
0055 SET_LABEL(outDump, ps);
0056 SET_LABEL(outHist, ps);
0057 if (outDump.empty())
0058 fPtr = &cout;
0059 else
0060 fPtr = new ofstream(outDump.c_str());
0061 }
0062
0063 TestBPHRecoDecay::~TestBPHRecoDecay() {}
0064
0065 void TestBPHRecoDecay::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0066 edm::ParameterSetDescription desc;
0067 desc.add<string>("patMuonLabel", "");
0068 desc.add<string>("ccCandsLabel", "");
0069 desc.add<string>("pfCandsLabel", "");
0070 desc.add<string>("pcCandsLabel", "");
0071 desc.add<string>("gpCandsLabel", "");
0072 desc.add<string>("outDump", "dump.txt");
0073 desc.add<string>("outHist", "hist.root");
0074 descriptions.add("testBPHRecoDecay", desc);
0075 return;
0076 }
0077
0078 void TestBPHRecoDecay::beginJob() {
0079 *fPtr << "TestBPHRecoDecay::beginJob" << endl;
0080 createHisto("massJPsi", 50, 2.85, 3.35);
0081 createHisto("mcstJPsi", 50, 2.85, 3.35);
0082 createHisto("massPhi", 50, 0.995, 1.045);
0083 createHisto("massBu", 20, 5.0, 5.5);
0084 createHisto("mcstBu", 20, 5.0, 5.5);
0085 createHisto("massBs", 20, 5.1, 5.6);
0086 createHisto("mcstBs", 20, 5.1, 5.6);
0087 createHisto("massBsPhi", 50, 0.995, 1.045);
0088 return;
0089 }
0090
0091 void TestBPHRecoDecay::analyze(const edm::Event& ev, const edm::EventSetup& es) {
0092 ostream& outF = *fPtr;
0093 outF << "--------- event " << ev.id().run() << " / " << ev.id().event() << " ---------" << endl;
0094
0095
0096
0097
0098
0099 int nrc = 0;
0100
0101
0102 edm::Handle<vector<reco::PFCandidate>> pfCands;
0103 if (usePF) {
0104 pfCandsToken.get(ev, pfCands);
0105 nrc = pfCands->size();
0106 if (pfCands.isValid())
0107 outF << nrc << " pfCands found" << endl;
0108 else
0109 outF << "no pfCands" << endl;
0110 }
0111
0112
0113
0114
0115
0116 edm::Handle<vector<BPHTrackReference::candidate>> pcCands;
0117 if (usePC) {
0118 pcCandsToken.get(ev, pcCands);
0119 nrc = pcCands->size();
0120 if (pcCands.isValid())
0121 outF << nrc << " pcCands found" << endl;
0122 else
0123 outF << "no pcCands" << endl;
0124 }
0125
0126
0127 edm::Handle<vector<pat::GenericParticle>> gpCands;
0128 if (useGP) {
0129 gpCandsToken.get(ev, gpCands);
0130 nrc = gpCands->size();
0131 if (gpCands.isValid())
0132 outF << nrc << " gpCands found" << endl;
0133 else
0134 outF << "no gpCands" << endl;
0135 }
0136
0137
0138 edm::Handle<pat::MuonCollection> patMuon;
0139 if (usePM) {
0140 patMuonToken.get(ev, patMuon);
0141 int n = patMuon->size();
0142 if (patMuon.isValid())
0143 outF << n << " muons found" << endl;
0144 else
0145 outF << "no muons" << endl;
0146 }
0147
0148
0149
0150 vector<const reco::Candidate*> muDaugs;
0151 set<const pat::Muon*> muonSet;
0152 if (useCC) {
0153 edm::Handle<vector<pat::CompositeCandidate>> ccCands;
0154 ccCandsToken.get(ev, ccCands);
0155 int n = ccCands->size();
0156 if (ccCands.isValid())
0157 outF << n << " ccCands found" << endl;
0158 else
0159 outF << "no ccCands" << endl;
0160 muDaugs.clear();
0161 muDaugs.reserve(n);
0162 muonSet.clear();
0163 set<const pat::Muon*>::const_iterator iter;
0164 set<const pat::Muon*>::const_iterator iend;
0165 int i;
0166 for (i = 0; i < n; ++i) {
0167 const pat::CompositeCandidate& cc = ccCands->at(i);
0168 int j;
0169 int m = cc.numberOfDaughters();
0170 for (j = 0; j < m; ++j) {
0171 const reco::Candidate* dp = cc.daughter(j);
0172 const pat::Muon* mp = dynamic_cast<const pat::Muon*>(dp);
0173 iter = muonSet.begin();
0174 iend = muonSet.end();
0175 bool add = (mp != nullptr) && (muonSet.find(mp) == iend);
0176 while (add && (iter != iend)) {
0177 if (BPHRecoBuilder::sameTrack(mp, *iter++, 1.0e-5))
0178 add = false;
0179 }
0180 if (add)
0181 muonSet.insert(mp);
0182 }
0183 }
0184 iter = muonSet.begin();
0185 iend = muonSet.end();
0186 while (iter != iend)
0187 muDaugs.push_back(*iter++);
0188 if ((n = muDaugs.size()))
0189 outF << n << " muons found" << endl;
0190 }
0191
0192
0193
0194
0195
0196
0197 class MuonChargeSelect : public BPHRecoSelect {
0198 public:
0199 MuonChargeSelect(int c) : charge(c) {}
0200 ~MuonChargeSelect() override {}
0201 bool accept(const reco::Candidate& cand) const override {
0202 const pat::Muon* p = dynamic_cast<const pat::Muon*>(&cand);
0203 if (p == nullptr)
0204 return false;
0205 return ((charge * cand.charge()) > 0);
0206 }
0207
0208 private:
0209 int charge;
0210 };
0211
0212
0213 class MuonPtSelect : public BPHRecoSelect {
0214 public:
0215 MuonPtSelect(float pt) : ptCut(pt) {}
0216 ~MuonPtSelect() override {}
0217 bool accept(const reco::Candidate& cand) const override {
0218 const pat::Muon* p = dynamic_cast<const pat::Muon*>(&cand);
0219 if (p == nullptr)
0220 return false;
0221 return (p->p4().pt() > ptCut);
0222 }
0223
0224 private:
0225 float ptCut;
0226 };
0227
0228
0229 class MuonEtaSelect : public BPHRecoSelect {
0230 public:
0231 MuonEtaSelect(float eta) : etaCut(eta) {}
0232 ~MuonEtaSelect() override {}
0233 bool accept(const reco::Candidate& cand) const override {
0234 const pat::Muon* p = dynamic_cast<const pat::Muon*>(&cand);
0235 if (p == nullptr)
0236 return false;
0237 return (fabs(p->p4().eta()) < etaCut);
0238 }
0239
0240 private:
0241 float etaCut;
0242 };
0243
0244
0245 class KaonChargeSelect : public BPHRecoSelect {
0246 public:
0247 KaonChargeSelect(int c) : charge(c) {}
0248 ~KaonChargeSelect() override {}
0249 bool accept(const reco::Candidate& cand) const override { return ((charge * cand.charge()) > 0); }
0250
0251 private:
0252 int charge;
0253 };
0254
0255 class KaonNeutralVeto : public BPHRecoSelect {
0256 public:
0257 KaonNeutralVeto() {}
0258 ~KaonNeutralVeto() override {}
0259 bool accept(const reco::Candidate& cand) const override { return lround(fabs(cand.charge())); }
0260 };
0261
0262
0263 class KaonPtSelect : public BPHRecoSelect {
0264 public:
0265 KaonPtSelect(float pt) : ptCut(pt) {}
0266 ~KaonPtSelect() override {}
0267 bool accept(const reco::Candidate& cand) const override { return (cand.p4().pt() > ptCut); }
0268
0269 private:
0270 float ptCut;
0271 };
0272
0273
0274 class KaonEtaSelect : public BPHRecoSelect {
0275 public:
0276 KaonEtaSelect(float eta) : etaCut(eta) {}
0277 ~KaonEtaSelect() override {}
0278 bool accept(const reco::Candidate& cand) const override { return (fabs(cand.p4().eta()) < etaCut); }
0279
0280 private:
0281 float etaCut;
0282 };
0283
0284
0285
0286
0287
0288
0289 class MassSelect : public BPHMomentumSelect {
0290 public:
0291 MassSelect(double minMass, double maxMass) : mMin(minMass), mMax(maxMass) {}
0292 ~MassSelect() override {}
0293 bool accept(const BPHDecayMomentum& cand) const override {
0294 double mass = cand.composite().mass();
0295 return ((mass > mMin) && (mass < mMax));
0296 }
0297
0298 private:
0299 double mMin;
0300 double mMax;
0301 };
0302
0303
0304 class Chi2Select : public BPHVertexSelect {
0305 public:
0306 Chi2Select(double minProb) : mProb(minProb) {}
0307 ~Chi2Select() override {}
0308 bool accept(const BPHDecayVertex& cand) const override {
0309 const reco::Vertex& v = cand.vertex();
0310 if (v.isFake())
0311 return false;
0312 if (!v.isValid())
0313 return false;
0314 return (TMath::Prob(v.chi2(), lround(v.ndof())) > mProb);
0315 }
0316
0317 private:
0318 double mProb;
0319 };
0320
0321
0322
0323 outF << "build and dump JPsi" << endl;
0324 MuonPtSelect muPt(4.0);
0325 MuonEtaSelect muEta(2.1);
0326 string muPos = "MuPos";
0327 string muNeg = "MuNeg";
0328 BPHRecoBuilder bJPsi(es);
0329 if (usePM) {
0330 bJPsi.add(muPos, BPHRecoBuilder::createCollection(patMuon, "cfmig"), 0.105658);
0331 bJPsi.add(muNeg, BPHRecoBuilder::createCollection(patMuon, "cfmig"), 0.105658);
0332 } else if (useCC) {
0333 bJPsi.add(muPos, BPHRecoBuilder::createCollection(muDaugs, "cfmig"), 0.105658);
0334 bJPsi.add(muNeg, BPHRecoBuilder::createCollection(muDaugs, "cfmig"), 0.105658);
0335 }
0336 bJPsi.filter(muPos, muPt);
0337 bJPsi.filter(muNeg, muPt);
0338 bJPsi.filter(muPos, muEta);
0339 bJPsi.filter(muNeg, muEta);
0340
0341 MassSelect massJPsi(2.5, 3.7);
0342 Chi2Select chi2Valid(0.0);
0343 bJPsi.filter(massJPsi);
0344 bJPsi.filter(chi2Valid);
0345 vector<BPHPlusMinusConstCandPtr> lJPsi = BPHPlusMinusCandidate::build(bJPsi, muPos, muNeg, 3.096916, 0.00004);
0346
0347
0348
0349
0350
0351
0352
0353
0354
0355 int iJPsi;
0356 int nJPsi = lJPsi.size();
0357 outF << nJPsi << " JPsi cand found" << endl;
0358 for (iJPsi = 0; iJPsi < nJPsi; ++iJPsi)
0359 dumpRecoCand("JPsi", lJPsi[iJPsi].get());
0360
0361
0362
0363 outF << "build and dump Phi" << endl;
0364 BPHRecoBuilder bPhi(es);
0365 KaonChargeSelect tkPos(+1);
0366 KaonChargeSelect tkNeg(-1);
0367 KaonPtSelect tkPt(0.7);
0368 string kPos = "KPos";
0369 string kNeg = "KNeg";
0370 if (usePF) {
0371 bPhi.add(kPos, BPHRecoBuilder::createCollection(pfCands), 0.493677);
0372 bPhi.add(kNeg, BPHRecoBuilder::createCollection(pfCands), 0.493677);
0373 } else if (usePC) {
0374 bPhi.add(kPos, BPHRecoBuilder::createCollection(pcCands), 0.493677);
0375 bPhi.add(kNeg, BPHRecoBuilder::createCollection(pcCands), 0.493677);
0376 } else if (useGP) {
0377 bPhi.add(kPos, BPHRecoBuilder::createCollection(gpCands), 0.493677);
0378 bPhi.add(kNeg, BPHRecoBuilder::createCollection(gpCands), 0.493677);
0379 }
0380 bPhi.filter(kPos, tkPos);
0381 bPhi.filter(kNeg, tkNeg);
0382 bPhi.filter(kPos, tkPt);
0383 bPhi.filter(kNeg, tkPt);
0384
0385 MassSelect massPhi(1.00, 1.04);
0386 bPhi.filter(massPhi);
0387 vector<BPHRecoConstCandPtr> lPhi = BPHRecoCandidate::build(bPhi);
0388
0389
0390
0391
0392
0393
0394 int iPhi;
0395 int nPhi = lPhi.size();
0396 outF << nPhi << " Phi cand found" << endl;
0397 for (iPhi = 0; iPhi < nPhi; ++iPhi)
0398 dumpRecoCand("Phi", lPhi[iPhi].get());
0399
0400
0401
0402 if (nJPsi && nPhi) {
0403 outF << "build and dump Bs" << endl;
0404 BPHRecoBuilder bBs(es);
0405 bBs.setMinPDiffererence(1.0e-5);
0406 bBs.add("JPsi", lJPsi);
0407 bBs.add("Phi", lPhi);
0408 MassSelect mJPsi(2.946916, 3.246916);
0409 MassSelect mPhi(1.009461, 1.029461);
0410 bBs.filter("JPsi", mJPsi);
0411 bBs.filter("Phi", mPhi);
0412 Chi2Select chi2Bs(0.02);
0413 bBs.filter(chi2Bs);
0414 vector<BPHRecoConstCandPtr> lBs = BPHRecoCandidate::build(bBs);
0415 int iBs;
0416 int nBs = lBs.size();
0417 outF << nBs << " Bs cand found" << endl;
0418
0419 for (iBs = 0; iBs < nBs; ++iBs) {
0420
0421 const BPHRecoCandidate* cptr = lBs[iBs].get();
0422 cptr->kinematicTree("JPsi", 3.096916, 0.000040);
0423 }
0424 for (iBs = 0; iBs < nBs; ++iBs)
0425 dumpRecoCand("Bs", lBs[iBs].get());
0426 }
0427
0428
0429
0430 if (nJPsi && nrc) {
0431 outF << "build and dump Bu" << endl;
0432 BPHRecoBuilder bBu(es);
0433 bBu.setMinPDiffererence(1.0e-5);
0434 bBu.add("JPsi", lJPsi);
0435 if (usePF) {
0436 bBu.add("Kaon", BPHRecoBuilder::createCollection(pfCands), 0.493677);
0437 } else if (usePC) {
0438 bBu.add("Kaon", BPHRecoBuilder::createCollection(pcCands), 0.493677);
0439 } else if (useGP) {
0440 bBu.add("Kaon", BPHRecoBuilder::createCollection(gpCands), 0.493677);
0441 }
0442 MassSelect mJPsi(2.946916, 3.246916);
0443 KaonNeutralVeto knv;
0444 bBu.filter("JPsi", mJPsi);
0445 bBu.filter("Kaon", tkPt);
0446 bBu.filter("Kaon", knv);
0447 Chi2Select chi2Bu(0.02);
0448 bBu.filter(chi2Bu);
0449 vector<BPHRecoConstCandPtr> lBu = BPHRecoCandidate::build(bBu);
0450 int iBu;
0451 int nBu = lBu.size();
0452 outF << nBu << " Bu cand found" << endl;
0453
0454 for (iBu = 0; iBu < nBu; ++iBu) {
0455
0456 const BPHRecoCandidate* cptr = lBu[iBu].get();
0457 cptr->kinematicTree("JPsi", 3.096916, 0.000040);
0458 }
0459 for (iBu = 0; iBu < nBu; ++iBu)
0460 dumpRecoCand("Bu", lBu[iBu].get());
0461 }
0462
0463 return;
0464 }
0465
0466 void TestBPHRecoDecay::endJob() {
0467 *fPtr << "TestBPHRecoDecay::endJob" << endl;
0468 TDirectory* currentDir = gDirectory;
0469 TFile file(outHist.c_str(), "RECREATE");
0470 map<string, TH1F*>::iterator iter = histoMap.begin();
0471 map<string, TH1F*>::iterator iend = histoMap.end();
0472 while (iter != iend)
0473 iter++->second->Write();
0474 currentDir->cd();
0475 return;
0476 }
0477
0478 void TestBPHRecoDecay::dumpRecoCand(const string& name, const BPHRecoCandidate* cand) {
0479 fillHisto(name, cand);
0480
0481 ostream& outF = *fPtr;
0482
0483 static string cType = " cowboy";
0484 static string sType = " sailor";
0485 static string dType = "";
0486 string* type;
0487 const BPHPlusMinusCandidate* pmCand = dynamic_cast<const BPHPlusMinusCandidate*>(cand);
0488 if (pmCand != nullptr) {
0489 if (pmCand->isCowboy())
0490 type = &cType;
0491 else
0492 type = &sType;
0493 } else
0494 type = &dType;
0495
0496 outF << "****** " << name << " cand mass: " << cand->composite().mass() << " momentum " << cand->composite().px()
0497 << " " << cand->composite().py() << " " << cand->composite().pz() << *type << endl;
0498
0499 bool validFit = cand->isValidFit();
0500 const RefCountedKinematicTree kt = cand->kinematicTree();
0501 const RefCountedKinematicParticle kp = cand->currentParticle();
0502 if (validFit) {
0503 outF << "****** " << name << " constr mass: " << cand->p4().mass() << " momentum " << cand->p4().px() << " "
0504 << cand->p4().py() << " " << cand->p4().pz() << endl;
0505 }
0506
0507 const reco::Vertex& vx = cand->vertex();
0508 const reco::Vertex::Point& vp = vx.position();
0509 double chi2 = vx.chi2();
0510 int ndof = lround(vx.ndof());
0511 double prob = TMath::Prob(chi2, ndof);
0512 string tdca = "";
0513 if (pmCand != nullptr) {
0514 stringstream sstr;
0515 sstr << " - " << pmCand->cAppInRPhi().distance();
0516 tdca = sstr.str();
0517 }
0518 outF << "****** " << name << " vertex: " << vx.isFake() << " " << vx.isValid() << " - " << chi2 << " " << ndof << " "
0519 << prob << " - " << vp.X() << " " << vp.Y() << " " << vp.Z() << tdca << endl;
0520
0521 const vector<string>& dl = cand->daugNames();
0522 int i;
0523 int n = dl.size();
0524 for (i = 0; i < n; ++i) {
0525 const string& name = dl[i];
0526 const reco::Candidate* dp = cand->getDaug(name);
0527 GlobalPoint gp(vp.X(), vp.Y(), vp.Z());
0528 GlobalVector dm(0.0, 0.0, 0.0);
0529 const reco::TransientTrack* tt = cand->getTransientTrack(dp);
0530 if (tt != nullptr) {
0531 TrajectoryStateClosestToPoint tscp = tt->trajectoryStateClosestToPoint(gp);
0532 dm = tscp.momentum();
0533
0534
0535 }
0536 outF << "daughter " << i << " " << name << " " << (dp->charge() > 0 ? '+' : '-') << " momentum: " << dp->px() << " "
0537 << dp->py() << " " << dp->pz() << " - at vertex: " << dm.x() << " " << dm.y() << " " << dm.z() << endl;
0538 }
0539 const vector<string>& dc = cand->compNames();
0540 int j;
0541 int m = dc.size();
0542 for (j = 0; j < m; ++j) {
0543 const string& name = dc[j];
0544 const BPHRecoCandidate* dp = cand->getComp(name).get();
0545 outF << "composite daughter " << j << " " << name << " momentum: " << dp->composite().px() << " "
0546 << dp->composite().py() << " " << dp->composite().pz() << endl;
0547 }
0548
0549 if (validFit) {
0550 const RefCountedKinematicVertex kv = cand->currentDecayVertex();
0551 GlobalPoint gp = kv->position();
0552 outF << " kin fit vertex: " << gp.x() << " " << gp.y() << " " << gp.z() << endl;
0553 vector<RefCountedKinematicParticle> dk = kt->finalStateParticles();
0554 int k;
0555 int l = dk.size();
0556 for (k = 0; k < l; ++k) {
0557 const reco::TransientTrack& tt = dk[k]->refittedTransientTrack();
0558 TrajectoryStateClosestToPoint tscp = tt.trajectoryStateClosestToPoint(gp);
0559 GlobalVector dm = tscp.momentum();
0560
0561
0562 outF << "daughter " << k << " refitted: " << dm.x() << " " << dm.y() << " " << dm.z() << endl;
0563 }
0564 }
0565
0566 return;
0567 }
0568
0569 void TestBPHRecoDecay::fillHisto(const string& name, const BPHRecoCandidate* cand) {
0570 string mass = "mass";
0571 string mcst = "mcst";
0572 fillHisto(mass + name, cand->composite().mass());
0573 if (cand->isValidFit())
0574 fillHisto(mcst + name, cand->p4().mass());
0575
0576 const vector<string>& dc = cand->compNames();
0577 int i;
0578 int n = dc.size();
0579 for (i = 0; i < n; ++i) {
0580 const string& daug = dc[i];
0581 const BPHRecoCandidate* dptr = cand->getComp(daug).get();
0582 fillHisto(mass + name + daug, dptr->composite().mass());
0583 }
0584 return;
0585 }
0586
0587 void TestBPHRecoDecay::fillHisto(const string& name, float x) {
0588 map<string, TH1F*>::iterator iter = histoMap.find(name);
0589 map<string, TH1F*>::iterator iend = histoMap.end();
0590 if (iter == iend)
0591 return;
0592 iter->second->Fill(x);
0593 return;
0594 }
0595
0596 void TestBPHRecoDecay::createHisto(const string& name, int nbin, float hmin, float hmax) {
0597 histoMap[name] = new TH1F(name.c_str(), name.c_str(), nbin, hmin, hmax);
0598 return;
0599 }
0600
0601 #include "FWCore/Framework/interface/MakerMacros.h"
0602
0603 DEFINE_FWK_MODULE(TestBPHRecoDecay);