ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/Development/Plotting/Modules/DustyModule.C
Revision: 1.2
Committed: Wed Sep 5 20:17:20 2012 UTC (12 years, 8 months ago) by buchmann
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Changes since 1.1: +195 -0 lines
Log Message:
Added experimental met z prediction to dusty module; the refined version will be uploaded to MetPlotting in a couple of seconds

File Contents

# Content
1 /****
2
3 Off peak status (RestrictToMassPeak) :
4
5 THIS MODULE WILL NOT BE ADAPTED AS IT ONLY CONTAINS OUTDATED FUNCTIONS!
6
7
8 . Necessary adaptations identified
9 . Started working on necessary adaptations
10 . Necessary adaptations implemented
11 . Necessary adaptations tested
12
13
14
15 ****/
16
17 /*
18
19 Contains functions that were used in the past but are out now ...
20
21 */
22 #include <iostream>
23 #include <vector>
24 #include <sys/stat.h>
25
26 #include <TCut.h>
27 #include <TROOT.h>
28 #include <TCanvas.h>
29 #include <TMath.h>
30 #include <TColor.h>
31 #include <TPaveText.h>
32 #include <TRandom.h>
33 #include <TH1.h>
34 #include <TH2.h>
35 #include <TF1.h>
36 #include <TSQLResult.h>
37 #include <TProfile.h>
38
39 string centralupdownname(int num) {
40 if (num==0) return "central";
41 if (num==1) return "down";
42 if (num==2) return "up";
43 return "centralupdownnameERROR";
44 }
45
46 string dataormc(int isdata) {
47 if(isdata) return "Data";
48 else return "MC";
49 }
50
51 int histocounter=0;
52 string give_histo_number(string basename="", int isdata=-1, int centralcounter=-1) {
53 histocounter++;
54 stringstream histocounternum;
55 if(isdata==-1) {
56 if(basename=="") histocounternum << histocounter;
57 else histocounternum << basename << "_" << histocounter;
58 }
59 else {
60 histocounternum << basename << "_" << dataormc(isdata) << "_" << centralupdownname(centralcounter);
61 }
62 return histocounternum.str();
63 }
64
65 int central=0;
66 int down=1;
67 int up=2;
68
69 void crunch_the_numbers(TH1F *RcorrJZBeemmop[3][2],TH1F *RcorrJZBeeop[3][2],TH1F *RcorrJZBmmop[3][2],TH1F *LcorrJZBeemmop[3][2],TH1F *LcorrJZBeeop[3][2],TH1F *LcorrJZBmmop[3][2],TH1F *RcorrJZBemop[3][2],TH1F *LcorrJZBemop[3][2],float zjetsestimate[3][2][20],float ttbarestimate[3][2][2][20],float predicted[3][2][20],float observed[3][2][20],int isdata,vector<float> jzbcuts,float eemmmumu[2][3][20]) {
70 // dout << "Crunching the numbers for is_data=" << isdata << endl;
71 for(int k=0;k<20;k++) {for(int i=0;i<2;i++){for(int j=0;j<3;j++) {eemmmumu[i][j][k]=0;}}}
72 for(int icut=0;icut<jzbcuts.size()-1;icut++) {//do calculation for each JZB cut
73 for(int ibin=1;ibin<jzbcuts.size();ibin++) {//and to do that, we need to consider each bin.
74 if(RcorrJZBeemmop[0][isdata]->GetBinLowEdge(ibin)<jzbcuts[icut]) continue;
75 // if(icut==0) dout << "Predicted+=" << LcorrJZBeemmop[central][isdata]->GetBinContent(ibin) << " eemm,L" << endl;
76 predicted[central][isdata][icut]+=LcorrJZBeemmop[central][isdata]->GetBinContent(ibin);
77 predicted[up][isdata][icut]+=LcorrJZBeemmop[up][isdata]->GetBinContent(ibin);
78 predicted[down][isdata][icut]+=LcorrJZBeemmop[down][isdata]->GetBinContent(ibin);
79
80 // if(icut==0) dout << "Predicted+=" << RcorrJZBemop[central][isdata]->GetBinContent(ibin) << " em, R" << endl;
81 predicted[central][isdata][icut]+=RcorrJZBemop[central][isdata]->GetBinContent(ibin);
82 predicted[up][isdata][icut]+=RcorrJZBemop[up][isdata]->GetBinContent(ibin);
83 predicted[down][isdata][icut]+=RcorrJZBemop[down][isdata]->GetBinContent(ibin);
84
85 // if(icut==0) dout << "Predicted-=" << LcorrJZBemop[central][isdata]->GetBinContent(ibin) <<" em, L" << endl;
86 predicted[central][isdata][icut]-=LcorrJZBemop[central][isdata]->GetBinContent(ibin);
87 predicted[up][isdata][icut]-=LcorrJZBemop[up][isdata]->GetBinContent(ibin);
88 predicted[down][isdata][icut]-=LcorrJZBemop[down][isdata]->GetBinContent(ibin);
89
90 // if(icut==0) dout << "Observed=" << RcorrJZBeemmop[central][isdata]->GetBinContent(ibin) << "eemm, R" << endl;
91 observed[central][isdata][icut]+=RcorrJZBeemmop[central][isdata]->GetBinContent(ibin);
92 observed[up][isdata][icut]+=RcorrJZBeemmop[up][isdata]->GetBinContent(ibin);
93 observed[down][isdata][icut]+=RcorrJZBeemmop[down][isdata]->GetBinContent(ibin);
94
95 zjetsestimate[central][isdata][icut]+=LcorrJZBeemmop[central][isdata]->GetBinContent(ibin);
96 ttbarestimate[central][isdata][0][icut]+=RcorrJZBemop[central][isdata]->GetBinContent(ibin);
97 ttbarestimate[central][isdata][1][icut]+=LcorrJZBemop[central][isdata]->GetBinContent(ibin);
98
99 eemmmumu[0][0][icut]+=LcorrJZBeeop[central][isdata]->GetBinContent(ibin);
100 eemmmumu[0][1][icut]+=LcorrJZBmmop[central][isdata]->GetBinContent(ibin);
101 eemmmumu[0][2][icut]+=LcorrJZBemop[central][isdata]->GetBinContent(ibin);
102
103 eemmmumu[1][0][icut]+=RcorrJZBeeop[central][isdata]->GetBinContent(ibin);
104 eemmmumu[1][1][icut]+=RcorrJZBmmop[central][isdata]->GetBinContent(ibin);
105 eemmmumu[1][2][icut]+=RcorrJZBemop[central][isdata]->GetBinContent(ibin);
106
107 }//end of ibin loop
108 }//end of icut loop
109 }
110
111 void make_table(float eemmmumu[2][3][20],int icut,float JZBcut) {
112 vector< vector <string> > entries;
113 vector <string> line;
114 line.push_back("");
115 line.push_back("eemm (ee/mm)");
116 line.push_back("em");
117 entries.push_back(line);
118 line.clear();
119 line.push_back("JZB>"+any2string(JZBcut));
120 line.push_back(any2string(eemmmumu[1][0][icut]+eemmmumu[1][1][icut])+"("+any2string(eemmmumu[1][0][icut])+"/"+any2string(eemmmumu[1][1][icut])+")");
121 line.push_back(any2string(eemmmumu[1][2][icut]));
122 entries.push_back(line);
123 line.clear();
124 line.push_back("JZB<-"+any2string(JZBcut));
125 line.push_back(any2string(eemmmumu[0][0][icut]+eemmmumu[0][1][icut])+"("+any2string(eemmmumu[0][0][icut])+"/"+any2string(eemmmumu[0][1][icut])+")");
126 line.push_back(any2string(eemmmumu[0][2][icut]));
127 entries.push_back(line);
128 make_nice_jzb_table(entries);
129 /*
130 dout << " \t\t | \t eemm (ee/mm) \t | \t em" << endl;
131 dout << "JZB>peak\t | \t " << eemmmumu[1][0]+eemmmumu[1][1] << "(" << eemmmumu[1][0] << "/" << eemmmumu[1][1] << ")\t | \t" << eemmmumu[1][2] << endl;
132 dout << "JZB<peak\t | \t " << eemmmumu[0][0]+eemmmumu[0][1] << "(" << eemmmumu[0][0] << "/" << eemmmumu[0][1] << ")\t | \t" << eemmmumu[0][2] << endl;
133 dout << endl;
134 */
135 }
136
137 void present_result(vector<float> &jzbcuts,float predicted[3][2][20],float observed[3][2][20],float zjetsestimate[3][2][20],float ttbarestimate[3][2][2][20],int icut, float eemmmumu[2][3][20]) {
138
139 //blublu
140 dout << endl << endl;
141 dout << "--------------------------------------------------------------------------------------" << endl;
142 dout << "DATA: " << endl;
143 dout << " For JZB>" << jzbcuts[icut];
144 float max,downvar,upvar;
145 stringstream printtitle;
146 printtitle << "JZB>" << jzbcuts[icut] << " (data)";
147 ComputePoissonError(zjetsestimate[central][data][icut],ttbarestimate[central][data][0][icut],ttbarestimate[central][data][1][icut],max,downvar,upvar,printtitle.str());
148 //dout << " Predicted: " << print_range(predicted[central][data][icut],predicted[up][data][icut],predicted[down][data][icut]) << " (stat) +" << 0.2*zjetsestimate[central][data][icut] << "-" << 0.4*zjetsestimate[central][data][icut] << " (sys) " << " using Poisson: " << print_range(max,max+upvar,max-downvar) << ")" << endl;
149 float pred=predicted[central][data][icut];
150 float sysP=abs(predicted[up][data][icut]-pred);
151 float sysN=abs(pred-predicted[down][data][icut]);
152 sysN = sysN + pred*(1-fitresultconstdata);//fitresultconst comes from the fit in the 0-30 GeV range!
153 dout << " Predicted: " << pred << "+" << statErrorP(pred) << "-" << statErrorN(pred) << " (stat) +" << sysP << " - " << sysN << " (sys) " << " using Poisson: " << print_range(max,max+upvar,max-downvar) << ")" << endl;
154 //0.2*zjetsestimate[central][data][icut] << "-" << 0.4*zjetsestimate[central][data][icut] << " (sys) " << " using Poisson: " << print_range(max,max+upvar,max-downvar) << ")" << endl;
155 dout << " Details: ZJetsEstimate= " << zjetsestimate[central][data][icut] << ", TTbar estimate=" << ttbarestimate[central][data][0][icut]-ttbarestimate[central][data][1][icut] << " [" << ttbarestimate[central][data][0][icut] << " R -" << ttbarestimate[central][data][1][icut] << " L ]" << endl;
156 dout << " For JZB>" << jzbcuts[icut] << " Observed: " << observed[central][data][icut] << endl;
157 dout << "TABLE:" << endl;
158 make_table(eemmmumu,icut,jzbcuts[icut]);
159 dout << "MC: " << endl;
160 printtitle.str("");
161 printtitle<<"JZB>"<<jzbcuts[icut]<<" (MC)";
162 ComputePoissonError(zjetsestimate[central][mc][icut],ttbarestimate[central][mc][0][icut],ttbarestimate[central][mc][1][icut],max,downvar,upvar,printtitle.str());
163 pred=predicted[central][mc][icut];
164 sysP=abs(predicted[up][mc][icut]-pred);
165 sysN=abs(pred-predicted[down][mc][icut]);
166 sysN = sysN + pred*(1-fitresultconstmc);//fitresultconst comes from the fit in the 0-30 GeV range!
167 dout << " Predicted: " << pred << "+" << statErrorP(pred) << "-" << statErrorN(pred) << " (stat) +" << sysP << " - " << sysN << " (sys) " << " using Poisson: " << print_range(max,max+upvar,max-downvar) << ")" << endl;
168 // dout << " Predicted: " << print_range(pred,statErrorP(pred),statErrorN(pred)) << " (stat) +" << sysP << " - " << sysN << " (sys) " << " using Poisson: " << print_range(max,max+upvar,max-downvar) << ")" << endl;
169 //dout << " Predicted: " << print_range(predicted[central][mc][icut],predicted[down][mc][icut],predicted[up][mc][icut]) << " (stat) +" << 0.2*zjetsestimate[central][mc][icut] << "-" << 0.4*zjetsestimate[central][mc][icut] << " (sys) " << " using Poisson: " << print_range(max,max+upvar,max-downvar) << ")" << endl;
170 //print_range(predicted[central][mc][icut],predicted[central][mc][icut]+upvar,predicted[central][mc][icut]-downvar) << ")" << endl;
171 dout << " Details: ZJetsEstimate= " << zjetsestimate[central][mc][icut] << ", TTbar estimate=" << ttbarestimate[central][mc][0][icut]-ttbarestimate[central][mc][1][icut] << " [" << ttbarestimate[central][mc][0][icut] << " R -" << ttbarestimate[central][mc][1][icut] << " L ] )" << endl;
172 dout << " For JZB>" << jzbcuts[icut] << " Observed: " << observed[central][mc][icut] << endl;
173 dout << endl << endl;
174 }
175
176 void calculate_predicted_and_observed_eemm(float MCPeak,float MCPeakError,float DataPeak,float DataPeakError,vector<float> jzbcuts) { // DO NOT CHANGE THIS TO ... &jzbcuts !!! we add one!
177 jzbcuts.push_back(14000);
178
179 /*
180 * We want the following numbers: For all JZB cuts, we want: ee,mm,em; observed, predicted. we want final results with errors from the peak and statistical error.
181 * How to accomplish this; we draw histos for MC&data, once with the peak (to extract obs/pred), once with peak+sigma, once with peak-sigma. this gives us the peak error.
182 * Statistical error: Comes from Poissonian errors ... for this we have a special little macro that will give us the value for three given parameters (cool, right?)
183 */
184
185 string jzbpeak[3][2];
186 jzbpeak[central][data]=give_jzb_expression(DataPeak,data);
187 jzbpeak[up][data]=give_jzb_expression(DataPeak+DataPeakError,data);
188 jzbpeak[down][data]=give_jzb_expression(DataPeak-DataPeakError,data);
189 jzbpeak[central][mc]=give_jzb_expression(MCPeak,mc);
190 jzbpeak[up][mc]=give_jzb_expression(MCPeak+MCPeakError,mc);
191 jzbpeak[down][mc]=give_jzb_expression(MCPeak-MCPeakError,mc);
192
193 TH1F *RcorrJZBeemmop[3][2];
194 TH1F *LcorrJZBeemmop[3][2];
195 TH1F *RcorrJZBeeop[3][2];
196 TH1F *LcorrJZBeeop[3][2];
197 TH1F *RcorrJZBmmop[3][2];
198 TH1F *LcorrJZBmmop[3][2];
199 TH1F *RcorrJZBemop[3][2];
200 TH1F *LcorrJZBemop[3][2];
201 float zjetsestimate[3][2][20];
202 float predicted[3][2][20];
203 float observed[3][2][20];
204 float ttbarestimate[3][2][2][20];
205 float eemmmumu[2][3][20];
206
207
208 for(int isdata=0;isdata<=1;isdata++) {
209 for(int centraldownup=0;centraldownup<=2;centraldownup++) {
210 RcorrJZBeemmop[centraldownup][isdata] = allsamples.Draw(give_histo_number("RcorrJZBeemmop",isdata,centraldownup),jzbpeak[centraldownup][isdata].c_str(),jzbcuts, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,isdata, luminosity);
211 LcorrJZBeemmop[centraldownup][isdata] = allsamples.Draw(give_histo_number("LcorrJZBeemmop",isdata,centraldownup),("-"+jzbpeak[centraldownup][isdata]).c_str(),jzbcuts, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,isdata, luminosity);
212
213 RcorrJZBemop[centraldownup][isdata] = allsamples.Draw(give_histo_number("RcorrJZBemop",isdata,centraldownup),jzbpeak[centraldownup][isdata].c_str(),jzbcuts, "JZB [GeV]", "events", cutmass&&cutOSOF&&cutnJets,isdata, luminosity);
214 LcorrJZBemop[centraldownup][isdata] = allsamples.Draw(give_histo_number("LcorrJZBemop",isdata,centraldownup),("-"+jzbpeak[centraldownup][isdata]).c_str(),jzbcuts, "JZB [GeV]", "events", cutmass&&cutOSOF&&cutnJets,isdata,luminosity);
215 RcorrJZBemop[centraldownup][isdata]->Add(allsamples.Draw(give_histo_number("RcorrJZBemopSB",isdata,centraldownup),jzbpeak[centraldownup][isdata].c_str(),jzbcuts, "JZB [GeV]", "events", sidebandcut&&cutOSOF&&cutnJets,isdata, luminosity));
216 LcorrJZBemop[centraldownup][isdata]->Add(allsamples.Draw(give_histo_number("LcorrJZBemopSB",isdata,centraldownup),("-"+jzbpeak[centraldownup][isdata]).c_str(),jzbcuts, "JZB [GeV]", "events", sidebandcut&&cutOSOF&&cutnJets,isdata,luminosity));
217 RcorrJZBemop[centraldownup][isdata]->Add(allsamples.Draw(give_histo_number("RcorrJZBeemmopSB",isdata,centraldownup),jzbpeak[centraldownup][isdata].c_str(),jzbcuts, "JZB [GeV]", "events", sidebandcut&&cutOSSF&&cutnJets,isdata, luminosity));
218 LcorrJZBemop[centraldownup][isdata]->Add(allsamples.Draw(give_histo_number("LcorrJZBeemmopSB",isdata,centraldownup),("-"+jzbpeak[centraldownup][isdata]).c_str(),jzbcuts, "JZB [GeV]", "events", sidebandcut&&cutOSSF&&cutnJets,isdata,luminosity));
219 LcorrJZBemop[centraldownup][isdata]->Scale(1.0/3);
220 RcorrJZBemop[centraldownup][isdata]->Scale(1.0/3);
221 RcorrJZBeeop[centraldownup][isdata] = allsamples.Draw(give_histo_number("RcorrJZBeeop",isdata,centraldownup),jzbpeak[centraldownup][isdata].c_str(),jzbcuts, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets&&"id1==0",isdata, luminosity);
222 LcorrJZBeeop[centraldownup][isdata] = allsamples.Draw(give_histo_number("LcorrJZBeeop",isdata,centraldownup),("-"+jzbpeak[centraldownup][isdata]).c_str(),jzbcuts, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets&&"id1==0",isdata, luminosity);
223 RcorrJZBmmop[centraldownup][isdata] = allsamples.Draw(give_histo_number("RcorrJZBmmop",isdata,centraldownup),jzbpeak[centraldownup][isdata].c_str(),jzbcuts, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets&&"id1==1",isdata, luminosity);
224 LcorrJZBmmop[centraldownup][isdata] = allsamples.Draw(give_histo_number("LcorrJZBmmop",isdata,centraldownup),("-"+jzbpeak[centraldownup][isdata]).c_str(),jzbcuts, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets&&"id1==1",isdata, luminosity);
225 for(int icut=0;icut<=jzbcuts.size();icut++) {
226 zjetsestimate[centraldownup][isdata][icut]=0;
227 predicted[centraldownup][isdata][icut]=0;
228 observed[centraldownup][isdata][icut]=0;
229 ttbarestimate[centraldownup][isdata][0][icut]=0;
230 ttbarestimate[centraldownup][isdata][1][icut]=0;
231 }
232 }//end of for central,up,down loop
233 crunch_the_numbers(RcorrJZBeemmop,RcorrJZBeeop,RcorrJZBmmop,LcorrJZBeemmop,LcorrJZBeeop,LcorrJZBmmop,RcorrJZBemop,LcorrJZBemop,zjetsestimate,ttbarestimate,predicted,observed,isdata,jzbcuts,eemmmumu);
234
235 // void crunch_the_numbers(TH1F *RcorrJZBeemmop[3][2],TH1F *RcorrJZBeeop[3][2],TH1F *RcorrJZBmmop[3][2],TH1F *LcorrJZBeemmop[3][2],TH1F *LcorrJZBeeop[3][2],TH1F *LcorrJZBmmop[3][2],TH1F *RcorrJZBemop[3][2],TH1F *LcorrJZBemop[3][2],float zjetsestimate[3][2][20],float ttbarestimate[3][2][2][20],float predicted[3][2][20],float observed[3][2][20],int isdata,vector<float> jzbcuts,float eemmmumu[2][3]) {
236
237 }//end of is data loop
238
239
240 dout << "We obtain the following results: " << endl;
241 for(int icut=0;icut<jzbcuts.size()-1;icut++) {
242 present_result(jzbcuts,predicted,observed,zjetsestimate,ttbarestimate,icut,eemmmumu);
243 }
244
245
246 }
247
248
249
250
251
252
253
254
255
256
257
258
259
260 void alpha_scan_susy_space(string mcjzb, string datajzb) {
261 TCanvas *c3 = new TCanvas("c3","c3");
262 vector<float> binning;
263 binning=allsamples.get_optimal_binsize(mcjzb,cutmass&&cutOSSF&&cutnJets,20,50,800);
264 float arrbinning[binning.size()];
265 for(int i=0;i<binning.size();i++) arrbinning[i]=binning[i];
266 TH1F *puredata = allsamples.Draw("puredata", datajzb,binning, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,data,luminosity);
267 puredata->SetMarkerSize(DataMarkerSize);
268 TH1F *allbgs = allsamples.Draw("allbgs", "-"+datajzb,binning, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,data,luminosity);
269 TH1F *allbgsb = allsamples.Draw("allbgsb", "-"+datajzb,binning, "JZB [GeV]", "events", cutmass&&cutOSOF&&cutnJets,data,luminosity);
270 TH1F *allbgsc = allsamples.Draw("allbgsc", datajzb,binning, "JZB [GeV]", "events", cutmass&&cutOSOF&&cutnJets,data,luminosity);
271 allbgs->Add(allbgsb,-1);
272 allbgs->Add(allbgsc);
273 int ndata=puredata->Integral();
274 ofstream myfile;
275 myfile.open ("susyscan_log.txt");
276 TFile *susyscanfile = new TFile("/scratch/fronga/SMS/T5z_SqSqToQZQZ_38xFall10.root");
277 TTree *suevents = (TTree*)susyscanfile->Get("events");
278 TH2F *exclusionmap = new TH2F("exclusionmap","",20,0,500,20,0,1000);
279 TH2F *exclusionmap1s = new TH2F("exclusionmap1s","",20,0,500,20,0,1000);
280 TH2F *exclusionmap2s = new TH2F("exclusionmap2s","",20,0,500,20,0,1000);
281 TH2F *exclusionmap3s = new TH2F("exclusionmap3s","",20,0,500,20,0,1000);
282
283 susy_scan_axis_labeling(exclusionmap);
284 susy_scan_axis_labeling(exclusionmap1s);
285 susy_scan_axis_labeling(exclusionmap2s);
286 susy_scan_axis_labeling(exclusionmap3s);
287
288 Int_t MyPalette[100];
289 Double_t r[] = {0., 0.0, 1.0, 1.0, 1.0};
290 Double_t g[] = {0., 0.0, 0.0, 1.0, 1.0};
291 Double_t b[] = {0., 1.0, 0.0, 0.0, 1.0};
292 Double_t stop[] = {0., .25, .50, .75, 1.0};
293 Int_t FI = TColor::CreateGradientColorTable(5, stop, r, g, b, 100);
294 for (int i=0;i<100;i++) MyPalette[i] = FI+i;
295
296 gStyle->SetPalette(100, MyPalette);
297
298 for(int m23=50;m23<500;m23+=25) {
299 for (int m0=(2*(m23-50)+150);m0<=1000;m0+=50)
300 {
301 c3->cd();
302 stringstream drawcondition;
303 drawcondition << "pfJetGoodNum>=3&&(TMath::Abs(masses[0]-"<<m0<<")<10&&TMath::Abs(masses[2]-masses[3]-"<<m23<<")<10)&&mll>5&&id1==id2";
304 TH1F *puresignal = new TH1F("puresignal","puresignal",binning.size()-1,arrbinning);
305 TH1F *puresignall= new TH1F("puresignall","puresignal",binning.size()-1,arrbinning);
306 stringstream drawvar,drawvar2;
307 drawvar<<mcjzb<<">>puresignal";
308 drawvar2<<"-"<<mcjzb<<">>puresignall";
309 suevents->Draw(drawvar.str().c_str(),drawcondition.str().c_str());
310 suevents->Draw(drawvar2.str().c_str(),drawcondition.str().c_str());
311 if(puresignal->Integral()<60) {
312 delete puresignal;
313 continue;
314 }
315 puresignal->Add(puresignall,-1);//we need to correct for the signal contamination - we effectively only see (JZB>0)-(JZB<0) !!
316 puresignal->Scale(ndata/(20*puresignal->Integral()));//normalizing it to 5% of the data
317 stringstream saveas;
318 saveas<<"Model_Scan/m0_"<<m0<<"__m23_"<<m23;
319 dout << "PLEASE KEEP IN MIND THAT SIGNAL CONTAMINATION IS NOT REALLY TAKEN CARE OF YET DUE TO LOW STATISTICS! SHOULD BE SOMETHING LIKE THIS : "<< endl;
320 // TH1F *signalpredlo = allsamples.Draw("signalpredlo", "-"+mcjzb, binning, "JZB [GeV]", "events", cutmass&&cutOSOF&&cutnJets,mc, luminosity,allsamples.FindSample("LM4"));
321 // TH1F *signalpredro = allsamples.Draw("signalpredro", mcjzb, binning, "JZB [GeV]", "events", cutmass&&cutOSOF&&cutnJets,mc, luminosity,allsamples.FindSample("LM4"));
322 // TH1F *puredata = allsamples.Draw("puredata", datajzb,binning, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,data,luminosity);
323 // signalpred->Add(signalpredlo,-1);
324 // signalpred->Add(signalpredro);
325 // puresignal->Add(signalpred,-1);//subtracting signal contamination
326 //---------------------
327 // dout << "(m0,m23)=("<<m0<<","<<m23<<") contains " << puresignal->Integral() << endl;
328 // TH1F *puresignal = allsamples.Draw("puresignal",mcjzb, binning, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,mc, luminosity,allsamples.FindSample("LM4"));
329 vector<float> results=establish_upper_limits(puredata,allbgs,puresignal,saveas.str(),myfile);
330 if(results.size()==0) {
331 delete puresignal;
332 continue;
333 }
334 exclusionmap->Fill(m23,m0,results[0]);
335 exclusionmap1s->Fill(m23,m0,results[1]);
336 exclusionmap2s->Fill(m23,m0,results[2]);
337 exclusionmap3s->Fill(m23,m0,results[3]);
338 delete puresignal;
339 dout << "(m0,m23)=("<<m0<<","<<m23<<") : 3 sigma at " << results[3] << endl;
340 }
341 }//end of model scan for loop
342
343 dout << "Exclusion Map contains" << exclusionmap->Integral() << " (integral) and entries: " << exclusionmap->GetEntries() << endl;
344 c3->cd();
345 exclusionmap->Draw("CONTZ");
346 CompleteSave(c3,"Model_Scan/CONT/Model_Scan_Mean_values");
347 exclusionmap->Draw("COLZ");
348 CompleteSave(c3,"Model_Scan/COL/Model_Scan_Mean_values");
349
350 exclusionmap1s->Draw("CONTZ");
351 CompleteSave(c3,"Model_Scan/CONT/Model_Scan_1sigma_values");
352 exclusionmap1s->Draw("COLZ");
353 CompleteSave(c3,"Model_Scan/COL/Model_Scan_1sigma_values");
354
355 exclusionmap2s->Draw("CONTZ");
356 CompleteSave(c3,"Model_Scan/CONT/Model_Scan_2sigma_values");
357 exclusionmap2s->Draw("COLZ");
358 CompleteSave(c3,"Model_Scan/COL/Model_Scan_2sigma_values");
359
360 exclusionmap3s->Draw("CONTZ");
361 CompleteSave(c3,"Model_Scan/CONT/Model_Scan_3sigma_values");
362 exclusionmap3s->Draw("COLZ");
363 CompleteSave(c3,"Model_Scan/COL/Model_Scan_3sigma_values");
364
365 TFile *exclusion_limits = new TFile("exclusion_limits.root","RECREATE");
366 exclusionmap->Write();
367 exclusionmap1s->Write();
368 exclusionmap2s->Write();
369 exclusionmap3s->Write();
370 exclusion_limits->Close();
371 susyscanfile->Close();
372
373 myfile.close();
374 }
375
376
377
378
379 //-------------------------------------------------------------
380
381 /*
382
383
384 not in use anymore:
385
386 -------------------------------------<Spring11 - used up to Jul 14>----------------------------------------------
387
388 // allsamples.AddSample("/scratch/buchmann/AllData_Jun10___DCS_648pb_MoreTriggers4_DCS.root","Data",1,1,true,false,0,kBlack);
389 // allsamples.AddSample("/scratch/buchmann/AllData_Jun10___486pb_MoreTriggers3.root","Data",1,1,true,false,0,kBlack);
390 //allsamples.AddSample("/scratch/buchmann/AllData_Jun17_160404-166861_real_complete_2.root","Data",1,1,true,false,0,kBlack);
391 allsamples.AddSample("/scratch/buchmann/ntuples2011/data/AllData_Jul6_1078pb.root","Data",1,1,true,false,0,kBlack);
392
393 // allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/TTJets_TuneZ2_7TeV-madgraph-tauola.root","t#bar{t}",totEventsTTbar,TTbarCrossSection,false,false,1,ttbar_color);//spring11
394 allsamples.AddSample("/scratch/fronga/MC_v1.29/TTJets_TuneZ2_7TeV-madgraph-tauola_Summer11_3.root","TTbar",totEventsTTbar,TTbarCrossSection,false,false,1,ttbar_color);//summer11
395 // allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/WJetsToLNu_TuneZ2_7TeV-madgraph-tauola.root","WJets",totEventsWJets,WJetsCrossSection,false,false,3,wjets_color); //spring11
396 allsamples.AddSample("/scratch/fronga/MC_v1.29/WJetsToLNu_TuneZ2_7TeV-madgraph-tauola.root","WJets",totEventsWJets,WJetsCrossSection,false,false,3,wjets_color);//summer
397 allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/VVJetsTo4L_TuneD6T_7TeV-madgraph-tauola.root","DiBosons",totEventsVVJets,VVJetsCrossSection,false,false,4,diboson_color);
398 allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/TToBLNu_TuneZ2_s-channel_7TeV-madgraph_2.root","SingleTop",totEventsSingleTopS,SingleTopSCrossSection,false,false,2,singletop_color);
399 allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/TToBLNu_TuneZ2_t-channel_7TeV-madgraph.root","SingleTop",totEventsSingleTopT,SingleTopTCrossSection,false,false,2,singletop_color);
400 allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/TToBLNu_TuneZ2_tW-channel_7TeV-madgraph.root","SingleTop",totEventsSingleTopU,SingleTopUCrossSection,false,false,2,singletop_color);
401 allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/ZinvisibleJets_7TeV-madgraph.root","Z nunu",totEventsZnunu,ZnunuCrossSection,false,false,4,diboson_color);
402 // allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/QCD_TuneD6T_HT-100To250_7TeV-madgraph.root","QCD",totEventsQCD100to250,QCD100to250CrossSection,false,false,5,qcd_color);
403 // allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/QCD_TuneD6T_HT-250To500_7TeV-madgraph.root","QCD",totEventsQCD250to500,QCD250to500CrossSection,false,false,5,qcd_color);
404 // allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/QCD_TuneD6T_HT-500To1000_7TeV-madgraph.root","QCD",totEventsQCD500to1000,QCD500to1000CrossSection,false,false,5,qcd_color);
405 // allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/QCD_TuneD6T_HT-1000_7TeV-madgraph.root","QCD",totEventsQCD1000toInf,QCD1000toInfCrossSection,false,false,5,qcd_color);
406 // allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/DYJetsToLL_TuneZ2_M-50_7TeV-madgraph-tauola.root","ZJets",totEventsZjets,ZJetsCrossSection,false,false,6,dy_color);//spring11
407 allsamples.AddSample("/scratch/fronga/MC_v1.29/DYJetsToLL_TuneZ2_M-50_7TeV-madgraph-tauola_Summer11.root","Z+Jets",totEventsZjets,ZJetsCrossSection,false,false,6,dy_color);//summer11
408 allsamples.AddSample("/scratch/buchmann/MCSpring2011PU__inclPFvarss__PURW_160404_166861/LM4_SUSY_sftsht_7TeV-pythia6.root","LM4",218536,1.879,false,true,7,lm_color);
409
410 -------------------------------------</Spring11 - used up to Jul 14>----------------------------------------------
411
412
413
414
415
416 //how this works: AddSample("FilePath,"Sample name",Nevents,XS,is_data,is_signal,groupindex
417 allsamples.AddSample("/scratch/buchmann/AllDataNov4.root","Data",1,1,true,false,0,kBlack);
418 allsamples.AddSample("/scratch/buchmann/CutAnalysis/WJetsToLNu_TuneZ2_7TeV-madgraph-tauola.root","WJets to LNu",14350756,31314,false,false,2,kGreen);
419 // allsamples.AddSample("/scratch/buchmann/CutAnalysis/DYJetsToLL_TuneD6T_M-10To50_7TeV-madgraph-tauola.root","DY,M10To50",186891,400,false,false,3);
420 allsamples.AddSample("/scratch/buchmann/CutAnalysis/DYJetsToLL_TuneZ2_M-50_7TeV-madgraph-tauola.root","DY,M50",2552567,3048,false,false,3,kYellow);
421 allsamples.AddSample("/scratch/buchmann/CutAnalysis/TTJets_TuneZ2_7TeV-madgraph-tauola.root","TTJets",1144028,157.5,false,false,4,kMagenta-1);
422 allsamples.AddSample("/scratch/buchmann/CutAnalysis/TToBLNu_TuneZ2_s-channel_7TeV-madgraph.root","T to blnu,s",489472,1.49,false,false,5,kBlue);
423 allsamples.AddSample("/scratch/buchmann/CutAnalysis/TToBLNu_TuneZ2_t-channel_7TeV-madgraph.root","T to blnu,t",477610,20.54,false,false,5,kBlue);
424 allsamples.AddSample("/scratch/buchmann/CutAnalysis/TToBLNu_TuneZ2_tW-channel_7TeV-madgraph.root","T to blnu,u",477599,10.6,false,false,5,kBlue);
425 allsamples.AddSample("/scratch/buchmann/CutAnalysis/VVJetsTo4L_TuneD6T_7TeV-madgraph-tauola.root","VVJets to 4L",509072,4.8,false,false,1,kGreen);
426 allsamples.AddSample("/scratch/buchmann/CutAnalysis/LM4_SUSY_sftsht_7TeV-pythia6.root","20xLM4",218536,1.879*100,false,true,6,kMagenta);
427
428
429
430 /*
431 Int_t fermiwjets = TColor::GetColor("#75f075");
432 Int_t ferminonW = TColor::GetColor("#7d99d1");
433 Int_t fermiwwwz = TColor::GetColor("#2b26ac");
434 Int_t fermizjets = TColor::GetColor("#0044ff");
435 Int_t fermitop = TColor::GetColor("#353535");
436 Int_t fermihiggs = TColor::GetColor("#b30000");
437
438
439 dy_color = fermiwjets;
440 wjets_color = ferminonW;
441 ttbar_color = fermizjets;
442 singletop_color = fermitop;
443 qcd_color = nice_orange;
444 diboson_color = fermiwwwz;
445 lm_color = fermihiggs;
446 */
447
448
449 /*
450 Int_t dy_color = TColor::GetColor("#ff1300");
451 Int_t wjets_color = TColor::GetColor("#ff7a00");
452 Int_t ttbar_color = TColor::GetColor("#03899c");
453 Int_t singletop_color = TColor::GetColor("#00c322");
454 Int_t qcd_color = TColor::GetColor("#4a11ae");
455 Int_t diboson_color = TColor::GetColor("#8f04a8");
456 Int_t lm_color = TColor::GetColor("#dcf900");
457
458 //fugly set of colors:
459 dy_color = TColor::GetColor("#ff1300");
460 wjets_color = TColor::GetColor("#ff7a00");
461 ttbar_color = TColor::GetColor("#03899c");
462 singletop_color = TColor::GetColor("#00c322");
463 qcd_color = TColor::GetColor("#4a11ae");
464 diboson_color = TColor::GetColor("#8f04a8");
465 lm_color = TColor::GetColor("#dcf900");
466
467 //twilight
468 dy_color = TColor::GetColor("#141414");
469 wjets_color = TColor::GetColor("#f8f8f8");
470 ttbar_color = TColor::GetColor("#cf6a4c");
471 singletop_color = TColor::GetColor("#cda869");
472 qcd_color = TColor::GetColor("#e9c062");
473 diboson_color = TColor::GetColor("#5f5a60");
474 lm_color = TColor::GetColor("#7587a6");
475 //http://www.colorschemer.com/schemes/viewscheme.php?id=9126
476
477 */
478
479 /* //how this works: AddSample("FilePath,"Sample name",Nevents,XS,is_data,is_signal,groupindex
480 allsamples.AddSample("/scratch/fronga/JSON_160404-163869_v1.29/AllCertified191_v1.29.root","Data",1,1,true,false,0,kBlack);
481 allsamples.AddSample("/scratch/fronga/MC_v1.29/DYJetsToLL_TuneZ2_M-50_7TeV-madgraph-tauola.root","ZJets",totEventsZjets,ZJetsCrossSection,false,false,3,kYellow);
482 allsamples.AddSample("/scratch/fronga/MC_v1.29/WJetsToLNu_TuneZ2_7TeV-madgraph-tauola.root","WJets",totEventsWJets,WJetsCrossSection,false,false,2,kGray);
483 allsamples.AddSample("/scratch/fronga/MC_v1.29/TTJets_TuneZ2_7TeV-madgraph-tauola.root","t#bar{t}",totEventsTTbar,TTbarCrossSection,false,false,4,kMagenta-1);
484 allsamples.AddSample("/scratch/fronga/MC_v1.29/TToBLNu_TuneZ2_s-channel_7TeV-madgraph_2.root","SingleTop",totEventsSingleTopS,SingleTopSCrossSection,false,false,5,kBlue);
485 allsamples.AddSample("/scratch/fronga/MC_v1.29/TToBLNu_TuneZ2_t-channel_7TeV-madgraph.root","SingleTop",totEventsSingleTopT,SingleTopTCrossSection,false,false,5,kBlue);
486 allsamples.AddSample("/scratch/fronga/MC_v1.29/TToBLNu_TuneZ2_tW-channel_7TeV-madgraph.root","SingleTop",totEventsSingleTopU,SingleTopUCrossSection,false,false,5,kBlue);
487 allsamples.AddSample("/scratch/fronga/MC_v1.29/VVJetsTo4L_TuneD6T_7TeV-madgraph-tauola.root","Dibosons",totEventsVVJets,VVJetsCrossSection,false,false,1,kGreen);
488 allsamples.AddSample("/scratch/fronga/MC_v1.29/ZinvisibleJets_7TeV-madgraph.root","Z nunu",totEventsZnunu,ZnunuCrossSection,false,false,1,kGreen);
489 allsamples.AddSample("/scratch/fronga/MC_v1.29/QCD_TuneD6T_HT-100To250_7TeV-madgraph.root","QCD",totEventsQCD100to250,QCD100to250CrossSection,false,false,6,kRed);
490 allsamples.AddSample("/scratch/fronga/MC_v1.29/QCD_TuneD6T_HT-250To500_7TeV-madgraph.root","QCD",totEventsQCD250to500,QCD250to500CrossSection,false,false,6,kRed);
491 // allsamples.AddSample("/scratch/fronga/MC_v1.29/QCD_TuneD6T_HT-500To1000_7TeV-madgraph.root","QCD",totEventsQCD500to1000,QCD500to1000CrossSection,false,false,6,kRed);
492 allsamples.AddSample("/scratch/fronga/MC_v1.29/QCD_TuneD6T_HT-1000_7TeV-madgraph.root","QCD",totEventsQCD1000toInf,QCD1000toInfCrossSection,false,false,6,kRed);
493 allsamples.AddSample("/scratch/fronga/MC_v1.29/LM4_SUSY_sftsht_7TeV-pythia6.root","20xLM4",218536,1.879*20,false,true,7,kMagenta);
494 */
495
496 //how this works: AddSample("FilePath,"Sample name",Nevents,XS,is_data,is_signal,groupindex
497 // allsamples.AddSample("/scratch/fronga/JSON_160404-163869_v1.29/AllCertified191_v1.29.root","Data",1,1,true,false,0,kBlack);
498 // allsamples.AddSample("/scratch/buchmann/May2011ReRecoAllMayJSON.root","Data",1,1,true,false,0,kBlack);
499
500 /*
501 allsamples.AddSample("/scratch/buchmann/DataMay2011ReReco_CertifiedReRecoJSON.root","Data",1,1,true,false,0,kBlack);
502 allsamples.AddSample("/scratch/fronga/MC_v1.29/TTJets_TuneZ2_7TeV-madgraph-tauola.root","t#bar{t}",totEventsTTbar,TTbarCrossSection,false,false,1,ttbar_color);
503 allsamples.AddSample("/scratch/fronga/MC_v1.29/WJetsToLNu_TuneZ2_7TeV-madgraph-tauola.root","WJets",totEventsWJets,WJetsCrossSection,false,false,3,wjets_color);
504 allsamples.AddSample("/scratch/fronga/MC_v1.29/VVJetsTo4L_TuneD6T_7TeV-madgraph-tauola.root","Dibosons",totEventsVVJets,VVJetsCrossSection,false,false,4,diboson_color);
505 allsamples.AddSample("/scratch/fronga/MC_v1.29/TToBLNu_TuneZ2_s-channel_7TeV-madgraph_2.root","SingleTop",totEventsSingleTopS,SingleTopSCrossSection,false,false,2,singletop_color);
506 allsamples.AddSample("/scratch/fronga/MC_v1.29/TToBLNu_TuneZ2_t-channel_7TeV-madgraph.root","SingleTop",totEventsSingleTopT,SingleTopTCrossSection,false,false,2,singletop_color);
507 allsamples.AddSample("/scratch/fronga/MC_v1.29/TToBLNu_TuneZ2_tW-channel_7TeV-madgraph.root","SingleTop",totEventsSingleTopU,SingleTopUCrossSection,false,false,2,singletop_color);
508 allsamples.AddSample("/scratch/fronga/MC_v1.29/ZinvisibleJets_7TeV-madgraph.root","Z nunu",totEventsZnunu,ZnunuCrossSection,false,false,4,diboson_color);
509 // allsamples.AddSample("/scratch/fronga/MC_v1.29/QCD_TuneD6T_HT-100To250_7TeV-madgraph.root","QCD",totEventsQCD100to250,QCD100to250CrossSection,false,false,5,qcd_color);
510 // allsamples.AddSample("/scratch/fronga/MC_v1.29/QCD_TuneD6T_HT-250To500_7TeV-madgraph.root","QCD",totEventsQCD250to500,QCD250to500CrossSection,false,false,5,qcd_color);
511 allsamples.AddSample("/scratch/fronga/MC_v1.29/QCD_TuneD6T_HT-500To1000_7TeV-madgraph.root","QCD",totEventsQCD500to1000,QCD500to1000CrossSection,false,false,5,qcd_color);
512 allsamples.AddSample("/scratch/fronga/MC_v1.29/QCD_TuneD6T_HT-1000_7TeV-madgraph.root","QCD",totEventsQCD1000toInf,QCD1000toInfCrossSection,false,false,5,qcd_color);
513 allsamples.AddSample("/scratch/fronga/MC_v1.29/DYJetsToLL_TuneZ2_M-50_7TeV-madgraph-tauola.root","ZJets",totEventsZjets,ZJetsCrossSection,false,false,6,dy_color);
514
515
516 allsamples.AddSample("/scratch/fronga/MC_v1.29/LM4_SUSY_sftsht_7TeV-pythia6.root","20xLM4",218536,1.879*20,false,true,7,lm_color);
517
518
519
520
521
522
523
524 // allsamples.AddSample("/scratch/buchmann/DataMay2011ReReco_CertifiedReRecoJSON.root","Data",1,1,true,false,0,kBlack);
525 // allsamples.AddSample("/scratch/buchmann/AllDataJun3_again.root","Data",1,1,true,false,0,kBlack);
526 allsamples.AddSample("/scratch/buchmann/AllData_Jun10___486pb_merged.root","Data",1,1,true,false,0,kBlack);
527 allsamples.AddSample("/scratch/buchmann/MC_Spring11_PU_PF/TTJets_TuneZ2_7TeV-madgraph-tauola.root","t#bar{t}",totEventsTTbar,TTbarCrossSection,false,false,1,ttbar_color);
528 allsamples.AddSample("/scratch/buchmann/MC_Spring11_PU_PF/WJetsToLNu_TuneZ2_7TeV-madgraph-tauola.root","WJets",totEventsWJets,WJetsCrossSection,false,false,3,wjets_color);
529 allsamples.AddSample("/scratch/buchmann/MC_Spring11_PU_PF/VVJetsTo4L_TuneD6T_7TeV-madgraph-tauola.root","Dibosons",totEventsVVJets,VVJetsCrossSection,false,false,4,diboson_color);
530 allsamples.AddSample("/scratch/buchmann/MC_Spring11_PU_PF/TToBLNu_TuneZ2_s-channel_7TeV-madgraph_2.root","SingleTop",totEventsSingleTopS,SingleTopSCrossSection,false,false,2,singletop_color);
531 allsamples.AddSample("/scratch/buchmann/MC_Spring11_PU_PF/TToBLNu_TuneZ2_t-channel_7TeV-madgraph.root","SingleTop",totEventsSingleTopT,SingleTopTCrossSection,false,false,2,singletop_color);
532 allsamples.AddSample("/scratch/buchmann/MC_Spring11_PU_PF/TToBLNu_TuneZ2_tW-channel_7TeV-madgraph.root","SingleTop",totEventsSingleTopU,SingleTopUCrossSection,false,false,2,singletop_color);
533 // allsamples.AddSample("/scratch/fronga/MC_v1.29/ZinvisibleJets_7TeV-madgraph.root","Z nunu",totEventsZnunu,ZnunuCrossSection,false,false,4,diboson_color);
534 // allsamples.AddSample("/scratch/fronga/MC_v1.29/QCD_TuneD6T_HT-100To250_7TeV-madgraph.root","QCD",totEventsQCD100to250,QCD100to250CrossSection,false,false,5,qcd_color);
535 // allsamples.AddSample("/scratch/fronga/MC_v1.29/QCD_TuneD6T_HT-250To500_7TeV-madgraph.root","QCD",totEventsQCD250to500,QCD250to500CrossSection,false,false,5,qcd_color);
536 // allsamples.AddSample("/scratch/fronga/MC_v1.29/QCD_TuneD6T_HT-500To1000_7TeV-madgraph.root","QCD",totEventsQCD500to1000,QCD500to1000CrossSection,false,false,5,qcd_color);
537 // allsamples.AddSample("/scratch/fronga/MC_v1.29/QCD_TuneD6T_HT-1000_7TeV-madgraph.root","QCD",totEventsQCD1000toInf,QCD1000toInfCrossSection,false,false,5,qcd_color);
538 allsamples.AddSample("/scratch/buchmann/MC_Spring11_PU_PF/DYJetsToLL_TuneZ2_M-50_7TeV-madgraph-tauola.root","ZJets",totEventsZjets,ZJetsCrossSection,false,false,6,dy_color);
539
540
541 // allsamples.AddSample("/scratch/fronga/MC_v1.29/LM4_SUSY_sftsht_7TeV-pythia6.root","20xLM4",218536,1.879*20,false,true,7,lm_color);
542 allsamples.AddSample("/scratch/buchmann/LM4_SUSY_sftsht_7TeV-pythia.root","20xLM4",218536,1.879*20,false,true,7,lm_color);
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558 */
559
560
561 void DoNewExperimentalMetPrediction() {
562 // Do it without Pt correction
563 float MCPeakNoPtCorr=0,MCPeakErrorNoPtCorr=0,DataPeakNoPtCorr=0,DataPeakErrorNoPtCorr=0,MCSigma=0,DataSigma=0;
564 stringstream resultsNoPtCorr;
565 stringstream NoPtCorrdatajzb;
566 stringstream NoPtCorrmcjzb;
567
568 string originalJZBvariableDATA = PlottingSetup::jzbvariabledata;
569 string originalJZBvariableMC = PlottingSetup::jzbvariablemc;
570
571 PlottingSetup::jzbvariabledata="jzb[1]";
572 PlottingSetup::jzbvariablemc="jzb[1]";
573
574 stringstream mcjzbNoPtCorr;
575 find_peaks(MCPeakNoPtCorr,MCPeakErrorNoPtCorr, DataPeakNoPtCorr, DataPeakErrorNoPtCorr,resultsNoPtCorr,true,NoPtCorrdatajzb,NoPtCorrmcjzb);
576 string mcjzb = NoPtCorrmcjzb.str();
577
578 PlottingSetup::jzbvariabledata=originalJZBvariableDATA;
579 PlottingSetup::jzbvariablemc=originalJZBvariableMC;
580
581 cout << "Have found a peak in mc at " << MCPeakNoPtCorr << " +/- " << MCPeakErrorNoPtCorr << " ; will use this result." << endl;
582
583 string CorrectedMet="met[4]+"+any2string(abs(2*(MCPeakNoPtCorr)));
584 if(2*(MCPeakNoPtCorr)<0) CorrectedMet="met[4]-"+any2string(2*(MCPeakNoPtCorr));
585 string CorrectedMetUp="met[4]+"+any2string(abs(2*(MCPeakNoPtCorr+MCPeakErrorNoPtCorr)));
586 if(2*(MCPeakNoPtCorr+MCPeakErrorNoPtCorr)<0) "met[4]-"+any2string(abs(2*(MCPeakNoPtCorr+MCPeakErrorNoPtCorr)));
587 string CorrectedMetDown="met[4]+"+any2string(abs(2*(MCPeakNoPtCorr-MCPeakErrorNoPtCorr)));
588 if(2*(MCPeakNoPtCorr-MCPeakErrorNoPtCorr)<0) "met[4]-"+any2string(abs(2*(MCPeakNoPtCorr+MCPeakErrorNoPtCorr)));
589
590 do_experimental_zjets_prediction("met[4]",CorrectedMet,"jzb[1]",MCPeakNoPtCorr,"CovertMapping");
591 do_experimental_zjets_prediction("met[4]",CorrectedMetDown,"jzb[1]",(MCPeakNoPtCorr-MCPeakErrorNoPtCorr),"CovertMapping_SysDown");
592 do_experimental_zjets_prediction("met[4]",CorrectedMetUp,"jzb[1]",(MCPeakNoPtCorr+MCPeakErrorNoPtCorr),"CovertMapping_SysUp");
593
594
595 }
596
597 void DoNewExperimentalMetPrediction_WithPtCorrection() {
598 // Do it with Pt correction
599 float MCPeakWithPtCorr=0,MCPeakErrorWithPtCorr=0,DataPeakWithPtCorr=0,DataPeakErrorWithPtCorr=0,MCSigma=0,DataSigma=0;
600 stringstream resultsWithPtCorr;
601 stringstream WithPtCorrdatajzb;
602 stringstream WithPtCorrmcjzb;
603
604 string originalJZBvariableDATA = PlottingSetup::jzbvariabledata;
605 string originalJZBvariableMC = PlottingSetup::jzbvariablemc;
606
607 stringstream mcjzbWithPtCorr;
608 // find_peaks(MCPeakWithPtCorr,MCPeakErrorWithPtCorr, DataPeakWithPtCorr, DataPeakErrorWithPtCorr,resultsWithPtCorr,true,WithPtCorrdatajzb,WithPtCorrmcjzb);
609 write_warning(__FUNCTION__,"Hardcoded peak position");MCPeakWithPtCorr=3.58273,MCPeakErrorWithPtCorr=0.208997;
610
611 string mcjzb = WithPtCorrmcjzb.str();
612
613 cout << "Have found a peak in mc at " << MCPeakWithPtCorr << " +/- " << MCPeakErrorWithPtCorr << " ; will use this result." << endl;
614
615 string CorrectedMet="met[4]+"+any2string(abs(2*(MCPeakWithPtCorr)));
616 if(2*(MCPeakWithPtCorr)<0) CorrectedMet="met[4]-"+any2string(2*(MCPeakWithPtCorr));
617 string CorrectedMetUp="met[4]+"+any2string(abs(2*(MCPeakWithPtCorr+MCPeakErrorWithPtCorr)));
618 if(2*(MCPeakWithPtCorr+MCPeakErrorWithPtCorr)<0) "met[4]-"+any2string(abs(2*(MCPeakWithPtCorr+MCPeakErrorWithPtCorr)));
619 string CorrectedMetDown="met[4]+"+any2string(abs(2*(MCPeakWithPtCorr-MCPeakErrorWithPtCorr)));
620 if(2*(MCPeakWithPtCorr-MCPeakErrorWithPtCorr)<0) "met[4]-"+any2string(abs(2*(MCPeakWithPtCorr+MCPeakErrorWithPtCorr)));
621
622 do_experimental_zjets_prediction("met[4]",CorrectedMet,originalJZBvariableMC,MCPeakWithPtCorr,"CovertMapping_WithPtCorr_AdaptedPredMet");
623 // do_experimental_zjets_prediction("met[4]",CorrectedMetDown,originalJZBvariableMC,(MCPeakWithPtCorr-MCPeakErrorWithPtCorr),"CovertMapping_WithPtCorr_AdaptedPredMet_SysDown");
624 // do_experimental_zjets_prediction("met[4]",CorrectedMetUp,originalJZBvariableMC,(MCPeakWithPtCorr+MCPeakErrorWithPtCorr),"CovertMapping_WithPtCorr_AdaptedPredMet_SysUp");
625
626 do_experimental_zjets_prediction("met[4]",CorrectedMet+"+0.036*pt",originalJZBvariableMC,MCPeakWithPtCorr,"CovertMapping_WithPtCorr_AdaptedPredMetWithPtDep");
627
628 // do_experimental_zjets_prediction("met[4]","met[4]",originalJZBvariableMC,MCPeakWithPtCorr,"CovertMapping_WithPtCorr");
629
630
631 }
632
633
634 void OUTDATED_ExperimentalMetPrediction() {
635
636 // compute_r_in_out();
637
638 // DoNewExperimentalMetPrediction();
639 DoNewExperimentalMetPrediction_WithPtCorrection();
640
641 cout << "About to switch off the weight : " << (const char*) cutWeight << endl;
642 TCut weightbackup=cutWeight;
643 cutWeight="1.0";
644 cout << "Have set cut weight to 1.0 (no PURW)" << endl;
645
646 // do_experimental_zjets_prediction("met[4]",PlottingSetup::jzbvariablemc,0,"PtCorrectedJZB");
647 /*
648 // Do it without PU re-weighting
649 float MCPeakNoPU=0,MCPeakErrorNoPU=0,DataPeakNoPU=0,DataPeakErrorNoPU=0,MCSigma=0,DataSigma=0;
650 stringstream resultsNoPU;
651 stringstream noPUdatajzb;
652 stringstream noPUmcjzb;
653
654 stringstream mcjzbnoPU;
655 find_peaks(MCPeakNoPU,MCPeakErrorNoPU, DataPeakNoPU, DataPeakErrorNoPU,resultsNoPU,true,noPUdatajzb,noPUmcjzb);
656 string mcjzb = noPUmcjzb.str();
657 */
658 // do_experimental_zjets_prediction("met[4]","jzb[1]",0,"RAWJZB");
659
660 // do_experimental_zjets_prediction("met[4]","((jzb[1]))",-0.614166,"PeakCorrectedJZB");/*
661 // do_experimental_zjets_prediction("met[4]",mcjzb,"FullyCorrectedJZB");
662 /*
663 do_experimental_zjets_prediction(tauevents.GetTree(),"met[4]","jzb[1]","tau_RAWJZB");
664 do_experimental_zjets_prediction(tauevents.GetTree(),"met[4]","((jzb[1])- 3.53399 )","tau_PeakCorrectedJZB");
665 do_experimental_zjets_prediction(tauevents.GetTree(),"met[4]","((jzb[1]+0.059979*pt))","tau_PtCorrectedJZB");
666 do_experimental_zjets_prediction(tauevents.GetTree(),"met[4]","((jzb[1]+0.059979*pt)- 3.53399 )","tau_FullyCorrectedJZB");
667 */
668 /**/
669
670 cutWeight=weightbackup;
671 cout << "Have reset cut weight (PURW reactivated)" << endl;
672 }
673
674 oid do_experimental_zjets_prediction(string metvariable, string predvariable, string jzbvariable, float epsilonZero, string SaveAsName) {
675 TCanvas *can = new TCanvas("can","can");
676
677 cout << "Going to use the variable " << metvariable << " to draw met " << endl;
678 cout << "Going to use the variable " << predvariable<< " to draw the met PREDICTION" << endl;
679 cout << "Going to use the variable " << jzbvariable << " to cut on JZB" << endl;
680 int nbins=50;
681 float xmin=0;
682 float xmax=100;
683
684 string basiccut="pfJetGoodNum40>=3&&abs(mll-91.2)<20&&pt1>20&&pt2>20&&(abs(jzb[1]-"+any2string(epsilonZero)+")>abs("+any2string(epsilonZero)+"))&&((";
685
686 cout << "The basic cut is " << basiccut << endl;
687 stringstream ObsCommand;
688 ObsCommand << basiccut << jzbvariable << ")>0)&&(id1==id2)&&abs("+metvariable+"-75)>5";
689 stringstream PredCommand;
690 PredCommand << basiccut << jzbvariable << ")<0)&&(id1==id2)&&abs("+predvariable+"-75)>5";
691 stringstream PredCommand3;
692 PredCommand3 << basiccut << jzbvariable << ")>0)&&(id1!=id2)&&abs("+predvariable+"-75)>5";
693 stringstream PredCommand4;
694 PredCommand4 << basiccut << jzbvariable << ")<0)&&(id1!=id2)&&abs("+predvariable+"-75)>5";
695
696 TH1F *obs = allsamples.Draw("obs",metvariable,nbins,xmin,xmax,"MET [GeV]","events",TCut(ObsCommand.str().c_str()),mc,PlottingSetup::luminosity,allsamples.FindSample("Z_em_DYJetsToLL"));
697 TH1F *pred = allsamples.Draw("pred", predvariable,nbins,xmin,xmax,"MET [GeV]","events",TCut(PredCommand.str().c_str()),mc,PlottingSetup::luminosity,allsamples.FindSample("Z_em_DYJetsToLL"));
698 TH1F *pre3 = allsamples.Draw("pre3",predvariable,nbins,xmin,xmax,"MET [GeV]","events",TCut(PredCommand3.str().c_str()),mc,PlottingSetup::luminosity,allsamples.FindSample("Z_em_DYJetsToLL"));
699 TH1F *pre4 = allsamples.Draw("pre4",predvariable,nbins,xmin,xmax,"MET [GeV]","events",TCut(PredCommand4.str().c_str()),mc,PlottingSetup::luminosity,allsamples.FindSample("Z_em_DYJetsToLL"));
700
701 obs->SetLineColor(kBlue);
702 obs->SetMarkerSize(0.3);
703 obs->SetMarkerColor(kBlue);
704
705 pred->Add(pre3);
706 pred->Add(pre4,-1);
707
708 pred->SetMarkerColor(kRed);
709 pred->SetLineColor(kRed);
710
711 can->SetLogy(0);
712 obs->Draw("e1");
713 pred->Draw("histo,same");
714 obs->Draw("same,e1");
715
716 TLegend *leg = new TLegend(0.7,0.7,0.89,0.89);
717 leg->SetFillColor(kWhite);
718 leg->AddEntry(obs,"Observed","lp");
719 leg->AddEntry(pred,"Predicted","lp");
720 leg->Draw("same");
721
722 stringstream SaveAsName1;
723 SaveAsName1 << "MetPredPlots/" << SaveAsName;
724 CompleteSave(can,SaveAsName1.str());
725
726 can->SetLogy(0);
727 TH1F *ratio = (TH1F*)obs->Clone("ratio");
728 ratio->SetTitle("obs/pred");
729 ratio->SetName("ratio");
730 ratio->Divide(pred);
731 ratio->GetYaxis()->SetTitle("obs/pred");
732 ratio->GetYaxis()->CenterTitle();
733 ratio->GetYaxis()->SetRangeUser(0,5);
734 TF1 *func = new TF1("func","[0]+[1]*x",10,100);
735 ratio->Fit("func","R");
736 ratio->Draw();
737 TLine *line = new TLine(0.0,1.0,100,1.0);
738 line->SetLineColor(kBlue);
739 line->SetLineWidth(2);
740 line->Draw("same");
741
742 stringstream SaveAsName2;
743 SaveAsName2 << "MetPredPlots/" << SaveAsName << "_ratio";
744 CompleteSave(can,SaveAsName2.str());
745
746 delete line;
747 delete ratio;
748 delete pred;
749 delete pre3;
750 delete pre4;
751 delete obs;
752 delete can;
753 }