ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/Selection/src/applySelection.cc
Revision: 1.9
Committed: Mon Oct 17 16:23:43 2011 UTC (13 years, 7 months ago) by khahn
Content type: text/plain
Branch: MAIN
Changes since 1.8: +25 -5 lines
Log Message:
*** empty log message ***

File Contents

# Content
1 //
2 // System headers
3 //
4 #include <vector> // STL vector class
5 #include <iostream> // standard I/O
6 #include <iomanip> // functions to format standard I/O
7 #include <fstream> // functions for file I/O
8 #include <string> // C++ string class
9 #include <sstream> // class for parsing strings
10 #include <assert.h>
11 #include <stdlib.h>
12 #include <getopt.h>
13 using namespace std;
14
15 //
16 // ROOT headers
17 //
18 #include <TROOT.h> // access to gROOT, entry point to ROOT system
19 #include <TSystem.h> // interface to OS
20 #include <TFile.h> // file handle class
21 #include <TNtuple.h>
22 #include <TTree.h> // class to access ntuples
23 #include <TChain.h> //
24 #include <TBranch.h> // class to access branches in TTree
25 #include <TClonesArray.h> // ROOT array class
26 #include <TCanvas.h> // class for drawing
27 #include <TH1F.h> // 1D histograms
28 #include <TBenchmark.h> // class to track macro running statistics
29 #include <TLorentzVector.h> // 4-vector class
30 #include <TVector3.h> // 3D vector class
31
32 //
33 // ntuple format headers
34 //
35 #include "HiggsAnaDefs.hh"
36 #include "TEventInfo.hh"
37 #include "TElectron.hh"
38 #include "TMuon.hh"
39 #include "TJet.hh"
40 #include "SiMVAElectronSelection.h"
41 #include "RunLumiRangeMap.h"
42
43 //
44 // utility headers
45 //
46 #include "ParseArgs.h"
47 #include "SampleWeight.h"
48 #include "Selection.h"
49 #include "HZZCiCElectronSelection.h"
50 #include "HZZLikelihoodElectronSelection.h"
51 #include "HZZBDTElectronSelection.h"
52 #include "SampleWeight.h"
53
54 //#define BDTFAIL_THEIR_EVENTS
55 //#define THEIR_EVENTS
56
57 //=== MAIN =================================================================================================
58 int main(int argc, char** argv)
59 {
60
61 vector<vector<string> > inputFiles;
62 inputFiles.push_back(vector<string>());
63
64 //
65 // args
66 //--------------------------------------------------------------------------------------------------------------
67 ControlFlags ctrl;
68 parse_args( argc, argv, ctrl );
69 if( ctrl.inputfile.empty() ||
70 ctrl.inputfile.empty() ) {
71 cerr << "usage: applySelection.exe <flags> " << endl;
72 cerr << "\tmandoatory flags : " << endl;
73 cerr << "\t--inputfile | file containing a list of ntuples to run over" << endl;
74 cerr << "\t--outputfile | file to store selected evet" << endl;
75 return 1;
76 }
77 ctrl.dump();
78
79 map<string,double> entrymap; // number of unskimmed entries in each file
80
81 //
82 // File I/O
83 //--------------------------------------------------------------------------------------------------------------
84 TChain * chain = new TChain("Events");
85 ifstream f(ctrl.inputfile.c_str());
86 string fname;
87 while (f >> fname) {
88 if( !(strncmp( fname.c_str(), "#", 1 ) ) ) continue; // skip commented lines
89 cout << "adding inputfile : " << fname.c_str() << endl;
90 entrymap[string(fname.c_str())] = unskimmedEntries(fname.c_str());
91 chain->AddFile(fname.c_str());
92 }
93
94 // table of cross section values
95 SimpleTable xstab("./data/xs.dat");
96
97 //
98 // Setup
99 //--------------------------------------------------------------------------------------------------------------
100 TH1F * h_evt = new TH1F("hevt", "hevt", 2, 0, 2 );
101 TH1F * h_evtfail = new TH1F("hevtfail", "hevtfail", 1024, 0, 1024 );
102 // TNtuple * passtuple = new TNtuple( "passtuple", "passtuple", "run:evt:lumi:channel:mZ1:mZ2:m4l:pt4l:w" );
103 TTree * passtuple = new TTree("passtuple", "passtuple" );
104 unsigned run, evt, lumi, channel;
105 float mZ1, mZ2, m4l, pt4l;
106 unsigned tZ1, tZ2;
107 double w;
108 passtuple->Branch("run", &run);
109 passtuple->Branch("evt", &evt);
110 passtuple->Branch("lumi", &lumi);
111 passtuple->Branch("mZ1", &mZ1);
112 passtuple->Branch("mZ2", &mZ2);
113 passtuple->Branch("tZ1", &tZ1);
114 passtuple->Branch("tZ2", &tZ2);
115 passtuple->Branch("m4l", &m4l);
116 passtuple->Branch("pt4l", &pt4l);
117 passtuple->Branch("w", &w);
118 initCiCSelection();
119 if(ctrl.eleSele=="lik") initLikSelection();
120 if(ctrl.eleSele=="bdt") initBDTSelection();
121 if(ctrl.eleSele=="si" ) initSiMVAElectronSelection();
122 initRunLumiRangeMap();
123
124 //
125 // Loop
126 //--------------------------------------------------------------------------------------------------------------
127
128 //
129 // Access samples and fill histograms
130 TFile *inputFile=0;
131 TTree *eventTree=0;
132
133 // Data structures to store info from TTrees
134 mithep::TEventInfo *info = new mithep::TEventInfo();
135 TClonesArray *electronArr = new TClonesArray("mithep::TElectron");
136 TClonesArray *muonArr = new TClonesArray("mithep::TMuon");
137
138
139 TBranch *infoBr;
140 TBranch *electronBr;
141 TBranch *muonBr;
142
143 // chain->SetBranchStatus("*",0); //disable all branches
144 // chain->SetBranchStatus("Info", 1);
145 // chain->SetBranchStatus("Muon", 1);
146 // chain->SetBranchStatus("Electron", 1);
147
148 chain->SetBranchAddress("Info", &info);
149 chain->SetBranchAddress("Electron", &electronArr);
150 chain->SetBranchAddress("Muon", &muonArr);
151
152
153 int count=0, pass=0;
154 cout << "nEntries: " << chain->GetEntries() << endl;
155 for(UInt_t ientry=0; ientry<chain->GetEntries(); ientry++) {
156
157 chain->GetEntry(ientry);
158 if(!(ientry%100000)) cout << "entry: " << ientry << endl;
159
160 // get event weight for this file
161 double eventweight=1;
162 if(ctrl.mc) eventweight = getWeight(xstab,entrymap,chain);
163
164 #ifdef THEIR_EVENTS
165 if( !( info->evtNum == 504867308 ||
166 info->evtNum == 368148849 ||
167 // info->evtNum == 129514273 ||
168 info->evtNum == 286336207 ||
169 info->evtNum == 344708580 ||
170 info->evtNum == 30998576 ||
171 info->evtNum == 155679852 ||
172 info->evtNum == 394010457 ||
173 info->evtNum == 917379387 ||
174 info->evtNum == 78213037 ||
175 info->evtNum == 862270386 ||
176 info->evtNum == 337493970 || // not baseline anymore?
177 info->evtNum == 1491724484 ||
178 info->evtNum == 480301165 ||
179 info->evtNum == 1038911933 ||
180 info->evtNum == 876658967 ||
181 info->evtNum == 966824024 ||
182 info->evtNum == 141954801 ||
183 info->evtNum == 160966858 ||
184 info->evtNum == 191231387 ||
185 info->evtNum == 66033190 ||
186 info->evtNum == 10347106 ||
187 info->evtNum == 107360878 ||
188 info->evtNum == 2554393033 ||
189 info->evtNum == 933807102 ||
190 info->evtNum == 1188043146 ||
191 info->evtNum == 559839432 ||
192 info->evtNum == 16706390 ||
193 info->evtNum == 65557571 ||
194 info->evtNum == 389185367 ||
195 info->evtNum == 2722114329 ) ) continue;
196 #endif
197
198
199 #ifdef BDTFAIL_THEIR_EVENTS
200 if( !( info->evtNum == 78213037 ||
201 info->evtNum == 876658967 ) ) continue;
202
203 #endif
204
205 unsigned evtfail = fails_HZZ4L_selection(ctrl, info, electronArr, muonArr, eventweight, passtuple);
206 h_evtfail->Fill( evtfail, eventweight );
207
208 count++;
209 if( !evtfail ) pass++;
210
211 } //end loop over data
212
213 if( ctrl.mc ) {
214 cout << "--------------" << endl;
215 cout << "denom: " << count << endl;
216 cout << "pass: " << pass << endl;
217 cout << "axe: " << (float)pass/count << endl;
218 cout << "--------------" << endl;
219 }
220
221 delete info;
222 delete electronArr;
223 delete muonArr;
224
225
226
227 //--------------------------------------------------------------------------------------------------------------
228 // Save Histograms;
229 //==============================================================================================================
230 const char * ofname;
231 if( strcmp( ctrl.outputfile.c_str(), "") ) {
232 ofname = ctrl.outputfile.c_str();
233 } else {
234 ofname = "tmp.root";
235 }
236
237 TFile *file = new TFile(ofname, "RECREATE");
238 h_evt->Write();
239 h_evtfail->Write();
240 passtuple->Write();
241 file->Close();
242 delete file;
243
244 // map<TString,TMVA::Reader*>::iterator it;
245 // for(it=readers.begin(); it!=readers.end(); it++) delete (*it).second;
246
247 cerr << "done!" << endl;
248 }
249
250
251