ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/csander/JetResolutionFromMC/MCResolutions/src/MCResolutions.cc
Revision: 1.10
Committed: Sat May 14 12:10:16 2011 UTC (13 years, 11 months ago) by csander
Content type: text/plain
Branch: MAIN
Changes since 1.9: +1006 -1153 lines
Log Message:
small improvements

File Contents

# User Rev Content
1 csander 1.1 // -*- C++ -*-
2     //
3     // Package: MCResolutions
4     // Class: MCResolutions
5 csander 1.2 //
6 csander 1.1 /**\class MCResolutions MCResolutions.cc JetResolutionFromMC/MCResolutions/src/MCResolutions.cc
7 csander 1.2
8 csander 1.1 Description: [one line class summary]
9 csander 1.2
10 csander 1.1 Implementation:
11 csander 1.10 [Notes on implementation]
12     */
13 csander 1.1 //
14     // Original Author: Christian Sander,,,
15     // Created: Wed Oct 6 18:22:21 CEST 2010
16 csander 1.10 // $Id: MCResolutions.cc,v 1.9 2011/01/21 08:02:26 csander Exp $
17 csander 1.1 //
18     //
19    
20 csander 1.2 #include "JetResolutionFromMC/MCResolutions/interface/MCResolutions.h"
21 csander 1.1
22     //
23     // constructors and destructor
24     //
25 csander 1.3 MCResolutions::MCResolutions(const edm::ParameterSet& iConfig) {
26 csander 1.10 //now do what ever initialization is needed
27     _jetTag = iConfig.getParameter<edm::InputTag> ("jetTag");
28     _muonTag = iConfig.getParameter<edm::InputTag> ("muonTag");
29     _genJetTag = iConfig.getParameter<edm::InputTag> ("genJetTag");
30     _weightName = iConfig.getParameter<edm::InputTag> ("weightName");
31     _bTag = iConfig.getParameter<edm::InputTag> ("bTag");
32     _EBRecHits = iConfig.getParameter<edm::InputTag> ("EBRecHits");
33     _EERecHits = iConfig.getParameter<edm::InputTag> ("EERecHits");
34     _jetMultPtCut = iConfig.getParameter<double> ("jetMultPtCut");
35     _jetMultEtaCut = iConfig.getParameter<double> ("jetMultEtaCut");
36     _deltaPhiDiJet = iConfig.getParameter<double> ("deltaPhiDiJet");
37     _absCut3rdJet = iConfig.getParameter<double> ("absCut3rdJet");
38     _relCut3rdJet = iConfig.getParameter<double> ("relCut3rdJet");
39     _deltaRMatch = iConfig.getParameter<double> ("deltaRMatch");
40     _deltaRMatchVeto = iConfig.getParameter<double> ("deltaRMatchVeto");
41     _absPtVeto = iConfig.getParameter<double> ("absPtVeto");
42     _relPtVeto = iConfig.getParameter<double> ("relPtVeto");
43     _deltaRDeadECal = iConfig.getParameter<double> ("deltaRDeadECal");
44     _GenJetPtCut = iConfig.getParameter<double> ("GenJetPtCut");
45     _bTagCut = iConfig.getParameter<double> ("bTagCut");
46     _bTagDeltaR = iConfig.getParameter<double> ("bTagDelatR");
47     _Bid = iConfig.getParameter<std::string> ("Bid");
48     _maskedEcalChannelStatusThreshold = iConfig.getParameter<int> ("maskedEcalChannelStatusThreshold");
49     _fileName = iConfig.getParameter<std::string> ("fileName");
50    
51     hfile = new TFile(_fileName.c_str(), "RECREATE", "Jet response in pT and eta bins");
52     //hfile->mkdir(_dirName.c_str(), _dirName.c_str());
53    
54     EBinEdges.push_back(0);
55     EBinEdges.push_back(20);
56     EBinEdges.push_back(30);
57     EBinEdges.push_back(50);
58     EBinEdges.push_back(80);
59     EBinEdges.push_back(120);
60     EBinEdges.push_back(170);
61     EBinEdges.push_back(230);
62     EBinEdges.push_back(300);
63     EBinEdges.push_back(380);
64     EBinEdges.push_back(470);
65     EBinEdges.push_back(570);
66     EBinEdges.push_back(680);
67     EBinEdges.push_back(800);
68     EBinEdges.push_back(1000);
69     EBinEdges.push_back(1300);
70     EBinEdges.push_back(1700);
71     EBinEdges.push_back(2200);
72    
73     PtBinEdges.push_back(0);
74     PtBinEdges.push_back(20);
75     PtBinEdges.push_back(30);
76     PtBinEdges.push_back(50);
77     PtBinEdges.push_back(80);
78     PtBinEdges.push_back(120);
79     PtBinEdges.push_back(170);
80     PtBinEdges.push_back(230);
81     PtBinEdges.push_back(300);
82     PtBinEdges.push_back(380);
83     PtBinEdges.push_back(470);
84     PtBinEdges.push_back(570);
85     PtBinEdges.push_back(680);
86     PtBinEdges.push_back(800);
87     PtBinEdges.push_back(1000);
88     PtBinEdges.push_back(1300);
89     PtBinEdges.push_back(1700);
90     PtBinEdges.push_back(2200);
91    
92     EtaBinEdges.push_back(0.0);
93     EtaBinEdges.push_back(0.5);
94     EtaBinEdges.push_back(1.1);
95     EtaBinEdges.push_back(1.7);
96     EtaBinEdges.push_back(2.3);
97     EtaBinEdges.push_back(3.2);
98     EtaBinEdges.push_back(5.0);
99    
100     //// Array of histograms for jet resolutions (all jet multiplicities)
101     h_tot_DiJet_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
102     h_b_DiJet_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
103     h_nob_DiJet_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
104     h_dead_DiJet_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
105     h_deadb_DiJet_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
106     for (std::vector<std::vector<TH1F*> >::iterator it = h_tot_DiJet_JetResPt_Pt.begin(); it
107     != h_tot_DiJet_JetResPt_Pt.end(); ++it) {
108     it->resize(PtBinEdges.size() - 1);
109     }
110     for (std::vector<std::vector<TH1F*> >::iterator it = h_b_DiJet_JetResPt_Pt.begin(); it
111     != h_b_DiJet_JetResPt_Pt.end(); ++it) {
112     it->resize(PtBinEdges.size() - 1);
113     }
114     for (std::vector<std::vector<TH1F*> >::iterator it = h_nob_DiJet_JetResPt_Pt.begin(); it
115     != h_nob_DiJet_JetResPt_Pt.end(); ++it) {
116     it->resize(PtBinEdges.size() - 1);
117     }
118     for (std::vector<std::vector<TH1F*> >::iterator it = h_dead_DiJet_JetResPt_Pt.begin(); it
119     != h_dead_DiJet_JetResPt_Pt.end(); ++it) {
120     it->resize(PtBinEdges.size() - 1);
121     }
122     for (std::vector<std::vector<TH1F*> >::iterator it = h_deadb_DiJet_JetResPt_Pt.begin(); it
123     != h_deadb_DiJet_JetResPt_Pt.end(); ++it) {
124     it->resize(PtBinEdges.size() - 1);
125     }
126     // h_tot_DiJet_JetResPt_E.resize(EtaBinEdges.size() - 1);
127     // h_b_DiJet_JetResPt_E.resize(EtaBinEdges.size() - 1);
128     // h_nob_DiJet_JetResPt_E.resize(EtaBinEdges.size() - 1);
129     // h_dead_DiJet_JetResPt_E.resize(EtaBinEdges.size() - 1);
130     // h_deadb_DiJet_JetResPt_E.resize(EtaBinEdges.size() - 1);
131     // for (std::vector<std::vector<TH1F*> >::iterator it = h_tot_DiJet_JetResPt_E.begin(); it != h_tot_DiJet_JetResPt_E.end(); ++it) {
132     // it->resize(EBinEdges.size() - 1);
133     // }
134     // for (std::vector<std::vector<TH1F*> >::iterator it = h_b_DiJet_JetResPt_E.begin(); it != h_b_DiJet_JetResPt_E.end(); ++it) {
135     // it->resize(EBinEdges.size() - 1);
136     // }
137     // for (std::vector<std::vector<TH1F*> >::iterator it = h_nob_DiJet_JetResPt_E.begin(); it != h_nob_DiJet_JetResPt_E.end(); ++it) {
138     // it->resize(EBinEdges.size() - 1);
139     // }
140     // for (std::vector<std::vector<TH1F*> >::iterator it = h_dead_DiJet_JetResPt_E.begin(); it != h_dead_DiJet_JetResPt_E.end(); ++it) {
141     // it->resize(EBinEdges.size() - 1);
142     // }
143     // for (std::vector<std::vector<TH1F*> >::iterator it = h_deadb_DiJet_JetResPt_E.begin(); it != h_deadb_DiJet_JetResPt_E.end(); ++it) {
144     // it->resize(EBinEdges.size() - 1);
145     // }
146    
147     //// Array of histograms for jet resolutions (all jet multiplicities)
148     h_tot_NJetAll_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
149     h_b_NJetAll_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
150     h_nob_NJetAll_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
151     h_dead_NJetAll_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
152     h_deadb_NJetAll_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
153     for (std::vector<std::vector<TH1F*> >::iterator it = h_tot_NJetAll_JetResPt_Pt.begin(); it
154     != h_tot_NJetAll_JetResPt_Pt.end(); ++it) {
155     it->resize(PtBinEdges.size() - 1);
156     }
157     for (std::vector<std::vector<TH1F*> >::iterator it = h_b_NJetAll_JetResPt_Pt.begin(); it
158     != h_b_NJetAll_JetResPt_Pt.end(); ++it) {
159     it->resize(PtBinEdges.size() - 1);
160     }
161     for (std::vector<std::vector<TH1F*> >::iterator it = h_nob_NJetAll_JetResPt_Pt.begin(); it
162     != h_nob_NJetAll_JetResPt_Pt.end(); ++it) {
163     it->resize(PtBinEdges.size() - 1);
164     }
165     for (std::vector<std::vector<TH1F*> >::iterator it = h_dead_NJetAll_JetResPt_Pt.begin(); it
166     != h_dead_NJetAll_JetResPt_Pt.end(); ++it) {
167     it->resize(PtBinEdges.size() - 1);
168     }
169     for (std::vector<std::vector<TH1F*> >::iterator it = h_deadb_NJetAll_JetResPt_Pt.begin(); it
170     != h_deadb_NJetAll_JetResPt_Pt.end(); ++it) {
171     it->resize(PtBinEdges.size() - 1);
172     }
173     // h_tot_NJetAll_JetResPt_E.resize(EtaBinEdges.size() - 1);
174     // h_b_NJetAll_JetResPt_E.resize(EtaBinEdges.size() - 1);
175     // h_nob_NJetAll_JetResPt_E.resize(EtaBinEdges.size() - 1);
176     // h_dead_NJetAll_JetResPt_E.resize(EtaBinEdges.size() - 1);
177     // h_deadb_NJetAll_JetResPt_E.resize(EtaBinEdges.size() - 1);
178     // for (std::vector<std::vector<TH1F*> >::iterator it = h_tot_NJetAll_JetResPt_E.begin(); it != h_tot_NJetAll_JetResPt_E.end(); ++it) {
179     // it->resize(EBinEdges.size() - 1);
180     // }
181     // for (std::vector<std::vector<TH1F*> >::iterator it = h_b_NJetAll_JetResPt_E.begin(); it != h_b_NJetAll_JetResPt_E.end(); ++it) {
182     // it->resize(EBinEdges.size() - 1);
183     // }
184     // for (std::vector<std::vector<TH1F*> >::iterator it = h_nob_NJetAll_JetResPt_E.begin(); it != h_nob_NJetAll_JetResPt_E.end(); ++it) {
185     // it->resize(EBinEdges.size() - 1);
186     // }
187     // for (std::vector<std::vector<TH1F*> >::iterator it = h_dead_NJetAll_JetResPt_E.begin(); it != h_dead_NJetAll_JetResPt_E.end(); ++it) {
188     // it->resize(EBinEdges.size() - 1);
189     // }
190     // for (std::vector<std::vector<TH1F*> >::iterator it = h_deadb_NJetAll_JetResPt_E.begin(); it != h_deadb_NJetAll_JetResPt_E.end(); ++it) {
191     // it->resize(EBinEdges.size() - 1);
192     // }
193    
194     //// Array of histograms for jet resolutions (jet multiplicities = 2)
195     // h_tot_NJet2_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
196     // h_b_NJet2_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
197     // h_nob_NJet2_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
198     // h_dead_NJet2_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
199     // h_deadb_NJet2_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
200     // for (std::vector<std::vector<TH1F*> >::iterator it = h_tot_NJet2_JetResPt_Pt.begin(); it != h_tot_NJet2_JetResPt_Pt.end(); ++it) {
201     // it->resize(PtBinEdges.size() - 1);
202     // }
203     // for (std::vector<std::vector<TH1F*> >::iterator it = h_b_NJet2_JetResPt_Pt.begin(); it != h_b_NJet2_JetResPt_Pt.end(); ++it) {
204     // it->resize(PtBinEdges.size() - 1);
205     // }
206     // for (std::vector<std::vector<TH1F*> >::iterator it = h_nob_NJet2_JetResPt_Pt.begin(); it != h_nob_NJet2_JetResPt_Pt.end(); ++it) {
207     // it->resize(PtBinEdges.size() - 1);
208     // }
209     // for (std::vector<std::vector<TH1F*> >::iterator it = h_dead_NJet2_JetResPt_Pt.begin(); it != h_dead_NJet2_JetResPt_Pt.end(); ++it) {
210     // it->resize(PtBinEdges.size() - 1);
211     // }
212     // for (std::vector<std::vector<TH1F*> >::iterator it = h_deadb_NJet2_JetResPt_Pt.begin(); it != h_deadb_NJet2_JetResPt_Pt.end(); ++it) {
213     // it->resize(PtBinEdges.size() - 1);
214     // }
215     // h_tot_NJet2_JetResPt_E.resize(EtaBinEdges.size() - 1);
216     // h_b_NJet2_JetResPt_E.resize(EtaBinEdges.size() - 1);
217     // h_nob_NJet2_JetResPt_E.resize(EtaBinEdges.size() - 1);
218     // h_dead_NJet2_JetResPt_E.resize(EtaBinEdges.size() - 1);
219     // h_deadb_NJet2_JetResPt_E.resize(EtaBinEdges.size() - 1);
220     // for (std::vector<std::vector<TH1F*> >::iterator it = h_tot_NJet2_JetResPt_E.begin(); it != h_tot_NJet2_JetResPt_E.end(); ++it) {
221     // it->resize(EBinEdges.size() - 1);
222     // }
223     // for (std::vector<std::vector<TH1F*> >::iterator it = h_b_NJet2_JetResPt_E.begin(); it != h_b_NJet2_JetResPt_E.end(); ++it) {
224     // it->resize(EBinEdges.size() - 1);
225     // }
226     // for (std::vector<std::vector<TH1F*> >::iterator it = h_nob_NJet2_JetResPt_E.begin(); it != h_nob_NJet2_JetResPt_E.end(); ++it) {
227     // it->resize(EBinEdges.size() - 1);
228     // }
229     // for (std::vector<std::vector<TH1F*> >::iterator it = h_dead_NJet2_JetResPt_E.begin(); it != h_dead_NJet2_JetResPt_E.end(); ++it) {
230     // it->resize(EBinEdges.size() - 1);
231     // }
232     // for (std::vector<std::vector<TH1F*> >::iterator it = h_deadb_NJet2_JetResPt_E.begin(); it != h_deadb_NJet2_JetResPt_E.end(); ++it) {
233     // it->resize(EBinEdges.size() - 1);
234     // }
235    
236     //// Array of histograms for jet resolutions (jet multiplicities = 3)
237     // h_tot_NJet3_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
238     // h_b_NJet3_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
239     // h_nob_NJet3_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
240     // h_dead_NJet3_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
241     // h_deadb_NJet3_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
242     // for (std::vector<std::vector<TH1F*> >::iterator it = h_tot_NJet3_JetResPt_Pt.begin(); it != h_tot_NJet3_JetResPt_Pt.end(); ++it) {
243     // it->resize(PtBinEdges.size() - 1);
244     // }
245     // for (std::vector<std::vector<TH1F*> >::iterator it = h_b_NJet3_JetResPt_Pt.begin(); it != h_b_NJet3_JetResPt_Pt.end(); ++it) {
246     // it->resize(PtBinEdges.size() - 1);
247     // }
248     // for (std::vector<std::vector<TH1F*> >::iterator it = h_nob_NJet3_JetResPt_Pt.begin(); it != h_nob_NJet3_JetResPt_Pt.end(); ++it) {
249     // it->resize(PtBinEdges.size() - 1);
250     // }
251     // for (std::vector<std::vector<TH1F*> >::iterator it = h_dead_NJet3_JetResPt_Pt.begin(); it != h_dead_NJet3_JetResPt_Pt.end(); ++it) {
252     // it->resize(PtBinEdges.size() - 1);
253     // }
254     // for (std::vector<std::vector<TH1F*> >::iterator it = h_deadb_NJet3_JetResPt_Pt.begin(); it != h_deadb_NJet3_JetResPt_Pt.end(); ++it) {
255     // it->resize(PtBinEdges.size() - 1);
256     // }
257     // h_tot_NJet3_JetResPt_E.resize(EtaBinEdges.size() - 1);
258     // h_b_NJet3_JetResPt_E.resize(EtaBinEdges.size() - 1);
259     // h_nob_NJet3_JetResPt_E.resize(EtaBinEdges.size() - 1);
260     // h_dead_NJet3_JetResPt_E.resize(EtaBinEdges.size() - 1);
261     // h_deadb_NJet3_JetResPt_E.resize(EtaBinEdges.size() - 1);
262     // for (std::vector<std::vector<TH1F*> >::iterator it = h_tot_NJet3_JetResPt_E.begin(); it != h_tot_NJet3_JetResPt_E.end(); ++it) {
263     // it->resize(EBinEdges.size() - 1);
264     // }
265     // for (std::vector<std::vector<TH1F*> >::iterator it = h_b_NJet3_JetResPt_E.begin(); it != h_b_NJet3_JetResPt_E.end(); ++it) {
266     // it->resize(EBinEdges.size() - 1);
267     // }
268     // for (std::vector<std::vector<TH1F*> >::iterator it = h_nob_NJet3_JetResPt_E.begin(); it != h_nob_NJet3_JetResPt_E.end(); ++it) {
269     // it->resize(EBinEdges.size() - 1);
270     // }
271     // for (std::vector<std::vector<TH1F*> >::iterator it = h_dead_NJet3_JetResPt_E.begin(); it != h_dead_NJet3_JetResPt_E.end(); ++it) {
272     // it->resize(EBinEdges.size() - 1);
273     // }
274     // for (std::vector<std::vector<TH1F*> >::iterator it = h_deadb_NJet3_JetResPt_E.begin(); it != h_deadb_NJet3_JetResPt_E.end(); ++it) {
275     // it->resize(EBinEdges.size() - 1);
276     // }
277    
278     //// Array of histograms for jet resolutions (jet multiplicities = 4)
279     // h_tot_NJet4_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
280     // h_b_NJet4_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
281     // h_nob_NJet4_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
282     // h_dead_NJet4_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
283     // h_deadb_NJet4_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
284     // for (std::vector<std::vector<TH1F*> >::iterator it = h_tot_NJet4_JetResPt_Pt.begin(); it != h_tot_NJet4_JetResPt_Pt.end(); ++it) {
285     // it->resize(PtBinEdges.size() - 1);
286     // }
287     // for (std::vector<std::vector<TH1F*> >::iterator it = h_b_NJet4_JetResPt_Pt.begin(); it != h_b_NJet4_JetResPt_Pt.end(); ++it) {
288     // it->resize(PtBinEdges.size() - 1);
289     // }
290     // for (std::vector<std::vector<TH1F*> >::iterator it = h_nob_NJet4_JetResPt_Pt.begin(); it != h_nob_NJet4_JetResPt_Pt.end(); ++it) {
291     // it->resize(PtBinEdges.size() - 1);
292     // }
293     // for (std::vector<std::vector<TH1F*> >::iterator it = h_dead_NJet4_JetResPt_Pt.begin(); it != h_dead_NJet4_JetResPt_Pt.end(); ++it) {
294     // it->resize(PtBinEdges.size() - 1);
295     // }
296     // for (std::vector<std::vector<TH1F*> >::iterator it = h_deadb_NJet4_JetResPt_Pt.begin(); it != h_deadb_NJet4_JetResPt_Pt.end(); ++it) {
297     // it->resize(PtBinEdges.size() - 1);
298     // }
299     // h_tot_NJet4_JetResPt_E.resize(EtaBinEdges.size() - 1);
300     // h_b_NJet4_JetResPt_E.resize(EtaBinEdges.size() - 1);
301     // h_nob_NJet4_JetResPt_E.resize(EtaBinEdges.size() - 1);
302     // h_dead_NJet4_JetResPt_E.resize(EtaBinEdges.size() - 1);
303     // h_deadb_NJet4_JetResPt_E.resize(EtaBinEdges.size() - 1);
304     // for (std::vector<std::vector<TH1F*> >::iterator it = h_tot_NJet4_JetResPt_E.begin(); it != h_tot_NJet4_JetResPt_E.end(); ++it) {
305     // it->resize(EBinEdges.size() - 1);
306     // }
307     // for (std::vector<std::vector<TH1F*> >::iterator it = h_b_NJet4_JetResPt_E.begin(); it != h_b_NJet4_JetResPt_E.end(); ++it) {
308     // it->resize(EBinEdges.size() - 1);
309     // }
310     // for (std::vector<std::vector<TH1F*> >::iterator it = h_nob_NJet4_JetResPt_E.begin(); it != h_nob_NJet4_JetResPt_E.end(); ++it) {
311     // it->resize(EBinEdges.size() - 1);
312     // }
313     // for (std::vector<std::vector<TH1F*> >::iterator it = h_dead_NJet4_JetResPt_E.begin(); it != h_dead_NJet4_JetResPt_E.end(); ++it) {
314     // it->resize(EBinEdges.size() - 1);
315     // }
316     // for (std::vector<std::vector<TH1F*> >::iterator it = h_deadb_NJet4_JetResPt_E.begin(); it != h_deadb_NJet4_JetResPt_E.end(); ++it) {
317     // it->resize(EBinEdges.size() - 1);
318     // }
319    
320     //// Array of histograms for jet resolutions (jet multiplicities >= 5)
321     // h_tot_NJet5p_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
322     // h_b_NJet5p_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
323     // h_nob_NJet5p_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
324     // h_dead_NJet5p_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
325     // h_deadb_NJet5p_JetResPt_Pt.resize(EtaBinEdges.size() - 1);
326     // for (std::vector<std::vector<TH1F*> >::iterator it = h_tot_NJet5p_JetResPt_Pt.begin(); it != h_tot_NJet5p_JetResPt_Pt.end(); ++it) {
327     // it->resize(PtBinEdges.size() - 1);
328     // }
329     // for (std::vector<std::vector<TH1F*> >::iterator it = h_b_NJet5p_JetResPt_Pt.begin(); it != h_b_NJet5p_JetResPt_Pt.end(); ++it) {
330     // it->resize(PtBinEdges.size() - 1);
331     // }
332     // for (std::vector<std::vector<TH1F*> >::iterator it = h_nob_NJet5p_JetResPt_Pt.begin(); it != h_nob_NJet5p_JetResPt_Pt.end(); ++it) {
333     // it->resize(PtBinEdges.size() - 1);
334     // }
335     // for (std::vector<std::vector<TH1F*> >::iterator it = h_dead_NJet5p_JetResPt_Pt.begin(); it != h_dead_NJet5p_JetResPt_Pt.end(); ++it) {
336     // it->resize(PtBinEdges.size() - 1);
337     // }
338     // for (std::vector<std::vector<TH1F*> >::iterator it = h_deadb_NJet5p_JetResPt_Pt.begin(); it != h_deadb_NJet5p_JetResPt_Pt.end(); ++it) {
339     // it->resize(PtBinEdges.size() - 1);
340     // }
341     // h_tot_NJet5p_JetResPt_E.resize(EtaBinEdges.size() - 1);
342     // h_b_NJet5p_JetResPt_E.resize(EtaBinEdges.size() - 1);
343     // h_nob_NJet5p_JetResPt_E.resize(EtaBinEdges.size() - 1);
344     // h_dead_NJet5p_JetResPt_E.resize(EtaBinEdges.size() - 1);
345     // h_deadb_NJet5p_JetResPt_E.resize(EtaBinEdges.size() - 1);
346     // for (std::vector<std::vector<TH1F*> >::iterator it = h_tot_NJet5p_JetResPt_E.begin(); it != h_tot_NJet5p_JetResPt_E.end(); ++it) {
347     // it->resize(EBinEdges.size() - 1);
348     // }
349     // for (std::vector<std::vector<TH1F*> >::iterator it = h_b_NJet5p_JetResPt_E.begin(); it != h_b_NJet5p_JetResPt_E.end(); ++it) {
350     // it->resize(EBinEdges.size() - 1);
351     // }
352     // for (std::vector<std::vector<TH1F*> >::iterator it = h_nob_NJet5p_JetResPt_E.begin(); it != h_nob_NJet5p_JetResPt_E.end(); ++it) {
353     // it->resize(EBinEdges.size() - 1);
354     // }
355     // for (std::vector<std::vector<TH1F*> >::iterator it = h_dead_NJet5p_JetResPt_E.begin(); it != h_dead_NJet5p_JetResPt_E.end(); ++it) {
356     // it->resize(EBinEdges.size() - 1);
357     // }
358     // for (std::vector<std::vector<TH1F*> >::iterator it = h_deadb_NJet5p_JetResPt_E.begin(); it != h_deadb_NJet5p_JetResPt_E.end(); ++it) {
359     // it->resize(EBinEdges.size() - 1);
360     // }
361 csander 1.1
362     }
363    
364 csander 1.3 MCResolutions::~MCResolutions() {
365 csander 1.2
366 csander 1.10 // do anything here that needs to be done at desctruction time
367     // (e.g. close files, deallocate resources etc.)
368 csander 1.1
369     }
370    
371     //
372     // member functions
373     //
374    
375     // ------------ method called to for each event ------------
376 csander 1.10 void MCResolutions::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
377     using namespace std;
378 csander 1.5
379 csander 1.10 // Event setup
380     envSet(iSetup);
381 csander 1.9
382 csander 1.10 getChannelStatusMaps();
383 csander 1.9
384 csander 1.10 //Weight
385     edm::Handle<double> event_weight;
386     bool findWeight = iEvent.getByLabel(_weightName, event_weight);
387     weight = (event_weight.isValid() ? (*event_weight) : 1.0);
388     if (!findWeight) {
389     cout << "Weight not found!" << endl;
390     }
391    
392     //GenJets
393     edm::Handle<edm::View<reco::GenJet> > Jets_gen;
394     iEvent.getByLabel(_genJetTag, Jets_gen);
395    
396     //RecoJets
397     edm::Handle<edm::View<pat::Jet> > Jets_rec;
398     iEvent.getByLabel(_jetTag, Jets_rec);
399    
400     //softMuonBJetTags
401     edm::Handle<reco::JetTagCollection> bTagHandle;
402     iEvent.getByLabel("trackCountingHighPurBJetTags", bTagHandle);
403     const reco::JetTagCollection & bTags = *(bTagHandle.product());
404    
405     edm::Handle<reco::GenParticleCollection> genParticles;
406     iEvent.getByLabel("genParticles", genParticles);
407    
408     //// Do DiJet selelction
409     bool isDiJet = false;
410     const reco::GenJet* gj1 = 0;
411     const reco::GenJet* gj2 = 0;
412     const reco::GenJet* gj3 = 0;
413     int k = 0;
414     for (edm::View<reco::GenJet>::const_iterator it = Jets_gen->begin(); it != Jets_gen->end(); ++it) {
415     ++k;
416     if (k == 1)
417     gj1 = &(*it);
418     if (k == 2)
419     gj2 = &(*it);
420     if (k == 3) {
421     gj3 = &(*it);
422     break;
423     }
424     }
425     if (k > 1) {
426     if (deltaPhi(gj1->phi(), gj2->phi()) > _deltaPhiDiJet) {
427     if (k == 2)
428     isDiJet = true;
429     if (k == 3) {
430     //if (gj3->pt()<(_relCut3rdJet*(gj1->pt()+gj2->pt())/2)){
431     if (gj3->pt() < (_relCut3rdJet * gj1->pt())) {
432     isDiJet = true;
433 csander 1.4 }
434 csander 1.10 }
435     }
436     }
437    
438     //// Get jet multiplicity
439     double NGenJet = 0;
440     for (edm::View<reco::GenJet>::const_iterator it = Jets_gen-> begin(); it != Jets_gen-> end(); ++it) {
441     if (it->pt() > _jetMultPtCut && abs(it->eta()) < _jetMultEtaCut)
442     ++NGenJet;
443     }
444    
445     int JetNo = 0;
446     for (edm::View<reco::GenJet>::const_iterator it = Jets_gen->begin(); it != Jets_gen->end(); ++it) {
447     ++JetNo;
448    
449     //use only two leading jets
450     if (JetNo > 2)
451     break;
452    
453     if (it->pt() < _GenJetPtCut)
454     continue;
455     //// First look if there is no significant GenJet near the tested GenJet
456     double dRgenjet = 999.;
457     double PtdRmin = 0;
458     for (edm::View<reco::GenJet>::const_iterator kt = Jets_gen-> begin(); kt != Jets_gen->end(); ++kt) {
459     if (&(*it) != &(*kt))
460 csander 1.3 continue;
461 csander 1.10 double dR = deltaR(*it, *kt);
462     if (dR < dRgenjet) {
463     dRgenjet = dR;
464     PtdRmin = kt->pt();
465     }
466     }
467     if (dRgenjet < _deltaRMatchVeto && PtdRmin / it->pt() < _relPtVeto)
468     continue;
469     const pat::Jet* matchedJet = 0;
470     const pat::Jet* nearestJet = 0;
471     double dRmatched = 999.;
472     double dRnearest = 999.;
473     for (edm::View<pat::Jet>::const_iterator jt = Jets_rec-> begin(); jt != Jets_rec->end(); ++jt) {
474     double dR = deltaR(*it, *jt);
475     if (dR < dRmatched) {
476     nearestJet = matchedJet;
477     dRnearest = dRmatched;
478     matchedJet = &(*jt);
479     dRmatched = dR;
480     } else if (dR < dRnearest) {
481     nearestJet = &(*jt);
482     dRnearest = dR;
483     }
484     }
485     //// look if there is no further significant CaloJet near the genJet
486     if (dRmatched < _deltaRMatch && (nearestJet == 0 || dRnearest > _deltaRMatchVeto || (nearestJet->pt()
487     < _absPtVeto && nearestJet->pt() / matchedJet->pt() < _relPtVeto))) {
488     //cout << dRmatched << endl;
489     double res = matchedJet->pt() / it->pt();
490     h_tot_NJetAll_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
491     //h_tot_NJetAll_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
492     if (isDiJet && JetNo < 3) {
493     h_tot_DiJet_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
494     //h_tot_DiJet_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
495     }
496     // if (NGenJet == 2) {
497     // h_tot_NJet2_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
498     // h_tot_NJet2_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
499     // }
500     // if (NGenJet == 3) {
501     // h_tot_NJet3_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
502     // h_tot_NJet3_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
503     // }
504     // if (NGenJet == 4) {
505     // h_tot_NJet4_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
506     // h_tot_NJet4_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
507     // }
508     // if (NGenJet >= 5) {
509     // h_tot_NJet5p_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
510     // h_tot_NJet5p_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
511     // }
512     //// look if matched jet can be matched to a b-tagged jet
513     bool bTag = false;
514     if (_Bid == "btag") {
515     //// Use matched parton for heavy flavour ID
516     for (size_t j = 0; j < genParticles->size(); ++j) {
517     const reco::GenParticle & p = (*genParticles)[j];
518     if (p.status() != 2)
519     continue;
520     int id = abs(p.pdgId());
521     if (id != 4 && id != 5)
522     continue;
523     float dR = deltaR(p, *it);
524     //// This is my own criterion ... should be not that bad but has to be validated/redifined
525     if (dR < 0.25 && p.pt() / it->pt() > 0.25) {
526     //cout << "Jet pT: " << it->pt() << endl;
527     //cout << "Parton pT: " << p.pt() << endl;
528     bTag = true;
529     break;
530     }
531 csander 1.3 }
532 csander 1.10 } else if (_Bid == "own") {
533     //// Use bTags for heavy flavour ID
534     for (int i = 0; i != (int) bTags.size(); ++i) {
535     double dR_btag = deltaR(*matchedJet, *(bTags[i].first));
536     if (dR_btag < _bTagDeltaR && bTags[i].second > _bTagCut) {
537     bTag = true;
538     //cout << "Is B_own" << endl;
539     break;
540     }
541 csander 1.5 }
542 csander 1.10 } else if (_Bid == "algo") {
543     //// Use algorithmic matching for heavy flavour ID
544     if (fabs(matchedJet->partonFlavour()) == 4 || fabs(matchedJet->partonFlavour()) == 5) {
545     bTag = true;
546     //cout << "Is B_algo" << endl;
547 csander 1.3 }
548 csander 1.10 }
549 csander 1.6
550 csander 1.10 if (bTag) {
551     //cout << "matched jet has b tag discriminator = "<<bTags[i].second << endl;
552     h_b_NJetAll_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
553     //h_b_NJetAll_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
554     if (isDiJet && JetNo < 3) {
555     h_b_DiJet_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
556     //h_b_DiJet_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
557     }
558     // if (NGenJet == 2) {
559     // h_b_NJet2_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
560     // h_b_NJet2_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
561     // }
562     // if (NGenJet == 3) {
563     // h_b_NJet3_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
564     // h_b_NJet3_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
565     // }
566     // if (NGenJet == 4) {
567     // h_b_NJet4_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
568     // h_b_NJet4_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
569     // }
570     // if (NGenJet >= 5) {
571     // h_b_NJet5p_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
572     // h_b_NJet5p_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
573     // }
574     } else {
575     h_nob_NJetAll_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
576     //h_nob_NJetAll_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
577     if (isDiJet && JetNo < 3) {
578     h_nob_DiJet_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
579     //h_nob_DiJet_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
580 csander 1.6 }
581 csander 1.10 // if (NGenJet == 2) {
582     // h_nob_NJet2_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
583     // h_nob_NJet2_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
584     // }
585     // if (NGenJet == 3) {
586     // h_nob_NJet3_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
587     // h_nob_NJet3_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
588     // }
589     // if (NGenJet == 4) {
590     // h_nob_NJet4_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
591     // h_nob_NJet4_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
592     // }
593     // if (NGenJet >= 5) {
594     // h_nob_NJet5p_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
595     // h_nob_NJet5p_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
596     // }
597     }
598    
599     //// look if matched jet points into direction of masked ECAL cluster
600     for (std::map<DetId, std::vector<double> >::iterator kt = EcalAllDeadChannelsValMap.begin(); kt
601     != EcalAllDeadChannelsValMap.end(); ++kt) {
602     math::PtEtaPhiMLorentzVectorD Evec(100., kt->second.at(0), kt->second.at(1), 0.);
603     double dR_dead = deltaR(*matchedJet, Evec);
604     if (dR_dead < _deltaRDeadECal) {
605     //cout << "jet is pointing to dead Ecal cluster (eta, phi)= "<< it->eta() << ", " << it->phi() << endl;
606     h_dead_NJetAll_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
607     //h_dead_NJetAll_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
608     if (bTag) {
609     h_deadb_NJetAll_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
610     //h_deadb_NJetAll_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
611     }
612     if (isDiJet && JetNo < 3) {
613     h_dead_DiJet_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
614     //h_dead_DiJet_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
615     if (bTag) {
616     h_deadb_DiJet_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
617     //h_deadb_DiJet_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
618     }
619     }
620     // if (NGenJet == 2) {
621     // h_dead_NJet2_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
622     // h_dead_NJet2_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
623     // if (bTag) {
624     // h_deadb_NJet2_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
625     // h_deadb_NJet2_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
626     // }
627     // }
628     // if (NGenJet == 3) {
629     // h_dead_NJet3_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
630     // h_dead_NJet3_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
631     // if (bTag) {
632     // h_deadb_NJet3_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
633     // h_deadb_NJet3_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
634     // }
635     // }
636     // if (NGenJet == 4) {
637     // h_dead_NJet4_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
638     // h_dead_NJet4_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
639     // if (bTag) {
640     // h_deadb_NJet4_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
641     // h_deadb_NJet4_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
642     // }
643     // }
644     // if (NGenJet >= 5) {
645     // h_dead_NJet5p_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
646     // h_dead_NJet5p_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
647     // if (bTag) {
648     // h_deadb_NJet5p_JetResPt_Pt.at(EtaBin(it->eta())).at(PtBin(it->pt()))->Fill(res, weight);
649     // h_deadb_NJet5p_JetResPt_E.at(EtaBin(it->eta())).at(EBin(it->energy()))->Fill(res, weight);
650     // }
651     // }
652     break;
653 csander 1.3 }
654 csander 1.10 }
655     }
656     }
657 csander 1.1
658     }
659    
660     // ------------ method called once each job just before starting event loop ------------
661 csander 1.10 void MCResolutions::beginJob() {
662 csander 1.3
663 csander 1.10 for (unsigned int i_pt = 0; i_pt < PtBinEdges.size() - 1; ++i_pt) {
664     for (unsigned int i_eta = 0; i_eta < EtaBinEdges.size() - 1; ++i_eta) {
665     char hname[100];
666     //// Book histograms (all jet multiplicities)
667     sprintf(hname, "h_tot_NJetAll_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
668     h_tot_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
669     h_tot_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
670     sprintf(hname, "h_b_NJetAll_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
671     h_b_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
672     h_b_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
673     sprintf(hname, "h_nob_NJetAll_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
674     h_nob_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
675     h_nob_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
676     sprintf(hname, "h_dead_NJetAll_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
677     h_dead_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
678     h_dead_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
679     sprintf(hname, "h_deadb_NJetAll_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
680     h_deadb_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
681     h_deadb_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
682     //// Book histograms (DiJets)
683     sprintf(hname, "h_tot_DiJet_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
684     h_tot_DiJet_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
685     h_tot_DiJet_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
686     sprintf(hname, "h_b_DiJet_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
687     h_b_DiJet_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
688     h_b_DiJet_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
689     sprintf(hname, "h_nob_DiJet_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
690     h_nob_DiJet_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
691     h_nob_DiJet_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
692     sprintf(hname, "h_dead_DiJet_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
693     h_dead_DiJet_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
694     h_dead_DiJet_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
695     sprintf(hname, "h_deadb_DiJet_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
696     h_deadb_DiJet_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
697     h_deadb_DiJet_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
698     //// Book histograms (jet multiplicities = 2)
699     // sprintf(hname, "h_tot_NJet2_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
700     // h_tot_NJet2_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
701     // h_tot_NJet2_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
702     // sprintf(hname, "h_b_NJet2_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
703     // h_b_NJet2_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
704     // h_b_NJet2_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
705     // sprintf(hname, "h_nob_NJet2_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
706     // h_nob_NJet2_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
707     // h_nob_NJet2_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
708     // sprintf(hname, "h_dead_NJet2_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
709     // h_dead_NJet2_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
710     // h_dead_NJet2_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
711     // sprintf(hname, "h_deadb_NJet2_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
712     // h_deadb_NJet2_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
713     // h_deadb_NJet2_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
714     //// Book histograms (jet multiplicities = 3)
715     // sprintf(hname, "h_tot_NJet3_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
716     // h_tot_NJet3_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
717     // h_tot_NJet3_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
718     // sprintf(hname, "h_b_NJet3_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
719     // h_b_NJet3_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
720     // h_b_NJet3_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
721     // sprintf(hname, "h_nob_NJet3_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
722     // h_nob_NJet3_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
723     // h_nob_NJet3_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
724     // sprintf(hname, "h_dead_NJet3_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
725     // h_dead_NJet3_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
726     // h_dead_NJet3_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
727     // sprintf(hname, "h_deadb_NJet3_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
728     // h_deadb_NJet3_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
729     // h_deadb_NJet3_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
730     //// Book histograms (jet multiplicities = 4)
731     // sprintf(hname, "h_tot_NJet4_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
732     // h_tot_NJet4_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
733     // h_tot_NJet4_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
734     // sprintf(hname, "h_b_NJet4_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
735     // h_b_NJet4_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
736     // h_b_NJet4_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
737     // sprintf(hname, "h_nob_NJet4_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
738     // h_nob_NJet4_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
739     // h_nob_NJet4_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
740     // sprintf(hname, "h_dead_NJet4_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
741     // h_dead_NJet4_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
742     // h_dead_NJet4_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
743     // sprintf(hname, "h_deadb_NJet4_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
744     // h_deadb_NJet4_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
745     // h_deadb_NJet4_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
746     //// Book histograms (jet multiplicities >= 5)
747     // sprintf(hname, "h_tot_NJet5p_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
748     // h_tot_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
749     // h_tot_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
750     // sprintf(hname, "h_b_NJet5p_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
751     // h_b_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
752     // h_b_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
753     // sprintf(hname, "h_nob_NJet5p_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
754     // h_nob_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
755     // h_nob_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
756     // sprintf(hname, "h_dead_NJet5p_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
757     // h_dead_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
758     // h_dead_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
759     // sprintf(hname, "h_deadb_NJet5p_ResponsePt_Pt%i_Eta%i", i_pt, i_eta);
760     // h_deadb_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt) = new TH1F(hname, hname, 300, 0., 3.);
761     // h_deadb_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt)->Sumw2();
762     }
763     }
764    
765     // for (unsigned int i_e = 0; i_e < EBinEdges.size() - 1; ++i_e) {
766     // for (unsigned int i_eta = 0; i_eta < EtaBinEdges.size() - 1; ++i_eta) {
767     // char hname[100];
768     // //// Book histograms (all jet multiplicities)
769     // sprintf(hname, "h_tot_NJetAll_ResponsePt_E%i_Eta%i", i_e, i_eta);
770     // h_tot_NJetAll_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
771     // h_tot_NJetAll_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
772     // sprintf(hname, "h_b_NJetAll_ResponsePt_E%i_Eta%i", i_e, i_eta);
773     // h_b_NJetAll_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
774     // h_b_NJetAll_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
775     // sprintf(hname, "h_nob_NJetAll_ResponsePt_E%i_Eta%i", i_e, i_eta);
776     // h_nob_NJetAll_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
777     // h_nob_NJetAll_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
778     // sprintf(hname, "h_dead_NJetAll_ResponsePt_E%i_Eta%i", i_e, i_eta);
779     // h_dead_NJetAll_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
780     // h_dead_NJetAll_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
781     // sprintf(hname, "h_deadb_NJetAll_ResponsePt_E%i_Eta%i", i_e, i_eta);
782     // h_deadb_NJetAll_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
783     // h_deadb_NJetAll_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
784     // //// Book histograms (DiJets)
785     // sprintf(hname, "h_tot_DiJet_ResponsePt_E%i_Eta%i", i_e, i_eta);
786     // h_tot_DiJet_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
787     // h_tot_DiJet_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
788     // sprintf(hname, "h_b_DiJet_ResponsePt_E%i_Eta%i", i_e, i_eta);
789     // h_b_DiJet_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
790     // h_b_DiJet_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
791     // sprintf(hname, "h_nob_DiJet_ResponsePt_E%i_Eta%i", i_e, i_eta);
792     // h_nob_DiJet_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
793     // h_nob_DiJet_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
794     // sprintf(hname, "h_dead_DiJet_ResponsePt_E%i_Eta%i", i_e, i_eta);
795     // h_dead_DiJet_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
796     // h_dead_DiJet_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
797     // sprintf(hname, "h_deadb_DiJet_ResponsePt_E%i_Eta%i", i_e, i_eta);
798     // h_deadb_DiJet_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
799     // h_deadb_DiJet_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
800     // //// Book histograms (jet multiplicities = 2)
801     // sprintf(hname, "h_tot_NJet2_ResponsePt_E%i_Eta%i", i_e, i_eta);
802     // h_tot_NJet2_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
803     // h_tot_NJet2_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
804     // sprintf(hname, "h_b_NJet2_ResponsePt_E%i_Eta%i", i_e, i_eta);
805     // h_b_NJet2_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
806     // h_b_NJet2_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
807     // sprintf(hname, "h_nob_NJet2_ResponsePt_E%i_Eta%i", i_e, i_eta);
808     // h_nob_NJet2_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
809     // h_nob_NJet2_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
810     // sprintf(hname, "h_dead_NJet2_ResponsePt_E%i_Eta%i", i_e, i_eta);
811     // h_dead_NJet2_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
812     // h_dead_NJet2_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
813     // sprintf(hname, "h_deadb_NJet2_ResponsePt_E%i_Eta%i", i_e, i_eta);
814     // h_deadb_NJet2_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
815     // h_deadb_NJet2_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
816     // //// Book histograms (jet multiplicities = 3)
817     // sprintf(hname, "h_tot_NJet3_ResponsePt_E%i_Eta%i", i_e, i_eta);
818     // h_tot_NJet3_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
819     // h_tot_NJet3_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
820     // sprintf(hname, "h_b_NJet3_ResponsePt_E%i_Eta%i", i_e, i_eta);
821     // h_b_NJet3_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
822     // h_b_NJet3_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
823     // sprintf(hname, "h_nob_NJet3_ResponsePt_E%i_Eta%i", i_e, i_eta);
824     // h_nob_NJet3_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
825     // h_nob_NJet3_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
826     // sprintf(hname, "h_dead_NJet3_ResponsePt_E%i_Eta%i", i_e, i_eta);
827     // h_dead_NJet3_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
828     // h_dead_NJet3_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
829     // sprintf(hname, "h_deadb_NJet3_ResponsePt_E%i_Eta%i", i_e, i_eta);
830     // h_deadb_NJet3_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
831     // h_deadb_NJet3_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
832     // //// Book histograms (jet multiplicities = 4)
833     // sprintf(hname, "h_tot_NJet4_ResponsePt_E%i_Eta%i", i_e, i_eta);
834     // h_tot_NJet4_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
835     // h_tot_NJet4_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
836     // sprintf(hname, "h_b_NJet4_ResponsePt_E%i_Eta%i", i_e, i_eta);
837     // h_b_NJet4_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
838     // h_b_NJet4_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
839     // sprintf(hname, "h_nob_NJet4_ResponsePt_E%i_Eta%i", i_e, i_eta);
840     // h_nob_NJet4_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
841     // h_nob_NJet4_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
842     // sprintf(hname, "h_dead_NJet4_ResponsePt_E%i_Eta%i", i_e, i_eta);
843     // h_dead_NJet4_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
844     // h_dead_NJet4_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
845     // sprintf(hname, "h_deadb_NJet4_ResponsePt_E%i_Eta%i", i_e, i_eta);
846     // h_deadb_NJet4_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
847     // h_deadb_NJet4_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
848     // //// Book histograms (jet multiplicities = 5p)
849     // sprintf(hname, "h_tot_NJet5p_ResponsePt_E%i_Eta%i", i_e, i_eta);
850     // h_tot_NJet5p_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
851     // h_tot_NJet5p_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
852     // sprintf(hname, "h_b_NJet5p_ResponsePt_E%i_Eta%i", i_e, i_eta);
853     // h_b_NJet5p_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
854     // h_b_NJet5p_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
855     // sprintf(hname, "h_nob_NJet5p_ResponsePt_E%i_Eta%i", i_e, i_eta);
856     // h_nob_NJet5p_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
857     // h_nob_NJet5p_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
858     // sprintf(hname, "h_dead_NJet5p_ResponsePt_E%i_Eta%i", i_e, i_eta);
859     // h_dead_NJet5p_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
860     // h_dead_NJet5p_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
861     // sprintf(hname, "h_deadb_NJet5p_ResponsePt_E%i_Eta%i", i_e, i_eta);
862     // h_deadb_NJet5p_JetResPt_E.at(i_eta).at(i_e) = new TH1F(hname, hname, 300, 0., 3.);
863     // h_deadb_NJet5p_JetResPt_E.at(i_eta).at(i_e)->Sumw2();
864     // }
865     // }
866 csander 1.3
867 csander 1.10 mapsReady = false;
868 csander 1.3 }
869 csander 1.1
870     // ------------ method called once each job just after ending the event loop ------------
871 csander 1.10 void MCResolutions::endJob() {
872 csander 1.5
873 csander 1.10 hfile->cd();
874     //hfile->cd(_dirName.c_str());
875     // Save all objects in this file
876     for (unsigned int i_pt = 0; i_pt < PtBinEdges.size() - 1; ++i_pt) {
877     for (unsigned int i_eta = 0; i_eta < EtaBinEdges.size() - 1; ++i_eta) {
878     // total
879     hfile->WriteTObject(h_tot_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt));
880     hfile->WriteTObject(h_tot_DiJet_JetResPt_Pt.at(i_eta).at(i_pt));
881     // hfile->WriteTObject(h_tot_NJet2_JetResPt_Pt.at(i_eta).at(i_pt));
882     // hfile->WriteTObject(h_tot_NJet3_JetResPt_Pt.at(i_eta).at(i_pt));
883     // hfile->WriteTObject(h_tot_NJet4_JetResPt_Pt.at(i_eta).at(i_pt));
884     // hfile->WriteTObject(h_tot_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt));
885     // with btag
886     hfile->WriteTObject(h_b_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt));
887     hfile->WriteTObject(h_b_DiJet_JetResPt_Pt.at(i_eta).at(i_pt));
888     // hfile->WriteTObject(h_b_NJet2_JetResPt_Pt.at(i_eta).at(i_pt));
889     // hfile->WriteTObject(h_b_NJet3_JetResPt_Pt.at(i_eta).at(i_pt));
890     // hfile->WriteTObject(h_b_NJet4_JetResPt_Pt.at(i_eta).at(i_pt));
891     // hfile->WriteTObject(h_b_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt));
892     // without btag
893     hfile->WriteTObject(h_nob_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt));
894     hfile->WriteTObject(h_nob_DiJet_JetResPt_Pt.at(i_eta).at(i_pt));
895     // hfile->WriteTObject(h_nob_NJet2_JetResPt_Pt.at(i_eta).at(i_pt));
896     // hfile->WriteTObject(h_nob_NJet3_JetResPt_Pt.at(i_eta).at(i_pt));
897     // hfile->WriteTObject(h_nob_NJet4_JetResPt_Pt.at(i_eta).at(i_pt));
898     // hfile->WriteTObject(h_nob_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt));
899     // in direction of dead ECAL cells
900     hfile->WriteTObject(h_dead_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt));
901     hfile->WriteTObject(h_dead_DiJet_JetResPt_Pt.at(i_eta).at(i_pt));
902     // hfile->WriteTObject(h_dead_NJet2_JetResPt_Pt.at(i_eta).at(i_pt));
903     // hfile->WriteTObject(h_dead_NJet3_JetResPt_Pt.at(i_eta).at(i_pt));
904     // hfile->WriteTObject(h_dead_NJet4_JetResPt_Pt.at(i_eta).at(i_pt));
905     // hfile->WriteTObject(h_dead_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt));
906     // in direction of dead ECAL cells and with b tag
907     hfile->WriteTObject(h_deadb_NJetAll_JetResPt_Pt.at(i_eta).at(i_pt));
908     hfile->WriteTObject(h_deadb_DiJet_JetResPt_Pt.at(i_eta).at(i_pt));
909     // hfile->WriteTObject(h_deadb_NJet2_JetResPt_Pt.at(i_eta).at(i_pt));
910     // hfile->WriteTObject(h_deadb_NJet3_JetResPt_Pt.at(i_eta).at(i_pt));
911     // hfile->WriteTObject(h_deadb_NJet4_JetResPt_Pt.at(i_eta).at(i_pt));
912     // hfile->WriteTObject(h_deadb_NJet5p_JetResPt_Pt.at(i_eta).at(i_pt));
913     }
914     }
915     // for (unsigned int i_e = 0; i_e < EBinEdges.size() - 1; ++i_e) {
916     // for (unsigned int i_eta = 0; i_eta < EtaBinEdges.size() - 1; ++i_eta) {
917     // // total
918     // hfile->WriteTObject(h_tot_NJetAll_JetResPt_E.at(i_eta).at(i_e));
919     // hfile->WriteTObject(h_tot_DiJet_JetResPt_E.at(i_eta).at(i_e));
920     // hfile->WriteTObject(h_tot_NJet2_JetResPt_E.at(i_eta).at(i_e));
921     // hfile->WriteTObject(h_tot_NJet3_JetResPt_E.at(i_eta).at(i_e));
922     // hfile->WriteTObject(h_tot_NJet4_JetResPt_E.at(i_eta).at(i_e));
923     // hfile->WriteTObject(h_tot_NJet5p_JetResPt_E.at(i_eta).at(i_e));
924     // // with btag
925     // hfile->WriteTObject(h_b_NJetAll_JetResPt_E.at(i_eta).at(i_e));
926     // hfile->WriteTObject(h_b_DiJet_JetResPt_E.at(i_eta).at(i_e));
927     // hfile->WriteTObject(h_b_NJet2_JetResPt_E.at(i_eta).at(i_e));
928     // hfile->WriteTObject(h_b_NJet3_JetResPt_E.at(i_eta).at(i_e));
929     // hfile->WriteTObject(h_b_NJet4_JetResPt_E.at(i_eta).at(i_e));
930     // hfile->WriteTObject(h_b_NJet5p_JetResPt_E.at(i_eta).at(i_e));
931     // // without btag
932     // hfile->WriteTObject(h_nob_NJetAll_JetResPt_E.at(i_eta).at(i_e));
933     // hfile->WriteTObject(h_nob_DiJet_JetResPt_E.at(i_eta).at(i_e));
934     // hfile->WriteTObject(h_nob_NJet2_JetResPt_E.at(i_eta).at(i_e));
935     // hfile->WriteTObject(h_nob_NJet3_JetResPt_E.at(i_eta).at(i_e));
936     // hfile->WriteTObject(h_nob_NJet4_JetResPt_E.at(i_eta).at(i_e));
937     // hfile->WriteTObject(h_nob_NJet5p_JetResPt_E.at(i_eta).at(i_e));
938     // // in direction of dead ECAL cells
939     // hfile->WriteTObject(h_dead_NJetAll_JetResPt_E.at(i_eta).at(i_e));
940     // hfile->WriteTObject(h_dead_DiJet_JetResPt_E.at(i_eta).at(i_e));
941     // hfile->WriteTObject(h_dead_NJet2_JetResPt_E.at(i_eta).at(i_e));
942     // hfile->WriteTObject(h_dead_NJet3_JetResPt_E.at(i_eta).at(i_e));
943     // hfile->WriteTObject(h_dead_NJet4_JetResPt_E.at(i_eta).at(i_e));
944     // hfile->WriteTObject(h_dead_NJet5p_JetResPt_E.at(i_eta).at(i_e));
945     // // in direction of dead ECAL cells and with b tag
946     // hfile->WriteTObject(h_deadb_NJetAll_JetResPt_E.at(i_eta).at(i_e));
947     // hfile->WriteTObject(h_deadb_DiJet_JetResPt_E.at(i_eta).at(i_e));
948     // hfile->WriteTObject(h_deadb_NJet2_JetResPt_E.at(i_eta).at(i_e));
949     // hfile->WriteTObject(h_deadb_NJet3_JetResPt_E.at(i_eta).at(i_e));
950     // hfile->WriteTObject(h_deadb_NJet4_JetResPt_E.at(i_eta).at(i_e));
951     // hfile->WriteTObject(h_deadb_NJet5p_JetResPt_E.at(i_eta).at(i_e));
952     // }
953     // }
954    
955     hfile->cd();
956     //hfile->WriteObject(&EBinEdges, "EBinEdges");
957     hfile->WriteObject(&PtBinEdges, "PtBinEdges");
958     hfile->WriteObject(&EtaBinEdges, "EtaBinEdges");
959     //hfile->ls();
960 csander 1.5
961 csander 1.10 // Close the file.
962     hfile->Close();
963 csander 1.5
964     }
965 csander 1.3
966     int MCResolutions::EBin(const double& e) {
967 csander 1.10 int i_e = -1;
968     for (std::vector<double>::const_iterator it = EBinEdges.begin(); it != EBinEdges.end(); ++it) {
969     if ((*it) > e)
970     break;
971     ++i_e;
972     }
973     if (i_e < 0)
974     i_e = 0;
975     if (i_e > (int) EBinEdges.size() - 2)
976     i_e = (int) EBinEdges.size() - 2;
977 csander 1.3
978 csander 1.10 return i_e;
979 csander 1.3 }
980    
981     int MCResolutions::PtBin(const double& pt) {
982 csander 1.10 int i_pt = -1;
983     for (std::vector<double>::const_iterator it = PtBinEdges.begin(); it != PtBinEdges.end(); ++it) {
984     if ((*it) > pt)
985     break;
986     ++i_pt;
987     }
988     if (i_pt < 0)
989     i_pt = 0;
990     if (i_pt > (int) PtBinEdges.size() - 2)
991     i_pt = (int) PtBinEdges.size() - 2;
992 csander 1.3
993 csander 1.10 return i_pt;
994 csander 1.3 }
995    
996     int MCResolutions::EtaBin(const double& eta) {
997 csander 1.10 int i_eta = -1;
998     for (std::vector<double>::const_iterator it = EtaBinEdges.begin(); it != EtaBinEdges.end(); ++it) {
999     if ((*it) > fabs(eta))
1000     break;
1001     ++i_eta;
1002     }
1003     if (i_eta < 0)
1004     i_eta = 0;
1005     if (i_eta > (int) EtaBinEdges.size() - 2)
1006     i_eta = (int) EtaBinEdges.size() - 2;
1007     return i_eta;
1008 csander 1.3 }
1009    
1010     void MCResolutions::envSet(const edm::EventSetup& iSetup) {
1011    
1012 csander 1.10 ecalScale.setEventSetup(iSetup);
1013 csander 1.3
1014 csander 1.10 iSetup.get<EcalChannelStatusRcd> ().get(ecalStatus);
1015     iSetup.get<CaloGeometryRecord> ().get(geometry);
1016 csander 1.3
1017 csander 1.10 if (!ecalStatus.isValid())
1018     throw "Failed to get ECAL channel status!";
1019     if (!geometry.isValid())
1020     throw "Failed to get the geometry!";
1021 csander 1.3
1022     }
1023    
1024     int MCResolutions::getChannelStatusMaps() {
1025    
1026 csander 1.10 if (mapsReady)
1027     return -1;
1028 csander 1.3
1029 csander 1.10 EcalAllDeadChannelsValMap.clear();
1030 csander 1.3
1031 csander 1.10 // Loop over EB ...
1032     for (int ieta = -85; ieta <= 85; ieta++) {
1033     for (int iphi = 0; iphi <= 360; iphi++) {
1034     if (!EBDetId::validDetId(ieta, iphi))
1035     continue;
1036    
1037     const EBDetId detid = EBDetId(ieta, iphi, EBDetId::ETAPHIMODE);
1038     EcalChannelStatus::const_iterator chit = ecalStatus->find(detid);
1039     // refer https://twiki.cern.ch/twiki/bin/viewauth/CMS/EcalChannelStatus
1040     int status = (chit != ecalStatus->end()) ? chit->getStatusCode() : -1;
1041     //std::cout << ieta << " " << iphi << " " << status << std:: endl;
1042     const CaloSubdetectorGeometry* subGeom = geometry->getSubdetectorGeometry(detid);
1043     const CaloCellGeometry* cellGeom = subGeom->getGeometry(detid);
1044     double eta = cellGeom->getPosition().eta();
1045     double phi = cellGeom->getPosition().phi();
1046     double theta = cellGeom->getPosition().theta();
1047    
1048     if (status >= _maskedEcalChannelStatusThreshold) {
1049     std::vector<double> valVec;
1050     valVec.push_back(eta);
1051     valVec.push_back(phi);
1052     valVec.push_back(theta);
1053     //std::cout << eta << " " << phi << std::endl;
1054     EcalAllDeadChannelsValMap.insert(std::make_pair(detid, valVec));
1055     }
1056     } // end loop iphi
1057     } // end loop ieta
1058    
1059     // Loop over EE detid
1060     for (int ix = 0; ix <= 100; ix++) {
1061     for (int iy = 0; iy <= 100; iy++) {
1062     for (int iz = -1; iz <= 1; iz++) {
1063     if (iz == 0)
1064     continue;
1065     if (!EEDetId::validDetId(ix, iy, iz))
1066     continue;
1067    
1068     const EEDetId detid = EEDetId(ix, iy, iz, EEDetId::XYMODE);
1069     EcalChannelStatus::const_iterator chit = ecalStatus->find(detid);
1070     int status = (chit != ecalStatus->end()) ? chit->getStatusCode() : -1;
1071    
1072     const CaloSubdetectorGeometry* subGeom = geometry->getSubdetectorGeometry(detid);
1073     const CaloCellGeometry* cellGeom = subGeom->getGeometry(detid);
1074     double eta = cellGeom->getPosition().eta();
1075     double phi = cellGeom->getPosition().phi();
1076 csander 1.3 double theta = cellGeom->getPosition().theta();
1077    
1078 csander 1.10 if (status >= _maskedEcalChannelStatusThreshold) {
1079     std::vector<double> valVec;
1080     valVec.push_back(eta);
1081     valVec.push_back(phi);
1082     valVec.push_back(theta);
1083     //std::cout << eta << " " << phi << std::endl;
1084     EcalAllDeadChannelsValMap.insert(std::make_pair(detid, valVec));
1085 csander 1.3 }
1086 csander 1.10 } // end loop iz
1087     } // end loop iy
1088     } // end loop ix
1089 csander 1.3
1090 csander 1.10 mapsReady = true;
1091     return 1;
1092 csander 1.3 }
1093 csander 1.1
1094     //define this as a plug-in
1095 csander 1.10 DEFINE_FWK_MODULE( MCResolutions);