JetMCAlgos/test/jetMatch.cc

0001 // system include files
0002 #include <memory>
0003 #include <string>
0004 #include <iostream>
0005 #include <algorithm>
0006 #include <vector>
0007
0008 // user include files
0009 #include "FWCore/Framework/interface/Frameworkfwd.h"
0010 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/EventSetup.h"
0013 #include "FWCore/Framework/interface/MakerMacros.h"
0014
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "DataFormats/Common/interface/Ref.h"
0017 #include "DataFormats/JetReco/interface/Jet.h"
0018 #include "DataFormats/JetReco/interface/CaloJet.h"
0019 #include "DataFormats/JetReco/interface/GenJet.h"
0020
0021 #include "DataFormats/Candidate/interface/CandMatchMap.h"
0022 #include "DataFormats/Candidate/interface/CandMatchMapMany.h"
0023
0024 #include "TFile.h"
0025 #include "TH1.h"
0026
0027 #include <Math/VectorUtil.h>
0028
0029 using namespace std;
0030 using namespace reco;
0031 using namespace edm;
0032 using namespace ROOT::Math::VectorUtil;
0033
0034 class jetMatch : public edm::one::EDAnalyzer<> {
0035 public:
0036   explicit jetMatch(const edm::ParameterSet&);
0037   void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) override;
0038   void beginJob() override;
0039   void endJob() override;
0040
0041 private:
0042   EDGetTokenT<CandidateCollection> sourceToken_;
0043   EDGetTokenT<CandidateCollection> matchedToken_;
0044   EDGetTokenT<CandViewMatchMap> matchedjetsOneToken_;
0045   EDGetTokenT<CandMatchMapMany> matchedjetsManyToken_;
0046   string fOutputFileName_;
0047
0048   Handle<CandidateCollection> source;
0049   Handle<CandidateCollection> matched;
0050   Handle<CandViewMatchMap> matchedjetsOne;
0051   Handle<CandMatchMapMany> matchedjetsMany;
0052
0053   TFile* hOutputFile;
0054   TH1D* hTotalLenght;
0055 };
0056
0057 jetMatch::jetMatch(const edm::ParameterSet& iConfig) {
0058   sourceToken_ = consumes<CandidateCollection>(iConfig.getParameter<InputTag>("src"));
0059   matchedToken_ = consumes<CandidateCollection>(iConfig.getParameter<InputTag>("matched"));
0060   matchedjetsOneToken_ = consumes<CandViewMatchMap>(iConfig.getParameter<InputTag>("matchMapOne"));
0061   matchedjetsManyToken_ = consumes<CandMatchMapMany>(iConfig.getParameter<InputTag>("matchMapMany"));
0062   fOutputFileName_ = iConfig.getUntrackedParameter<string>("HistOutFile", std::string("testMatch.root"));
0063 }
0064
0065 void jetMatch::beginJob() {
0066   hOutputFile = new TFile(fOutputFileName_.c_str(), "RECREATE");
0067   hTotalLenght = new TH1D("hTotalLenght", "Total Lenght", 100, 0., 5.);
0068   return;
0069 }
0070
0071 void jetMatch::endJob() {
0072   hOutputFile->Write();
0073   hOutputFile->Close();
0074   return;
0075 }
0076
0077 void jetMatch::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0078   cout << "[GenJetTest] analysing event " << iEvent.id() << endl;
0079
0080   try {
0081     iEvent.getByToken(sourceToken_, source);
0082     iEvent.getByToken(matchedToken_, matched);
0083     iEvent.getByToken(matchedjetsOneToken_, matchedjetsOne);
0084     iEvent.getByToken(matchedjetsManyToken_, matchedjetsMany);
0085   } catch (std::exception& ce) {
0086     cerr << "[GenJetTest] caught std::exception " << ce.what() << endl;
0087     return;
0088   }
0089
0090   //
0091   // Printout for OneToOne matcher
0092   //
0093   double dR = -1.;
0094   cout << "**********************" << endl;
0095   cout << "* OneToOne Printout  *" << endl;
0096   cout << "**********************" << endl;
0097   for (CandViewMatchMap::const_iterator f = matchedjetsOne->begin(); f != matchedjetsOne->end(); f++) {
0098     const Candidate* sourceRef = &*(f->key);
0099     const Candidate* matchRef = &*(f->val);
0100     dR = DeltaR(sourceRef->p4(), matchRef->p4());
0101
0102     printf("[GenJetTest] (pt,eta,phi) source = %6.2f %5.2f %5.2f matched = %6.2f %5.2f %5.2f dR=%5.3f\n",
0103            sourceRef->et(),
0104            sourceRef->eta(),
0105            sourceRef->phi(),
0106            matchRef->et(),
0107            matchRef->eta(),
0108            matchRef->phi(),
0109            dR);
0110   }
0111
0112   hTotalLenght->Fill(dR);
0113
0114   //
0115   // Printout for OneToMany matcher
0116   //
0117   cout << "**********************" << endl;
0118   cout << "* OneToMany Printout *" << endl;
0119   cout << "**********************" << endl;
0120   for (CandMatchMapMany::const_iterator f = matchedjetsMany->begin(); f != matchedjetsMany->end(); f++) {
0121     const Candidate* sourceRef = &*(f->key);
0122
0123     printf(
0124         "[GenJetTest] (pt,eta,phi) source = %6.2f %5.2f %5.2f\n", sourceRef->et(), sourceRef->eta(), sourceRef->phi());
0125
0126     const vector<pair<CandidateRef, double> > vectMatched = f->val;
0127     vector<pair<CandidateRef, double> >::const_iterator matchIT;
0128     for (matchIT = vectMatched.begin(); matchIT != vectMatched.end(); matchIT++) {
0129       const Candidate* matchedRef = &*((*matchIT).first);
0130       double deltaR = (*matchIT).second;
0131       printf("             (pt,eta,phi) matched = %6.2f %5.2f %5.2f - dR=%5.2f\n",
0132              matchedRef->et(),
0133              matchedRef->eta(),
0134              matchedRef->phi(),
0135              deltaR);
0136     }
0137   }
0138 }
0139
0140 DEFINE_FWK_MODULE(jetMatch);