File indexing completed on 2024-04-06 12:23:34
0001
0002 #include <memory>
0003 #include <string>
0004 #include <iostream>
0005 #include <algorithm>
0006 #include <vector>
0007
0008
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
0023 #include "TFile.h"
0024 #include "TH1.h"
0025
0026 #include <Math/VectorUtil.h>
0027
0028 using namespace std;
0029 using namespace reco;
0030 using namespace edm;
0031 using namespace ROOT::Math::VectorUtil;
0032
0033 class matchOneToOne : public edm::one::EDAnalyzer<> {
0034 public:
0035 explicit matchOneToOne(const edm::ParameterSet&);
0036 ~matchOneToOne() {}
0037 virtual void analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup);
0038 virtual void beginJob();
0039 virtual void endJob();
0040
0041 private:
0042 EDGetTokenT<CandidateCollection> sourceToken_;
0043 EDGetTokenT<CandidateCollection> matchedToken_;
0044 EDGetTokenT<CandViewMatchMap> matchedjetsOne1Token_;
0045 EDGetTokenT<CandViewMatchMap> matchedjetsOne2Token_;
0046 string fOutputFileName_;
0047
0048 Handle<CandidateCollection> source;
0049 Handle<CandidateCollection> matched;
0050 Handle<CandViewMatchMap> matchedjetsOne1;
0051 Handle<CandViewMatchMap> matchedjetsOne2;
0052
0053 TFile* hOutputFile;
0054 TH1D* hTotalLenght;
0055 TH1F* hDR;
0056 };
0057
0058 matchOneToOne::matchOneToOne(const edm::ParameterSet& iConfig) {
0059 sourceToken_ = consumes<CandidateCollection>(iConfig.getParameter<InputTag>("src"));
0060 matchedToken_ = consumes<CandidateCollection>(iConfig.getParameter<InputTag>("matched"));
0061 matchedjetsOne1Token_ = consumes<CandViewMatchMap>(iConfig.getParameter<InputTag>("matchMapOne1"));
0062 matchedjetsOne2Token_ = consumes<CandViewMatchMap>(iConfig.getParameter<InputTag>("matchMapOne2"));
0063 fOutputFileName_ = iConfig.getUntrackedParameter<string>("HistOutFile", std::string("testMatch.root"));
0064 }
0065
0066 void matchOneToOne::beginJob() {
0067 hOutputFile = new TFile(fOutputFileName_.c_str(), "RECREATE");
0068 hDR = new TH1F("hDR", "", 1000, 0., 10.);
0069 hTotalLenght = new TH1D("hTotalLenght", "Total Lenght", 1000, 0., 5.);
0070 return;
0071 }
0072
0073 void matchOneToOne::endJob() {
0074 hOutputFile->Write();
0075 hOutputFile->Close();
0076 return;
0077 }
0078
0079 void matchOneToOne::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0080 cout << "[matchOneToOne] analysing event " << iEvent.id() << endl;
0081
0082 try {
0083 iEvent.getByToken(sourceToken_, source);
0084 iEvent.getByToken(matchedToken_, matched);
0085 iEvent.getByToken(matchedjetsOne1Token_, matchedjetsOne1);
0086 iEvent.getByToken(matchedjetsOne2Token_, matchedjetsOne2);
0087 } catch (std::exception& ce) {
0088 cerr << "[matchOneToOne] caught std::exception " << ce.what() << endl;
0089 return;
0090 }
0091
0092
0093
0094
0095 double dR = -1.;
0096 float totalLenght = 0;
0097 cout << "**********************" << endl;
0098 cout << "* OneToOne Printout *" << endl;
0099 cout << "**********************" << endl;
0100 for (CandViewMatchMap::const_iterator f = matchedjetsOne1->begin(); f != matchedjetsOne1->end(); f++) {
0101 const Candidate* sourceRef = &*(f->key);
0102 const Candidate* matchRef = &*(f->val);
0103 dR = DeltaR(sourceRef->p4(), matchRef->p4());
0104 totalLenght += dR;
0105
0106 printf("[matchOneToOne src2mtc] (pt,eta,phi) source = %6.2f %5.2f %5.2f matched = %6.2f %5.2f %5.2f dR=%5.3f\n",
0107 sourceRef->et(),
0108 sourceRef->eta(),
0109 sourceRef->phi(),
0110 matchRef->et(),
0111 matchRef->eta(),
0112 matchRef->phi(),
0113 dR);
0114
0115 hDR->Fill(dR);
0116 }
0117
0118 cout << "-----------------" << endl;
0119
0120 for (CandViewMatchMap::const_iterator f = matchedjetsOne2->begin(); f != matchedjetsOne2->end(); f++) {
0121 const Candidate* sourceRef = &*(f->key);
0122 const Candidate* matchRef = &*(f->val);
0123 dR = DeltaR(sourceRef->p4(), matchRef->p4());
0124 printf("[matchOneToOne mtc2src] (pt,eta,phi) source = %6.2f %5.2f %5.2f matched = %6.2f %5.2f %5.2f dR=%5.3f\n",
0125 sourceRef->et(),
0126 sourceRef->eta(),
0127 sourceRef->phi(),
0128 matchRef->et(),
0129 matchRef->eta(),
0130 matchRef->phi(),
0131 dR);
0132 }
0133
0134 hTotalLenght->Fill(totalLenght);
0135 }
0136
0137 DEFINE_FWK_MODULE(matchOneToOne);