1 |
schiefer |
1.1 |
////////////////////////////////////////////////////////////////////////////////
|
2 |
|
|
//
|
3 |
|
|
// mcresponse_x
|
4 |
|
|
// ------------
|
5 |
|
|
//
|
6 |
|
|
// 06/16/2007 Philipp Schieferdecker <philipp.schieferdecker@cern.ch>
|
7 |
|
|
////////////////////////////////////////////////////////////////////////////////
|
8 |
|
|
|
9 |
|
|
|
10 |
|
|
#include "utils/cmdline.h"
|
11 |
schiefer |
1.2 |
#include "utils/hist.h"
|
12 |
schiefer |
1.1 |
|
13 |
|
|
#include <TROOT.h>
|
14 |
|
|
#include <TRint.h>
|
15 |
|
|
#include <TFile.h>
|
16 |
|
|
#include <TTree.h>
|
17 |
|
|
#include <TEventList.h>
|
18 |
|
|
#include <TH1F.h>
|
19 |
|
|
#include <TF1.h>
|
20 |
|
|
|
21 |
|
|
#include <iostream>
|
22 |
|
|
#include <iomanip>
|
23 |
|
|
#include <cmath>
|
24 |
|
|
#include <string>
|
25 |
|
|
#include <vector>
|
26 |
|
|
|
27 |
|
|
|
28 |
|
|
using namespace std;
|
29 |
|
|
|
30 |
|
|
|
31 |
|
|
////////////////////////////////////////////////////////////////////////////////
|
32 |
|
|
// defaults
|
33 |
|
|
////////////////////////////////////////////////////////////////////////////////
|
34 |
schiefer |
1.4 |
string dflt_etbins ="8,12,15,20,25,30,45,65,90,120,180";
|
35 |
schiefer |
1.3 |
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";
|
36 |
schiefer |
1.1 |
string dflt_phibins="-2.75,-2.25,-1.75.,-1.25,-.75,-.25,.25,.75,1.25,1.75,2.25,2.75";
|
37 |
schiefer |
1.4 |
string dflt_emfbins="0.2,0.4,0.6,0.8";
|
38 |
schiefer |
1.1 |
|
39 |
|
|
|
40 |
|
|
////////////////////////////////////////////////////////////////////////////////
|
41 |
|
|
// declare global functions
|
42 |
|
|
////////////////////////////////////////////////////////////////////////////////
|
43 |
schiefer |
1.4 |
void replaceHistos(int nbins,
|
44 |
|
|
TH1F**& hVar,TH1F**& hEtCorr,TH1F**& hAbsRsp,TH1F**& hRelRsp,
|
45 |
|
|
TTree**& tVar);
|
46 |
|
|
void replaceHistos(int nbins1,int nbins2,
|
47 |
|
|
TH1F***&hVar,TH1F***&hEtCorr,TH1F***&hAbsRsp,TH1F***&hRelRsp,
|
48 |
|
|
TTree***& tVar);
|
49 |
|
|
void replaceHistos(int nbins1,int nbins2,int nbins3,
|
50 |
|
|
TH1F****&hVar,TH1F****&hEtCorr,
|
51 |
|
|
TH1F****&hAbsRsp,TH1F****&hRelRsp,
|
52 |
|
|
TTree****& tVar);
|
53 |
schiefer |
1.5 |
void fitHistos (TH1F** histos,int n,double nsigma);
|
54 |
|
|
void fitHistos (TH1F*** histos,int n1,int n2,double nsigma);
|
55 |
|
|
void fitHistos (TH1F**** histos,int n1,int n2,int n3,double nsigma);
|
56 |
|
|
void saveHistos (TH1F**& histos,int n);
|
57 |
|
|
void saveHistos (TH1F***& histos,int n1,int n2);
|
58 |
|
|
void saveHistos (TH1F****& histos,int n1,int n2,int n3);
|
59 |
schiefer |
1.1 |
|
60 |
|
|
|
61 |
|
|
|
62 |
|
|
////////////////////////////////////////////////////////////////////////////////
|
63 |
|
|
// main
|
64 |
|
|
////////////////////////////////////////////////////////////////////////////////
|
65 |
|
|
|
66 |
|
|
//______________________________________________________________________________
|
67 |
|
|
int main(int argc,char**argv)
|
68 |
|
|
{
|
69 |
|
|
cmdline cl;
|
70 |
|
|
cl.parse(argc,argv);
|
71 |
schiefer |
1.2 |
|
72 |
schiefer |
1.1 |
vector<string> input = cl.get_vector<string>("input");
|
73 |
schiefer |
1.4 |
string output = cl.get_value<string> ("output", "mcrsp.root");
|
74 |
|
|
string selection = cl.get_value<string> ("selection", "njt>0");
|
75 |
|
|
float drmax = cl.get_value<float> ("drmax", 0.3);
|
76 |
|
|
short applyjes = cl.get_value<short> ("applyjes", 0);
|
77 |
|
|
string treename = cl.get_value<string> ("treename", "t");
|
78 |
|
|
float nsigma = cl.get_value<float> ("nsigma", 3.0);
|
79 |
|
|
vector<float> etbins = cl.get_vector<float> ("etbins", dflt_etbins);
|
80 |
|
|
vector<float> etabins = cl.get_vector<float> ("etabins",dflt_etabins);
|
81 |
|
|
vector<float> phibins = cl.get_vector<float> ("phibins",dflt_phibins);
|
82 |
|
|
vector<float> emfbins = cl.get_vector<float> ("emfbins",dflt_emfbins);
|
83 |
|
|
bool abseta = cl.get_value<bool> ("abseta", true);
|
84 |
schiefer |
1.1 |
|
85 |
|
|
if (!cl.check()) return 0;
|
86 |
|
|
cl.print();
|
87 |
|
|
|
88 |
|
|
|
89 |
|
|
// etabins
|
90 |
|
|
if (!abseta) {
|
91 |
|
|
int neta=(int)etabins.size();
|
92 |
|
|
std::reverse(etabins.begin(),etabins.end());
|
93 |
|
|
for (int ieta=neta-1;ieta>=0;ieta--) etabins.push_back(etabins[ieta]);
|
94 |
|
|
for (int ieta=0;ieta<neta;ieta++) etabins[ieta]*=-1;
|
95 |
|
|
}
|
96 |
|
|
|
97 |
|
|
|
98 |
|
|
// jetalgs
|
99 |
|
|
vector<string> jetalgs;
|
100 |
|
|
for (vector<string>::const_iterator it=input.begin();it!=input.end();++it) {
|
101 |
|
|
string jetalg = *it;
|
102 |
|
|
string::size_type pos;
|
103 |
|
|
while ((pos=jetalg.find('/'))!=string::npos) jetalg = jetalg.substr(pos+1);
|
104 |
|
|
jetalg = jetalg.substr(0,jetalg.find(".root"));
|
105 |
|
|
jetalgs.push_back(jetalg);
|
106 |
|
|
}
|
107 |
|
|
|
108 |
|
|
// prepare branch variables and histograms / graphs
|
109 |
|
|
float weight;
|
110 |
|
|
char njt;
|
111 |
schiefer |
1.3 |
float jtet[100];
|
112 |
schiefer |
1.1 |
float jtpt[100];
|
113 |
|
|
float jteta[100];
|
114 |
|
|
float jtphi[100];
|
115 |
schiefer |
1.4 |
float jtemf[100];
|
116 |
schiefer |
1.3 |
float jtgenet[100];
|
117 |
schiefer |
1.1 |
float jtgenpt[100];
|
118 |
|
|
float jtgendr[100];
|
119 |
|
|
float jtjes[100][3];
|
120 |
|
|
|
121 |
schiefer |
1.3 |
int netbins = etbins.size()+1;
|
122 |
schiefer |
1.1 |
int netabins = etabins.size()+1;
|
123 |
|
|
int nphibins = phibins.size()+1;
|
124 |
schiefer |
1.4 |
int nemfbins = emfbins.size()+1;
|
125 |
schiefer |
1.1 |
|
126 |
schiefer |
1.3 |
argc=3; argv[1]="-l"; argv[2]="-b";
|
127 |
schiefer |
1.4 |
TRint* app = new TRint(argv[0],&argc,argv); app->Argc();
|
128 |
schiefer |
1.1 |
|
129 |
|
|
|
130 |
|
|
//
|
131 |
|
|
// loop over all input files
|
132 |
|
|
//
|
133 |
schiefer |
1.3 |
TFile* plotfile = new TFile(output.c_str(),"RECREATE");
|
134 |
schiefer |
1.2 |
|
135 |
schiefer |
1.1 |
for (unsigned int i=0;i<input.size();++i) {
|
136 |
schiefer |
1.3 |
TFile* f = new TFile(input[i].c_str(),"UPDATE"); if (!f->IsOpen()) return 0;
|
137 |
schiefer |
1.1 |
TTree* t = (TTree*)f->Get(treename.c_str());
|
138 |
|
|
|
139 |
|
|
TEventList* el = new TEventList("sel","sel");
|
140 |
|
|
t->Draw(">>sel",selection.c_str());
|
141 |
|
|
|
142 |
|
|
t->SetBranchAddress("weight", &weight);
|
143 |
|
|
t->SetBranchAddress("njt", &njt);
|
144 |
schiefer |
1.3 |
t->SetBranchAddress("jtet", jtet);
|
145 |
schiefer |
1.1 |
t->SetBranchAddress("jtpt", jtpt);
|
146 |
|
|
t->SetBranchAddress("jteta", jteta);
|
147 |
|
|
t->SetBranchAddress("jtphi", jtphi);
|
148 |
schiefer |
1.4 |
t->SetBranchAddress("jtemf", jtemf);
|
149 |
schiefer |
1.3 |
t->SetBranchAddress("jtgenet",jtgenet);
|
150 |
schiefer |
1.1 |
t->SetBranchAddress("jtgenpt",jtgenpt);
|
151 |
|
|
t->SetBranchAddress("jtgendr",jtgendr);
|
152 |
|
|
|
153 |
|
|
if (t->FindBranch("jtjes")) t->SetBranchAddress("jtjes", jtjes);
|
154 |
|
|
else
|
155 |
|
|
for (int i1=0;i1<100;i1++)
|
156 |
|
|
for (int i2=0;i2<3;i2++) jtjes[i1][i2]=1.0;
|
157 |
|
|
|
158 |
schiefer |
1.3 |
string absrsp_xtitle = "E_{T}^{gen}-E_{T} [GeV]";
|
159 |
|
|
string relrsp_xtitle = "E_{T}^{rec}/E_{T}";
|
160 |
schiefer |
1.5 |
|
161 |
|
|
TH1F** hGenEt = hist::initHistos("GenEt",jetalgs[i],100,
|
162 |
|
|
"jet E_{T}^{gen} [GeV]","GenEt",etbins);
|
163 |
|
|
TH1F** hEtCorrGenEt = hist::initHistos("EtCorrGenEt",jetalgs[i],100,
|
164 |
|
|
"jet E_{T}^{gen} [GeV]","GenEt",etbins);
|
165 |
|
|
TH1F** hAbsRspGenEt = hist::initHistos("AbsRspGenEt",jetalgs[i],100,-50,150,
|
166 |
|
|
absrsp_xtitle,"GenEt",etbins);
|
167 |
|
|
TH1F** hRelRspGenEt = hist::initHistos("RelRspGenEt",jetalgs[i],100,0,2,
|
168 |
|
|
relrsp_xtitle,"GenEt",etbins);
|
169 |
|
|
|
170 |
|
|
TH1F** hEt = hist::initHistos("Et",jetalgs[i],100,
|
171 |
|
|
"jet E_{T} [GeV]","Et",etbins);
|
172 |
|
|
TH1F** hEtCorrEt = hist::initHistos("EtCorrEt",jetalgs[i],100,
|
173 |
|
|
"jet E_{T} [GeV]","Et",etbins);
|
174 |
|
|
TH1F** hAbsRspEt = hist::initHistos("AbsRspEt",jetalgs[i],100,-50,150,
|
175 |
|
|
absrsp_xtitle,"Et",etbins);
|
176 |
|
|
TH1F** hRelRspEt = hist::initHistos("RelRspEt",jetalgs[i],100,0,2,
|
177 |
|
|
relrsp_xtitle,"Et",etbins);
|
178 |
|
|
|
179 |
|
|
TH1F** hEta = hist::initHistos("Eta",jetalgs[i],100,
|
180 |
|
|
"jet #eta","Eta",etabins);
|
181 |
|
|
TH1F** hEtCorrEta = hist::initHistos("EtCorrEta",jetalgs[i],100,
|
182 |
|
|
"jet E_{T} [GeV]","Eta",etabins);
|
183 |
|
|
TH1F** hAbsRspEta = hist::initHistos("AbsRspEta",jetalgs[i],100,-50,150,
|
184 |
|
|
absrsp_xtitle,"Eta",etabins);
|
185 |
|
|
TH1F** hRelRspEta = hist::initHistos("RelRspEta",jetalgs[i],100,0,2,
|
186 |
|
|
relrsp_xtitle,"Eta",etabins);
|
187 |
|
|
|
188 |
|
|
TH1F** hPhi = hist::initHistos("Phi",jetalgs[i],100,
|
189 |
|
|
"jet #phi","Phi",phibins);
|
190 |
|
|
TH1F** hEtCorrPhi = hist::initHistos("EtCorrPhi",jetalgs[i],100,
|
191 |
|
|
"jet E_{T} [Gev]","Phi",phibins);
|
192 |
|
|
TH1F** hAbsRspPhi = hist::initHistos("AbsRspPhi",jetalgs[i],100,-50,150,
|
193 |
|
|
absrsp_xtitle,"Phi",phibins);
|
194 |
|
|
TH1F** hRelRspPhi = hist::initHistos("RelRspPhi",jetalgs[i],100,0,2,
|
195 |
|
|
relrsp_xtitle,"Phi",phibins);
|
196 |
|
|
|
197 |
|
|
TH1F** hEmf = hist::initHistos("Emf",jetalgs[i],100,
|
198 |
|
|
"jet emf","Emf",emfbins);
|
199 |
|
|
TH1F** hEtCorrEmf = hist::initHistos("EtCorrEmf",jetalgs[i],100,
|
200 |
|
|
"jet E_{T} [Gev]","Emf",emfbins);
|
201 |
|
|
TH1F** hAbsRspEmf = hist::initHistos("AbsRspEmf",jetalgs[i],100,-50,150,
|
202 |
|
|
absrsp_xtitle,"Emf",emfbins);
|
203 |
|
|
TH1F** hRelRspEmf = hist::initHistos("RelRspEmf",jetalgs[i],100,0,2,
|
204 |
|
|
relrsp_xtitle,"Emf",emfbins);
|
205 |
|
|
|
206 |
|
|
TH1F*** hEtEtaEt = hist::initHistos("Et",jetalgs[i],100,"jet E_{T} [GeV]",
|
207 |
|
|
"Eta",etabins,"Et",etbins);
|
208 |
|
|
TH1F*** hEtCorrEtaEt = hist::initHistos("EtCorrEtaEt",jetalgs[i],100,
|
209 |
schiefer |
1.4 |
"jet E_{T} [GeV]",
|
210 |
schiefer |
1.5 |
"Eta",etabins,"Et",etbins);
|
211 |
|
|
TH1F*** hAbsRspEtaEt = hist::initHistos("AbsRspEtaEt",jetalgs[i],100,-50,150,
|
212 |
schiefer |
1.4 |
absrsp_xtitle,
|
213 |
schiefer |
1.5 |
"Eta",etabins,"Et",etbins);
|
214 |
|
|
TH1F*** hRelRspEtaEt = hist::initHistos("RelRspEtaEt",jetalgs[i],100,0,2,
|
215 |
schiefer |
1.4 |
relrsp_xtitle,
|
216 |
schiefer |
1.5 |
"Eta",etabins,"Et",etbins);
|
217 |
schiefer |
1.1 |
|
218 |
schiefer |
1.5 |
TH1F*** hEtEmfEt = hist::initHistos("Et",jetalgs[i],100,"jet E_{T} [GeV]",
|
219 |
|
|
"Emf",emfbins,"Et",etbins);
|
220 |
|
|
TH1F*** hEtCorrEmfEt = hist::initHistos("EtCorrEmfEt",jetalgs[i],100,
|
221 |
schiefer |
1.4 |
"jet E_{T} [GeV]",
|
222 |
schiefer |
1.5 |
"Emf",emfbins,"Et",etbins);
|
223 |
|
|
TH1F*** hAbsRspEmfEt = hist::initHistos("AbsRspEmfEt",jetalgs[i],100,-50,150,
|
224 |
schiefer |
1.4 |
absrsp_xtitle,
|
225 |
schiefer |
1.5 |
"Emf",emfbins,"Et",etbins);
|
226 |
|
|
TH1F*** hRelRspEmfEt = hist::initHistos("RelRspEmfEt",jetalgs[i],100,0,2,
|
227 |
schiefer |
1.4 |
relrsp_xtitle,
|
228 |
schiefer |
1.5 |
"Emf",emfbins,"Et",etbins);
|
229 |
schiefer |
1.4 |
|
230 |
schiefer |
1.5 |
TH1F**** hEtEmfEtaEt = hist::initHistos("Et",jetalgs[i],100,"jet E_{T} [GeV]",
|
231 |
|
|
"Emf",emfbins,
|
232 |
|
|
"Eta",etabins,
|
233 |
|
|
"Et", etbins);
|
234 |
|
|
TH1F**** hEtCorrEmfEtaEt = hist::initHistos("EtCorrEmfEtaEt",jetalgs[i],100,
|
235 |
|
|
"jet E_{T} [GeV]",
|
236 |
|
|
"Emf",emfbins,
|
237 |
|
|
"Eta",etabins,
|
238 |
|
|
"Et", etbins);
|
239 |
|
|
TH1F**** hAbsRspEmfEtaEt = hist::initHistos("AbsRspEmfEtaEt",jetalgs[i],
|
240 |
|
|
100,-50,150,absrsp_xtitle,
|
241 |
|
|
"Emf",emfbins,
|
242 |
|
|
"Eta",etabins,
|
243 |
|
|
"Et", etbins);
|
244 |
|
|
TH1F**** hRelRspEmfEtaEt = hist::initHistos("RelRspEmfEtaEt",jetalgs[i],
|
245 |
|
|
100,0,2,relrsp_xtitle,
|
246 |
|
|
"Emf",emfbins,
|
247 |
|
|
"Eta",etabins,
|
248 |
|
|
"Et", etbins);
|
249 |
schiefer |
1.4 |
|
250 |
schiefer |
1.1 |
|
251 |
|
|
int nevts = el->GetN();
|
252 |
|
|
cout<<jetalgs[i]<<": "<<nevts<<" events selected."<<endl;
|
253 |
schiefer |
1.4 |
|
254 |
schiefer |
1.3 |
// prepare intermediate trees
|
255 |
schiefer |
1.4 |
float genet,et,eta,phi,emf,etcorr,absrsp,relrsp;
|
256 |
|
|
TTree** tGenEt = new TTree*[netbins];
|
257 |
|
|
TTree** tEt = new TTree*[netbins];
|
258 |
|
|
TTree** tEta = new TTree*[netabins];
|
259 |
|
|
TTree** tPhi = new TTree*[nphibins];
|
260 |
|
|
TTree** tEmf = new TTree*[nemfbins];
|
261 |
|
|
TTree*** tEtEtaEt = new TTree**[netabins];
|
262 |
|
|
TTree*** tEtEmfEt = new TTree**[nemfbins];
|
263 |
|
|
TTree**** tEtEmfEtaEt = new TTree***[nemfbins];
|
264 |
|
|
|
265 |
schiefer |
1.3 |
for (int iet=0;iet<netbins;iet++) {
|
266 |
schiefer |
1.4 |
|
267 |
schiefer |
1.3 |
tGenEt[iet]=new TTree("tGenEt","tGenEt");
|
268 |
|
|
tGenEt[iet]->Branch("val", &genet, "val/F");
|
269 |
|
|
tGenEt[iet]->Branch("etcorr",&etcorr,"etcorr/F");
|
270 |
schiefer |
1.4 |
tGenEt[iet]->Branch("absrsp",&absrsp,"absrsp/F");
|
271 |
|
|
tGenEt[iet]->Branch("relrsp",&relrsp,"relrsp/F");
|
272 |
schiefer |
1.3 |
tGenEt[iet]->Branch("weight",&weight,"weight/F");
|
273 |
schiefer |
1.4 |
|
274 |
schiefer |
1.3 |
tEt[iet]=new TTree("tEt","tEt");
|
275 |
|
|
tEt[iet]->Branch("val", &et, "val/F");
|
276 |
|
|
tEt[iet]->Branch("etcorr",&etcorr,"etcorr/F");
|
277 |
schiefer |
1.4 |
tEt[iet]->Branch("absrsp",&absrsp,"absrsp/F");
|
278 |
|
|
tEt[iet]->Branch("relrsp",&relrsp,"relrsp/F");
|
279 |
schiefer |
1.3 |
tEt[iet]->Branch("weight",&weight,"weight/F");
|
280 |
|
|
}
|
281 |
|
|
for (int ieta=0;ieta<netabins;ieta++) {
|
282 |
|
|
tEta[ieta]=new TTree("tEta","tEta");
|
283 |
|
|
tEta[ieta]->Branch("val", &eta, "val/F");
|
284 |
|
|
tEta[ieta]->Branch("etcorr",&etcorr,"etcorr/F");
|
285 |
schiefer |
1.4 |
tEta[ieta]->Branch("absrsp",&absrsp,"absrsp/F");
|
286 |
|
|
tEta[ieta]->Branch("relrsp",&relrsp,"relrsp/F");
|
287 |
schiefer |
1.3 |
tEta[ieta]->Branch("weight",&weight,"weight/F");
|
288 |
|
|
}
|
289 |
|
|
for (int iphi=0;iphi<nphibins;iphi++) {
|
290 |
|
|
tPhi[iphi]=new TTree("tPhi","tPhi");
|
291 |
|
|
tPhi[iphi]->Branch("val", &phi, "val/F");
|
292 |
|
|
tPhi[iphi]->Branch("etcorr",&etcorr,"etcorr/F");
|
293 |
schiefer |
1.4 |
tPhi[iphi]->Branch("absrsp",&absrsp,"absrsp/F");
|
294 |
|
|
tPhi[iphi]->Branch("relrsp",&relrsp,"relrsp/F");
|
295 |
schiefer |
1.3 |
tPhi[iphi]->Branch("weight",&weight,"weight/F");
|
296 |
|
|
}
|
297 |
schiefer |
1.4 |
for (int iemf=0;iemf<nemfbins;iemf++) {
|
298 |
|
|
tEmf[iemf]=new TTree("tEmf","tEmf");
|
299 |
|
|
tEmf[iemf]->Branch("val", &emf, "val/F");
|
300 |
|
|
tEmf[iemf]->Branch("etcorr",&etcorr,"etcorr/F");
|
301 |
|
|
tEmf[iemf]->Branch("absrsp",&absrsp,"absrsp/F");
|
302 |
|
|
tEmf[iemf]->Branch("relrsp",&relrsp,"relrsp/F");
|
303 |
|
|
tEmf[iemf]->Branch("weight",&weight,"weight/F");
|
304 |
|
|
}
|
305 |
schiefer |
1.3 |
for (int ieta=0;ieta<netabins;ieta++) {
|
306 |
|
|
tEtEtaEt[ieta] = new TTree*[netbins];
|
307 |
|
|
for (int iet=0;iet<netbins;iet++) {
|
308 |
|
|
tEtEtaEt[ieta][iet] = new TTree("tEtEtaEt","tEtEtaEt");
|
309 |
schiefer |
1.4 |
tEtEtaEt[ieta][iet]->Branch("val", &et, "val/F");
|
310 |
schiefer |
1.3 |
tEtEtaEt[ieta][iet]->Branch("etcorr",&etcorr,"etcorr/F");
|
311 |
schiefer |
1.4 |
tEtEtaEt[ieta][iet]->Branch("absrsp",&absrsp,"absrsp/F");
|
312 |
|
|
tEtEtaEt[ieta][iet]->Branch("relrsp",&relrsp,"relrsp/F");
|
313 |
schiefer |
1.3 |
tEtEtaEt[ieta][iet]->Branch("weight",&weight,"weight/F");
|
314 |
|
|
}
|
315 |
|
|
}
|
316 |
schiefer |
1.4 |
for (int iemf=0;iemf<nemfbins;iemf++) {
|
317 |
|
|
tEtEmfEt[iemf] = new TTree*[netbins];
|
318 |
|
|
for (int iet=0;iet<netbins;iet++) {
|
319 |
|
|
tEtEmfEt[iemf][iet] = new TTree("tEtEmfEt","tEtEmfEt");
|
320 |
|
|
tEtEmfEt[iemf][iet]->Branch("val", &et, "val/F");
|
321 |
|
|
tEtEmfEt[iemf][iet]->Branch("etcorr",&etcorr,"etcorr/F");
|
322 |
|
|
tEtEmfEt[iemf][iet]->Branch("absrsp",&absrsp,"absrsp/F");
|
323 |
|
|
tEtEmfEt[iemf][iet]->Branch("relrsp",&relrsp,"relrsp/F");
|
324 |
|
|
tEtEmfEt[iemf][iet]->Branch("weight",&weight,"weight/F");
|
325 |
|
|
}
|
326 |
|
|
}
|
327 |
|
|
for (int iemf=0;iemf<nemfbins;iemf++) {
|
328 |
|
|
tEtEmfEtaEt[iemf] = new TTree**[netabins];
|
329 |
|
|
for (int ieta=0;ieta<netabins;ieta++) {
|
330 |
|
|
tEtEmfEtaEt[iemf][ieta] = new TTree*[netbins];
|
331 |
|
|
for (int iet=0;iet<netbins;iet++) {
|
332 |
|
|
tEtEmfEtaEt[iemf][ieta][iet] = new TTree("tEtEmfEtaEt","tEtEmfEtaEt");
|
333 |
|
|
tEtEmfEtaEt[iemf][ieta][iet]->Branch("val", &et, "val/F");
|
334 |
|
|
tEtEmfEtaEt[iemf][ieta][iet]->Branch("etcorr",&etcorr,"etcorr/F");
|
335 |
|
|
tEtEmfEtaEt[iemf][ieta][iet]->Branch("absrsp",&absrsp,"absrsp/F");
|
336 |
|
|
tEtEmfEtaEt[iemf][ieta][iet]->Branch("relrsp",&relrsp,"relrsp/F");
|
337 |
|
|
tEtEmfEtaEt[iemf][ieta][iet]->Branch("weight",&weight,"weight/F");
|
338 |
|
|
}
|
339 |
|
|
}
|
340 |
|
|
}
|
341 |
schiefer |
1.1 |
|
342 |
schiefer |
1.5 |
cout<<jetalgs[i]<<": trees created."<<endl;
|
343 |
|
|
|
344 |
schiefer |
1.1 |
|
345 |
|
|
// loop over all events
|
346 |
|
|
for (int ievt=0;ievt<nevts;ievt++) {
|
347 |
|
|
|
348 |
|
|
t->GetEntry(el->GetEntry(ievt));
|
349 |
|
|
|
350 |
|
|
// loop over all jets in the event
|
351 |
|
|
for (int ijt=0;ijt<njt;ijt++) {
|
352 |
|
|
|
353 |
|
|
if (jtgendr[ijt]<0.0||jtgendr[ijt]>drmax) continue;
|
354 |
|
|
|
355 |
|
|
if (ijt>0) continue;
|
356 |
|
|
|
357 |
schiefer |
1.3 |
et = jtet[ijt];
|
358 |
|
|
genet = jtgenet[ijt];
|
359 |
|
|
eta = (abseta) ? std::abs(jteta[ijt]) : jteta[ijt];
|
360 |
|
|
phi = jtphi[ijt];
|
361 |
schiefer |
1.4 |
emf = jtemf[ijt];
|
362 |
schiefer |
1.3 |
etcorr = (applyjes==0) ? et : et*jtjes[ijt][applyjes-1];
|
363 |
schiefer |
1.1 |
|
364 |
schiefer |
1.4 |
absrsp = genet - etcorr;
|
365 |
|
|
relrsp = etcorr/genet;
|
366 |
schiefer |
1.2 |
|
367 |
schiefer |
1.3 |
// genet
|
368 |
|
|
int igenet = hist::get_ibin(genet,etbins);
|
369 |
|
|
tGenEt[igenet]->Fill();
|
370 |
schiefer |
1.1 |
|
371 |
schiefer |
1.3 |
// et
|
372 |
|
|
int iet = hist::get_ibin(et,etbins);
|
373 |
|
|
tEt[iet]->Fill();
|
374 |
schiefer |
1.1 |
|
375 |
|
|
// eta
|
376 |
schiefer |
1.2 |
int ieta = hist::get_ibin(eta,etabins);
|
377 |
schiefer |
1.3 |
tEta[ieta]->Fill();
|
378 |
|
|
|
379 |
schiefer |
1.1 |
// phi
|
380 |
schiefer |
1.2 |
int iphi = hist::get_ibin(phi,phibins);
|
381 |
schiefer |
1.3 |
tPhi[iphi]->Fill();
|
382 |
|
|
|
383 |
schiefer |
1.4 |
// emf
|
384 |
|
|
int iemf = hist::get_ibin(emf,emfbins);
|
385 |
|
|
tEmf[iemf]->Fill();
|
386 |
|
|
|
387 |
schiefer |
1.3 |
// et-eta
|
388 |
|
|
tEtEtaEt[ieta][iet]->Fill();
|
389 |
schiefer |
1.1 |
|
390 |
schiefer |
1.4 |
// et-emf
|
391 |
|
|
tEtEmfEt[iemf][iet]->Fill();
|
392 |
|
|
|
393 |
|
|
// et-eta-emf
|
394 |
|
|
tEtEmfEtaEt[iemf][ieta][iet]->Fill();
|
395 |
|
|
|
396 |
schiefer |
1.1 |
} // jets
|
397 |
|
|
} // evts
|
398 |
|
|
|
399 |
|
|
cout<<jetalgs[i]<<": all events processed / histograms filled."<<endl;
|
400 |
schiefer |
1.3 |
|
401 |
schiefer |
1.2 |
|
402 |
schiefer |
1.3 |
// replace histograms
|
403 |
schiefer |
1.4 |
|
404 |
|
|
replaceHistos(netbins,
|
405 |
schiefer |
1.5 |
hGenEt,hEtCorrGenEt,hAbsRspGenEt,hRelRspGenEt,tGenEt);
|
406 |
|
|
replaceHistos(netbins,hEt,hEtCorrEt,hAbsRspEt,hRelRspEt,tEt);
|
407 |
|
|
replaceHistos(netabins,hEta,hEtCorrEta,hAbsRspEta,hRelRspEta,tEta);
|
408 |
|
|
replaceHistos(nphibins,hPhi,hEtCorrPhi,hAbsRspPhi,hRelRspPhi,tPhi);
|
409 |
|
|
replaceHistos(nemfbins,hEmf,hEtCorrEmf,hAbsRspEmf,hRelRspEmf,tEmf);
|
410 |
schiefer |
1.4 |
replaceHistos(netabins,netbins,
|
411 |
schiefer |
1.5 |
hEtEtaEt,hEtCorrEtaEt,hAbsRspEtaEt,hRelRspEtaEt,
|
412 |
schiefer |
1.4 |
tEtEtaEt);
|
413 |
|
|
replaceHistos(nemfbins,netbins,
|
414 |
schiefer |
1.5 |
hEtEmfEt,hEtCorrEmfEt,hAbsRspEmfEt,hRelRspEmfEt,
|
415 |
schiefer |
1.4 |
tEtEmfEt);
|
416 |
|
|
replaceHistos(nemfbins,netabins,netbins,
|
417 |
schiefer |
1.5 |
hEtEmfEtaEt,hEtCorrEmfEtaEt,
|
418 |
|
|
hAbsRspEmfEtaEt,hRelRspEmfEtaEt,
|
419 |
schiefer |
1.4 |
tEtEmfEtaEt);
|
420 |
schiefer |
1.2 |
|
421 |
schiefer |
1.5 |
cout<<jetalgs[i]<<": all histograms replaced, and tree deleted."<<endl;
|
422 |
|
|
|
423 |
schiefer |
1.2 |
TDirectory* d = plotfile->mkdir(jetalgs[i].c_str());
|
424 |
|
|
d->cd();
|
425 |
schiefer |
1.1 |
|
426 |
schiefer |
1.3 |
// fit response histos
|
427 |
schiefer |
1.5 |
fitHistos(hAbsRspGenEt,netbins,nsigma);
|
428 |
|
|
fitHistos(hRelRspGenEt,netbins,nsigma);
|
429 |
schiefer |
1.3 |
|
430 |
schiefer |
1.5 |
fitHistos(hAbsRspEt, netbins,nsigma);
|
431 |
|
|
fitHistos(hRelRspEt, netbins,nsigma);
|
432 |
schiefer |
1.3 |
|
433 |
schiefer |
1.5 |
fitHistos(hAbsRspEta, netabins,nsigma);
|
434 |
|
|
fitHistos(hRelRspEta, netabins,nsigma);
|
435 |
schiefer |
1.3 |
|
436 |
schiefer |
1.5 |
fitHistos(hAbsRspPhi, nphibins,nsigma);
|
437 |
|
|
fitHistos(hRelRspPhi, nphibins,nsigma);
|
438 |
schiefer |
1.3 |
|
439 |
schiefer |
1.5 |
fitHistos(hAbsRspEmf, nemfbins,nsigma);
|
440 |
|
|
fitHistos(hRelRspEmf, nemfbins,nsigma);
|
441 |
schiefer |
1.4 |
|
442 |
schiefer |
1.5 |
fitHistos(hAbsRspEtaEt,netabins,netbins,nsigma);
|
443 |
|
|
fitHistos(hRelRspEtaEt,netabins,netbins,nsigma);
|
444 |
schiefer |
1.3 |
|
445 |
schiefer |
1.5 |
fitHistos(hAbsRspEmfEt,nemfbins,netbins,nsigma);
|
446 |
|
|
fitHistos(hRelRspEmfEt,nemfbins,netbins,nsigma);
|
447 |
schiefer |
1.4 |
|
448 |
schiefer |
1.5 |
fitHistos(hAbsRspEmfEtaEt,nemfbins,netabins,netbins,nsigma);
|
449 |
|
|
fitHistos(hRelRspEmfEtaEt,nemfbins,netabins,netbins,nsigma);
|
450 |
|
|
|
451 |
|
|
cout<<jetalgs[i]<<": fits done."<<endl;
|
452 |
schiefer |
1.4 |
|
453 |
schiefer |
1.3 |
|
454 |
|
|
// save histograms
|
455 |
schiefer |
1.5 |
saveHistos(hGenEt, netbins);
|
456 |
|
|
saveHistos(hEtCorrGenEt,netbins);
|
457 |
|
|
saveHistos(hAbsRspGenEt,netbins);
|
458 |
|
|
saveHistos(hRelRspGenEt,netbins);
|
459 |
|
|
|
460 |
|
|
saveHistos(hEt, netbins);
|
461 |
|
|
saveHistos(hEtCorrEt, netbins);
|
462 |
|
|
saveHistos(hAbsRspEt, netbins);
|
463 |
|
|
saveHistos(hRelRspEt, netbins);
|
464 |
|
|
|
465 |
|
|
saveHistos(hEta, netabins);
|
466 |
|
|
saveHistos(hEtCorrEta, netabins);
|
467 |
|
|
saveHistos(hAbsRspEta, netabins);
|
468 |
|
|
saveHistos(hRelRspEta, netabins);
|
469 |
|
|
|
470 |
|
|
saveHistos(hPhi, nphibins);
|
471 |
|
|
saveHistos(hEtCorrPhi, nphibins);
|
472 |
|
|
saveHistos(hAbsRspPhi, nphibins);
|
473 |
|
|
saveHistos(hRelRspPhi, nphibins);
|
474 |
|
|
|
475 |
|
|
saveHistos(hEmf, nemfbins);
|
476 |
|
|
saveHistos(hEtCorrEmf, nemfbins);
|
477 |
|
|
saveHistos(hAbsRspEmf, nemfbins);
|
478 |
|
|
saveHistos(hRelRspEmf, nemfbins);
|
479 |
|
|
|
480 |
|
|
saveHistos(hEtEtaEt, netabins,netbins);
|
481 |
|
|
saveHistos(hEtCorrEtaEt,netabins,netbins);
|
482 |
|
|
saveHistos(hAbsRspEtaEt,netabins,netbins);
|
483 |
|
|
saveHistos(hRelRspEtaEt,netabins,netbins);
|
484 |
|
|
|
485 |
|
|
saveHistos(hEtEmfEt, nemfbins,netbins);
|
486 |
|
|
saveHistos(hEtCorrEmfEt,nemfbins,netbins);
|
487 |
|
|
saveHistos(hAbsRspEmfEt,nemfbins,netbins);
|
488 |
|
|
saveHistos(hRelRspEmfEt,nemfbins,netbins);
|
489 |
|
|
|
490 |
|
|
saveHistos(hEtEmfEtaEt, nemfbins,netabins,netbins);
|
491 |
|
|
saveHistos(hEtCorrEmfEtaEt,nemfbins,netabins,netbins);
|
492 |
|
|
saveHistos(hAbsRspEmfEtaEt,nemfbins,netabins,netbins);
|
493 |
|
|
saveHistos(hRelRspEmfEtaEt,nemfbins,netabins,netbins);
|
494 |
|
|
|
495 |
|
|
plotfile->Write();
|
496 |
|
|
|
497 |
|
|
cout<<jetalgs[i]<<": histograms saved."<<endl;
|
498 |
schiefer |
1.1 |
}
|
499 |
|
|
|
500 |
schiefer |
1.2 |
plotfile->Write();
|
501 |
|
|
plotfile->Close();
|
502 |
|
|
delete plotfile;
|
503 |
schiefer |
1.1 |
|
504 |
|
|
return 0;
|
505 |
|
|
}
|
506 |
|
|
|
507 |
|
|
|
508 |
|
|
|
509 |
|
|
////////////////////////////////////////////////////////////////////////////////
|
510 |
|
|
// implementation of global functions
|
511 |
|
|
////////////////////////////////////////////////////////////////////////////////
|
512 |
|
|
|
513 |
|
|
//______________________________________________________________________________
|
514 |
schiefer |
1.4 |
void replaceHistos(int nbins,
|
515 |
|
|
TH1F**& hVar,TH1F**& hEtCorr,TH1F**& hAbsRsp,TH1F**& hRelRsp,
|
516 |
|
|
TTree**& tVar)
|
517 |
schiefer |
1.1 |
{
|
518 |
schiefer |
1.3 |
TH1::SetDefaultSumw2();
|
519 |
schiefer |
1.1 |
for (int ibin=0;ibin<nbins;ibin++) {
|
520 |
schiefer |
1.4 |
|
521 |
schiefer |
1.3 |
TH1F* h(0);
|
522 |
schiefer |
1.4 |
|
523 |
|
|
if (tVar[ibin]->GetEntries()==0) continue;
|
524 |
schiefer |
1.1 |
|
525 |
schiefer |
1.3 |
tVar[ibin]->Project("h(50)","val","weight*(1)");
|
526 |
|
|
h = (TH1F*)gROOT->FindObject("h");
|
527 |
|
|
h->SetNameTitle(hVar[ibin]->GetName(),hVar[ibin]->GetTitle());
|
528 |
|
|
h->SetXTitle(hVar[ibin]->GetXaxis()->GetTitle());
|
529 |
|
|
h->SetYTitle(hVar[ibin]->GetYaxis()->GetTitle());
|
530 |
|
|
delete hVar[ibin];
|
531 |
|
|
hVar[ibin] = h;
|
532 |
|
|
|
533 |
|
|
tVar[ibin]->Project("h(50)","etcorr","weight*(1)");
|
534 |
|
|
h = (TH1F*)gROOT->FindObject("h");
|
535 |
|
|
h->SetNameTitle(hEtCorr[ibin]->GetName(),hEtCorr[ibin]->GetTitle());
|
536 |
|
|
h->SetXTitle(hEtCorr[ibin]->GetXaxis()->GetTitle());
|
537 |
|
|
h->SetYTitle(hEtCorr[ibin]->GetYaxis()->GetTitle());
|
538 |
|
|
delete hEtCorr[ibin];
|
539 |
|
|
hEtCorr[ibin] = h;
|
540 |
schiefer |
1.4 |
|
541 |
|
|
tVar[ibin]->Project("h(50)","absrsp","weight*(1)");
|
542 |
|
|
h = (TH1F*)gROOT->FindObject("h");
|
543 |
|
|
h->SetNameTitle(hAbsRsp[ibin]->GetName(),hAbsRsp[ibin]->GetTitle());
|
544 |
|
|
h->SetXTitle(hAbsRsp[ibin]->GetXaxis()->GetTitle());
|
545 |
|
|
h->SetYTitle(hAbsRsp[ibin]->GetYaxis()->GetTitle());
|
546 |
|
|
delete hAbsRsp[ibin];
|
547 |
|
|
hAbsRsp[ibin] = h;
|
548 |
|
|
|
549 |
|
|
tVar[ibin]->Project("h(50)","relrsp","weight*(1)");
|
550 |
|
|
h = (TH1F*)gROOT->FindObject("h");
|
551 |
|
|
h->SetNameTitle(hRelRsp[ibin]->GetName(),hRelRsp[ibin]->GetTitle());
|
552 |
|
|
h->SetXTitle(hRelRsp[ibin]->GetXaxis()->GetTitle());
|
553 |
|
|
h->SetYTitle(hRelRsp[ibin]->GetYaxis()->GetTitle());
|
554 |
|
|
delete hRelRsp[ibin];
|
555 |
|
|
hRelRsp[ibin] = h;
|
556 |
|
|
|
557 |
schiefer |
1.3 |
delete tVar[ibin];
|
558 |
schiefer |
1.1 |
}
|
559 |
|
|
|
560 |
schiefer |
1.3 |
delete [] tVar;
|
561 |
schiefer |
1.1 |
}
|
562 |
|
|
|
563 |
schiefer |
1.4 |
|
564 |
|
|
//______________________________________________________________________________
|
565 |
|
|
void replaceHistos(int nbins1,int nbins2,
|
566 |
|
|
TH1F***& hVar,TH1F***& hEtCorr,TH1F***& hAbsRsp,TH1F***& hRelRsp,
|
567 |
|
|
TTree***& tVar)
|
568 |
|
|
{
|
569 |
|
|
for (int i1=0;i1<nbins1;i1++)
|
570 |
|
|
replaceHistos(nbins2,hVar[i1],hEtCorr[i1],hAbsRsp[i1],hRelRsp[i1],tVar[i1]);
|
571 |
schiefer |
1.5 |
delete [] tVar;
|
572 |
schiefer |
1.4 |
}
|
573 |
|
|
|
574 |
|
|
|
575 |
|
|
//______________________________________________________________________________
|
576 |
|
|
void replaceHistos(int nbins1,int nbins2,int nbins3,
|
577 |
|
|
TH1F****& hVar,TH1F****& hEtCorr,TH1F****&hAbsRsp,TH1F****&hRelRsp,
|
578 |
|
|
TTree****& tVar)
|
579 |
|
|
{
|
580 |
|
|
for (int i1=0;i1<nbins1;i1++)
|
581 |
|
|
replaceHistos(nbins2,nbins3,
|
582 |
|
|
hVar[i1],hEtCorr[i1],hAbsRsp[i1],hRelRsp[i1],tVar[i1]);
|
583 |
schiefer |
1.5 |
delete [] tVar;
|
584 |
schiefer |
1.4 |
}
|
585 |
|
|
|
586 |
|
|
|
587 |
schiefer |
1.1 |
//______________________________________________________________________________
|
588 |
schiefer |
1.3 |
void fitHistos(TH1F** histos,int n,double nsigma)
|
589 |
schiefer |
1.1 |
{
|
590 |
schiefer |
1.3 |
for (int i=0;i<n;i++) {
|
591 |
|
|
float nrm = histos[i]->Integral();
|
592 |
|
|
float mean = histos[i]->GetMean();
|
593 |
|
|
float rms = histos[i]->GetRMS();
|
594 |
|
|
TF1* f = new TF1("gaus","gaus",mean-nsigma*rms,mean+nsigma*rms);
|
595 |
|
|
f->SetParameter(0,nrm);
|
596 |
|
|
f->SetParameter(1,mean);
|
597 |
|
|
f->SetParameter(2,rms);
|
598 |
|
|
f->SetLineWidth(1);
|
599 |
|
|
f->SetLineColor(kRed);
|
600 |
|
|
histos[i]->Fit(f,"QR");
|
601 |
schiefer |
1.1 |
}
|
602 |
|
|
}
|
603 |
|
|
|
604 |
|
|
|
605 |
|
|
//______________________________________________________________________________
|
606 |
schiefer |
1.3 |
void fitHistos(TH1F*** histos,int n1,int n2,double nsigma)
|
607 |
schiefer |
1.1 |
{
|
608 |
schiefer |
1.3 |
for (int i1=0;i1<n1;i1++) fitHistos(histos[i1],n2,nsigma);
|
609 |
schiefer |
1.1 |
}
|
610 |
|
|
|
611 |
|
|
|
612 |
|
|
//______________________________________________________________________________
|
613 |
schiefer |
1.4 |
void fitHistos(TH1F**** histos,int n1,int n2,int n3,double nsigma)
|
614 |
|
|
{
|
615 |
|
|
for (int i1=0;i1<n1;i1++) fitHistos(histos[i1],n2,n3,nsigma);
|
616 |
|
|
}
|
617 |
|
|
|
618 |
|
|
|
619 |
|
|
//______________________________________________________________________________
|
620 |
schiefer |
1.5 |
void saveHistos(TH1F**& histos,int n)
|
621 |
schiefer |
1.1 |
{
|
622 |
schiefer |
1.5 |
for (int i=0;i<n;i++) {
|
623 |
|
|
histos[i]->Write();
|
624 |
|
|
delete histos[i];
|
625 |
|
|
}
|
626 |
|
|
delete [] histos;
|
627 |
schiefer |
1.1 |
}
|
628 |
|
|
|
629 |
|
|
|
630 |
|
|
//______________________________________________________________________________
|
631 |
schiefer |
1.5 |
void saveHistos(TH1F***& histos,int n1,int n2)
|
632 |
schiefer |
1.1 |
{
|
633 |
schiefer |
1.3 |
for (int i1=0;i1<n1;i1++) saveHistos(histos[i1],n2);
|
634 |
schiefer |
1.5 |
delete [] histos;
|
635 |
schiefer |
1.1 |
}
|
636 |
schiefer |
1.4 |
|
637 |
|
|
|
638 |
|
|
//______________________________________________________________________________
|
639 |
schiefer |
1.5 |
void saveHistos(TH1F****& histos,int n1,int n2,int n3)
|
640 |
schiefer |
1.4 |
{
|
641 |
|
|
for (int i1=0;i1<n1;i1++) saveHistos(histos[i1],n2,n3);
|
642 |
schiefer |
1.5 |
delete [] histos;
|
643 |
schiefer |
1.4 |
}
|