ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/SchieferD/jetcalib/mcresponse_x.cpp
(Generate patch)

Comparing UserCode/SchieferD/jetcalib/mcresponse_x.cpp (file contents):
Revision 1.6 by schiefer, Fri Aug 31 13:56:04 2007 UTC vs.
Revision 1.11 by schiefer, Mon Sep 3 08:14:28 2007 UTC

# Line 10 | Line 10
10   #include "utils/cmdline.h"
11   #include "utils/hist.h"
12  
13 + #include <TSystem.h>
14   #include <TROOT.h>
15   #include <TRint.h>
16   #include <TFile.h>
# Line 30 | Line 31 | using namespace std;
31  
32  
33   ////////////////////////////////////////////////////////////////////////////////
33 – // defaults
34 – ////////////////////////////////////////////////////////////////////////////////
35 – string dflt_etbins ="8,12,15,20,25,30,45,65,90,120,180";
36 – string dflt_etabins=".2,.4,.6,.8,1.,1.2,1.4,1.6,1.8,2.,2.3,2.5,2.7,3.0,3.3";
37 – string dflt_phibins="-2.75,-2.25,-1.75.,-1.25,-.75,-.25,.25,.75,1.25,1.75,2.25,2.75";
38 – string dflt_emfbins="0.2,0.4,0.6,0.8";
39 –
40 –
41 – ////////////////////////////////////////////////////////////////////////////////
34   // declare global functions
35   ////////////////////////////////////////////////////////////////////////////////
36   void     replaceHistos(int nbins,
# Line 77 | Line 69 | int main(int argc,char**argv)
69    short          applyjes  = cl.get_value<short>  ("applyjes",          0);
70    string         treename  = cl.get_value<string> ("treename",        "t");
71    float          nsigma    = cl.get_value<float>  ("nsigma",          3.0);
72 <  vector<float>  etbins    = cl.get_vector<float> ("etbins",  dflt_etbins);
73 <  vector<float>  etabins   = cl.get_vector<float> ("etabins",dflt_etabins);
74 <  vector<float>  phibins   = cl.get_vector<float> ("phibins",dflt_phibins);
75 <  vector<float>  emfbins   = cl.get_vector<float> ("emfbins",dflt_emfbins);
72 >  vector<float>  etbins    = cl.get_vector<float> ("etbins",           "");
73 >  vector<float>  etabins   = cl.get_vector<float> ("etabins",          "");
74 >  vector<float>  phibins   = cl.get_vector<float> ("phibins",          "");
75 >  vector<float>  emfbins   = cl.get_vector<float> ("emfbins",          "");
76    bool           abseta    = cl.get_value<bool>   ("abseta",         true);
77    
78    if (!cl.check()) return 0;
79    cl.print();
80    
81 <
81 >  int netbins  = (etbins.size() >0) ? etbins.size() +1 : 0;
82 >  int netabins = (etabins.size()>0) ? etabins.size()+1 : 0;
83 >  int nphibins = (phibins.size()>0) ? phibins.size()+1 : 0;
84 >  int nemfbins = (emfbins.size()>0) ? emfbins.size()+1 : 0;
85 >  
86 >  if (netbins+netabins+nphibins+nemfbins==0) {
87 >    cout<<"Must specify bins: etbins, etabins, phibins, AND/OR emfbins!"<<endl;
88 >    return 0;
89 >  }
90 >  
91    // etabins
92 <  if (!abseta) {
92 >  if (netabins>0&&!abseta) {
93      int neta=(int)etabins.size();
94      std::reverse(etabins.begin(),etabins.end());
95      for (int ieta=neta-1;ieta>=0;ieta--) etabins.push_back(etabins[ieta]);
# Line 117 | Line 118 | int main(int argc,char**argv)
118    float jtgendr[100];
119    float jtjes[100][3];
120    
120 –  int netbins  = etbins.size()+1;
121 –  int netabins = etabins.size()+1;
122 –  int nphibins = phibins.size()+1;
123 –  int nemfbins = emfbins.size()+1;
124 –  
121    argc=3; argv[1]="-l"; argv[2]="-b";
122    TRint* app = new TRint(argv[0],&argc,argv); app->Argc();
123    
# Line 132 | Line 128 | int main(int argc,char**argv)
128    TFile* plotfile = new TFile(output.c_str(),"UPDATE");
129    
130    for (unsigned int i=0;i<input.size();++i) {
131 <    TFile* f = new TFile(input[i].c_str(),"UPDATE"); if (!f->IsOpen()) return 0;
131 >    TFile* f = new TFile(input[i].c_str(),"READ"); if (!f->IsOpen()) return 0;
132      TTree* t = (TTree*)f->Get(treename.c_str());
133      
134      TEventList* el = new TEventList("sel","sel");
# Line 146 | Line 142 | int main(int argc,char**argv)
142      t->SetBranchAddress("jtemf",  jtemf);
143      t->SetBranchAddress("jtgenet",jtgenet);
144      t->SetBranchAddress("jtgendr",jtgendr);
145 <
146 <    if (t->FindBranch("jtjes")) t->SetBranchAddress("jtjes", jtjes);
147 <    else
145 >    
146 >    vector<TBranch*> branches;
147 >    branches.push_back(t->GetBranch("weight"));
148 >    branches.push_back(t->GetBranch("njt"));
149 >    branches.push_back(t->GetBranch("jtet"));
150 >    branches.push_back(t->GetBranch("jteta"));
151 >    branches.push_back(t->GetBranch("jtphi"));
152 >    branches.push_back(t->GetBranch("jtemf"));
153 >    branches.push_back(t->GetBranch("jtgenet"));
154 >    branches.push_back(t->GetBranch("jtgendr"));
155 >    
156 >    if (t->FindBranch("jtjes")) {
157 >      t->SetBranchAddress("jtjes", jtjes);
158 >      branches.push_back(t->GetBranch("jtjes"));
159 >    }
160 >    else {
161        for (int i1=0;i1<100;i1++)
162          for (int i2=0;i2<3;i2++) jtjes[i1][i2]=1.0;
163 <    
163 >    }
164      string absrsp_xtitle = "E_{T}^{gen}-E_{T} [GeV]";
165      string relrsp_xtitle = "E_{T}^{rec}/E_{T}";
166  
158 –    TH1F** hGenEt       = hist::initHistos("GenEt",jetalgs[i],100,
159 –                                           "jet E_{T}^{gen} [GeV]","GenEt",etbins);
160 –    TH1F** hEtCorrGenEt = hist::initHistos("EtCorrGenEt",jetalgs[i],100,
161 –                                           "jet E_{T}^{gen} [GeV]","GenEt",etbins);
162 –    TH1F** hAbsRspGenEt = hist::initHistos("AbsRspGenEt",jetalgs[i],100,-50,150,
163 –                                           absrsp_xtitle,"GenEt",etbins);
164 –    TH1F** hRelRspGenEt = hist::initHistos("RelRspGenEt",jetalgs[i],100,0,2,
165 –                                           relrsp_xtitle,"GenEt",etbins);
166 –    
167 –    TH1F** hEt       = hist::initHistos("Et",jetalgs[i],100,
168 –                                        "jet E_{T} [GeV]","Et",etbins);
169 –    TH1F** hEtCorrEt = hist::initHistos("EtCorrEt",jetalgs[i],100,
170 –                                        "jet E_{T} [GeV]","Et",etbins);
171 –    TH1F** hAbsRspEt = hist::initHistos("AbsRspEt",jetalgs[i],100,-50,150,
172 –                                        absrsp_xtitle,"Et",etbins);
173 –    TH1F** hRelRspEt = hist::initHistos("RelRspEt",jetalgs[i],100,0,2,
174 –                                        relrsp_xtitle,"Et",etbins);
175 –    
176 –    TH1F** hEta       = hist::initHistos("Eta",jetalgs[i],100,
177 –                                         "jet #eta","Eta",etabins);
178 –    TH1F** hEtCorrEta = hist::initHistos("EtCorrEta",jetalgs[i],100,
179 –                                         "jet E_{T} [GeV]","Eta",etabins);
180 –    TH1F** hAbsRspEta = hist::initHistos("AbsRspEta",jetalgs[i],100,-50,150,
181 –                                         absrsp_xtitle,"Eta",etabins);
182 –    TH1F** hRelRspEta = hist::initHistos("RelRspEta",jetalgs[i],100,0,2,
183 –                                         relrsp_xtitle,"Eta",etabins);
184 –    
185 –    TH1F** hPhi       = hist::initHistos("Phi",jetalgs[i],100,
186 –                                         "jet #phi","Phi",phibins);
187 –    TH1F** hEtCorrPhi = hist::initHistos("EtCorrPhi",jetalgs[i],100,
188 –                                         "jet E_{T} [Gev]","Phi",phibins);
189 –    TH1F** hAbsRspPhi = hist::initHistos("AbsRspPhi",jetalgs[i],100,-50,150,
190 –                                         absrsp_xtitle,"Phi",phibins);
191 –    TH1F** hRelRspPhi = hist::initHistos("RelRspPhi",jetalgs[i],100,0,2,
192 –                                         relrsp_xtitle,"Phi",phibins);
193 –    
194 –    TH1F** hEmf       = hist::initHistos("Emf",jetalgs[i],100,
195 –                                         "jet emf","Emf",emfbins);
196 –    TH1F** hEtCorrEmf = hist::initHistos("EtCorrEmf",jetalgs[i],100,
197 –                                         "jet E_{T} [Gev]","Emf",emfbins);
198 –    TH1F** hAbsRspEmf = hist::initHistos("AbsRspEmf",jetalgs[i],100,-50,150,
199 –                                         absrsp_xtitle,"Emf",emfbins);
200 –    TH1F** hRelRspEmf = hist::initHistos("RelRspEmf",jetalgs[i],100,0,2,
201 –                                         relrsp_xtitle,"Emf",emfbins);
202 –    
203 –    TH1F*** hEtEtaEt     = hist::initHistos("Et",jetalgs[i],100,"jet E_{T} [GeV]",
204 –                                            "Eta",etabins,"Et",etbins);
205 –    TH1F*** hEtCorrEtaEt = hist::initHistos("EtCorrEtaEt",jetalgs[i],100,
206 –                                            "jet E_{T} [GeV]",
207 –                                            "Eta",etabins,"Et",etbins);
208 –    TH1F*** hAbsRspEtaEt = hist::initHistos("AbsRspEtaEt",jetalgs[i],100,-50,150,
209 –                                            absrsp_xtitle,
210 –                                            "Eta",etabins,"Et",etbins);
211 –    TH1F*** hRelRspEtaEt = hist::initHistos("RelRspEtaEt",jetalgs[i],100,0,2,
212 –                                            relrsp_xtitle,
213 –                                            "Eta",etabins,"Et",etbins);
214 –    
215 –    TH1F*** hEtEmfEt     = hist::initHistos("Et",jetalgs[i],100,"jet E_{T} [GeV]",
216 –                                            "Emf",emfbins,"Et",etbins);
217 –    TH1F*** hEtCorrEmfEt = hist::initHistos("EtCorrEmfEt",jetalgs[i],100,
218 –                                            "jet E_{T} [GeV]",
219 –                                            "Emf",emfbins,"Et",etbins);
220 –    TH1F*** hAbsRspEmfEt = hist::initHistos("AbsRspEmfEt",jetalgs[i],100,-50,150,
221 –                                            absrsp_xtitle,
222 –                                            "Emf",emfbins,"Et",etbins);
223 –    TH1F*** hRelRspEmfEt = hist::initHistos("RelRspEmfEt",jetalgs[i],100,0,2,
224 –                                            relrsp_xtitle,
225 –                                            "Emf",emfbins,"Et",etbins);
226 –    
227 –    TH1F**** hEtEmfEtaEt     = hist::initHistos("Et",jetalgs[i],100,"jet E_{T} [GeV]",
228 –                                                "Emf",emfbins,
229 –                                                "Eta",etabins,
230 –                                                "Et", etbins);
231 –    TH1F**** hEtCorrEmfEtaEt = hist::initHistos("EtCorrEmfEtaEt",jetalgs[i],100,
232 –                                                "jet E_{T} [GeV]",
233 –                                                "Emf",emfbins,
234 –                                                "Eta",etabins,
235 –                                                "Et", etbins);
236 –    TH1F**** hAbsRspEmfEtaEt = hist::initHistos("AbsRspEmfEtaEt",jetalgs[i],
237 –                                                100,-50,150,absrsp_xtitle,
238 –                                                "Emf",emfbins,
239 –                                                "Eta",etabins,
240 –                                                "Et", etbins);
241 –    TH1F**** hRelRspEmfEtaEt = hist::initHistos("RelRspEmfEtaEt",jetalgs[i],
242 –                                                100,0,2,relrsp_xtitle,
243 –                                                "Emf",emfbins,
244 –                                                "Eta",etabins,
245 –                                                "Et", etbins);
167      
168 +    TH1F** hGenEt(0);
169 +    TH1F** hEtCorrGenEt(0);
170 +    TH1F** hAbsRspGenEt(0);
171 +    TH1F** hRelRspGenEt(0);
172 +    
173 +    if (netbins>0) {
174 +      hGenEt       = hist::initHistos("GenEt",jetalgs[i],100,
175 +                                      "jet E_{T}^{gen} [GeV]",
176 +                                      "GenEt",etbins);
177 +      hEtCorrGenEt = hist::initHistos("EtCorrGenEt",jetalgs[i],100,
178 +                                      "jet E_{T}^{gen} [GeV]",
179 +                                      "GenEt",etbins);
180 +      hAbsRspGenEt = hist::initHistos("AbsRspGenEt",jetalgs[i],
181 +                                      100,-50,150,
182 +                                      absrsp_xtitle,
183 +                                      "GenEt",etbins);
184 +      hRelRspGenEt = hist::initHistos("RelRspGenEt",jetalgs[i],
185 +                                      100,0,2,
186 +                                      relrsp_xtitle,
187 +                                      "GenEt",etbins);
188 +    }
189 +    
190 +
191 +    TH1F** hEt(0);
192 +    TH1F** hEtCorrEt(0);
193 +    TH1F** hAbsRspEt(0);
194 +    TH1F** hRelRspEt(0);
195 +    
196 +    if (netbins>0) {
197 +      hEt       = hist::initHistos("Et",jetalgs[i],100,
198 +                                   "jet E_{T} [GeV]","Et",etbins);
199 +      hEtCorrEt = hist::initHistos("EtCorrEt",jetalgs[i],100,
200 +                                   "jet E_{T} [GeV]","Et",etbins);
201 +      hAbsRspEt = hist::initHistos("AbsRspEt",jetalgs[i],100,-50,150,
202 +                                   absrsp_xtitle,"Et",etbins);
203 +      hRelRspEt = hist::initHistos("RelRspEt",jetalgs[i],100,0,2,
204 +                                   relrsp_xtitle,"Et",etbins);
205 +    }
206 +    
207 +
208 +    TH1F** hEta(0);
209 +    TH1F** hEtCorrEta(0);
210 +    TH1F** hAbsRspEta(0);
211 +    TH1F** hRelRspEta(0);
212 +    
213 +    if (netabins>0) {
214 +      hEta       = hist::initHistos("Eta",jetalgs[i],100,
215 +                                    "jet #eta","Eta",etabins);
216 +      hEtCorrEta = hist::initHistos("EtCorrEta",jetalgs[i],100,
217 +                                    "jet E_{T} [GeV]","Eta",etabins);
218 +      hAbsRspEta = hist::initHistos("AbsRspEta",jetalgs[i],100,-50,150,
219 +                                    absrsp_xtitle,"Eta",etabins);
220 +      hRelRspEta = hist::initHistos("RelRspEta",jetalgs[i],100,0,2,
221 +                                    relrsp_xtitle,"Eta",etabins);
222 +    }
223 +    
224 +
225 +    TH1F** hPhi(0);
226 +    TH1F** hEtCorrPhi(0);
227 +    TH1F** hAbsRspPhi(0);
228 +    TH1F** hRelRspPhi(0);
229 +    
230 +    if (nphibins>0) {
231 +      hPhi       = hist::initHistos("Phi",jetalgs[i],100,
232 +                                    "jet #phi","Phi",phibins);
233 +      hEtCorrPhi = hist::initHistos("EtCorrPhi",jetalgs[i],100,
234 +                                    "jet E_{T} [Gev]","Phi",phibins);
235 +      hAbsRspPhi = hist::initHistos("AbsRspPhi",jetalgs[i],100,-50,150,
236 +                                    absrsp_xtitle,"Phi",phibins);
237 +      hRelRspPhi = hist::initHistos("RelRspPhi",jetalgs[i],100,0,2,
238 +                                    relrsp_xtitle,"Phi",phibins);
239 +    }
240 +    
241 +
242 +    TH1F** hEmf(0);
243 +    TH1F** hEtCorrEmf(0);
244 +    TH1F** hAbsRspEmf(0);
245 +    TH1F** hRelRspEmf(0);
246 +    
247 +    if (nemfbins>0) {
248 +      hEmf       = hist::initHistos("Emf",jetalgs[i],100,
249 +                                    "jet emf","Emf",emfbins);
250 +      hEtCorrEmf = hist::initHistos("EtCorrEmf",jetalgs[i],100,
251 +                                    "jet E_{T} [Gev]","Emf",emfbins);
252 +      hAbsRspEmf = hist::initHistos("AbsRspEmf",jetalgs[i],100,-50,150,
253 +                                    absrsp_xtitle,"Emf",emfbins);
254 +      hRelRspEmf = hist::initHistos("RelRspEmf",jetalgs[i],100,0,2,
255 +                                    relrsp_xtitle,"Emf",emfbins);
256 +    }
257 +    
258 +    
259 +    TH1F*** hEtEtaEt(0);
260 +    TH1F*** hEtCorrEtaEt(0);
261 +    TH1F*** hAbsRspEtaEt(0);
262 +    TH1F*** hRelRspEtaEt(0);
263 +    
264 +    if (netbins>0&&netabins>0) {
265 +      hEtEtaEt     = hist::initHistos("Et",jetalgs[i],100,"jet E_{T} [GeV]",
266 +                                      "Eta",etabins,"Et",etbins);
267 +      hEtCorrEtaEt = hist::initHistos("EtCorrEtaEt",jetalgs[i],100,
268 +                                      "jet E_{T} [GeV]",
269 +                                      "Eta",etabins,"Et",etbins);
270 +      hAbsRspEtaEt = hist::initHistos("AbsRspEtaEt",jetalgs[i],100,-50,150,
271 +                                      absrsp_xtitle,
272 +                                      "Eta",etabins,"Et",etbins);
273 +      hRelRspEtaEt = hist::initHistos("RelRspEtaEt",jetalgs[i],100,0,2,
274 +                                      relrsp_xtitle,
275 +                                      "Eta",etabins,"Et",etbins);
276 +    }
277 +
278 +
279 +    TH1F*** hEtEmfEt(0);
280 +    TH1F*** hEtCorrEmfEt(0);
281 +    TH1F*** hAbsRspEmfEt(0);
282 +    TH1F*** hRelRspEmfEt(0);
283 +    
284 +    if (netbins>0&&nemfbins>0) {
285 +      hEtEmfEt     = hist::initHistos("Et",jetalgs[i],100,"jet E_{T} [GeV]",
286 +                                      "Emf",emfbins,"Et",etbins);
287 +      hEtCorrEmfEt = hist::initHistos("EtCorrEmfEt",jetalgs[i],100,
288 +                                      "jet E_{T} [GeV]",
289 +                                      "Emf",emfbins,"Et",etbins);
290 +      hAbsRspEmfEt = hist::initHistos("AbsRspEmfEt",jetalgs[i],100,-50,150,
291 +                                      absrsp_xtitle,
292 +                                      "Emf",emfbins,"Et",etbins);
293 +      hRelRspEmfEt = hist::initHistos("RelRspEmfEt",jetalgs[i],100,0,2,
294 +                                      relrsp_xtitle,
295 +                                      "Emf",emfbins,"Et",etbins);
296 +    }
297 +    
298 +
299 +    TH1F**** hEtEmfEtaEt(0);
300 +    TH1F**** hEtCorrEmfEtaEt(0);
301 +    TH1F**** hAbsRspEmfEtaEt(0);
302 +    TH1F**** hRelRspEmfEtaEt(0);
303 +    
304 +    if (netbins>0&&netabins>0&&nemfbins>0) {
305 +      hEtEmfEtaEt     = hist::initHistos("Et",jetalgs[i],100,"jet E_{T} [GeV]",
306 +                                         "Emf",emfbins,
307 +                                         "Eta",etabins,
308 +                                         "Et", etbins);
309 +      hEtCorrEmfEtaEt = hist::initHistos("EtCorrEmfEtaEt",jetalgs[i],100,
310 +                                         "jet E_{T} [GeV]",
311 +                                         "Emf",emfbins,
312 +                                         "Eta",etabins,
313 +                                         "Et", etbins);
314 +      hAbsRspEmfEtaEt = hist::initHistos("AbsRspEmfEtaEt",jetalgs[i],
315 +                                         100,-50,150,absrsp_xtitle,
316 +                                         "Emf",emfbins,
317 +                                         "Eta",etabins,
318 +                                         "Et", etbins);
319 +      hRelRspEmfEtaEt = hist::initHistos("RelRspEmfEtaEt",jetalgs[i],
320 +                                         100,0,2,relrsp_xtitle,
321 +                                         "Emf",emfbins,
322 +                                         "Eta",etabins,
323 +                                         "Et", etbins);
324 +    }
325      
326      int nevts = el->GetN();
327      cout<<jetalgs[i]<<": "<<nevts<<" events selected."<<endl;
328      
329      // prepare intermediate trees
330      float genet,et,eta,phi,emf,etcorr,absrsp,relrsp;
331 <    TTree**   tGenEt      = new TTree*[netbins];
332 <    TTree**   tEt         = new TTree*[netbins];
333 <    TTree**   tEta        = new TTree*[netabins];
334 <    TTree**   tPhi        = new TTree*[nphibins];
335 <    TTree**   tEmf        = new TTree*[nemfbins];
336 <    TTree***  tEtEtaEt    = new TTree**[netabins];
337 <    TTree***  tEtEmfEt    = new TTree**[nemfbins];
338 <    TTree**** tEtEmfEtaEt = new TTree***[nemfbins];
331 >
332 >    TTree**   tGenEt(0);
333 >    TTree**   tEt(0);
334 >    TTree**   tEta(0);
335 >    TTree**   tPhi(0);
336 >    TTree**   tEmf(0);
337 >    TTree***  tEtEtaEt(0);
338 >    TTree***  tEtEmfEt(0);
339 >    TTree**** tEtEmfEtaEt(0);
340 >
341 >    gROOT->cd();
342 >    
343 >    if (netbins>0)  tGenEt = new TTree*[netbins];
344 >    if (netbins>0)  tEt    = new TTree*[netbins];
345 >    if (netabins>0) tEta   = new TTree*[netabins];
346 >    if (nphibins>0) tPhi   = new TTree*[nphibins];
347 >    if (nemfbins>0) tEmf   = new TTree*[nemfbins];
348 >    if (netbins>0&&netabins>0) tEtEtaEt = new TTree**[netabins];
349 >    if (netbins>0&&nemfbins>0) tEtEmfEt = new TTree**[nemfbins];
350 >    if (netbins>0&&netabins>0&&nemfbins>0) tEtEmfEtaEt = new TTree***[nemfbins];
351      
352      for (int iet=0;iet<netbins;iet++) {
353        
# Line 299 | Line 389 | int main(int argc,char**argv)
389        tEmf[iemf]->Branch("relrsp",&relrsp,"relrsp/F");
390        tEmf[iemf]->Branch("weight",&weight,"weight/F");
391      }
392 <    for (int ieta=0;ieta<netabins;ieta++) {
393 <      tEtEtaEt[ieta] = new TTree*[netbins];
304 <      for (int iet=0;iet<netbins;iet++) {
305 <        tEtEtaEt[ieta][iet] = new TTree("tEtEtaEt","tEtEtaEt");
306 <        tEtEtaEt[ieta][iet]->Branch("val",   &et,    "val/F");
307 <        tEtEtaEt[ieta][iet]->Branch("etcorr",&etcorr,"etcorr/F");
308 <        tEtEtaEt[ieta][iet]->Branch("absrsp",&absrsp,"absrsp/F");
309 <        tEtEtaEt[ieta][iet]->Branch("relrsp",&relrsp,"relrsp/F");
310 <        tEtEtaEt[ieta][iet]->Branch("weight",&weight,"weight/F");
311 <      }
312 <    }
313 <    for (int iemf=0;iemf<nemfbins;iemf++) {
314 <      tEtEmfEt[iemf] = new TTree*[netbins];
315 <      for (int iet=0;iet<netbins;iet++) {
316 <        tEtEmfEt[iemf][iet] = new TTree("tEtEmfEt","tEtEmfEt");
317 <        tEtEmfEt[iemf][iet]->Branch("val",   &et,    "val/F");
318 <        tEtEmfEt[iemf][iet]->Branch("etcorr",&etcorr,"etcorr/F");
319 <        tEtEmfEt[iemf][iet]->Branch("absrsp",&absrsp,"absrsp/F");
320 <        tEtEmfEt[iemf][iet]->Branch("relrsp",&relrsp,"relrsp/F");
321 <        tEtEmfEt[iemf][iet]->Branch("weight",&weight,"weight/F");
322 <      }
323 <    }
324 <    for (int iemf=0;iemf<nemfbins;iemf++) {
325 <      tEtEmfEtaEt[iemf] = new TTree**[netabins];
392 >
393 >    if (netbins>0) {
394        for (int ieta=0;ieta<netabins;ieta++) {
395 <        tEtEmfEtaEt[iemf][ieta] = new TTree*[netbins];
395 >        tEtEtaEt[ieta] = new TTree*[netbins];
396 >        for (int iet=0;iet<netbins;iet++) {
397 >          tEtEtaEt[ieta][iet] = new TTree("tEtEtaEt","tEtEtaEt");
398 >          tEtEtaEt[ieta][iet]->Branch("val",   &et,    "val/F");
399 >          tEtEtaEt[ieta][iet]->Branch("etcorr",&etcorr,"etcorr/F");
400 >          tEtEtaEt[ieta][iet]->Branch("absrsp",&absrsp,"absrsp/F");
401 >          tEtEtaEt[ieta][iet]->Branch("relrsp",&relrsp,"relrsp/F");
402 >          tEtEtaEt[ieta][iet]->Branch("weight",&weight,"weight/F");
403 >        }
404 >      }
405 >      for (int iemf=0;iemf<nemfbins;iemf++) {
406 >        tEtEmfEt[iemf] = new TTree*[netbins];
407          for (int iet=0;iet<netbins;iet++) {
408 <          tEtEmfEtaEt[iemf][ieta][iet] = new TTree("tEtEmfEtaEt","tEtEmfEtaEt");
409 <          tEtEmfEtaEt[iemf][ieta][iet]->Branch("val",   &et,    "val/F");
410 <          tEtEmfEtaEt[iemf][ieta][iet]->Branch("etcorr",&etcorr,"etcorr/F");
411 <          tEtEmfEtaEt[iemf][ieta][iet]->Branch("absrsp",&absrsp,"absrsp/F");
412 <          tEtEmfEtaEt[iemf][ieta][iet]->Branch("relrsp",&relrsp,"relrsp/F");
413 <          tEtEmfEtaEt[iemf][ieta][iet]->Branch("weight",&weight,"weight/F");
408 >          tEtEmfEt[iemf][iet] = new TTree("tEtEmfEt","tEtEmfEt");
409 >          tEtEmfEt[iemf][iet]->Branch("val",   &et,    "val/F");
410 >          tEtEmfEt[iemf][iet]->Branch("etcorr",&etcorr,"etcorr/F");
411 >          tEtEmfEt[iemf][iet]->Branch("absrsp",&absrsp,"absrsp/F");
412 >          tEtEmfEt[iemf][iet]->Branch("relrsp",&relrsp,"relrsp/F");
413 >          tEtEmfEt[iemf][iet]->Branch("weight",&weight,"weight/F");
414 >        }
415 >      }
416 >      if (netabins>0) {
417 >        for (int iemf=0;iemf<nemfbins;iemf++) {
418 >          tEtEmfEtaEt[iemf] = new TTree**[netabins];
419 >          for (int ieta=0;ieta<netabins;ieta++) {
420 >            tEtEmfEtaEt[iemf][ieta] = new TTree*[netbins];
421 >            for (int iet=0;iet<netbins;iet++) {
422 >              tEtEmfEtaEt[iemf][ieta][iet] = new TTree("tEtEmfEtaEt",
423 >                                                       "tEtEmfEtaEt");
424 >              tEtEmfEtaEt[iemf][ieta][iet]->Branch("val",   &et,    "val/F");
425 >              tEtEmfEtaEt[iemf][ieta][iet]->Branch("etcorr",&etcorr,"etcorr/F");
426 >              tEtEmfEtaEt[iemf][ieta][iet]->Branch("absrsp",&absrsp,"absrsp/F");
427 >              tEtEmfEtaEt[iemf][ieta][iet]->Branch("relrsp",&relrsp,"relrsp/F");
428 >              tEtEmfEtaEt[iemf][ieta][iet]->Branch("weight",&weight,"weight/F");
429 >            }
430 >          }
431          }
432        }
433      }
434      
435      cout<<jetalgs[i]<<": trees created."<<endl;
436  
341 –    
437      // loop over all events
438      for (int ievt=0;ievt<nevts;ievt++) {
439        
440 <      t->GetEntry(el->GetEntry(ievt));
441 <      
440 >      //t->GetEntry(el->GetEntry(ievt));
441 >      for (int ibrnch=0;ibrnch<(int)branches.size();ibrnch++)
442 >        branches[ibrnch]->GetEntry(el->GetEntry(ievt));
443 >
444        // loop over all jets in the event
445        for (int ijt=0;ijt<njt;ijt++) {
446          
447          if (jtgendr[ijt]<0.0||jtgendr[ijt]>drmax) continue;
448  
449 <        // if (ijt>0) continue;
450 <
449 >        //if (ijt>ijtmax) continue;
450 >        
451          et      = jtet[ijt];
452          genet   = jtgenet[ijt];
453          eta     = (abseta) ? std::abs(jteta[ijt]) : jteta[ijt];
# Line 364 | Line 461 | int main(int argc,char**argv)
461  
462          if (TMath::IsNaN(emf)) emf = 0.0;
463          
464 <        
465 <        // genet
466 <        int igenet = hist::get_ibin(genet,etbins);
467 <        tGenEt[igenet]->Fill();
468 <        
469 <        // et
470 <        int iet = hist::get_ibin(et,etbins);
471 <        tEt[iet]->Fill();
472 <        
473 <        // eta
474 <        int ieta = hist::get_ibin(eta,etabins);
475 <        tEta[ieta]->Fill();
476 <        
477 <        // phi
478 <        int iphi = hist::get_ibin(phi,phibins);
479 <        tPhi[iphi]->Fill();
480 <        
384 <        // emf
385 <        int iemf = hist::get_ibin(emf,emfbins);
386 <        tEmf[iemf]->Fill();
387 <        
388 <        // et-eta
389 <        tEtEtaEt[ieta][iet]->Fill();
390 <        
391 <        // et-emf
392 <        tEtEmfEt[iemf][iet]->Fill();
393 <        
394 <        // et-eta-emf
395 <        tEtEmfEtaEt[iemf][ieta][iet]->Fill();
464 >        int igenet = (netbins>0)  ? hist::get_ibin(genet,etbins) : -1;
465 >        int iet    = (netbins>0)  ? hist::get_ibin(et,etbins)    : -1;
466 >        int ieta   = (netabins>0) ? hist::get_ibin(eta,etabins)  : -1;
467 >        int iphi   = (nphibins>0) ? hist::get_ibin(phi,phibins)  : -1;
468 >        int iemf   = (nemfbins>0) ? hist::get_ibin(emf,emfbins)  : -1;
469 >        
470 >        if (netbins>0)  tGenEt[igenet]->Fill();
471 >        if (netbins>0)  tEt[iet]->Fill();
472 >        if (netabins>0) tEta[ieta]->Fill();
473 >        if (nphibins>0) tPhi[iphi]->Fill();
474 >        if (nemfbins>0) tEmf[iemf]->Fill();
475 >
476 >        if (netbins>0&&netabins>0) tEtEtaEt[ieta][iet]->Fill();
477 >        if (netbins>0&&nemfbins>0) tEtEmfEt[iemf][iet]->Fill();
478 >
479 >        if (netbins>0&&netabins>0&&nemfbins>0)
480 >          tEtEmfEtaEt[iemf][ieta][iet]->Fill();
481          
482        } // jets
483      } // evts
# Line 422 | Line 507 | int main(int argc,char**argv)
507      cout<<jetalgs[i]<<": all histograms replaced, and tree deleted."<<endl;
508  
509      
510 +    plotfile->cd();
511      TDirectory* d = plotfile->GetDirectory(jetalgs[i].c_str());
512      if (0!=d) plotfile->rmdir(jetalgs[i].c_str());
513      d = plotfile->mkdir(jetalgs[i].c_str());
# Line 430 | Line 516 | int main(int argc,char**argv)
516      // fit response histos
517      fitHistos(hAbsRspGenEt,netbins,nsigma);
518      fitHistos(hRelRspGenEt,netbins,nsigma);
519 <
519 >    
520      fitHistos(hAbsRspEt,   netbins,nsigma);
521      fitHistos(hRelRspEt,   netbins,nsigma);
522  
# Line 439 | Line 525 | int main(int argc,char**argv)
525  
526      fitHistos(hAbsRspPhi,  nphibins,nsigma);
527      fitHistos(hRelRspPhi,  nphibins,nsigma);
528 <    
528 >
529      fitHistos(hAbsRspEmf,  nemfbins,nsigma);
530      fitHistos(hRelRspEmf,  nemfbins,nsigma);
531 <    
531 >
532      fitHistos(hAbsRspEtaEt,netabins,netbins,nsigma);
533      fitHistos(hRelRspEtaEt,netabins,netbins,nsigma);
534  
# Line 496 | Line 582 | int main(int argc,char**argv)
582      saveHistos(hAbsRspEmfEtaEt,nemfbins,netabins,netbins);
583      saveHistos(hRelRspEmfEtaEt,nemfbins,netabins,netbins);
584  
499 –    plotfile->Write();
500 –  
585      cout<<jetalgs[i]<<": histograms saved."<<endl;
586    }
587    
588    plotfile->Write();
589    plotfile->Close();
590    delete plotfile;
591 +
592 +  app->Run();
593    
594    return 0;
595   }
# Line 519 | Line 605 | void replaceHistos(int nbins,
605                     TH1F**& hVar,TH1F**& hEtCorr,TH1F**& hAbsRsp,TH1F**& hRelRsp,
606                     TTree**& tVar)
607   {
608 +  if (nbins==0) return;
609 +  
610    TH1::SetDefaultSumw2();
611    for (int ibin=0;ibin<nbins;ibin++) {
612      
# Line 570 | Line 658 | void replaceHistos(int nbins1,int nbins2
658                     TH1F***& hVar,TH1F***& hEtCorr,TH1F***& hAbsRsp,TH1F***& hRelRsp,
659                     TTree***& tVar)
660   {
661 +  if (nbins1==0) return;
662 +  
663    for (int i1=0;i1<nbins1;i1++)
664      replaceHistos(nbins2,hVar[i1],hEtCorr[i1],hAbsRsp[i1],hRelRsp[i1],tVar[i1]);
665    delete [] tVar;
# Line 581 | Line 671 | void replaceHistos(int nbins1,int nbins2
671                     TH1F****& hVar,TH1F****& hEtCorr,TH1F****&hAbsRsp,TH1F****&hRelRsp,
672                     TTree****& tVar)
673   {
674 +  if (nbins1==0) return;
675    for (int i1=0;i1<nbins1;i1++)
676      replaceHistos(nbins2,nbins3,
677                    hVar[i1],hEtCorr[i1],hAbsRsp[i1],hRelRsp[i1],tVar[i1]);
# Line 601 | Line 692 | void fitHistos(TH1F** histos,int n,doubl
692      f->SetParameter(2,rms);
693      f->SetLineWidth(1);
694      f->SetLineColor(kRed);
695 <    histos[i]->Fit(f,"QR");
695 >    histos[i]->Fit(f,"QR"); // Avoid TCanvas?
696    }
697   }
698  
# Line 609 | Line 700 | void fitHistos(TH1F** histos,int n,doubl
700   //______________________________________________________________________________
701   void fitHistos(TH1F*** histos,int n1,int n2,double nsigma)
702   {
703 +  if (n1==0||n2==0) return;
704    for (int i1=0;i1<n1;i1++) fitHistos(histos[i1],n2,nsigma);
705   }
706  
# Line 616 | Line 708 | void fitHistos(TH1F*** histos,int n1,int
708   //______________________________________________________________________________
709   void fitHistos(TH1F**** histos,int n1,int n2,int n3,double nsigma)
710   {
711 +  if (n1==0||n2==0||n3==0) return;
712    for (int i1=0;i1<n1;i1++) fitHistos(histos[i1],n2,n3,nsigma);
713   }
714  
# Line 623 | Line 716 | void fitHistos(TH1F**** histos,int n1,in
716   //______________________________________________________________________________
717   void saveHistos(TH1F**& histos,int n)
718   {
719 +  if (n==0) return;
720    for (int i=0;i<n;i++) {
721      histos[i]->Write();
722      delete histos[i];
# Line 634 | Line 728 | void saveHistos(TH1F**& histos,int n)
728   //______________________________________________________________________________
729   void saveHistos(TH1F***& histos,int n1,int n2)
730   {
731 +  if (n1==0) return;
732    for (int i1=0;i1<n1;i1++) saveHistos(histos[i1],n2);
733    delete [] histos;
734   }
# Line 642 | Line 737 | void saveHistos(TH1F***& histos,int n1,i
737   //______________________________________________________________________________
738   void saveHistos(TH1F****& histos,int n1,int n2,int n3)
739   {
740 +  if (n1==0) return;
741    for (int i1=0;i1<n1;i1++) saveHistos(histos[i1],n2,n3);
742    delete [] histos;
743   }

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines