ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/AnalysisFramework/Plotting/Modules/SampleClass.C
Revision: 1.17
Committed: Thu Sep 8 15:20:56 2011 UTC (13 years, 8 months ago) by buchmann
Content type: text/plain
Branch: MAIN
Changes since 1.16: +9 -9 lines
Log Message:
Fixed over/underflow bug which lead to the first/last bin becoming the under/overflow bin instead of actually containing it as well

File Contents

# Content
1 #include <iostream>
2 #include <sstream>
3 #include <vector>
4 #include <utility>
5
6 #include <TFile.h>
7 #include <TTree.h>
8 #include <TH1.h>
9 #include <TH2.h>
10 #include <TCut.h>
11 #include <THStack.h>
12 #include <TColor.h>
13 #include <TCanvas.h>
14 #include <TError.h>
15 #include <TText.h>
16 #include <TLegend.h>
17 #include <TError.h>
18 #include <TTreeFormula.h>
19
20 #define SampleClassLoaded
21
22 #ifndef Verbosity
23 #define Verbosity 0
24 #endif
25 #ifndef HUSH
26 #define HUSH 1
27 #endif
28 #ifndef GeneralToolBoxLoaded
29 #include "GeneralToolBox.C"
30 #endif
31 using namespace std;
32
33 bool doesROOTFileExist(string filename);
34 bool addoverunderflowbins=false;
35 string treename="PFevents";
36
37 class sample
38 {
39 public:
40 string filename;
41 string samplename;
42 int Nentries;
43 float xs;//cross section
44 float weight;
45 bool is_data;
46 bool is_signal;
47 bool is_active;
48 int groupindex;
49 Color_t samplecolor;
50
51 TFile *tfile;
52 TTree *events;
53
54 sample(string m_filename,string m_samplename,int m_Nentries, float m_xs,bool is_data, bool m_is_signal, int m_groupindex, Color_t color);
55 void closeFile();
56 };
57
58 string write_mc_or_data(bool is_data)
59 {
60 if(is_data) return "data";
61 return "MC";
62 }
63
64 sample::sample(string m_filename, string m_samplename, int m_Nentries, float m_xs,bool m_is_data, bool m_is_signal, int m_groupindex, Color_t mycolor)
65 {
66 this->filename=m_filename;
67 this->samplename=m_samplename;
68 this->Nentries=m_Nentries;
69 this->xs=m_xs;
70 this->is_data=m_is_data;
71 this->is_signal=m_is_signal;
72 this->weight=(xs/(float)Nentries);
73 this->groupindex=m_groupindex;
74 this->is_active=true;
75 this->samplecolor=mycolor;
76 if(!doesROOTFileExist(this->filename)) {
77 stringstream message;
78 message << "The " << write_mc_or_data(is_data) << " sample " << this->samplename << " is invalid because the associated file path, " << this->filename << " is incorrect.";
79 write_error("Sample::Sample",message.str());
80 this->is_active=false;
81 }
82
83 if(doesROOTFileExist(this->filename)) {
84 //suppressing stupid 64/32 errors here (Warning in <TFile::ReadStreamerInfo>: /scratch/buchmann/MC_Spring11_PU_PF/TToBLNu_TuneZ2_t-channel_7TeV-madgraph.root: not a TStreamerInfo object)
85 Int_t currlevel=gErrorIgnoreLevel;
86 gErrorIgnoreLevel=5000;
87 this->tfile = new TFile(m_filename.c_str());
88 gErrorIgnoreLevel=currlevel;
89 this->events=(TTree*)(this->tfile)->Get(treename.c_str());
90 if(Verbosity>0) dout << "The " << write_mc_or_data(is_data) << " file " << this->filename << " has been added successfully to the list of samples. " << endl;
91 }
92 else {
93 this->is_active=false;
94 }
95
96 }
97
98 void sample::closeFile()
99 {
100 if(doesROOTFileExist(this->filename)) {
101 (this->tfile)->Close();
102 }
103 else {
104 dout << "SAMPLE " << this->samplename << " cannot be closed as the underlying file (" << this->filename << ") does not exist!" << endl;
105 this->is_active=false;
106 }
107 }
108
109 bool doesROOTFileExist(string filename)
110 {
111 //suppressing stupid 64/32 errors here (Warning in <TFile::ReadStreamerInfo>: /scratch/buchmann/MC_Spring11_PU_PF/TToBLNu_TuneZ2_t-channel_7TeV-madgraph.root: not a TStreamerInfo object)
112 Int_t currlevel=gErrorIgnoreLevel;
113 gErrorIgnoreLevel=5000;
114 TFile *f = new TFile(filename.c_str());
115
116 if (f->IsZombie()) {
117 dout << "Error opening file" << filename << endl;
118 return 0;
119 }
120 f->Close();
121 gErrorIgnoreLevel=currlevel;
122 return 1;
123 }
124 //********************************************************
125
126 TCut essentialcut("mll>0");
127 // This will reweight all the events: use "weight" to correct MC for pileUP, "1.0" otherwise
128 TCut cutWeight("1.0");
129
130 void setessentialcut(TCut ess) {
131 essentialcut=ess;
132 }
133
134 class samplecollection
135 {
136 public:
137 string name;
138 vector<sample> collection;
139 int nsamples;
140 int ndatasamples;
141 int nmcsamples;
142 int ngroups;
143
144 samplecollection(string m_name);
145 void AddSample(string m_filename,string m_samplename,int m_Nentries, float m_xs,bool is_data, bool m_is_signal, int groupindex, Color_t color);
146
147 //vector quantities
148 TH1F* Draw(string m_histoname,string m_var, vector<float> binning, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, vector<int> onlyindex, bool drawsignal);
149 TH1F* Draw(string m_histoname,string m_var, vector<float> binning, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, bool drawsignal);
150
151 //minx,maxs quantities
152 TH1F* Draw(string m_histoname,string m_var, int m_nbins, float m_minx, float m_maxx, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, bool drawsignal);
153 TH1F* Draw(string m_histoname,string m_var, int m_nbins, float m_minx, float m_maxx, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, vector<int> onlyindex, bool drawsignal);
154
155 TH2F* Draw(string m_histoname,string m_var, vector<float> binningx, vector<float> binningy, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, vector<int> onlyindex, bool drawsignal);
156
157 vector<float> get_optimal_binsize(string variable, TCut cut,int nbins, float low, float high);
158
159 THStack DrawStack(string m_histoname,string m_var, int m_nbins, float m_minx, float m_maxx, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, bool drawsignal);
160 THStack DrawStack(string m_histoname,string m_var, vector<float> binning, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, bool drawsignal);
161 vector<int> FindSample(string what);
162 void ListSamples();
163 bool do_sample(int thissample, vector<int> &selected_samples);
164 Color_t GetColor(string filename);
165 Color_t GetColor(int sampleindex);
166
167 TLegend* allbglegend(string title,float x, float y);
168 TLegend* allbglegend(string title, TH1F *data, float x, float y);
169
170 void PickUpFromThisFile(int isamp, string cut, vector<string> &output, vector<string> &pickupfile);
171 void PickUpEvents(string cut);
172 string find_units(string&);
173 };
174
175 samplecollection::samplecollection(string m_name)
176 {
177 this->name=m_name;
178 if(Verbosity>0) dout << "Initiated sample collection " << this->name << " " << endl;
179 }
180
181 void samplecollection::ListSamples()
182 {
183 dout << "---------------------------------------------------------------------------------------------------" << endl;
184 dout << "Listing all " << this->nsamples << " sample(s) of the sample collection called " << this->name << " : " << endl;
185 if(this->ndatasamples>0) {
186 dout << "Data sample(s): " << endl;
187 for (int isamp=0;isamp<this->collection.size();isamp++)
188 {
189 if((this->collection)[isamp].is_data) dout << " - " << (this->collection)[isamp].samplename << " from " << (this->collection)[isamp].filename << endl;
190 }
191 }
192 if(this->nmcsamples>0) {
193 dout << "MC sample(s): " << endl;
194 for (int isamp=0;isamp<this->collection.size();isamp++)
195 {
196 if(!(this->collection)[isamp].is_data) dout << " - " << (this->collection)[isamp].samplename << " from " << (this->collection)[isamp].filename << endl;
197 }
198 }
199 dout << "---------------------------------------------------------------------------------------------------" << endl;
200 }
201
202 void samplecollection::AddSample(string m_filename,string m_samplename,int m_Nentries, float m_xs,bool is_data, bool m_is_signal, int groupindex, Color_t newcolor)
203 {
204 sample NewSample(m_filename,m_samplename,m_Nentries,m_xs,is_data,m_is_signal,groupindex,newcolor);
205 (this->collection).push_back(NewSample);
206 if((this->collection).size()==1) {
207 this->nmcsamples=0;
208 this->ndatasamples=0;
209 this->nsamples=0;
210 this->ngroups=0;
211 }
212 if(groupindex>this->ngroups) this->ngroups=groupindex;
213 if(is_data) this->ndatasamples++;
214 else this->nmcsamples++;
215 this->nsamples=this->collection.size();
216 }
217
218 bool samplecollection::do_sample(int thissample, vector<int> &selected_samples)
219 {
220 bool drawit=false;
221 for(int isel=0;isel<selected_samples.size();isel++)
222 {
223 if(selected_samples[isel]==thissample) {
224 drawit=true;
225 break;
226 }
227 }
228 return drawit;
229 }
230
231 TH1F* samplecollection::Draw(string m_histoname,string m_var, int m_nbins, float m_minx, float m_maxx, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, bool drawsignal=false)
232 {
233 vector<int> emptyvector;
234 TH1F *histo = this->Draw(m_histoname,m_var, m_nbins, m_minx, m_maxx, m_xlabel, m_ylabel, Cut, m_is_data, luminosity,emptyvector,drawsignal);
235 return histo;
236 }
237
238 TH1F* samplecollection::Draw(string m_histoname,string m_var, vector<float> binning, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, bool drawsignal=false)
239 {
240 vector<int> emptyvector;
241 TH1F *histo = this->Draw(m_histoname,m_var, binning, m_xlabel, m_ylabel, Cut, m_is_data, luminosity,emptyvector,drawsignal);
242 return histo;
243 }
244
245 TH1F* samplecollection::Draw(string m_histoname,string m_var, vector<float> binning, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, vector<int> onlyindex, bool drawsignal=false) {
246 if(Verbosity>0) dout << endl << endl;
247 if(Verbosity>0) dout << "-------------------------------------------------------------------------------------" << endl;
248 if(Verbosity>0) dout << "histoname : " << m_histoname << " , m_var = " << m_var << ", m_xlabel=" << m_xlabel << " m_ylabel=" << m_ylabel << " m_is_data = " << m_is_data << " lumi = " << luminosity << " onlyindex size: " << onlyindex.size() << endl;
249 if(HUSH==0) dout << "Drawing histo called " << m_histoname << "... " << endl;
250 bool do_only_selected_samples=false;
251 if(onlyindex.size()>0&&onlyindex[0]!=-1) {
252 if(Verbosity>0) {dout << "Requested to only draw sample corresponding to the following sample(s) : " << endl;}
253 for(int is=0;is<onlyindex.size();is++) {
254 if(Verbosity>0) dout << " - " << (this->collection)[onlyindex[is]].filename << " (sample index: " << onlyindex[is] << ")" << endl;
255 }
256 do_only_selected_samples=true;
257 }
258 if(onlyindex.size()==1&&onlyindex[0]==-1) {
259 do_only_selected_samples=true; // this is the special case when we request a non-existing sample - we shouldn't draw anything then.
260 onlyindex.clear();
261 }
262 stringstream h_histoname;
263 h_histoname<<"h_"<<m_histoname;
264 float binningarray[binning.size()+1];
265 for(int i=0;i<binning.size();i++) {binningarray[i]=binning[i];}
266 TH1F *histo = new TH1F(m_histoname.c_str(),"",binning.size()-1,binningarray);
267 histo->Sumw2();
268
269 stringstream drawthis;
270 drawthis<<m_var<<">>tempdrawhisto";
271
272 for (unsigned int isample=0;isample<(this->collection).size();isample++) {
273 bool use_this_sample=false;
274 if(!(this->collection)[isample].is_active) continue; // skip inactive samples right away
275 if(((this->collection)[isample].is_data==m_is_data)&&(!do_only_selected_samples)) {
276 if(drawsignal==false&&(this->collection)[isample].is_signal==true) continue;
277 use_this_sample=true;
278 }
279 if(do_only_selected_samples&&this->do_sample(isample,onlyindex)) {
280 use_this_sample=true;
281 }
282 TH1F *tempdrawhisto = new TH1F("tempdrawhisto","tempdrawhisto",binning.size()-1,binningarray);
283 tempdrawhisto->Sumw2();
284 if(use_this_sample) {
285 if(Verbosity>0) dout << "[samplecollection::Draw] : Added contribution from sample " << (this->collection)[isample].samplename << endl;
286 (this->collection)[isample].events->Draw(drawthis.str().c_str(),(essentialcut&&Cut)*cutWeight);//this weight is based on PU etc. not XS
287 if(addoverunderflowbins) {
288 //now also adding the overflow & underflow bins:
289 tempdrawhisto->SetBinContent(tempdrawhisto->GetNbinsX(),tempdrawhisto->GetBinContent(tempdrawhisto->GetNbinsX()+1)+tempdrawhisto->GetBinContent(tempdrawhisto->GetNbinsX()));
290 tempdrawhisto->SetBinError(tempdrawhisto->GetNbinsX(),TMath::Sqrt(tempdrawhisto->GetBinContent(tempdrawhisto->GetNbinsX())));
291 tempdrawhisto->SetBinContent(1,tempdrawhisto->GetBinContent(0)+tempdrawhisto->GetBinContent(1));
292 tempdrawhisto->SetBinError(1,TMath::Sqrt(tempdrawhisto->GetBinContent(1)));
293 }
294
295 if(!(this->collection)[isample].is_data) tempdrawhisto->Scale(luminosity*((this->collection)[isample].weight));//weight applied here is XS & N(entries)
296 histo->Add(tempdrawhisto);
297 }
298 tempdrawhisto->Delete();
299 }//end of loop over isample
300 if(Verbosity>0) dout << "Histo has been filled and now contains " << histo->Integral() << " points (integral)" << endl;
301 histo->GetXaxis()->SetTitle(m_xlabel.c_str());
302 // Try to add bin width information: look for units in m_xlabel
303 string units = find_units(m_xlabel);
304 if ( units.length()>0 ) {
305 stringstream ylabel;
306 ylabel << m_ylabel << " / " << histo->GetBinWidth(1) << " " << units;
307 histo->GetYaxis()->SetTitle( ylabel.str().c_str() );
308 } else {
309 histo->GetYaxis()->SetTitle(m_ylabel.c_str());
310 }
311 histo->GetXaxis()->CenterTitle();
312 histo->GetYaxis()->CenterTitle();
313 if(do_only_selected_samples) histo->SetLineColor(this->GetColor(onlyindex[0]));
314 return histo;
315 }
316
317 TH2F* samplecollection::Draw(string m_histoname,string m_var, vector<float> binningx, vector<float> binningy, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, vector<int> onlyindex, bool drawsignal=false) {
318 if(Verbosity>0) dout << endl << endl;
319 if(Verbosity>0) dout << "-------------------------------------------------------------------------------------" << endl;
320 if(Verbosity>0) dout << "histoname : " << m_histoname << " , m_var = " << m_var << ", m_xlabel=" << m_xlabel << " m_ylabel=" << m_ylabel << " m_is_data = " << m_is_data << " lumi = " << luminosity << " onlyindex size: " << onlyindex.size() << endl;
321 if(HUSH==0) dout << "Drawing histo called " << m_histoname << "... " << endl;
322 bool do_only_selected_samples=false;
323 if(onlyindex.size()>0&&onlyindex[0]!=-1) {
324 if(Verbosity>0) {dout << "Requested to only draw sample corresponding to the following sample(s) : " << endl;}
325 for(int is=0;is<onlyindex.size();is++) {
326 if(Verbosity>0) dout << " - " << (this->collection)[onlyindex[is]].filename << " (sample index: " << onlyindex[is] << ")" << endl;
327 }
328 do_only_selected_samples=true;
329 }
330 if(onlyindex.size()==1&&onlyindex[0]==-1) {
331 do_only_selected_samples=true; // this is the special case when we request a non-existing sample - we shouldn't draw anything then.
332 onlyindex.clear();
333 }
334 stringstream h_histoname;
335 h_histoname<<"h_"<<m_histoname;
336 float binningarray[binningx.size()+1];
337 float binningyarray[binningy.size()+1];
338 for(int i=0;i<binningx.size();i++) {binningarray[i]=binningx[i];binningyarray[i]=binningy[i];}
339 TH2F *histo = new TH2F(m_histoname.c_str(),"",binningx.size()-1,binningarray,binningy.size(),binningyarray);
340 histo->Sumw2();
341
342 stringstream drawthis;
343 drawthis<<m_var<<">>tempdrawhisto";
344
345 for (unsigned int isample=0;isample<(this->collection).size();isample++) {
346 bool use_this_sample=false;
347 if(!(this->collection)[isample].is_active) continue; // skip inactive samples right away
348 if(((this->collection)[isample].is_data==m_is_data)&&(!do_only_selected_samples)) {
349 if(drawsignal==false&&(this->collection)[isample].is_signal==true) continue;
350 use_this_sample=true;
351 }
352 if(do_only_selected_samples&&this->do_sample(isample,onlyindex)) {
353 use_this_sample=true;
354 }
355 TH2F *tempdrawhisto = new TH2F("tempdrawhisto","tempdrawhisto",binningx.size()-1,binningarray,binningy.size()-1,binningyarray);
356 tempdrawhisto->Sumw2();
357 if(use_this_sample) {
358 if(Verbosity>0) dout << "[samplecollection::Draw] : Added contribution from sample " << (this->collection)[isample].samplename << endl;
359 (this->collection)[isample].events->Draw(drawthis.str().c_str(),(essentialcut&&Cut)*cutWeight);//this weight is based on PU etc. not XS
360 if(addoverunderflowbins) {
361 //now also adding the overflow & underflow bins:
362 tempdrawhisto->SetBinContent(tempdrawhisto->GetNbinsX(),tempdrawhisto->GetBinContent(tempdrawhisto->GetNbinsX()+1)+tempdrawhisto->GetBinContent(tempdrawhisto->GetNbinsX()));
363 tempdrawhisto->SetBinError(tempdrawhisto->GetNbinsX(),TMath::Sqrt(tempdrawhisto->GetBinContent(tempdrawhisto->GetNbinsX())));
364 tempdrawhisto->SetBinContent(1,tempdrawhisto->GetBinContent(0)+tempdrawhisto->GetBinContent(1));
365 tempdrawhisto->SetBinError(1,TMath::Sqrt(tempdrawhisto->GetBinContent(1)));
366 }
367
368 if(!(this->collection)[isample].is_data) tempdrawhisto->Scale(luminosity*((this->collection)[isample].weight));//weight applied here is XS & N(entries)
369 histo->Add(tempdrawhisto);
370 }
371 tempdrawhisto->Delete();
372 }//end of loop over isample
373 if(Verbosity>0) dout << "Histo has been filled and now contains " << histo->Integral() << " points (integral)" << endl;
374 histo->GetXaxis()->SetTitle(m_xlabel.c_str());
375 histo->GetYaxis()->SetTitle(m_ylabel.c_str());
376 histo->GetXaxis()->CenterTitle();
377 histo->GetYaxis()->CenterTitle();
378 if(do_only_selected_samples) histo->SetLineColor(this->GetColor(onlyindex[0]));
379 return histo;
380 }
381
382
383 TH1F* samplecollection::Draw(string m_histoname,string m_var, int m_nbins, float m_minx, float m_maxx, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, vector<int> onlyindex, bool drawsignal=false)
384 {
385 vector<float> binning;
386 for(int i=0;i<=m_nbins;i++)
387 {
388 binning.push_back(((float)(m_maxx-m_minx)/((float)m_nbins))*i+m_minx);
389 }
390
391 TH1F *histo = this->Draw(m_histoname,m_var, binning, m_xlabel, m_ylabel, Cut, m_is_data, luminosity, onlyindex,drawsignal);
392 return histo;
393 }
394
395 THStack samplecollection::DrawStack(string m_histoname,string m_var, vector<float> binning, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, bool drawsignal=false) {
396 stringstream h_histoname;
397 h_histoname<<"h_"<<m_histoname;
398 THStack thestack("thestack",m_histoname.c_str());
399
400 stringstream drawthis;
401 drawthis<<m_var<<">>tempdrawhisto";
402
403 TH1F *histogroups[this->ngroups+1];
404 int bookedhistos[this->ngroups+1];
405 for(int ih=0;ih<=this->ngroups;ih++) bookedhistos[ih]=0;
406
407 float binningarray[binning.size()+1];
408 for(int i=0;i<binning.size();i++) {binningarray[i]=binning[i];}
409
410 for (unsigned int isample=0;isample<(this->collection).size();isample++) {
411 if(!drawsignal&&(this->collection)[isample].is_signal) continue;
412 if((this->collection)[isample].is_active&&(this->collection)[isample].is_data==m_is_data) {//fills mc if we want mc, else fills data.
413 TH1F *tempdrawhisto = new TH1F("tempdrawhisto","",binning.size()-1,binningarray);
414 (this->collection)[isample].events->Draw(drawthis.str().c_str(),(essentialcut&&Cut)*cutWeight);
415
416 if(addoverunderflowbins) {
417 //now also adding the overflow & underflow bins:
418 tempdrawhisto->SetBinContent(tempdrawhisto->GetNbinsX(),tempdrawhisto->GetBinContent(tempdrawhisto->GetNbinsX()+1)+tempdrawhisto->GetBinContent(tempdrawhisto->GetNbinsX()));
419 // tempdrawhisto->SetBinError(tempdrawhisto->GetNbinsX(),TMath::Sqrt(tempdrawhisto->GetBinContent(tempdrawhisto->GetNbinsX()))); // errors not necessary for stack!
420 tempdrawhisto->SetBinContent(1,tempdrawhisto->GetBinContent(0)+tempdrawhisto->GetBinContent(1));
421 // tempdrawhisto->SetBinError(1,TMath::Sqrt(tempdrawhisto->GetBinContent(1))); // errors not necessary for stack!
422 }
423
424 if(!(this->collection)[isample].is_data) tempdrawhisto->Scale(luminosity*((this->collection)[isample].weight));
425 tempdrawhisto->SetFillColor(this->GetColor(isample));
426 if(bookedhistos[(this->collection)[isample].groupindex]==0) {
427 tempdrawhisto->SetName(GetNumericHistoName().c_str());
428 histogroups[(this->collection)[isample].groupindex]=(TH1F*)tempdrawhisto->Clone();
429 bookedhistos[(this->collection)[isample].groupindex]=1;
430 }
431 else
432 {
433 histogroups[(this->collection)[isample].groupindex]->Add(tempdrawhisto);
434 bookedhistos[(this->collection)[isample].groupindex]++;
435 }
436 tempdrawhisto->Delete();
437 }
438 }
439 vector<int> ordered_indices;
440 vector<float> ordered_integrals;
441
442 for(int ig=0;ig<=this->ngroups;ig++) {
443 if(bookedhistos[ig]>0) {
444 ordered_indices.push_back(ig);
445 ordered_integrals.push_back(histogroups[ig]->Integral());
446 }
447 }
448
449 /*
450 bubbleSort(ordered_integrals,ordered_indices);
451 // for(int index=ordered_indices.size()-1;index>=0;index--) {
452 for(int index=0;index<ordered_indices.size();index++) {
453 thestack.Add(histogroups[ordered_indices[index]]);
454 }
455 */
456 for(int index=0;index<ordered_indices.size();index++) {
457 thestack.Add(histogroups[ordered_indices[index]]);
458 }
459 return thestack;
460 }
461
462 THStack samplecollection::DrawStack(string m_histoname,string m_var, int m_nbins, float m_minx, float m_maxx, string m_xlabel, string m_ylabel, TCut Cut, bool m_is_data, float luminosity, bool
463 drawsignal=false)
464 {
465 vector<float> binning;
466 for(int i=0;i<=m_nbins;i++)
467 {
468 binning.push_back(((float)(m_maxx-m_minx)/((float)m_nbins))*i+m_minx);
469 }
470 return this->DrawStack(m_histoname,m_var, binning, m_xlabel, m_ylabel, Cut, m_is_data, luminosity, drawsignal);
471 }
472
473
474 vector<int> samplecollection::FindSample(string what)
475 {
476 vector<int> hitcollection;
477 for(unsigned int isam=0;isam<(this->collection).size();isam++)
478 {
479 if(((this->collection)[isam].filename).find(what)!=string::npos) {
480 hitcollection.push_back(isam);
481 } else {
482 }
483 }
484 if(hitcollection.size()==0) {
485 hitcollection.push_back(-1);
486 write_warning(__FUNCTION__,"Couldn't find sample "+string(what)+" using sample collection \""+string(this->name)+"\"");
487 }
488 return hitcollection;
489 }
490
491 Color_t samplecollection::GetColor(string filename) {
492 vector<int> corresponding_samples = this->FindSample(filename);
493 if(corresponding_samples.size()>0) return this->GetColor(corresponding_samples[0]);
494 return kBlack;
495 }
496
497 Color_t samplecollection::GetColor(int sampleindex) {
498 return this->collection[sampleindex].samplecolor;
499 }
500
501 TLegend* samplecollection::allbglegend(string title, TH1F *data,float posx=0.65, float posy=0.60) {
502 // TLegend *leg = new TLegend(0.65,0.60,0.89,0.77);
503 TLegend *leg = new TLegend(posx,posy,0.89,0.89);
504 if(title!="") leg->SetHeader(title.c_str());
505 if(data->GetName()!="nothing") leg->AddEntry(data,"Data","p");
506 leg->SetFillColor(kWhite);
507 leg->SetBorderSize(0);
508 leg->SetLineColor(kWhite);
509
510 TH1F *fakehistos[(this->collection).size()];
511 bool donealready[(this->collection).size()];
512 for(int i=0;i<(this->collection).size();i++) donealready[i]=false;
513 for(int isample=0;isample<(this->collection).size();isample++) {
514 if((this->collection)[isample].is_data||(this->collection)[isample].is_signal) continue;
515
516 if(!donealready[(this->collection)[isample].groupindex]) {
517 fakehistos[(this->collection)[isample].groupindex] = new TH1F(GetNumericHistoName().c_str(),"",1,0,1);
518 fakehistos[(this->collection)[isample].groupindex]->Fill(1);
519 fakehistos[(this->collection)[isample].groupindex]->SetFillColor(this->GetColor(isample));
520 leg->AddEntry(fakehistos[(this->collection)[isample].groupindex],((this->collection)[isample].samplename).c_str(),"f");
521 donealready[(this->collection)[isample].groupindex]=true;
522 }
523 }
524 DrawPrelim();
525 leg->SetTextFont(42);
526 return leg;
527 }
528
529 TLegend* samplecollection::allbglegend(string title="",float posx=0.65, float posy=0.60) {
530 Int_t currlevel=gErrorIgnoreLevel;
531 gErrorIgnoreLevel=5000;
532 TH1F *blub = new TH1F("nothing","nothing",1,0,1);
533 gErrorIgnoreLevel=currlevel;//we know this possibly replaces a previous histo, but we don't care since it's fake anyway.
534 return this->allbglegend(title,blub,posx,posy);
535 }
536
537
538 void set_treename(string treen="events") {
539 dout << "Treename has been set to " << treen << endl;
540 if(treen=="PFevents"||treen=="events") treename=treen;
541 else write_error("sample::set_treename","Setting the treename failed as you provided an invalid tree name.");
542 }
543
544 vector<float> samplecollection::get_optimal_binsize(string variable, TCut cut,int nbins, float low, float high) {
545 TH1F *histo = this->Draw("histo",variable,5000,low,high, "finding_optimal_binsize", "events", cut,0,1000);
546 float eventsperbin=histo->Integral()/nbins;
547 vector<float> binning;
548 binning.push_back(high);
549 float runningsum=0;
550 for(int i=histo->GetNbinsX();i>0;i--) {
551 runningsum+=histo->GetBinContent(i);
552 if(TMath::Abs(runningsum-eventsperbin)<0.05*eventsperbin||runningsum>eventsperbin) {
553 binning.push_back(histo->GetBinLowEdge(i));
554 runningsum=0;
555 }
556 }
557 if(runningsum>0) binning.push_back(low);
558 for(int i=0;i<int(binning.size()/2);i++) {
559 float temp=binning[i];
560 binning[i]=binning[binning.size()-1-i];
561 binning[binning.size()-1-i]=temp;
562 }
563
564 return binning;
565
566 }
567
568 void samplecollection::PickUpFromThisFile(int isamp, string cut, vector<string> &output, vector<string> &pickupfile) {
569 int lumi,eventNum,runNum;
570 float jzb[30];
571 (this->collection)[isamp].events->SetBranchAddress("jzb",&jzb);
572 (this->collection)[isamp].events->SetBranchAddress("lumi",&runNum);
573 (this->collection)[isamp].events->SetBranchAddress("eventNum",&eventNum);
574 (this->collection)[isamp].events->SetBranchAddress("lumi",&lumi);
575 (this->collection)[isamp].events->SetBranchAddress("runNum",&runNum);
576
577
578 TTreeFormula *select = new TTreeFormula("select", cut.c_str()&&essentialcut, (this->collection)[isamp].events);
579 int npickedup=0;
580 for (Int_t entry = 0 ; entry < (this->collection)[isamp].events->GetEntries() ; entry++) {
581 (this->collection)[isamp].events->LoadTree(entry);
582 if (select->EvalInstance()) {
583 (this->collection)[isamp].events->GetEntry(entry);
584 dout << runNum << ":" << lumi << ":" << eventNum << endl;
585 npickedup++;
586 }
587 }
588 dout << "Printed out a total of " << npickedup << " events" << endl;
589 }
590
591
592
593 void samplecollection::PickUpEvents(string cut) {
594 vector<string> output;
595 vector<string> pickupfile;
596 for (int isamp=0;isamp<this->collection.size();isamp++)
597 {
598 if((this->collection)[isamp].is_data) {
599 //we have a data sample !
600 this->PickUpFromThisFile(isamp,cut,output,pickupfile);
601 }
602 }
603 //do something with output and of course the pickup file!
604 }
605
606 //________________________________________________________________________________________
607 // Find units from histogram x-label (looks for '[...]')
608 string samplecollection::find_units(string& xlabel) {
609
610 string units;
611 size_t pos1 = xlabel.find("[");
612 if ( pos1 != string::npos ) {
613 size_t pos2 = xlabel.find("]");
614 units = xlabel.substr(pos1+1,pos2-pos1-1);
615 }// else {
616 // size_t pos1 = xlabel.find("(");
617 // if ( pos1 != string::npos ) {
618 // size_t pos2 = xlabel.find(")");
619 // units = xlabel.substr(pos1+1,pos2-pos1-1);
620 // }
621 // }
622 return units;
623
624 }
625
626 void switch_overunderflow(bool newpos=false) {
627 addoverunderflowbins=newpos;
628 }
629
630