ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/FGolf/ss/Tools/histtools.C
Revision: 1.2
Committed: Fri Sep 3 22:57:21 2010 UTC (14 years, 8 months ago) by fgolf
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Changes since 1.1: +0 -0 lines
State: FILE REMOVED
Error occurred while calculating annotation data.
Log Message:
remove this directory as these tools are now in a general Tools directory

File Contents

# Content
1 #ifndef HISTTOOLS_H
2 #define HISTTOOLS_H
3
4 #include "TList.h"
5 #include "TObjArray.h"
6 #include "TH1.h"
7 #include "TH2.h"
8 #include "TLegend.h"
9 #include "THStack.h"
10 #include "TCanvas.h"
11 #include "TFile.h"
12 #include "TRegexp.h"
13 #include "TKey.h"
14 #include <iostream>
15 #include <vector>
16
17 #include "histtools.h"
18 namespace hist {
19
20 //Add all histograms whose names match the given regular expression pattern
21 //or begin with the given prefix. If the final histogram named outHistName
22 //does not exist it is created.
23
24 void add(const char* outHistName, const char* patORpfx) {
25 TRegexp reg(patORpfx, kFALSE);
26
27 TList* list = gDirectory->GetList() ;
28 TIterator* iter = list->MakeIterator();
29
30 TObject* obj = 0;
31 TObject* hist = 0;
32 Bool_t makeOutHist = false;
33
34 hist = gDirectory->Get(outHistName);
35 //If out hist does not exist, remember to create it
36 if (! hist) makeOutHist = true;
37
38 while ((obj = iter->Next())) {
39 if (! obj->InheritsFrom(TH1::Class())) continue;
40
41 TString name = obj->GetName();
42 //Don't add out hist
43 if (name == TString(outHistName)) continue;
44
45 if (TString(patORpfx).MaybeRegexp()) {
46 if (TString(obj->GetName()).Index(reg) < 0 ) continue;
47 } else if (! name.BeginsWith(patORpfx)) continue;
48
49 if (makeOutHist) {
50 hist = obj->Clone(outHistName);
51
52 if (hist->InheritsFrom(TH2::Class()))
53 ((TH2*)hist)->Reset();
54 else
55 ((TH1*)hist)->Reset();
56
57 ((TH1*)hist)->SetTitle(outHistName);
58 ((TH1*)hist)->Sumw2();
59 makeOutHist = false;
60 }
61
62 ((TH1*)hist)->Add((TH1*)obj);
63 }
64 }
65
66 //Add all histograms whose names match one of ten possible regular expression
67 //patterns or begin with one of ten possible given prefixes. Feel free to
68 //mix and match regular expression patterns and prefixes. If the final hist-
69 //ogram named outHistName does not exist it is created.
70
71 void add(const char* outHistName, const char* patORpfx0, const char* patORpfx1, const char* patORpfx2,
72 const char* patORpfx3, const char* patORpfx4, const char* patORpfx5, const char* patORpfx6,
73 const char* patORpfx7, const char* patORpfx8, const char* patORpfx9)
74 {
75 add(outHistName, patORpfx0);
76 add(outHistName, patORpfx1);
77 if (patORpfx2) add(outHistName, patORpfx2);
78 if (patORpfx3) add(outHistName, patORpfx3);
79 if (patORpfx4) add(outHistName, patORpfx4);
80 if (patORpfx5) add(outHistName, patORpfx5);
81 if (patORpfx6) add(outHistName, patORpfx6);
82 if (patORpfx7) add(outHistName, patORpfx7);
83 if (patORpfx8) add(outHistName, patORpfx8);
84 if (patORpfx9) add(outHistName, patORpfx9);
85 }
86
87
88 //For all histograms whose names match the given regular expression pattern
89 //or begin with the given prefix, set the fill, line and marker colors to the
90 //given value.
91
92 void color(const char* patORpfx, Color_t color) {
93
94 TRegexp reg(patORpfx, kFALSE);
95
96 TList* list = gDirectory->GetList() ;
97
98 TIterator* iter = list->MakeIterator();
99
100 TObject* obj = 0;
101
102 while ((obj = iter->Next())) {
103 if (! obj->InheritsFrom(TH1::Class())) continue;
104
105 TString name = obj->GetName();
106
107 if (TString(patORpfx).MaybeRegexp()) {
108 if (TString(obj->GetName()).Index(reg) < 0 ) continue;
109 } else if (! name.BeginsWith(patORpfx)) continue;
110
111 ((TH1*)obj)->SetFillColor(color);
112 ((TH1*)obj)->SetLineColor(color);
113 ((TH1*)obj)->SetMarkerColor(color);
114 }
115 }
116
117
118 void pattern(const char* patORpfx, Color_t pattern) {
119 TRegexp reg(patORpfx, kFALSE);
120
121 TList* list = gDirectory->GetList() ;
122 TIterator* iter = list->MakeIterator();
123
124 TObject* obj = 0;
125
126 while ((obj = iter->Next())) {
127 if (! obj->InheritsFrom(TH1::Class())) continue;
128
129 TString name = obj->GetName();
130
131 if (TString(patORpfx).MaybeRegexp()) {
132 if (TString(obj->GetName()).Index(reg) < 0 ) continue;
133 } else if (! name.BeginsWith(patORpfx)) continue;
134
135 ((TH1*)obj)->SetFillStyle(pattern);
136 }
137 }
138
139
140
141
142 //Return a pointer to a TLegend with an entry for each histogram drawn on a
143 //given TCanvas. Display either the line, point or fill values. Optionally
144 //apply colors to all histograms. By default, entry labels are the names of
145 //their respective histograms. Optionally, if histogram names are of the
146 //form XX_YY_ZZ_WW, entry labels can be XX (token=0), YY (token=1), etc.
147
148 TLegend* legend(TCanvas* canvas, Option_t* option, Bool_t addColor, Int_t token,
149 Float_t xmin, Float_t ymin, Float_t xmax, Float_t ymax) {
150 if(! canvas) return 0;
151
152 TLegend* leg = new TLegend(xmin, ymin, xmax, ymax);
153 TList* list = canvas->GetListOfPrimitives();
154 // TIterator* iter = list->MakeIterator();
155
156 TObject* obj = 0;
157
158 //Hist color iterator
159 Int_t colorIt = 1;
160
161 //while (obj = iter->Next()) {
162 for(int i = list->GetSize() - 1 ; i >=0 ; i--) {
163 obj = list->At(i);
164 if (! obj->InheritsFrom(TH1::Class())) continue;
165
166 if (addColor) {
167 hist::color(obj->GetName(), colorIt);
168 ++colorIt;
169 }
170
171
172 if (token == -1)
173 leg->AddEntry(obj, obj->GetName(), option);
174 else {
175 TString name(obj->GetName());
176 TObjArray* a = name.Tokenize("_");
177 if (a->GetEntries() <= token)
178 leg->AddEntry(obj, obj->GetName(), option);
179 else
180 leg->AddEntry(obj, a->At(token)->GetName(), option);
181 }
182 }
183
184 return leg;
185 }
186
187 //Return a pointer to a TLegend with an entry for each histogram added to a
188 //given THStack. Display either the line, point or fill values. Optionally
189 //apply colors to all histograms. By default, entry labels are the names of
190 //their respective histograms. Optionally, if histogram names are of the
191 //form XX_YY_ZZ_WW, entry labels can be XX (token=0), YY (token=1), etc.
192
193 TLegend* legend(THStack* stack, Option_t* option, Bool_t addColor, Int_t token,
194 Float_t xmin, Float_t ymin, Float_t xmax, Float_t ymax,
195 TH1F *hdata) {
196
197 if(! stack) return 0;
198
199 TLegend* leg = new TLegend(xmin, ymin, xmax, ymax);
200 leg->SetFillColor(kWhite);
201 TList* list = stack->GetHists();
202 // TIterator* iter = list->MakeIterator();
203
204 TObject* obj = 0;
205
206
207 if(hdata != NULL)
208 leg->AddEntry(hdata, "Data", "P");
209
210
211 //Hist color iterator
212 Int_t colorIt = 1;
213
214 //while (obj = iter->Next()) {
215 for(int i = list->GetSize() - 1 ; i >=0 ; i--) {
216 obj = list->At(i);
217 if (! obj->InheritsFrom(TH1::Class())) continue;
218 TString name(obj->GetName());
219
220 if(name.Contains("QCD"))
221 continue;
222
223 if (addColor) {
224 hist::color(obj->GetName(), colorIt);
225 ++colorIt;
226 }
227
228 if (token == -1)
229 leg->AddEntry(obj, obj->GetName(), option);
230 else {
231 TString name(obj->GetName());
232 TObjArray* a = name.Tokenize("_");
233 TString entry(a->GetEntries() <= token ? obj->GetName() : a->At(token)->GetName());
234
235 if (entry == "t"){
236 entry = "Single top";
237 }else if (entry == "ttbar"){
238 entry = "#font[12]{t#bar{t}}";
239 } else if (entry == "ttdil"){
240 entry = "#font[12]{t#bar{t}} signal";
241 } else if (entry == "ttotr"){
242 entry = "#font[12]{t#bar{t}} other";
243 } else if (entry == "VV"){
244 //entry = "#font[12]{ZZ}, #font[12]{WZ}, #font[12]{ZZ}";
245 entry = "#font[12]{VV}";
246 } else if (entry == "DYtautau"){
247 //entry = "DY#rightarrow #tau#tau";
248 entry = "Z/#gamma*#rightarrow#tau^{+}#tau^{-}";
249 } else if (entry == "DYeemm"){
250 //entry = "DY#rightarrow #font[12]{ee}, #mu#mu";
251 entry = "Z/#gamma*#rightarrowl^{+}l^{-}";
252 } else if (entry == "wjets"){
253 //entry = "#font[12]{W}+jets";
254 entry = "#font[12]{W}#rightarrowl#nu";
255 }else {
256 entry.ReplaceAll("tautau", "#tau#tau");
257 entry.ReplaceAll("mm", "#mu#mu");
258 }
259 leg->AddEntry(obj, entry, option);
260 }
261 }
262
263 return leg;
264 }
265
266 //Normalize to one all histograms whose names match the given regular exp-
267 //ression pattern or begin with the given prefix.
268
269 void normalize(const char* patORpfx) {
270 TRegexp reg(patORpfx, kFALSE);
271
272 TList* list = gDirectory->GetList() ;
273 TIterator* iter = list->MakeIterator();
274
275 TObject* obj = 0;
276
277 while ((obj = iter->Next())) {
278 if (! obj->InheritsFrom(TH1::Class())) continue;
279
280 TString name = obj->GetName();
281
282 if (TString(patORpfx).MaybeRegexp()) {
283 if (TString(obj->GetName()).Index(reg) < 0 ) continue;
284 } else if (! name.BeginsWith(patORpfx)) continue;
285
286 Double_t integral = 0;
287
288 if (obj->InheritsFrom(TH2::Class()))
289 integral = ((TH2*)obj)->Integral();
290 else
291 integral = ((TH1*)obj)->Integral();
292
293 if (integral) {
294 ((TH1*)obj)->Sumw2();
295 ((TH1*)obj)->Scale(1./integral);
296 }
297 }
298 }
299
300 //Scale by the given value all histograms whose names match the given regular
301 //expression pattern or begin with the given prefix.
302
303 void scale(const char* patORpfx, Double_t scale) {
304 TRegexp reg(patORpfx, kFALSE);
305
306 TList* list = gDirectory->GetList() ;
307 TIterator* iter = list->MakeIterator();
308
309 TObject* obj = 0;
310
311 while ((obj = iter->Next())) {
312 if (! obj->InheritsFrom(TH1::Class())) continue;
313
314 TString name = obj->GetName();
315
316 if (TString(patORpfx).MaybeRegexp()) {
317 if (TString(obj->GetName()).Index(reg) < 0 ) continue;
318 } else if (! name.BeginsWith(patORpfx)) continue;
319
320
321 if ((((TH1*)obj)->GetNbinsX()+2) != ((TH1*)obj)->GetSumw2N())
322 ((TH1*)obj)->Sumw2();
323
324 ((TH1*)obj)->Scale(scale);
325 }
326 }
327
328 //Don't you hate it when you draw multiple histograms on the same canvas only
329 //to find that the bottom histogram's range does not encompass those of the
330 //histograms drawn on top? This method determines the maximum and minimum y
331 //range of all the histograms drawn on a given TCanvas and appropriately re-
332 //sizes the bottom histogram.
333
334 void setrangey(TCanvas* canvas) {
335 if(! canvas) return ;
336
337 TList* list = canvas->GetListOfPrimitives();
338 TIterator* iter = list->MakeIterator();
339
340 TObject* obj = 0;
341 TObject* top = 0;
342
343 //Extremes
344 Double_t maxy = -999999;
345 Double_t miny = 999999;
346
347 while ((obj = iter->Next())) {
348 if (! obj->InheritsFrom(TH1::Class())) continue;
349
350 if (! top) top = obj;
351
352 if (((TH1*)obj)->GetMaximum() > maxy) maxy = ((TH1*)obj)->GetMaximum();
353 if (((TH1*)obj)->GetMinimum() < miny) miny = ((TH1*)obj)->GetMinimum();
354 }
355
356 ((TH1*)top)->SetMaximum(maxy*1.3);
357 //Protect against log scale
358 if (canvas->GetLogy() && ! miny)
359 ((TH1*)top)->SetMinimum(1E-3);
360 else
361 ((TH1*)top)->SetMinimum(miny*0.7);
362 }
363
364 //Create a stacked histogram consisting of all histograms whose names match
365 //the given regular expression pattern or begin with the given prefix. If
366 //the THStack named stackHistName does not exist it is created. Optionally
367 //apply colors to all histograms. Set drawOption to "nostack" if you do not
368 //want to stack, to "hist" to display histograms without errors, to "histe"
369 //to display histograms with errors, etc.
370
371 void stack(const char* stackHistName, const char* patORpfx, Bool_t addColor, Option_t* drawOption, Int_t orderScheme,
372 const char* bsmName, bool doRefPats, std::string prefixToSkip){//, float yMinimum) {
373 TRegexp reg(patORpfx, kTRUE);
374
375 TList* list = gDirectory->GetList() ;
376
377 TObject* obj = 0;
378 TObject* stack = 0;
379 Bool_t makeStackHist = false;
380
381 stack = gDirectory->Get(stackHistName);
382 //If stack hist does not exist, remember to create it
383 if (! stack) makeStackHist = true;
384
385 //Hist color iterator
386 Int_t colorIt = 1;
387
388 std::vector<TString*> pats(0);
389 if (orderScheme == 0){
390 pats.push_back(new TString(""));
391 //pats.push_back(new TString(""));
392 //pats.push_back(new TString("ttdil_"));
393 //pats.push_back(new TString("ttotr_"));
394 //pats.push_back(new TString("DYeemm_"));
395 //pats.push_back(new TString("DYtautau_"));
396 //pats.push_back(new TString("wjets_"));
397 //pats.push_back(new TString("VV_"));
398 //pats.push_back(new TString("tw_"));
399 //pats.push_back(new TString("data_"));
400 } else if (orderScheme == 1) {
401 //this requires some care: only the histograms beginning with what's in the pattern will be processed
402 pats.push_back(new TString("ttdil_"));
403 pats.push_back(new TString("ttotr_"));
404 pats.push_back(new TString("t_"));
405 pats.push_back(new TString("VV_"));
406 pats.push_back(new TString("DYtautau_"));
407 pats.push_back(new TString("DYeemm_"));
408 pats.push_back(new TString("wjets_"));
409 pats.push_back(new TString("QCD_"));
410 // pats.push_back(new TString("Vgamma_"));
411
412 if (doRefPats){
413 pats.push_back(new TString("ref_ttdil_"));
414 pats.push_back(new TString("ref_ttotr_"));
415 pats.push_back(new TString("ref_t_"));
416 pats.push_back(new TString("ref_VV_"));
417 pats.push_back(new TString("ref_DYtautau_"));
418 pats.push_back(new TString("ref_DYeemm_"));
419 pats.push_back(new TString("ref_wjets_"));
420 pats.push_back(new TString("ref_QCD_"));
421 // pats.push_back(new TString("ref_Vgamma_"));
422 }
423 } else if (orderScheme == 2) {
424 //this requires some care: only the histograms beginning with what's in the pattern will be processed
425 pats.push_back(new TString("ttdil_"));
426 pats.push_back(new TString("ttotr_"));
427 pats.push_back(new TString("t_"));
428 pats.push_back(new TString("VV_"));
429 pats.push_back(new TString("DYtautau_"));
430 pats.push_back(new TString("DYeemm_"));
431 pats.push_back(new TString("wjets_"));
432 pats.push_back(new TString("QCD_"));
433 pats.push_back(new TString("Vgamma_"));
434 if (bsmName!=0) pats.push_back(new TString(Form("%s_",bsmName)));
435
436 if (doRefPats){
437 pats.push_back(new TString("ref_ttdil_"));
438 pats.push_back(new TString("ref_ttotr_"));
439 pats.push_back(new TString("ref_t_"));
440 pats.push_back(new TString("ref_VV_"));
441 pats.push_back(new TString("ref_DYtautau_"));
442 pats.push_back(new TString("ref_DYeemm_"));
443 pats.push_back(new TString("ref_wjets_"));
444 pats.push_back(new TString("ref_QCD_"));
445 pats.push_back(new TString("ref_Vgamma_"));
446 if(bsmName!=0) pats.push_back(new TString(Form("ref_%s_",bsmName)));
447 }
448 } else {
449 pats.push_back(new TString(""));
450 }
451
452 for (unsigned int iPat = 0; iPat < pats.size(); ++iPat) {
453 TIterator* iter = list->MakeIterator();
454 while ((obj = iter->Next())) {
455 if (! obj->InheritsFrom(TH1::Class())) continue;
456 //((TH1F*)obj)->SetMinimum(yMinimum);
457 TString name = obj->GetName();
458 if(prefixToSkip != "" && name.BeginsWith(prefixToSkip.c_str())) continue;
459 //if(name.Contains("data"))
460 if (!( name.BeginsWith(pats[iPat]->Data()))) continue;
461
462 if (TString(patORpfx).MaybeRegexp()) {
463 if (TString(obj->GetName()).Index(reg) < 0 ) continue;
464 } else if (! name.BeginsWith(patORpfx)) continue;
465
466 if (makeStackHist) {
467 stack = new THStack(stackHistName, stackHistName);
468 makeStackHist = false;
469 }
470
471 if (addColor) {
472 hist::color(obj->GetName(), colorIt);
473 ++colorIt;
474 }
475 ((THStack*)stack)->Add((TH1*)obj, drawOption);
476
477 }
478 // delete pats[iPat];//cleanup
479 }
480
481
482 // Currently breaks .ls
483 //gDirectory->Append(stack);
484 }
485
486 //Set the x-axis title of all histograms whose names match the given regular
487 //expression pattern or begin with the given prefix.
488
489 void xaxis(const char* patORpfx, const char* title) {
490 TRegexp reg(patORpfx, kFALSE);
491
492 TList* list = gDirectory->GetList() ;
493 TIterator* iter = list->MakeIterator();
494
495 TObject* obj = 0;
496
497 while ((obj = iter->Next())) {
498 if (! (obj->InheritsFrom(TH1::Class()) || obj->InheritsFrom(THStack::Class()))) continue;
499
500 TString name = obj->GetName();
501
502 if (TString(patORpfx).MaybeRegexp()) {
503 if (TString(obj->GetName()).Index(reg) < 0 ) continue;
504 } else if (! name.BeginsWith(patORpfx)) continue;
505
506 if (obj->InheritsFrom(TH1::Class()))
507 ((TH1*)obj)->GetXaxis()->SetTitle(title);
508 if (obj->InheritsFrom(THStack::Class())) {
509 ((THStack*)obj)->Draw();
510 ((THStack*)obj)->GetXaxis()->SetTitle(title);
511 }
512 }
513 }
514
515 //Set the y-axis title of all histograms whose names match the given regular
516 //expression pattern or begin with the given prefix.
517
518 void yaxis(const char* patORpfx, const char* title) {
519 TRegexp reg(patORpfx, kFALSE);
520
521 TList* list = gDirectory->GetList() ;
522 TIterator* iter = list->MakeIterator();
523
524 TObject* obj = 0;
525
526 while ((obj = iter->Next())) {
527 if (! (obj->InheritsFrom(TH1::Class()) || obj->InheritsFrom(THStack::Class()))) continue;
528
529 TString name = obj->GetName();
530
531 if (TString(patORpfx).MaybeRegexp()) {
532 if (TString(obj->GetName()).Index(reg) < 0 ) continue;
533 } else if (! name.BeginsWith(patORpfx)) continue;
534
535 if (obj->InheritsFrom(TH1::Class()))
536 ((TH1*)obj)->GetYaxis()->SetTitle(title);
537 if (obj->InheritsFrom(THStack::Class())) {
538 ((THStack*)obj)->Draw();
539 ((THStack*)obj)->GetYaxis()->SetTitle(title);
540 }
541 }
542 }
543
544 // Input: 2 histogram
545 // Output: one histogram which is the efficiency:
546 // h1 : TOTAL NUMBER OF EVENTS
547 // h2 : NUMBER OF EVENTS THAT PASS
548
549 #include "TH1.h"
550
551 // Method by pointer
552 TH1F* eff(TH1F* h1, TH1F* h2, const char* name){
553
554 // first, verify that all histograms have same binning
555 // nx is the number of visible bins
556 // nxtot = nx+2 includes underflow and overflow
557 Int_t nx = h1->GetNbinsX();
558 if (h2->GetNbinsX() != nx) {
559 cout << "Histograms must have same number of bins" << endl;
560 return 0;
561 }
562
563 // get the new histogram
564 TH1F* temp = (TH1F*) h1->Clone(name);
565 temp->SetTitle(name);
566 temp->Reset();
567 temp->Sumw2();
568
569 // Do the calculation
570 temp->Divide(h2,h1,1.,1.,"B");
571
572 // Done
573 return temp;
574 }
575
576
577 // Method by name
578 TH1F* eff(const char* name1, const char* name2, const char* name){
579
580 // Get a list of object and their iterator
581 TList* list = gDirectory->GetList() ;
582 TIterator* iter = list->MakeIterator();
583
584 // Loop over objects, set the pointers
585 TObject* obj;
586 TH1F* h1=0;
587 TH1F* h2=0;
588 TString str1 = Form("%s",name1);
589 TString str2 = Form("%s",name2);
590 while((obj=iter->Next())) {
591 TString objName = obj->GetName();
592 if (objName == str1) h1 = (TH1F*) obj;
593 if (objName == str2) h2 = (TH1F*) obj;
594 }
595
596 // quit if not found
597 if (h1 == 0) {
598 cout << "Histogram " << name1 << " not found" << endl;
599 return 0;
600 }
601 if (h2 == 0) {
602 cout << "Histogram " << name2 << " not found" << endl;
603 return 0;
604 }
605
606 // Call the method by pointer
607 TH1F* temp = eff(h1, h2, name);
608 return temp;
609 }
610 // Input: 4 histogram
611 // Output: one histogram which is the BG subtracted efficiency:
612 // h1 : TOTAL NUMBER OF EVENTS, SIGNAL REGION
613 // h2 : NUMBER OF EVENTS THAT PASS, SIGNAL REGION
614 // h3 : TOTAL NUMBER OF EVENTS, SIDE BAND
615 // h4 : NUMBER OF EVENTS THAT PASS, SIDE BAND
616
617 #include "TH1.h"
618
619
620 TH1F* eff_bg(TH1F* h1, TH1F* h2, TH1F* h3, TH1F* h4, const char* name){
621
622 // first, verify that all histograms have same binning
623 // nx is the number of visible bins
624 // nxtot = nx+2 includes underflow and overflow
625 Int_t nx = h1->GetNbinsX();
626 Int_t nxtot = nx + 2;
627 if (h2->GetNbinsX() != nx) {
628 cout << "Histograms must have same number of bins" << endl;
629 return 0;
630 }
631 if (h3->GetNbinsX() != nx) {
632 cout << "Histograms must have same number of bins" << endl;
633 return 0;
634 }
635 if (h3->GetNbinsX() != nx) {
636 cout << "Histograms must have same number of bins" << endl;
637 return 0;
638 }
639
640 // get the new histogram
641 TH1F* temp = (TH1F*) h1->Clone(name);
642 temp->SetTitle(name);
643 temp->Reset();
644 temp->Sumw2();
645
646 // Loop over bins, calculate efficiency and error, put it in histogram
647 for (Int_t i=0; i<nxtot; i++) {
648 Double_t x1 = h1->GetBinContent(i);
649 Double_t x2 = h2->GetBinContent(i);
650 Double_t x3 = h3->GetBinContent(i);
651 Double_t x4 = h4->GetBinContent(i);
652 Double_t denom = x1 - x3;
653 Double_t eff;
654 if (denom == 0.) {
655 eff = 0;
656 } else {
657 eff = (x2-x4)/denom;
658 }
659 Double_t failSig = x1 - x2;
660 Double_t failBg = x3 - x4;
661 Double_t blah = (1-eff)*(1-eff)*(x2+x4) + eff*eff*(failSig+failBg);
662 if (blah <= 0.) blah=0.0;
663 Double_t err;
664 if (denom == 0) {
665 err = 0.;
666 } else {
667 err = sqrt(blah)/denom;
668 }
669 temp->SetBinContent(i,eff);
670 temp->SetBinError(i,err);
671 }
672
673 // Done
674 return temp;
675 }
676
677 #include <TList.h>
678 #include <TIterator.h>
679
680 void deleteHistos() {
681 // Delete all existing histograms in memory
682 TObject* obj;
683 TList* list = gDirectory->GetList() ;
684 TIterator* iter = list->MakeIterator();
685 while ((obj=iter->Next())) {
686 if (obj->IsA()->InheritsFrom(TH1::Class()) ||
687 obj->IsA()->InheritsFrom(TH2::Class()) ) {delete obj;}
688 }
689 }
690 void histio()
691 {
692 }
693
694 void saveHist(const char* filename, const char* pat)
695 {
696 TList* list = gDirectory->GetList() ;
697 TIterator* iter = list->MakeIterator();
698
699 TRegexp re(pat,kTRUE) ;
700
701 TFile outf(filename,"RECREATE") ;
702 while(TObject* obj=iter->Next()) {
703 if (TString(obj->GetName()).Index(re)>=0) {
704 obj->Write() ;
705 cout << "." ;
706 cout.flush() ;
707 }
708 }
709 cout << endl ;
710 outf.Close() ;
711
712 delete iter ;
713 }
714
715
716 void loadHist(const char* filename, const char* pfx, const char* pat, Bool_t doAdd)
717 {
718 TFile inf(filename) ;
719 //inf.ReadAll() ;
720 TList* list = inf.GetListOfKeys() ;
721 TIterator* iter = list->MakeIterator();
722
723 TRegexp re(pat,kTRUE) ;
724 if (pat!=0) cout << "pat = " << pat << endl ;
725 else cout<<"no pattern: read all"<<endl;
726
727 gDirectory->cd("Rint:") ;
728
729 TObject* obj ;
730 TKey* key ;
731 cout << "doAdd = " << (doAdd?"T":"F") << endl ;
732 cout << "loadHist: reading." ;
733 while((key=(TKey*)iter->Next())) {
734
735 // Int_t ridx = TString(key->GetName()).Index(re) ;
736 // cout<<"reading "<<key->GetName()<<std::endl;
737 if (pat!=0){
738 // cout<<"Check patt "<<pat<<endl;
739 if (TString(key->GetName()).Index(re)==-1) {
740 continue ;
741 }
742 }
743
744 obj = inf.Get(key->GetName()) ;
745 TObject* clone ;
746 if (pfx) {
747
748 // Find existing TH1-derived objects
749 TObject* oldObj = 0 ;
750 if (doAdd){
751 oldObj = gDirectory->Get(Form("%s_%s",pfx,obj->GetName())) ;
752 if (oldObj && !oldObj->IsA()->InheritsFrom(TH1::Class())) {
753 oldObj = 0 ;
754 }
755 }
756 if (oldObj) {
757 clone = oldObj ;
758 ((TH1*)clone)->Add((TH1*)obj) ;
759 } else {
760 clone = obj->Clone(Form("%s_%s",pfx,obj->GetName())) ;
761 }
762
763
764 } else {
765
766 // Find existing TH1-derived objects
767 TObject* oldObj = 0 ;
768 if (doAdd){
769 oldObj = gDirectory->Get(key->GetName()) ;
770 if (oldObj && !oldObj->IsA()->InheritsFrom(TH1::Class())) {
771 oldObj = 0 ;
772 }
773 }
774
775 if (oldObj) {
776 clone = oldObj ;
777 ((TH1*)clone)->Add((TH1*)obj) ;
778 } else {
779 clone = obj->Clone() ;
780 }
781 }
782 if (!gDirectory->GetList()->FindObject(clone)) {
783 gDirectory->Append(clone) ;
784 }
785 //cout << "." ;
786 cout.flush() ;
787 }
788 cout << endl;
789 inf.Close() ;
790 delete iter ;
791 }
792
793
794 //takes in a vector of prefixes, and a output prefix. Combines all histograms
795 //whose prefixes matches those in the input vector of TStrings into a new histogram
796 //with the the outPrefix
797 void combineHists(const std::vector<TString> v_prfxsToCombine, const TString outPrefix) {
798
799
800 if(v_prfxsToCombine.size() == 0) {
801 cout << "Input vector size is 0" << endl;
802 }
803 //histogram of common suffixes
804 vector<TString> v_suffixes;
805 //get the list of histograms that match the first entry in the prefix
806 TString temp = v_prfxsToCombine.at(0);
807 TList *rootlist = gDirectory->GetList();
808
809 for(int i = 0; i < rootlist->GetSize(); i++) {
810
811 TObject *obj = rootlist->At(i);
812 TString name = obj->GetName();
813 if(!obj->InheritsFrom(TH1::Class()))
814 continue;
815 if(!name.BeginsWith(temp+"_"))
816 continue;
817
818 name.ReplaceAll(temp+"_", "_");
819
820
821 if(obj->InheritsFrom(TH2::Class())) {
822 TH2F *h = dynamic_cast<TH2F*>(obj->Clone());
823 h->SetName((outPrefix+name).Data());
824 h->SetTitle((outPrefix+name).Data());
825 for(unsigned int j = 1; j < v_prfxsToCombine.size(); j++) {
826 TH2F *htemp = dynamic_cast<TH2F*>(gDirectory->Get((v_prfxsToCombine.at(j) + name).Data()));
827 h->Add(htemp);
828 }
829 } else if(obj->InheritsFrom(TH1::Class())) {
830 TH1F *h = dynamic_cast<TH1F*>(obj->Clone());
831 h->SetName((outPrefix+name).Data());
832 h->SetTitle((outPrefix+name).Data());
833 for(unsigned int j = 1; j < v_prfxsToCombine.size(); j++) {
834 TH1F *htemp = dynamic_cast<TH1F*>(gDirectory->Get((v_prfxsToCombine.at(j) + name).Data()));
835 //cout << "TH1F: " << v_prfxsToCombine.at(j) + name << endl;
836 if(htemp == NULL)
837 cout << "htemp is NULL" << endl;
838 h->Add(htemp);
839 }
840 }
841 }//rootlist loop
842
843
844 //now delete the histos that match the prfxs to combine
845 for(unsigned int i = 0; i < v_prfxsToCombine.size(); i++ ) {
846
847 // Delete all existing histograms in memory
848 TObject* obj;
849 TList* list = gDirectory->GetList() ;
850 TIterator* iter = list->MakeIterator();
851 while ((obj=iter->Next())) {
852 TString name = obj->GetName();
853 if(name.BeginsWith(outPrefix+"_"))
854 continue;
855 if(TString(obj->GetName()).BeginsWith(v_prfxsToCombine.at(i).Data())) {
856 if (obj->IsA()->InheritsFrom(TH1::Class()) ||
857 obj->IsA()->InheritsFrom(TH2::Class()) )
858 delete obj;
859 }
860 }
861 }//loop over prefixes
862
863 }//fnc declaration
864
865 }
866 #endif