ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/AnalysisFramework/Plotting/Modules/ResultModule.C
Revision: 1.15
Committed: Mon Oct 24 15:05:37 2011 UTC (13 years, 6 months ago) by buchmann
Content type: text/plain
Branch: MAIN
Changes since 1.14: +86 -34 lines
Log Message:
Upgrade from HoneyPot to IceCreamBowl: merged in offpeak stuff, different warning color for frederic, saving to rootfile, only 3 attempts when computing limits and much, much more

File Contents

# Content
1 #include <iostream>
2 #include <vector>
3 #include <sys/stat.h>
4
5 #include <TCut.h>
6 #include <TROOT.h>
7 #include <TCanvas.h>
8 #include <TMath.h>
9 #include <TColor.h>
10 #include <TPaveText.h>
11 #include <TRandom.h>
12 #include <TH1.h>
13 #include <TH2.h>
14 #include <TF1.h>
15 #include <TSQLResult.h>
16
17 //#include "TTbar_stuff.C"
18 using namespace std;
19
20 using namespace PlottingSetup;
21
22 void fill_result_histos(float &zossfp, float &zossfperr, float &zosofp, float &zossfn, float &zosofn, float &sbossfp, float &sbosofp, float &sbossfn, float &sbosofn,string datajzb,float cut, float cuthigh, int mcordata, float &result, vector<int> sel, samplecollection &sampleC, string addcut="") {
23 string xlabel="JZB [GeV] -- for algoritm internal use only!";
24 bool dosignal=false;
25 if(mcordata==mcwithsignal) {
26 dosignal=true;
27 mcordata=0;
28 }
29 TCut basiccutplus;
30 if(addcut=="") basiccutplus=basiccut;
31 else basiccutplus=basiccut&&addcut.c_str();
32 TH1F *ZOSSFP;
33 TH1F *ZOSOFP;
34 TH1F *ZOSSFN;
35 TH1F *ZOSOFN;
36
37 TH1F *SBOSSFP;
38 TH1F *SBOSOFP;
39 TH1F *SBOSSFN;
40 TH1F *SBOSOFN;
41
42 if(mcordata==mc||mcordata==data||mcordata==mcwithsignal) {
43 ZOSSFP = sampleC.Draw("ZOSSFP",datajzb,1,cut,cuthigh, xlabel, "events",cutmass&&cutOSSF&&cutnJets&&basiccutplus,mcordata,luminosity,dosignal);
44 ZOSOFP = sampleC.Draw("ZOSOFP",datajzb,1,cut,cuthigh, xlabel, "events",cutmass&&cutOSOF&&cutnJets&&basiccutplus,mcordata,luminosity,dosignal);
45 ZOSSFN = sampleC.Draw("ZOSSFN","-"+datajzb,1,cut,cuthigh, xlabel, "events",cutmass&&cutOSSF&&cutnJets&&basiccutplus,mcordata,luminosity,dosignal);
46 ZOSOFN = sampleC.Draw("ZOSOFN","-"+datajzb,1,cut,cuthigh, xlabel, "events",cutmass&&cutOSOF&&cutnJets&&basiccutplus,mcordata,luminosity,dosignal);
47
48 flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
49 if(PlottingSetup::RestrictToMassPeak) {
50 SBOSSFP = sampleC.Draw("SBOSSFP",datajzb,1,cut,cuthigh, xlabel, "events",cutOSSF&&cutnJets&&basiccutplus&&sidebandcut,mcordata,luminosity,dosignal);
51 SBOSOFP = sampleC.Draw("SBOSOFP",datajzb,1,cut,cuthigh, xlabel, "events",cutOSOF&&cutnJets&&basiccutplus&&sidebandcut,mcordata,luminosity,dosignal);
52 SBOSSFN = sampleC.Draw("SBOSSFN","-"+datajzb,1,cut,cuthigh, xlabel, "events",cutOSSF&&cutnJets&&basiccutplus&&sidebandcut,mcordata,luminosity,dosignal);
53 SBOSOFN = sampleC.Draw("SBOSOFN","-"+datajzb,1,cut,cuthigh, xlabel, "events",cutOSOF&&cutnJets&&basiccutplus&&sidebandcut,mcordata,luminosity,dosignal);
54 }
55 } else {
56 //doing signal only!
57 ZOSSFP = sampleC.Draw("ZOSSFP",datajzb,1,cut,cuthigh, xlabel, "events",cutmass&&cutOSSF&&cutnJets&&basiccutplus,mcordata,luminosity,sel);
58 ZOSOFP = sampleC.Draw("ZOSOFP",datajzb,1,cut,cuthigh, xlabel, "events",cutmass&&cutOSOF&&cutnJets&&basiccutplus,mcordata,luminosity,sel);
59 ZOSSFN = sampleC.Draw("ZOSSFN","-"+datajzb,1,cut,cuthigh, xlabel, "events",cutmass&&cutOSSF&&cutnJets&&basiccutplus,mcordata,luminosity,sel);
60 ZOSOFN = sampleC.Draw("ZOSOFN","-"+datajzb,1,cut,cuthigh, xlabel, "events",cutmass&&cutOSOF&&cutnJets&&basiccutplus,mcordata,luminosity,sel);
61
62 flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
63 if(PlottingSetup::RestrictToMassPeak) {
64 SBOSSFP = sampleC.Draw("SBOSSFP",datajzb,1,cut,cuthigh, xlabel, "events",cutOSSF&&cutnJets&&basiccutplus&&sidebandcut,mcordata,luminosity,sel);
65 SBOSOFP = sampleC.Draw("SBOSOFP",datajzb,1,cut,cuthigh, xlabel, "events",cutOSOF&&cutnJets&&basiccutplus&&sidebandcut,mcordata,luminosity,sel);
66 SBOSSFN = sampleC.Draw("SBOSSFN","-"+datajzb,1,cut,cuthigh, xlabel, "events",cutOSSF&&cutnJets&&basiccutplus&&sidebandcut,mcordata,luminosity,sel);
67 SBOSOFN = sampleC.Draw("SBOSOFN","-"+datajzb,1,cut,cuthigh, xlabel, "events",cutOSOF&&cutnJets&&basiccutplus&&sidebandcut,mcordata,luminosity,sel);
68 }
69 }
70
71 double err; // UGH!
72 zossfp=IntegralAndError(ZOSSFP,1,ZOSSFP->GetNbinsX(),err,"");//making this compatible with my computer which has an outdated version of root.
73 zossfperr = err;
74 zosofp=ZOSOFP->Integral();
75 zossfn=ZOSSFN->Integral();
76 zosofn=ZOSOFN->Integral();
77
78 flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
79 if(PlottingSetup::RestrictToMassPeak) {
80 sbossfp=SBOSSFP->Integral();
81 sbosofp=SBOSOFP->Integral();
82 sbossfn=SBOSSFN->Integral();
83 sbosofn=SBOSOFN->Integral();
84 } else {
85 sbossfp=0;
86 sbosofp=0;
87 sbossfn=0;
88 sbosofn=0;
89 }
90
91 delete ZOSSFP;
92 delete ZOSOFP;
93 delete ZOSSFN;
94 delete ZOSOFN;
95
96 if(PlottingSetup::RestrictToMassPeak) {
97 delete SBOSSFP;
98 delete SBOSOFP;
99 delete SBOSSFN;
100 delete SBOSOFN;
101 }
102
103 flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
104 if(PlottingSetup::RestrictToMassPeak) result = zossfn + (1.0/3.)*(zosofp-zosofn)+ (1.0/3.)*(sbossfp-sbossfn)+ (1.0/3.)*(sbosofp-sbosofn);
105 else result = zossfn + (zosofp-zosofn);
106
107 }
108
109 void fill_result_histos(float &zossfp, float &zossfperr, float &zosofp, float &zossfn, float &zosofn, float &sbossfp, float &sbosofp, float &sbossfn, float &sbosofn,string datajzb,float cut, float cuthigh, int mcordata, float &result, string addcut="") {
110 vector<int> emptyvector;
111 fill_result_histos(zossfp, zossfperr,zosofp, zossfn, zosofn, sbossfp, sbosofp, sbossfn, sbosofn,datajzb,cut,cuthigh,mcordata,result,emptyvector,allsamples,addcut);
112 }
113
114 /*
115 void get_result_above_one_fixed_jzb_value(float cut ,string mcjzb,string datajzb, int mcordata,float jzbpeakerrorMC, float jzbpeakerrorData, TCanvas *rescan, bool chatty=false, bool dopoisson=false) {
116 rescan->cd();
117 if(mcordata==data) dout << "***\nCrunching numbers for JZB>" << cut << " : " << endl;
118
119 float zossfp, zossfperr, zosofp, zossfn, zosofn, sbossfp, sbosofp, sbossfn, sbosofn,result;
120 if(mcordata==mc) fill_result_histos(zossfp, zossfperr, zosofp, zossfn, zosofn, sbossfp, sbosofp, sbossfn, sbosofn,mcjzb,cut,mcordata,result);
121 if(mcordata==data) fill_result_histos(zossfp, zossfperr, zosofp, zossfn, zosofn, sbossfp, sbosofp, sbossfn, sbosofn,datajzb,cut,mcordata,result);
122 if(mcordata==mcwithsignal) fill_result_histos(zossfp, zossfperr, zosofp, zossfn, zosofn, sbossfp, sbosofp, sbossfn, sbosofn,mcjzb,cut,mcordata,result);
123
124 float zossfpee, zossfpeeerr, zosofpee, zossfnee, zosofnee, sbossfpee, sbosofpee, sbossfnee, sbosofnee,resultee;
125 float zossfpmm, zossfpmmerr, zosofpmm, zossfnmm, zosofnmm, sbossfpmm, sbosofpmm, sbossfnmm, sbosofnmm,resultmm;
126 if(mcordata==data) {
127 fill_result_histos(zossfpee, zossfpeeerr, zosofpee, zossfnee, zosofnee, sbossfpee, sbosofpee, sbossfnee, sbosofnee,datajzb,cut,mcordata,resultee,"id1==0");
128 fill_result_histos(zossfpmm, zossfpmmerr, zosofpmm, zossfnmm, zosofnmm, sbossfpmm, sbosofpmm, sbossfnmm, sbosofnmm,datajzb,cut,mcordata,resultmm,"id1==1");
129 }
130
131 float ppzossfp, ppzossfperr, ppzosofp, ppzossfn, ppzosofn, ppsbossfp, ppsbosofp, ppsbossfn, ppsbosofn, ppresult;
132 if(mcordata==mc) fill_result_histos(ppzossfp, ppzossfperr,ppzosofp, ppzossfn, ppzosofn, ppsbossfp, ppsbosofp, ppsbossfn, ppsbosofn,newjzbexpression(mcjzb,jzbpeakerrorMC),cut,mcordata,ppresult);
133 if(mcordata==data) fill_result_histos(ppzossfp, ppzossfperr,ppzosofp, ppzossfn, ppzosofn, ppsbossfp, ppsbosofp, ppsbossfn, ppsbosofn,newjzbexpression(datajzb,jzbpeakerrorData),cut,mcordata,ppresult);
134 if(mcordata==mcwithsignal) fill_result_histos(ppzossfp, ppzossfperr,ppzosofp, ppzossfn, ppzosofn, ppsbossfp, ppsbosofp, ppsbossfn, ppsbosofn,newjzbexpression(mcjzb,jzbpeakerrorMC),cut,mcordata,ppresult);
135
136 float pnzossfp, pnzossfperr, pnzosofp, pnzossfn, pnzosofn, pnsbossfp, pnsbosofp, pnsbossfn, pnsbosofn, pnresult;
137 if(mcordata==mc) fill_result_histos(pnzossfp, pnzossfperr,pnzosofp, pnzossfn, pnzosofn, pnsbossfp, pnsbosofp, pnsbossfn, pnsbosofn,newjzbexpression(mcjzb,-jzbpeakerrorMC),cut,mcordata,pnresult);
138 if(mcordata==data) fill_result_histos(pnzossfp, pnzossfperr,pnzosofp, pnzossfn, pnzosofn, pnsbossfp, pnsbosofp, pnsbossfn, pnsbosofn,newjzbexpression(datajzb,-jzbpeakerrorData),cut,mcordata,pnresult);
139 if(mcordata==mcwithsignal) fill_result_histos(pnzossfp, pnzossfperr,pnzosofp, pnzossfn, pnzosofn, pnsbossfp, pnsbosofp, pnsbossfn, pnsbosofn,newjzbexpression(mcjzb,-jzbpeakerrorMC),cut,mcordata,pnresult);
140
141 float syserr=0;
142 float peakerr=0;
143 float staterr=0;
144 float poissonstaterrup=-999;
145 float poissonstaterrdown=-999;
146 if(fabs(result-pnresult)>fabs(result-ppresult)) peakerr=fabs(result-pnresult); else peakerr=fabs(result-ppresult);
147
148 float zjetsestimateuncert=0.5;
149 float emuncert=0.5;
150 float emsidebanduncert=0.5;
151 float eemmsidebanduncert=0.5;
152 syserr = (zjetsestimateuncert*zossfn)*(zjetsestimateuncert*zossfn);//first term
153 syserr+= ((zosofp)*(zosofp) + (zosofn)*(zosofn))*(1.0/9)*emuncert*emuncert;//sys err from emu method
154 syserr+= ((sbossfp)*(sbossfp)+(sbossfn)*(sbossfn))*(1.0/9)*eemmsidebanduncert*eemmsidebanduncert; // sys err from eemm sidebands
155 syserr+= ((sbosofp)*(sbosofp)+(sbosofn)*(sbosofn))*(1.0/9)*emsidebanduncert*emsidebanduncert; // sys err from emu sidebands
156 syserr=TMath::Sqrt(syserr);
157
158 staterr=TMath::Sqrt(zossfn + (1.0/9)*(zosofp+zosofn)+ (1.0/9)*(sbossfp+sbossfn)+ (1.0/9)*(sbosofp+sbosofn));
159 if(dopoisson) advanced_poisson(zossfn,zosofp,zosofn,sbossfp,sbossfn,sbosofp,sbosofn,poissonstaterrdown,poissonstaterrup);
160
161 float e_to_emu=0.5;
162 float m_to_emu=0.5;
163
164 if(mcordata==mc) dout << " MC :: ";
165 if(mcordata==mcwithsignal) dout << " MC with S :: ";
166 if(mcordata==data) dout << " ";
167 dout << "Observed : " << zossfp << "+-" << zossfperr << endl;
168 if(mcordata==data) dout << " Composition: " << zossfpee << " (ee), " << zossfpmm << " (mm) " << endl;
169
170 if(mcordata==mc) dout << " MC :: ";
171 if(mcordata==mcwithsignal) dout << " MC with S :: ";
172 if(mcordata==data) dout << " ";
173 dout << "Predicted: " << zossfn << " + (1/3)*(" << zosofp << "-" << zosofn<<") [e&mu]+ (1/3)*(" << sbossfp << "-" << sbossfn<<") [SF,SB]+ (1/3)*(" << sbosofp << "-" << sbosofn<<") [OF,SB] = ";
174 dout << zossfn + (1.0/3)*(zosofp-zosofn)+ (1.0/3)*(sbossfp-sbossfn)+ (1.0/3)*(sbosofp-sbosofn)<< " +/- " << peakerr << " (peak ) " << " +/- " << syserr << " (sys) +/- " << staterr << " (stat)" << endl;
175 if(dopoisson) dout << " stat error with Poisson : +" << poissonstaterrup << " - " << poissonstaterrdown << endl;
176
177 if(mcordata==data) {
178 dout << " Composition: (ee)" << zossfnee << " + (1/3)*(" << e_to_emu*zosofp << "-" << e_to_emu*zosofn<<") ["<<e_to_emu<<"*e&mu]+ (1/3)*(" << sbossfpee << "-" << sbossfnee<<") [SF,SB]+ (1/3)*(" << sbosofp*e_to_emu << "-" << sbosofn*e_to_emu<<") ["<<e_to_emu<<"*OF,SB] = ";
179 dout << zossfnee + (1.0/3)*(e_to_emu*zosofp-e_to_emu*zosofn)+ (1.0/3)*(sbossfpee-sbossfnee) + (1.0/3)*(sbosofp*e_to_emu - sbosofn*e_to_emu) << endl;
180 dout << " (mm)" << zossfnmm << " + (1/3)*(" << m_to_emu*zosofp << "-" << m_to_emu*zosofn<<") ["<<m_to_emu<<"*e&mu]+ (1/3)*(" << sbossfpmm << "-" << sbossfnmm<<") [SF,SB]+ (1/3)*(" << sbosofp*m_to_emu << "-" << sbosofn*m_to_emu<<") ["<<m_to_emu<<"*OF,SB] = ";
181 dout << zossfnmm + (1.0/3)*(m_to_emu*zosofp - m_to_emu*zosofn) + (1.0/3)*(sbossfpmm -sbossfnmm)+ (1.0/3)*(sbosofp*m_to_emu-sbosofn*e_to_emu) << endl;
182 }
183
184 if(chatty) {
185 dout << " Pred(ZJets ) \t " << zossfn << endl;
186 dout << " Pred(e&mu;]) \t " << zosofp << "-" << zosofn << " = " << zosofp-zosofn<<endl;
187 dout << " Pred(ossf,sb]) \t " << sbossfp << "-" << sbossfn<<" = "<<sbossfp-sbossfn<<endl;
188 dout << " Pred(osof,sb]) \t " << sbosofp << "-" << sbosofn<<" = "<<sbosofp-sbosofn<<endl;
189 }
190
191 if(mcordata==data) {
192 //store the result!
193 Nobs.push_back(zossfp);
194 Npred.push_back(zossfn + (1.0/3)*(zosofp-zosofn)+ (1.0/3)*(sbossfp-sbossfn)+ (1.0/3)*(sbosofp-sbosofn));
195 float totprederr=0;
196 totprederr+=peakerr*peakerr;
197 totprederr+=syserr*syserr;
198 totprederr+=staterr*staterr;
199 totprederr=TMath::Sqrt(totprederr);
200 Nprederr.push_back(totprederr);
201 }
202 }
203 */
204 vector<float> get_result_between_two_fixed_jzb_values(float cut , float cuthigh, string mcjzb,string datajzb, int mcordata,float jzbpeakerrorMC, float jzbpeakerrorData, TCanvas *rescan, bool chatty=false, bool dopoisson=false, bool writeanything=true) {
205 /*return vector of floats
206 [0] Bpred [1] Bpred uncert [2] Observed [3] Observed error
207 // if we use this for the ratio plot we don't want to see any of the results (hence writeanything=false)
208 */
209 rescan->cd();
210 if(writeanything&&mcordata==data) dout << "***\nCrunching numbers for JZB>" << cut << " : " << endl;
211
212 float zossfp, zossfperr, zosofp, zossfn, zosofn, sbossfp, sbosofp, sbossfn, sbosofn,result;
213 if(mcordata==mc) fill_result_histos(zossfp, zossfperr, zosofp, zossfn, zosofn, sbossfp, sbosofp, sbossfn, sbosofn,mcjzb,cut,cuthigh,mcordata,result);
214 if(mcordata==data) fill_result_histos(zossfp, zossfperr, zosofp, zossfn, zosofn, sbossfp, sbosofp, sbossfn, sbosofn,datajzb,cut,cuthigh,mcordata,result);
215 if(mcordata==mcwithsignal) fill_result_histos(zossfp, zossfperr, zosofp, zossfn, zosofn, sbossfp, sbosofp, sbossfn, sbosofn,mcjzb,cut,cuthigh,mcordata,result);
216
217 float zossfpee, zossfpeeerr, zosofpee, zossfnee, zosofnee, sbossfpee, sbosofpee, sbossfnee, sbosofnee,resultee;
218 float zossfpmm, zossfpmmerr, zosofpmm, zossfnmm, zosofnmm, sbossfpmm, sbosofpmm, sbossfnmm, sbosofnmm,resultmm;
219 if(mcordata==data) {
220 fill_result_histos(zossfpee, zossfpeeerr, zosofpee, zossfnee, zosofnee, sbossfpee, sbosofpee, sbossfnee, sbosofnee,datajzb,cut,cuthigh,mcordata,resultee,"id1==0");
221 fill_result_histos(zossfpmm, zossfpmmerr, zosofpmm, zossfnmm, zosofnmm, sbossfpmm, sbosofpmm, sbossfnmm, sbosofnmm,datajzb,cut,cuthigh,mcordata,resultmm,"id1==1");
222 }
223
224 // float ppzossfp, ppzossfperr, ppzosofp, ppzossfn, ppzosofn, ppsbossfp, ppsbosofp, ppsbossfn, ppsbosofn, ppresult;
225 // if(mcordata==mc) fill_result_histos(ppzossfp, ppzossfperr,ppzosofp, ppzossfn, ppzosofn, ppsbossfp, ppsbosofp, ppsbossfn, ppsbosofn,newjzbexpression(mcjzb,jzbpeakerrorMC),cut,cuthigh,mcordata,ppresult);
226 // if(mcordata==data) fill_result_histos(ppzossfp, ppzossfperr,ppzosofp, ppzossfn, ppzosofn, ppsbossfp, ppsbosofp, ppsbossfn, ppsbosofn,newjzbexpression(datajzb,jzbpeakerrorData),cut,cuthigh,mcordata,ppresult);
227 // if(mcordata==mcwithsignal) fill_result_histos(ppzossfp, ppzossfperr,ppzosofp, ppzossfn, ppzosofn, ppsbossfp, ppsbosofp, ppsbossfn, ppsbosofn,newjzbexpression(mcjzb,jzbpeakerrorMC),cut,cuthigh,mcordata,ppresult);
228
229 float pnzossfp, pnzossfperr, pnzosofp, pnzossfn, pnzosofn, pnsbossfp, pnsbosofp, pnsbossfn, pnsbosofn, pnresult;
230 if(mcordata==mc) fill_result_histos(pnzossfp, pnzossfperr,pnzosofp, pnzossfn, pnzosofn, pnsbossfp, pnsbosofp, pnsbossfn, pnsbosofn,newjzbexpression(mcjzb,-jzbpeakerrorMC),cut,cuthigh,mcordata,pnresult);
231 if(mcordata==data) fill_result_histos(pnzossfp, pnzossfperr,pnzosofp, pnzossfn, pnzosofn, pnsbossfp, pnsbosofp, pnsbossfn, pnsbosofn,newjzbexpression(datajzb,-jzbpeakerrorData),cut,cuthigh,mcordata,pnresult);
232 if(mcordata==mcwithsignal) fill_result_histos(pnzossfp, pnzossfperr,pnzosofp, pnzossfn, pnzosofn, pnsbossfp, pnsbosofp, pnsbossfn, pnsbosofn,newjzbexpression(mcjzb,-jzbpeakerrorMC),cut,cuthigh,mcordata,pnresult);
233
234 float syserr=0;
235 float peakerr=0;
236 float staterr=0;
237 float poissonstaterrup=-999;
238 float poissonstaterrdown=-999;
239 // if(fabs(result-pnresult)>fabs(result-ppresult)) peakerr=fabs(result-pnresult); else peakerr=fabs(result-ppresult);
240
241 float zjetsestimateuncert=0.25;
242 float emuncert=0.5;
243 float emsidebanduncert=0.5;
244 float eemmsidebanduncert=0.5;
245 syserr = (zjetsestimateuncert*zossfn)*(zjetsestimateuncert*zossfn);//first term
246 flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
247 if(PlottingSetup::RestrictToMassPeak) {
248 syserr+= ((zosofp)*(zosofp) + (zosofn)*(zosofn))*(1.0/9)*emuncert*emuncert;//sys err from emu method
249 syserr+= ((sbossfp)*(sbossfp)+(sbossfn)*(sbossfn))*(1.0/9)*eemmsidebanduncert*eemmsidebanduncert; // sys err from eemm sidebands
250 syserr+= ((sbosofp)*(sbosofp)+(sbosofn)*(sbosofn))*(1.0/9)*emsidebanduncert*emsidebanduncert; // sys err from emu sidebands
251 } else {
252 syserr+= ((zosofp)*(zosofp) + (zosofn)*(zosofn))*emuncert*emuncert;//sys err from emu method
253 }
254
255 syserr=TMath::Sqrt(syserr);
256
257 flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
258 if(PlottingSetup::RestrictToMassPeak) {
259 staterr=TMath::Sqrt(zossfn + (1.0/9)*(zosofp+zosofn)+ (1.0/9)*(sbossfp+sbossfn)+ (1.0/9)*(sbosofp+sbosofn));
260 } else {
261 staterr=TMath::Sqrt(zossfn + (zosofp+zosofn));
262 }
263 if(dopoisson) advanced_poisson(zossfn,zosofp,zosofn,sbossfp,sbossfn,sbosofp,sbosofn,poissonstaterrdown,poissonstaterrup);
264
265 float e_to_emu=0.5;
266 float m_to_emu=0.5;
267
268 if(writeanything) {
269 if(mcordata==mc) dout << " MC :: ";
270 if(mcordata==mcwithsignal) dout << " MC with S :: ";
271 if(mcordata==data) dout << " ";
272 dout << "Observed : " << zossfp << "+-" << zossfperr << endl;
273 if(mcordata==data) dout << " Composition: " << zossfpee << " (ee), " << zossfpmm << " (mm) " << endl;
274
275 if(mcordata==mc) dout << " MC :: ";
276 if(mcordata==mcwithsignal) dout << " MC with S :: ";
277 if(mcordata==data) dout << " ";
278 flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
279 if(PlottingSetup::RestrictToMassPeak) {
280 dout << "Predicted: " << zossfn << " + (1/3)*(" << zosofp << "-" << zosofn<<") [e&mu]+ (1/3)*(" << sbossfp << "-" << sbossfn<<") [SF,SB]+ (1/3)*(" << sbosofp << "-" << sbosofn<<") [OF,SB] = ";
281 dout << zossfn + (1.0/3)*(zosofp-zosofn)+ (1.0/3)*(sbossfp-sbossfn)+ (1.0/3)*(sbosofp-sbosofn)<< " +/- " << peakerr << " (peak ) " << " +/- " << syserr << " (sys) +/- " << staterr << " (stat)" << endl;
282 if(dopoisson) dout << " stat error with Poisson : +" << poissonstaterrup << " - " << poissonstaterrdown << endl;
283 } else {
284 dout << "Predicted: " << zossfn << " + (" << zosofp << "-" << zosofn<<") [e&mu]= ";
285 dout << zossfn + (zosofp-zosofn) << " +/- " << peakerr << " (peak ) " << " +/- " << syserr << " (sys) +/- " << staterr << " (stat)" << endl;
286 if(dopoisson) dout << " stat error with Poisson : +" << poissonstaterrup << " - " << poissonstaterrdown << endl;
287 }
288
289
290 if(mcordata==data) {
291 flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
292 if(PlottingSetup::RestrictToMassPeak) {
293 dout << " Composition: (ee)" << zossfnee << " + (1/3)*(" << e_to_emu*zosofp << "-" << e_to_emu*zosofn<<") ["<<e_to_emu<<"*e&mu]+ (1/3)*(" << sbossfpee << "-" << sbossfnee<<") [SF,SB]+ (1/3)*(" << sbosofp*e_to_emu << "-" << sbosofn*e_to_emu<<") ["<<e_to_emu<<"*OF,SB] = ";
294 dout << zossfnee + (1.0/3)*(e_to_emu*zosofp-e_to_emu*zosofn)+ (1.0/3)*(sbossfpee-sbossfnee) + (1.0/3)*(sbosofp*e_to_emu - sbosofn*e_to_emu) << endl;
295 dout << " (mm)" << zossfnmm << " + (1/3)*(" << m_to_emu*zosofp << "-" << m_to_emu*zosofn<<") ["<<m_to_emu<<"*e&mu]+ (1/3)*(" << sbossfpmm << "-" << sbossfnmm<<") [SF,SB]+ (1/3)*(" << sbosofp*m_to_emu << "-" << sbosofn*m_to_emu<<") ["<<m_to_emu<<"*OF,SB] = ";
296 dout << zossfnmm + (1.0/3)*(m_to_emu*zosofp - m_to_emu*zosofn) + (1.0/3)*(sbossfpmm -sbossfnmm)+ (1.0/3)*(sbosofp*m_to_emu-sbosofn*e_to_emu) << endl;
297 } else {
298 dout << " Composition: (ee)" << zossfnee << " + (" << e_to_emu*zosofp << "-" << e_to_emu*zosofn<<") ["<<e_to_emu<<"*e&mu] = ";
299 dout << zossfnee + (e_to_emu*zosofp-e_to_emu*zosofn) << endl;
300 dout << " (mm)" << zossfnmm << " + (" << m_to_emu*zosofp << "-" << m_to_emu*zosofn<<") ["<<m_to_emu<<"*e&mu] = ";
301 dout << zossfnmm + (m_to_emu*zosofp - m_to_emu*zosofn) << endl;
302 }
303 }
304
305 if(chatty) {
306 dout << " Pred(ZJets ) \t " << zossfn << endl;
307 dout << " Pred(e&mu;]) \t " << zosofp << "-" << zosofn << " = " << zosofp-zosofn<<endl;
308 flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
309 if(PlottingSetup::RestrictToMassPeak) {
310 dout << " Pred(ossf,sb]) \t " << sbossfp << "-" << sbossfn<<" = "<<sbossfp-sbossfn<<endl;
311 dout << " Pred(osof,sb]) \t " << sbosofp << "-" << sbosofn<<" = "<<sbosofp-sbosofn<<endl;
312 }
313 }
314 }//end of writeanything
315
316 if(mcordata==data) {
317 //store the result!
318 Nobs.push_back(zossfp);
319 flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
320 if(PlottingSetup::RestrictToMassPeak) Npred.push_back(zossfn + (1.0/3)*(zosofp-zosofn)+ (1.0/3)*(sbossfp-sbossfn)+ (1.0/3)*(sbosofp-sbosofn));
321 else Npred.push_back(zossfn + (zosofp-zosofn));
322 float totprederr=0;
323 totprederr+=peakerr*peakerr;
324 totprederr+=syserr*syserr;
325 totprederr+=staterr*staterr;
326 totprederr=TMath::Sqrt(totprederr);
327 Nprederr.push_back(totprederr);
328 }
329 vector<float> resultvector;
330 flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak
331 if(PlottingSetup::RestrictToMassPeak) resultvector.push_back(zossfn + (1.0/3)*(zosofp-zosofn)+ (1.0/3)*(sbossfp-sbossfn)+ (1.0/3)*(sbosofp-sbosofn));
332 else resultvector.push_back(zossfn + (zosofp-zosofn));
333
334 resultvector.push_back(TMath::Sqrt(peakerr*peakerr+syserr*syserr+staterr*staterr));
335 resultvector.push_back(zossfp);
336 resultvector.push_back(zossfperr);
337
338 return resultvector;
339 }
340
341 void get_result_above_one_fixed_jzb_value(float cut ,string mcjzb,string datajzb, int mcordata,float jzbpeakerrorMC, float jzbpeakerrorData, TCanvas *rescan, bool chatty=false, bool dopoisson=false) {
342 get_result_between_two_fixed_jzb_values(cut ,15000, mcjzb,datajzb, mcordata,jzbpeakerrorMC, jzbpeakerrorData, rescan, chatty, dopoisson);
343 }
344
345
346 void get_result(string mcjzb, string datajzb, float jzbpeakerrordata, float jzbpeakerrormc, vector<float> jzbcuts, bool chatty=false, bool dopoisson=false,bool doquick=false) {
347 TCanvas *rescan = new TCanvas("rescan","Result Canvas");
348 for(int icut=0;icut<jzbcuts.size();icut++) {
349 get_result_above_one_fixed_jzb_value(jzbcuts[icut],mcjzb,datajzb,data,jzbpeakerrormc,jzbpeakerrordata,rescan,chatty,dopoisson);
350 if(!doquick) {
351 get_result_above_one_fixed_jzb_value(jzbcuts[icut],mcjzb,datajzb,mc,jzbpeakerrormc,jzbpeakerrordata,rescan,chatty,dopoisson);
352 get_result_above_one_fixed_jzb_value(jzbcuts[icut],mcjzb,datajzb,mcwithsignal,jzbpeakerrormc,jzbpeakerrordata,rescan,chatty,dopoisson);
353 }
354 }
355 delete rescan;
356 }