ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/AnalysisFramework/DistributedModelCalculations/ShapeLimit/calc_limits.C
Revision: 1.2
Committed: Wed Aug 31 08:09:54 2011 UTC (13 years, 8 months ago) by buchmann
Content type: text/plain
Branch: MAIN
CVS Tags: cbaf_4_98ifb_paper, cbaf_4p7ifb, Honeypot, cbaf_2p1ifb, HEAD
Changes since 1.1: +46 -12 lines
Log Message:
Rewrote the retrieval of limits from limit tree

File Contents

# Content
1 #include <iostream>
2 #include <vector>
3 #include <sys/stat.h>
4 #include <sstream>
5 #include <algorithm>
6
7 #include <TCut.h>
8 #include <TLatex.h>
9 #include <TROOT.h>
10 #include <TCanvas.h>
11 #include <TFile.h>
12 #include <TTree.h>
13 #include <TString.h>
14 #include <TMath.h>
15 #include <stdio.h>
16 #include <stdlib.h>
17 #include <TH2.h>
18 #include <TColor.h>
19 #include <TStyle.h>
20 #include <TText.h>
21 using namespace std;
22
23 double Round(double num, unsigned int dig)
24 {
25 num *= pow(10, dig);
26 if (num >= 0)
27 num = floor(num + 0.5);
28 else
29 num = ceil(num - 0.5);
30 num/= pow(10, dig);
31 return num;
32 }
33
34 int main (int narg, char * args[]) {
35 if(narg<2||narg>2) {
36 cout << "USAGE : " << endl;
37 cout << args[0] << " fusefile.root" << endl;
38 cout << " where fusefile.root is the file that results from the fusion of all partial shape limit files" << endl;
39 return -1;
40 }
41 TFile *infile = new TFile(args[1]);
42 TTree *limtree = (TTree*)infile->Get("limit");
43 Double_t lim;
44 limtree->SetBranchAddress("limit",&lim);
45 vector<float> vlims;
46 float limsum=0;
47 for (int i=0;i<limtree->GetEntries();i++) {
48 limtree->GetEntry(i);
49 vlims.push_back(lim);
50 }
51
52 float mean=TMath::Mean(limtree->GetEntries(),&vlims[0]);
53 float median=TMath::Median(limtree->GetEntries(),&vlims[0]);
54 float meanerr=TMath::RMS(limtree->GetEntries(),&vlims[0])/TMath::Sqrt(limtree->GetEntries()-1);
55 int nentries=limtree->GetEntries();
56
57 sort(vlims.begin(),vlims.end());
58 unsigned neg95 = unsigned(Round(vlims.size() * (1 - 0.954499736104) * 0.5, 0));
59 unsigned neg68 = unsigned(Round(vlims.size() * (1 - 0.682689492137) * 0.5, 0));
60 unsigned neg99 = unsigned(Round(vlims.size() * (1 - 0.997300203937) * 0.5, 0));
61 unsigned central = unsigned(Round(vlims.size() * 0.5, 0));
62 unsigned pos68 = unsigned(Round(vlims.size() * (1 + 0.682689492137) * 0.5, 0));
63 unsigned pos95 = unsigned(Round(vlims.size() * (1 + 0.954499736104) * 0.5, 0));
64 unsigned pos99 = unsigned(Round(vlims.size() * (1 + 0.997300203937) * 0.5, 0));
65
66 if(pos68>=vlims.size()) pos68=vlims.size()-1;
67 if(pos95>=vlims.size()) pos95=vlims.size()-1;
68 if(pos99>=vlims.size()) pos99=vlims.size()-1;
69
70 cout << "_________________________________________" << endl;
71 cout << endl;
72 cout << " Marco's summary: " << endl;
73 cout << "mean expected limit: r < "<<mean<<" +/- "<<meanerr<<" @ 95%CL ("<<nentries<<" toyMC)" << endl;
74 cout << "median expected limit: r < "<<median<<" @ 95%CL ("<<nentries<<" toyMC)" << endl;
75 cout << "68% expected band : " << vlims[neg68] <<" < r < " << vlims[pos68] << endl;
76 cout << "95% expected band : " << vlims[neg95] <<" < r < " << vlims[pos95] << endl;
77 cout << "99% expected band : " << vlims[neg99] <<" < r < " << vlims[pos99] << endl;
78
79
80
81 return 0;
82 }