ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/cbrown/Development/Plotting/RunMetJZBAnalysis.C
Revision: 1.1
Committed: Thu May 30 08:01:03 2013 UTC (11 years, 11 months ago) by buchmann
Content type: text/plain
Branch: MAIN
Log Message:
Dedicated file for MetJZB analysis

File Contents

# User Rev Content
1 buchmann 1.1 #include <iostream>
2     #include <vector>
3     #include <sys/stat.h>
4     #include <getopt.h>
5     #include <stdio.h>
6     #include <stdlib.h>
7    
8     #ifndef Verbosity
9     #define Verbosity 0
10     #endif
11     #ifndef HUSH
12     #define HUSH 1
13     #endif
14    
15     #include "Modules/GeneralToolBox.C"
16     #include "Modules/SampleClass.C"
17     #include "Modules/setTDRStyle.C"
18     #include "Modules/Setup.C"
19     #include "Modules/Poisson_Calculator.C"
20     #include "Modules/JSON/JSONSampleLoader.C"
21     #include "Modules/ActiveSamples.C"
22     #include "Modules/PeakFinder.C"
23     #include "Modules/UpperLimitsWithShape.C"
24     #include "Modules/Plotting_Functions.C"
25     #include "Modules/LimitCalculation.C"
26     #include "Modules/ResultModule.C"
27     #include "Modules/CrossSectionReader.C"
28     #include "Modules/Systematics.C"
29     #include "Modules/SugarCoating.C"
30     #include "Modules/ExclusionPlot.C"
31     #include "Modules/SUSYScan.C"
32     #include "Modules/AachenCompatibility.C"
33     #include "Modules/MetPlotting.C"
34     #include "Modules/WZStudy.C"
35     //#include "Modules/FSRStudy.C"
36     #include "Modules/TTbarAnalysis.C"
37     #include "Modules/OverlapStudy.C"
38    
39    
40     #include <TCut.h>
41     #include <TROOT.h>
42     #include <TCanvas.h>
43     #include <TMath.h>
44     #include <TColor.h>
45     #include <TPaveText.h>
46     #include <TRandom.h>
47     #include <TH1.h>
48     #include <TH2.h>
49     #include <TF1.h>
50     #include <TSQLResult.h>
51    
52     using namespace PlottingSetup;
53    
54     void usage(int passed=0 ) {
55     std::cout << "USAGE : " << std::endl;
56     std::cout << "You can use different options when running this program : " << std::endl;
57    
58     std::cout << " ********* SIGNAL REGION OPTIONS *********" << std::endl;
59     std::cout << "\033[1;34m cata\033[0m \t\t category a only" << std::endl;
60     std::cout << "\033[1;34m catab\033[0m \t\t category a+b" << std::endl;
61     std::cout << "\033[1;34m eth\033[0m \t\t ETH selection (make sure you specify a or ab)" << std::endl;
62     std::cout << "\033[1;34m aachen\033[0m \t Aachen selection (make sure you specify a or ab)" << std::endl;
63     std::cout << std::endl << std::endl;
64     std::cout << " ********* ANALYSIS OPTIONS *********" << std::endl;
65     std::cout << "\033[1;34m yields\033[0m \t\t Produce all yields (just yields!)" << std::endl;
66     std::cout << "\033[1;34m plots \033[0m \t\t Produce all plots" << std::endl;
67     std::cout << "\033[1;34m comp \033[0m \t\t Produce comparison plots (i.e. SF vs OF plots)" << std::endl;
68     std::cout << "\033[1;34m pred \033[0m \t\t Get JZB-based DY prediction" << std::endl;
69     exit(-1);
70     }
71    
72     int main (int argc, char ** argv)
73     {
74     int do_a = false;
75     int do_ab = false;
76     int do_aachen = false;
77     int do_eth = false;
78    
79    
80     int do_yields = false;
81     int do_plots = false;
82     int do_comp = false;
83     int do_pred = false;
84    
85     int savepdf=true;
86     int saveC=true;
87     int saveRoot=true;
88     int savepng=true;
89     int saveeps=false;
90    
91     std::string directory="";
92     int option_iterator;
93     int option_counter=0;
94     bool moreoptions=true;
95    
96     int SwitchOffSidebands=0;
97    
98     string jzbcuts_string="";
99    
100     while(moreoptions) {
101     static struct option long_options[] =
102     {
103     /* These options set a flag. */
104     {"cata", no_argument, &do_a, 1},
105     {"catab", no_argument, &do_ab, 1},
106     {"eth", no_argument, &do_eth, 1},
107     {"aachen", no_argument, &do_aachen, 1},
108     {"yields", no_argument, &do_yields, 1},
109     {"plots", no_argument, &do_plots, 1},
110     {"comp", no_argument, &do_comp, 1},
111     {"pred", no_argument, &do_pred, 1},
112     {"png", no_argument, &savepng,1},
113     {"eps", no_argument, &saveeps,1},
114     {"pdf", no_argument, &savepdf,1},
115     {"root", no_argument, &saveRoot,1},
116     {"C", no_argument, &saveC,1},
117     {"dir", required_argument, 0, 'd'},
118     {0, 0, 0, 0}
119     };
120     int option_index = 0;
121     option_iterator = getopt_long(argc, argv, "d:",long_options, &option_index);
122     if(option_iterator == -1) moreoptions=false;
123     else {
124     option_counter++;
125     switch (option_iterator)
126     {
127     case 0:
128     if (long_options[option_index].flag != 0)
129     break;
130     printf ("option %s", long_options[option_index].name);
131     if (optarg)
132     printf (" with arg %s", optarg);
133     printf ("\n");
134     break;
135     case 'd':
136     directory=(std::string)optarg;
137     std::cout<<"Option directory was passed with argument " << optarg << std::endl;
138     break;
139     case '?':
140     usage(option_iterator);
141     break;
142     default:
143     usage(option_iterator);
144     }
145     }
146     }
147    
148     if(directory!="") PlottingSetup::directoryname=directory;
149     if(option_counter==0) usage();
150    
151     ///----------------------------------- BELOW THIS LINE: NO MORE OPTIONS BUT ACTUAL FUNCTION CALLS! ---------------------------------------------------------
152     gROOT->SetStyle("Plain");
153     bool do_fat_line=false; // if you want to have HistLineWidth=1 and FuncWidth=1 as it was before instead of 2
154     setTDRStyle(do_fat_line);
155     gStyle->SetTextFont(42);
156     bool showList=true;
157     set_directory(directoryname);//Indicate the directory name where you'd like to save the output files in Setup.C
158     set_treename("events");//you can set the treename here to be used; options are "events" (for reco) for "PFevents" (for particle flow)
159    
160     if(!do_aachen && !do_eth) {
161     write_error(__FUNCTION__,"You need to specify which analysis you want to run - either Aachen or ETH!");
162     return -1;
163     }
164     if(do_aachen && do_eth) {
165     write_error(__FUNCTION__,"You need to specify which analysis you want to run - either Aachen or ETH, but not both!");
166     return -1;
167     }
168    
169     if(!do_a && !do_ab) {
170     write_error(__FUNCTION__,"You need to specify which category (either a or ab)");
171     return -1;
172     }
173    
174     if(do_a && do_ab) {
175     write_error(__FUNCTION__,"You need to specify which category (either a or ab)");
176     return -1;
177     }
178    
179     dout << endl << endl;
180     dout << "Goint to initialize analysis using the following setup: " << endl;
181     dout << "Category \t" << (do_a? " A": " A+B") << endl;
182     dout << "Selection \t" << (do_aachen? " Aachen": " ETH") << endl;
183    
184     dout << "Adapting all cuts according to your selection now: " << endl;
185     if(do_a) {
186     string backup_basicqualitycut = (const char*) basicqualitycut;
187     string backup_essentialcut = (const char*) essentialcut;
188     string backup_basiccut = (const char*) basiccut;
189     string backup_leptoncut = (const char*) leptoncut;
190     string backup_essential = (const char*) essential;
191     string backup_cutNJets = (const char*) cutnJets;
192    
193     string Sbasicqualitycut = backup_basicqualitycut;
194     Sbasicqualitycut = ReplaceAll(Sbasicqualitycut,")<2.4",")<1.4");
195     basicqualitycut=TCut(Sbasicqualitycut.c_str());
196    
197     string Sbasiccut = backup_basiccut;
198     Sbasiccut = ReplaceAll(Sbasiccut,")<2.4",")<1.4");
199     basiccut=TCut(Sbasiccut.c_str());
200    
201     string Sessentialcut = backup_essentialcut;
202     Sessentialcut = ReplaceAll(Sessentialcut,")<2.4",")<1.4");
203     essentialcut=TCut(Sessentialcut.c_str());
204    
205     string Sleptoncut = backup_leptoncut;
206     Sleptoncut = ReplaceAll(Sleptoncut,")<2.4",")<1.4");
207     leptoncut=TCut(Sleptoncut.c_str());
208    
209     string Sessential = backup_essential;
210     Sessential = ReplaceAll(Sessential,")<2.4",")<1.4");
211     essential = TCut(Sessential.c_str());
212    
213     string ScutnJets = backup_cutNJets;
214     ScutnJets = ReplaceAll(ScutnJets,")<2.4",")<1.4");
215     cutnJets=TCut(ScutnJets.c_str());
216     }
217     if(do_ab) {
218     string backup_basicqualitycut = (const char*) basicqualitycut;
219     string backup_essentialcut = (const char*) essentialcut;
220     string backup_basiccut = (const char*) basiccut;
221     string backup_leptoncut = (const char*) leptoncut;
222     string backup_essential = (const char*) essential;
223     string backup_cutNJets = (const char*) cutnJets;
224    
225     string Sbasicqualitycut = backup_basicqualitycut;
226     Sbasicqualitycut = ReplaceAll(Sbasicqualitycut,")<1.4",")<2.4");
227     basicqualitycut=TCut(Sbasicqualitycut.c_str());
228    
229     string Sbasiccut = backup_basiccut;
230     Sbasiccut = ReplaceAll(Sbasiccut,")<1.4",")<2.4");
231     basiccut=TCut(Sbasiccut.c_str());
232    
233     string Sessentialcut = backup_essentialcut;
234     Sessentialcut = ReplaceAll(Sessentialcut,")<1.4",")<2.4");
235     essentialcut=TCut(Sessentialcut.c_str());
236    
237     string Sleptoncut = backup_leptoncut;
238     Sleptoncut = ReplaceAll(Sleptoncut,")<1.4",")<2.4");
239     leptoncut=TCut(Sleptoncut.c_str());
240    
241     string Sessential = backup_essential;
242     Sessential = ReplaceAll(Sessential,")<1.4",")<2.4");
243     essential = TCut(Sessential.c_str());
244    
245     string ScutnJets = backup_cutNJets;
246     ScutnJets = ReplaceAll(ScutnJets,")<1.4",")<2.4");
247     ScutnJets = ReplaceAll(ScutnJets,"pfJetGoodNum40>=3","pfJetGoodNum40>=2");
248     cutnJets=TCut(ScutnJets.c_str());
249     }
250    
251     cout << "Basic cut : " << (const char*) basiccut << endl;
252     cout << "Essential cut : " << (const char*) essentialcut << endl;
253     cout << "Essential cut ('): " << (const char*) essential << endl;
254     cout << "Lepton cut : " << (const char*) leptoncut << endl;
255    
256     define_samples(showList,allsamples,signalsamples,scansample,raresample,systsamples,qcdsamples,comparesamples);
257    
258     setlumi(luminosity);
259     setessentialcut(essential);//this sets the essential cut; this one is used in the draw command so it is AUTOMATICALLY applied everywhere. IMPORTANT: Do NOT store weights here!
260     stringstream resultsummary;
261    
262     //write_analysis_type(PlottingSetup::RestrictToMassPeak,PlottingSetup::DoBTag);
263     do_png(savepng);
264     do_pdf(savepdf);
265     do_eps(saveeps);
266     do_C(saveC);
267     do_root(saveRoot);
268    
269     //**** part 1 : peak finding
270     float MCPeak=0,MCPeakError=0,DataPeak=0,DataPeakError=0,MCSigma=10,DataSigma=10;
271     method=Kostasmethod;//Kostasmethod;//dogaus3sigma;// options: dogaus,doKM,dogaus2sigma,dogaus3sigma
272     stringstream datajzb;
273     stringstream mcjzb;
274    
275    
276     if ( 0 ) {
277     find_peaks(MCPeak,MCPeakError, DataPeak, DataPeakError,resultsummary,true,datajzb,mcjzb);
278    
279     if(datajzb.str().size()<1) datajzb<<"("<<jzbvariabledata<<")";
280     if(mcjzb.str().size()<1) mcjzb<<"("<<jzbvariablemc<<")";
281     } else {
282     write_warning(__FUNCTION__,"NOT DOING ANY PEAK FINDING ATM");
283     datajzb << "((jzb[1]+0.059979*pt)- 4.18335 )";
284     mcjzb << "((jzb[1]+0.034665*pt)- 3.58273 )";
285     }
286    
287     dout << "With peak correction, we get : " << endl;
288     dout << " Data : " << datajzb.str() << endl;
289     dout << " MC : " << mcjzb.str() << endl;
290    
291    
292     if(do_comp) ProduceOFSFPlots(mcjzb.str(),datajzb.str(),do_aachen);
293    
294     if(do_yields) DoMetPlots(datajzb.str(),mcjzb.str(),do_aachen,do_eth);
295    
296     if(do_plots) DoMetPlots(datajzb.str(),mcjzb.str(),do_aachen,do_eth);
297    
298     if(do_pred) ExperimentalMetPrediction(false, do_aachen, false);
299    
300     // if(do_kinematic_variables||(do_all&&!PlottingSetup::Approved)) do_kinematic_plots(mcjzb.str(),datajzb.str());
301    
302    
303     /*
304    
305     if(do_pick_up_events) {
306     TCut essentialcutbkp = essentialcut;
307     essentialcut=TCut("1.0");
308    
309     cout << "LOW MASS, SF, ETH " << endl;
310     pick_up_events("((((((passed_triggers ||!is_data)))&&((abs(eta1)<1.4 && abs(eta2)<1.4 && pt1>30 &&
311     pt2>20)))&&(((met[4]>100)&&((id1==id2)&&(ch1*ch2<0)))&&((((pfJetGoodNum40>=2&&pfJetGoodID[0]!=0)&&(pfJetGoodNum40>=2&&pfJetGoodID[1]!=0))&&((mll>2)&&((abs(eta1)<1.4 && abs(eta2)<1.4 && pt1>30 && pt2>20))))&&(pfJetGoodNum40>=3))&&mll>20&&mll<70&&(abs(eta1)<1.4 && abs(eta2)<1.4 && pt1>30 && pt2>20)&&pt1>20&&pt2>20)))*((weight*(weight<1000)*(is_data+(!is_data)*((id1==id2)*0.95+(id1!=id2)*0.94))))","ETH_SF.txt");
312    
313     cout << "LOW MASS, OF, ETH " << endl;
314     pick_up_events("((((((passed_triggers ||!is_data)))&&((abs(eta1)<1.4 && abs(eta2)<1.4 && pt1>30 && pt2>20)))&&(((met[4]>100)&&((id1!=id2)&&(ch1*ch2<0)))&&((((pfJetGoodNum40>=2&&pfJetGoodID[0]!=0)&&(pfJetGoodNum40>=2&&pfJetGoodID[1]!=0))&&((mll>2)&&((abs(eta1)<1.4 && abs(eta2)<1.4 && pt1>30 && pt2>20))))&&(pfJetGoodNum40>=3))&&mll>20&&mll<70&&(abs(eta1)<1.4 && abs(eta2)<1.4 && pt1>30 && pt2>20)&&pt1>20&&pt2>20)))*((weight*(weight<1000)*(is_data+(!is_data)*((id1==id2)*0.95+(id1!=id2)*0.94))))","ETH_OF.txt");
315    
316     cout << "LOW MASS, SF, Aachen " << endl;
317     pick_up_events("((((((passed_triggers ||!is_data)))&&((abs(eta1)<2.4 && abs(eta2)<2.4 && pt1>30 && pt2>10)))&&(((met[4]>150)&&((id1==id2)&&(ch1*ch2<0)))&&((((pfJetGoodNum40>=2&&pfJetGoodID[0]!=0)&&(pfJetGoodNum40>=2&&pfJetGoodID[1]!=0))&&((mll>2)&&((abs(eta1)<2.4 && abs(eta2)<2.4 && pt1>30 && pt2>10))))&&(pfJetGoodNum40>=2&&pfTightHT>100))&&mll>20&&mll<70&&(abs(eta1)<2.4 && abs(eta2)<2.4 && pt1>30 && pt2>10)&&pt1>20&&pt2>10)))*((weight*(weight<1000)*(is_data+(!is_data)*((id1==id2)*0.95+(id1!=id2)*0.94))))","Aachen_SF.txt");
318    
319     cout << "LOW MASS, OF, Aachen " << endl;
320     pick_up_events("((((((passed_triggers ||!is_data)))&&((abs(eta1)<2.4 && abs(eta2)<2.4 && pt1>30 && pt2>10)))&&(((met[4]>150)&&((id1!=id2)&&(ch1*ch2<0)))&&((((pfJetGoodNum40>=2&&pfJetGoodID[0]!=0)&&(pfJetGoodNum40>=2&&pfJetGoodID[1]!=0))&&((mll>2)&&((abs(eta1)<2.4 && abs(eta2)<2.4 && pt1>30 && pt2>10))))&&(pfJetGoodNum40>=2&&pfTightHT>100))&&mll>20&&mll<70&&(abs(eta1)<2.4 && abs(eta2)<2.4 && pt1>30 && pt2>10)&&pt1>20&&pt2>10)))*((weight*(weight<1000)*(is_data+(!is_data)*((id1==id2)*0.95+(id1!=id2)*0.94))))","Aachen_OF.txt");
321    
322     essentialcut=essentialcutbkp;
323     }
324    
325     if(do_metjzb) met_jzb_cut(datajzb.str(),mcjzb.str(),jzb_cut);
326    
327     if(do_metpred) ExperimentalMetPrediction(false,false);
328    
329     if(do_metplots) {
330     DoMetPlots(datajzb.str(),mcjzb.str());
331     }
332    
333     if(do_wz) WZstudy();
334    
335     if(do_rmue) ComputeRMuE();
336    
337     if(do_overlap) overlapanalysis(datajzb.str(),mcjzb.str());
338    
339     */
340    
341     return 0;
342     }
343