ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/csander/JetResolutionFromMC/MCResolutions/src/MCResolutions.cc
Revision: 1.7
Committed: Mon Jan 17 15:33:49 2011 UTC (14 years, 3 months ago) by csander
Content type: text/plain
Branch: MAIN
Changes since 1.6: +139 -78 lines
Log Message: 
new version

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