ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/LJMet/MultivariateAnalysis/root/toyMC.C
(Generate patch)

Comparing UserCode/LJMet/MultivariateAnalysis/root/toyMC.C (file contents):
Revision 1.3 by kukartse, Fri May 1 07:54:10 2009 UTC vs.
Revision 1.5 by kukartse, Mon May 11 16:11:54 2009 UTC

# Line 6 | Line 6
6   #include "TTree.h"
7   #include "TRandom3.h"
8   #include "TCanvas.h"
9 + #include "TSlider.h"
10 +
11 + //gROOT->SetStyle("Plain");
12 + //gStyle->SetOptStat(0);
13  
14   using namespace std;
15  
# Line 14 | Line 18 | class toyMC
18   public:
19    
20    toyMC();
21 +  toyMC(int nbins);
22    ~toyMC();
23    
24    void load_templates();
# Line 29 | Line 34 | public:
34    TH1F * wzjets;
35    TH1F * qcd;
36  
37 +  int n_bins;
38 +
39    TObjArray * mc;
40    TFractionFitter * fit;
41  
42 +  int n_gen_sig,n_gen_bkg1,n_gen_bkg2;
43 +
44    TH1F * toy_data;
45  
46   private:
# Line 47 | Line 56 | toyMC::toyMC(){
56    mc = 0;
57    fit = 0;
58    toy_data = 0;
59 +  n_bins = 50;
60 + }
61 +
62 +
63 + toyMC::toyMC(int nbins){
64 +  data = 0;
65 +  ttbar = 0;
66 +  wzjets = 0;
67 +  qcd = 0;
68 +  mc = 0;
69 +  fit = 0;
70 +  toy_data = 0;
71 +  n_bins = nbins;
72   }
73  
74  
# Line 62 | Line 84 | toyMC::~toyMC(){
84  
85  
86   void toyMC::load_templates(){
87 <  TFile * data_file = new TFile("./TMVApp-data-27apr2009.root","READ");
87 >  TFile * data_file = new TFile("./TMVApp-data-06may2009.root","READ");
88    TFile * ttbar_phys_template_file = new TFile("./TMVA.root","READ");
89 <  TFile * qcd_template_file = new TFile("./TMVApp-qcd-27apr2009.root","READ");
89 >  TFile * qcd_template_file = new TFile("./TMVApp-qcd-06may2009.root","READ");
90    
91    TTree * t_data       = (TTree *)data_file->Get("classifier");
92    TTree * t_ttbar_phys = (TTree *)ttbar_phys_template_file->Get("TestTree");
# Line 79 | Line 101 | void toyMC::load_templates(){
101    delete ttbar;
102    delete wzjets;
103    delete qcd;
104 <  data   = new TH1F("data"  ,"data"  ,30, -0.8, 0.8);
105 <  ttbar  = new TH1F("ttbar" ,"ttbar" ,30, -0.8, 0.8);
106 <  wzjets = new TH1F("wzjets","wzjets",30, -0.8, 0.8);
107 <  qcd    = new TH1F("qcd"   ,"qcd"   ,30, -0.8, 0.8);
104 >  data   = new TH1F("data"  ,"data"  ,n_bins, -0.8, 0.8);
105 >  ttbar  = new TH1F("ttbar" ,"ttbar" ,n_bins, -0.8, 0.8);
106 >  wzjets = new TH1F("wzjets","wzjets",n_bins, -0.8, 0.8);
107 >  qcd    = new TH1F("qcd"   ,"qcd"   ,n_bins, -0.8, 0.8);
108 >
109 >  //
110 >  data->SetMarkerStyle(21);
111 >  data->SetMarkerSize(0.7);
112 >  //data->SetFillStyle(4050);
113 >  ttbar->SetFillColor(42);
114 >  //ttbar->SetFillStyle(4050);
115 >  wzjets->SetFillColor(46);
116 >  //wzjets->SetFillStyle(4050);
117 >  qcd->SetFillColor(48);
118 >  //TSlider *slider = 0;
119 >  //
120  
121    TCanvas * c1 = new TCanvas("canvas1","canvas1",800,600);
122    c1->Divide(2,2);
123 <  c1->cd(1);
123 >  TVirtualPad * pad = c1->cd(1);
124    t_data   -> Draw("MVA_BDT>>data");
125 <  c1->cd(2);
125 >  pad->SaveAs("mu_data.eps");
126 >  pad = c1->cd(2);
127    t_ttbar_phys  -> Draw("MVA_BDT>>ttbar","type==1");
128 <  c1->cd(3);
128 >  pad->SaveAs("mu_ttbar_template.eps");
129 >  pad = c1->cd(3);
130    t_ttbar_phys -> Draw("MVA_BDT>>wzjets","type==0");
131 <  c1->cd(4);
131 >  pad->SaveAs("mu_wzjets_template.eps");
132 >  pad = c1->cd(4);
133    t_qcd    -> Draw("MVA_BDT>>qcd");
134 +  pad->SaveAs("mu_qcd_template.eps");
135  
136    set_overflow_bins(data);
137    set_overflow_bins(ttbar);
# Line 113 | Line 151 | void toyMC::load_templates(){
151  
152  
153   void toyMC::load_ejets_templates(){
154 <  //TFile * data_file = new TFile("./TMVApp-data-electron-01may2009.root","READ");
155 <  TFile * data_file = new TFile("./TMVApp-data-wzfastsim-electron-01may2009.root","READ");
156 <  TFile * ttbar_phys_template_file = new TFile("./TMVA-el.root","READ");
157 <  TFile * qcd_template_file = new TFile("./TMVApp-qcd-electron-01may2009.root","READ");
154 >  TFile * data_file = new TFile("./TMVApp-data-electron-08may2009.root","READ");
155 >  //TFile * data_file = new TFile("./TMVApp-data-electron_wzfastsim-08may2009.root","READ");
156 >  TFile * ttbar_phys_template_file = new TFile("./TMVA.root","READ");
157 >  TFile * qcd_template_file = new TFile("./TMVApp-qcd-electron-08may2009.root","READ");
158    
159    TTree * t_data       = (TTree *)data_file->Get("classifier");
160    TTree * t_ttbar_phys = (TTree *)ttbar_phys_template_file->Get("TestTree");
# Line 127 | Line 165 | void toyMC::load_ejets_templates(){
165    delete wzjets;
166    delete qcd;
167  
168 <  data   = new TH1F("data"  ,"data"  ,50, -0.9, 0.9);
169 <  ttbar  = new TH1F("ttbar" ,"ttbar" ,50, -0.9, 0.9);
170 <  wzjets = new TH1F("wzjets","wzjets",50, -0.9, 0.9);
171 <  qcd    = new TH1F("qcd"   ,"qcd"   ,50, -0.9, 0.9);
168 >  data   = new TH1F("data"  ,"data"  ,n_bins, -0.9, 0.9);
169 >  ttbar  = new TH1F("ttbar" ,"ttbar" ,n_bins, -0.9, 0.9);
170 >  wzjets = new TH1F("wzjets","wzjets",n_bins, -0.9, 0.9);
171 >  qcd    = new TH1F("qcd"   ,"qcd"   ,n_bins, -0.9, 0.9);
172 >
173 >  //
174 >  data->SetMarkerStyle(21);
175 >  data->SetMarkerSize(0.7);
176 >  //data->SetFillStyle(4050);
177 >  ttbar->SetFillColor(42);
178 >  //ttbar->SetFillStyle(4050);
179 >  wzjets->SetFillColor(46);
180 >  //wzjets->SetFillStyle(4050);
181 >  qcd->SetFillColor(48);
182 >  //TSlider *slider = 0;
183 >  //
184  
185 +  TCanvas * c1 = new TCanvas("canvas1","canvas1",800,600);
186 +  c1->Divide(2,2);
187 +  TVirtualPad * pad = c1->cd(1);
188    t_data   -> Draw("MVA_BDT>>data");
189 +  pad->SaveAs("e_data.eps");
190 +  pad = c1->cd(2);
191    t_ttbar_phys  -> Draw("MVA_BDT>>ttbar","type==1");
192 +  pad->SaveAs("e_ttbar_template.eps");
193 +  pad = c1->cd(3);
194    t_ttbar_phys -> Draw("MVA_BDT>>wzjets","type==0");
195 +  pad->SaveAs("e_wzjets_template.eps");
196 +  pad = c1->cd(4);
197    t_qcd    -> Draw("MVA_BDT>>qcd");
198 +  pad->SaveAs("e_qcd_template.eps");
199  
200    //data_file->Close();
201    //ttbar_phys_template_file->Close();
# Line 152 | Line 212 | void toyMC::load_ejets_templates(){
212   int toyMC::do_fit(TH1F * _data){
213    delete fit;
214    fit = new TFractionFitter(_data, mc); // initialise
215 +
216    fit->Constrain(1,0.0,1.0);               // constrain fraction 1 to be between 0 and 1
217    fit->Constrain(2,0.0,1.0);               // constrain fraction 1 to be between 0 and 1
218 <  fit->Constrain(3,0.0,1.0);               // constrain fraction 1 to be between 0 and 1
219 <  //fit->Constrain(3,0.07298,0.07299);               // constrain fraction 1 to be between 0 and 1
220 <  //fit->Constrain(3,0.2981,0.2982);               // constrain fraction 1 to be between 0 and 1
221 <  //fit->SetRangeX(1,30);                    // use only the first 15 bins in the fit
218 >  //fit->Constrain(3,0.0,1.0);               // constrain fraction 1 to be between 0 and 1
219 >
220 >  //fit->Constrain(1,0.15,0.85);               // constrain fraction 1 to be between 0 and 1
221 >  //fit->Constrain(2,0.1651865,0.1651866);               // constrain fraction 1 to be between 0 and 1
222 >  //fit->Constrain(3,0.15,0.85);               // constrain fraction 1 to be between 0 and 1
223 >  //fit->Constrain(3,0.024999,0.025001);               // constrain fraction 1 to be between 0 and 1
224 >  fit->Constrain(3,0.4422735,0.4422736);               // constrain fraction 1 to be between 0 and 1
225 >  //double qcd_frac = (double)n_gen_bkg2/(double)(n_gen_sig+n_gen_bkg1+n_gen_bkg2);
226 >  //fit->Constrain(3,qcd_frac-0.000001,qcd_frac+0.000001);               // constrain fraction 1 to be between 0 and 1
227 >  fit->SetRangeX(1,n_bins);                    // use only some bins to fit
228    Int_t status = fit->Fit();               // perform the fit
229    if (status!=0) status = fit->Fit();
230    cout << "fit status: " << status << endl;
231    return status;
232   }
233  
234 + /*
235 + int toyMC::fit1(TH1F * _data, double bkg2_low = 0.0, bkg2_high = 1.0){
236 +  delete fit;
237 +  fit = new TFractionFitter(_data, mc); // initialise
238 +
239 +  fit->Constrain(1,data_low,data_high); // constrain fraction 1 to be between 0 and 1
240 +  fit->Constrain(2,bkg1_low,bkg1_high); // constrain fraction 1 to be between 0 and 1
241 +  fit->Constrain(3,bkg2_low,bkg2_high); // constrain fraction 1 to be between 0 and 1
242 +  fit->SetRangeX(first_bin,last_bin);   // use only some bins to fit
243 +  Int_t status = fit->Fit();            // perform the fit
244 +  if (status!=0) status = fit->Fit();   // make one more attempt to fit if the first one fails
245 +  cout << "fit status: " << status << endl;
246 +  return status;
247 + }
248 + */
249  
250   int toyMC::generate_toy(int n_sig, int n_bg, int n_qcd){
251    int n_events=0;
252    delete toy_data;
253  
254 <  toy_data = new TH1F("toy_data"  ,"toy_data"  ,50, -0.9, 0.9);
254 >  toy_data = new TH1F("toy_data"  ,"toy_data"  ,n_bins, -0.9, 0.9);
255  
256    int _nsig = _random.Poisson(n_sig);
257    int _nbg  = _random.Poisson(n_bg);
258    int _nqcd = _random.Poisson(n_qcd);
259    //int _nqcd = n_qcd;
260 +
261 +  n_gen_sig = _nsig;
262 +  n_gen_bkg1 = _nbg;
263 +  n_gen_bkg2 = _nqcd;  
264    n_events = _nsig+_nbg+_nqcd;
265  
266    for (int i=0;i<_nsig;i++){
# Line 189 | Line 275 | int toyMC::generate_toy(int n_sig, int n
275    //toy_data->Draw();
276    set_overflow_bins(toy_data);  
277    return n_events;
278 +  //return _nsig;
279   }
280  
281   void toyMC::set_overflow_bins(TH1F * h){
# Line 216 | Line 303 | void toyMC::debug_hist(TH1F * h){
303  
304   int toyMC::fit_data(){
305    //load_templates();
306 <  load_ejets_templates();
306 >  //load_ejets_templates();
307    int _res = do_fit(data);
308    int n_data = data->GetEntries();
309    TH1F * result = (TH1F*) fit->GetPlot();
310 +  result->SetFillColor(16);
311 +  //result->SetFillStyle(4050);
312    data->SetMinimum(0);
313    TCanvas * c = new TCanvas("canvas2","canvas2",800,600);
314    data->Draw("Ep");
# Line 230 | Line 319 | int toyMC::fit_data(){
319    ttbar->Scale(n_data*value/ttbar->Integral());
320    ttbar->Draw("same");
321    //
233  fit->GetResult(1,value,error);
234  wzjets->Scale(n_data*value/wzjets->Integral());
235  wzjets->Draw("same");
236  //
322    fit->GetResult(2,value,error);
323    qcd->Scale(n_data*value/qcd->Integral());
324    qcd->Draw("same");
325    data->Draw("Ep:same");
326 +  //
327 +  fit->GetResult(1,value,error);
328 +  wzjets->Scale(n_data*value/wzjets->Integral());
329 +  wzjets->Draw("same");
330 +  //
331 +  c->SaveAs("fit.eps");
332    return _res;
333   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines