CMS 3D CMS Logo

JetCorrectorOnTheFly.cc
Go to the documentation of this file.
3 
17 
18 #include "TH1F.h"
19 
20 
21 using namespace edm;
22 using namespace reco;
23 using namespace std;
24 
25 //
26 // class declaration
27 //
28 template<class Jet>
30 public:
31  explicit JetCorrectorOnTheFly(const edm::ParameterSet&);
32  ~JetCorrectorOnTheFly() override;
33 
34 private:
35  typedef std::vector<Jet> JetCollection;
36  void beginJob() override ;
37  void analyze(const edm::Event&, const edm::EventSetup&) override;
38  void endJob() override ;
39 
44  double mMinRawJetPt;
45  bool mDebug;
46  TH1F *mRawPt,*mCorPt;
47 };
48 //
49 //----------- Class Implementation ------------------------------------------
50 //
51 //---------------------------------------------------------------------------
52 template<class Jet>
54 {
55  mJetCorrector = consumes<reco::JetCorrector>(iConfig.getParameter<edm::InputTag>("JetCorrector"));
56  mJetToken = consumes<JetCollection>(edm::InputTag(iConfig.getParameter<edm::InputTag> ("JetCollectionName")));
57  mVertexToken = consumes<reco::VertexCollection>(edm::InputTag("offlinePrimaryVertices"));
58  mMinRawJetPt = iConfig.getParameter<double> ("MinRawJetPt");
59  mDebug = iConfig.getParameter<bool> ("Debug");
60 }
61 //---------------------------------------------------------------------------
62 template<class Jet>
64 {
65 
66 }
67 //---------------------------------------------------------------------------
68 template<class Jet>
70 {
72  iEvent.getByToken(mJetCorrector, corrector);
74  iEvent.getByToken(mJetToken,jets);
75 
77  iEvent.getByToken(mVertexToken,recVtxs);
78  int NPV(0);
79  for(unsigned int ind=0;ind<recVtxs->size();ind++) {
80  if (!((*recVtxs)[ind].isFake())) {
81  NPV++;
82  }
83  }
84  typename JetCollection::const_iterator i_jet;
86  for(i_jet = jets->begin(); i_jet != jets->end(); i_jet++) {
87  if (i_jet->pt() < mMinRawJetPt) continue;
88  //double scale = corrector->correction(i_jet->p4());
89  double scale = corrector->correction(*i_jet);
90  if (mDebug) {
91  std::cout<<"energy = "<<i_jet->energy()<<", "
92  <<"eta = "<<i_jet->eta()<<", "
93  <<"raw pt = "<<i_jet->pt()<<", "
94  <<"NPV = "<<NPV<<", "
95  <<"correction = "<<scale<<", "
96  <<"cor pt = "<<scale*i_jet->pt()<<endl;
97  }
98  mRawPt->Fill(i_jet->pt());
99  mCorPt->Fill(scale*i_jet->pt());
100  }
101 
102 
103 }
104 //---------------------------------------------------------------------------
105 template<class Jet>
107 {
108  mRawPt = fs->make<TH1F>("RawJetPt","RawJetPt",1000,0,1000);
109  mCorPt = fs->make<TH1F>("CorJetPt","CorJetPt",1000,0,1000);
110 }
111 //---------------------------------------------------------------------------
112 template<class Jet>
114 {
115 
116 }
117 //---------------------------------------------------------------------------
120 
123 
126 
127 
128 
129 
130 
131 
132 
133 
134 
135 
136 
137 
138 
139 
140 
141 
142 
143 
144 
145 
146 
147 
148 
T getParameter(std::string const &) const
edm::EDGetTokenT< JetCollection > mJetToken
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
std::vector< Jet > JetCollection
double correction(const LorentzVector &fJet) const
get correction using Jet information only
Definition: JetCorrector.h:49
JetCorrectorOnTheFly(const edm::ParameterSet &)
edm::EDGetTokenT< reco::VertexCollection > mVertexToken
example_stream void analyze(const edm::Event &, const edm::EventSetup &) override
void beginJob()
Definition: Breakpoints.cc:14
edm::EDGetTokenT< reco::JetCorrector > mJetCorrector
int iEvent
Definition: GenABIO.cc:224
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
JetCorrectorOnTheFly< PFJet > PFJetCorrectorOnTheFly
edm::Service< TFileService > fs
vector< PseudoJet > jets
JetCorrectorOnTheFly< JPTJet > JPTJetCorrectorOnTheFly
void analyze(const edm::Event &, const edm::EventSetup &) override
JetCorrectorOnTheFly< CaloJet > CaloJetCorrectorOnTheFly
fixed size matrix
HLT enums.