Selection/src/applySelection.cc

//
// System headers
//
#include <vector>                   // STL vector class
#include <iostream>                 // standard I/O
#include <iomanip>                  // functions to format standard I/O
#include <fstream>                  // functions for file I/O
#include <string>                   // C++ string class
#include <sstream>                  // class for parsing strings
#include <assert.h>
#include <stdlib.h>   
#include <getopt.h>
using namespace std;

//
// ROOT headers
//
#include <TROOT.h>                  // access to gROOT, entry point to ROOT system
#include <TSystem.h>                // interface to OS
#include <TFile.h>                  // file handle class
#include <TNtuple.h> 
#include <TTree.h>                  // class to access ntuples
#include <TChain.h>                 // 
#include <TBranch.h>                // class to access branches in TTree
#include <TClonesArray.h>           // ROOT array class
#include <TCanvas.h>                // class for drawing
#include <TH1F.h>                   // 1D histograms
#include <TBenchmark.h>             // class to track macro running statistics
#include <TLorentzVector.h>         // 4-vector class
#include <TVector3.h>               // 3D vector class

//
// ntuple format headers
//
#include "HiggsAnaDefs.hh"
#include "TEventInfo.hh"
#include "TElectron.hh"
#include "TMuon.hh"
#include "TJet.hh"
#include "RunLumiRangeMap.h"

//
// utility headers
//
#include "ParseArgs.h"
#include "SampleWeight.h"
#include "Selection.h"
#include "HZZCiCElectronSelection.h"
#include "HZZLikelihoodElectronSelection.h"
#include "HZZBDTElectronSelection.h"
#include "SampleWeight.h"

//#define THEIR_EVENTS

//=== MAIN =================================================================================================
int main(int argc, char** argv) 
{
  
  vector<vector<string> > inputFiles;
  inputFiles.push_back(vector<string>());

  //
  // args
  //--------------------------------------------------------------------------------------------------------------
  ControlFlags ctrl;
  parse_args( argc, argv, ctrl );
  if( ctrl.inputfile.empty() || 
      ctrl.inputfile.empty() ) {
    cerr << "usage: applySelection.exe <flags> " << endl;
    cerr << "\tmandoatory flags : " << endl;
    cerr << "\t--inputfile |  file containing a list of ntuples to run over" << endl;
    cerr << "\t--outputfile | file to store  selected evet" << endl;
    return 1;
  }
  ctrl.dump();

  map<string,double> entrymap; // number of unskimmed entries in each file

  //
  // File I/O
  //--------------------------------------------------------------------------------------------------------------
  TChain * chain = new TChain("Events");
  ifstream f(ctrl.inputfile.c_str());
  string fname;
  while (f >> fname) { 
    if( !(strncmp( fname.c_str(), "#", 1 ) ) ) continue; // skip commented lines
    cout << "adding inputfile : " << fname.c_str() << endl;
    entrymap[string(fname.c_str())] = unskimmedEntries(fname.c_str());
    chain->AddFile(fname.c_str());
  }

  // table of cross section values   
  SimpleTable xstab("./data/xs.dat");

  //
  // Setup
  //--------------------------------------------------------------------------------------------------------------
  TH1F * h_evt = new TH1F("hevt", "hevt", 2, 0, 2 ); 
  TH1F * h_evtfail = new TH1F("hevtfail", "hevtfail", 1024, 0, 1024 );
  //  TNtuple * passtuple = new TNtuple( "passtuple", "passtuple", "run:evt:lumi:channel:mZ1:mZ2:m4l:pt4l:w" );
  TTree * passtuple = new TTree("passtuple", "passtuple" );
  unsigned run, evt, lumi, channel;
  float mZ1, mZ2, m4l, pt4l;
  double w;
  passtuple->Branch("run", &run);
  passtuple->Branch("evt", &evt);
  passtuple->Branch("lumi", &lumi);
  passtuple->Branch("mZ1", &mZ1);
  passtuple->Branch("mZ2", &mZ2);
  passtuple->Branch("m4l", &m4l);
  passtuple->Branch("pt4l", &pt4l);
  passtuple->Branch("w", &w);
  initCiCSelection();
  if(ctrl.eleSele=="lik") initLikSelection();
  if(ctrl.eleSele=="bdt") initBDTSelection();
  initRunLumiRangeMap();

  //
  // Loop
  //--------------------------------------------------------------------------------------------------------------
  
  //
  // Access samples and fill histograms
  TFile *inputFile=0;
  TTree *eventTree=0;  
   
  // Data structures to store info from TTrees
  mithep::TEventInfo *info    = new mithep::TEventInfo();
  TClonesArray *electronArr   = new TClonesArray("mithep::TElectron");
  TClonesArray *muonArr       = new TClonesArray("mithep::TMuon");
  

  TBranch *infoBr;
  TBranch *electronBr;
  TBranch *muonBr;

//   chain->SetBranchStatus("*",0); //disable all branches
//   chain->SetBranchStatus("Info", 1);
//   chain->SetBranchStatus("Muon", 1);
//   chain->SetBranchStatus("Electron", 1);

  chain->SetBranchAddress("Info",       &info);
  chain->SetBranchAddress("Electron", &electronArr);
  chain->SetBranchAddress("Muon", &muonArr);


  cout << "nEntries: " << chain->GetEntries() << endl;
  for(UInt_t ientry=0; ientry<chain->GetEntries(); ientry++) {

    chain->GetEntry(ientry);
    if(!(ientry%100000)) cout << "entry: " << ientry << endl;

    // get event weight for this file                                                     
    double eventweight=1;
    if(ctrl.mc) eventweight = getWeight(xstab,entrymap,chain);

#ifdef THEIR_EVENTS
    if( !( info->evtNum == 504867308  ||
           info->evtNum == 368148849  ||
           info->evtNum == 129514273  ||
           info->evtNum == 286336207  ||
           info->evtNum == 344708580  ||
           info->evtNum == 30998576   ||
           info->evtNum == 155679852  ||
           info->evtNum == 394010457  ||
           info->evtNum == 917379387  ||
           info->evtNum == 78213037   ||
           info->evtNum == 337493970  ||
           info->evtNum == 1491724484 ||
           info->evtNum == 480301165  ||
           info->evtNum == 1038911933 ||
           info->evtNum == 876658967  ||
           info->evtNum == 966824024  ||
           info->evtNum == 141954801  ||
           info->evtNum == 160966858  ||
           info->evtNum == 191231387  ||
           info->evtNum == 66033190   ||
           info->evtNum == 10347106   ||
           info->evtNum == 107360878  ) ) continue;
#endif

    unsigned evtfail = fails_HZZ4L_selection(ctrl, info, electronArr, muonArr, eventweight, passtuple);    
    h_evtfail->Fill( evtfail, eventweight );
    if( ctrl.debug ) cout << endl << endl;  
          
  } //end loop over data     


  delete info;
  delete electronArr;
  delete muonArr;


  //--------------------------------------------------------------------------------------------------------------
  // Save Histograms;
  //============================================================================================================== 
  const char * ofname;
  if( strcmp( ctrl.outputfile.c_str(), "") ) {
    ofname = ctrl.outputfile.c_str();
  } else {
    ofname = "tmp.root";
  }

  TFile *file = new TFile(ofname, "RECREATE");
  h_evt->Write();
  h_evtfail->Write();
  passtuple->Write();
  file->Close();
  delete file;

  // map<TString,TMVA::Reader*>::iterator it;
  // for(it=readers.begin(); it!=readers.end(); it++) delete (*it).second;

  cerr << "done!" << endl;
} 


Revision:	1.7
Committed:	Thu Oct 13 14:27:53 2011 UTC (13 years, 7 months ago) by khahn
Content type:	text/plain
Branch:	MAIN
Changes since 1.6:	+2 -1 lines
Log Message:	evtweight is double
#	User	Rev	Content
1	khahn	1.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	dkralph	1.2	#include "HiggsAnaDefs.hh"
36	khahn	1.1	#include "TEventInfo.hh"
37			#include "TElectron.hh"
38			#include "TMuon.hh"
39			#include "TJet.hh"
40			#include "RunLumiRangeMap.h"
41
42			//
43			// utility headers
44			//
45			#include "ParseArgs.h"
46	dkralph	1.2	#include "SampleWeight.h"
47	khahn	1.1	#include "Selection.h"
48			#include "HZZCiCElectronSelection.h"
49	dkralph	1.4	#include "HZZLikelihoodElectronSelection.h"
50	dkralph	1.5	#include "HZZBDTElectronSelection.h"
51	dkralph	1.2	#include "SampleWeight.h"
52	khahn	1.1
53			//#define THEIR_EVENTS
54
55			//=== MAIN =================================================================================================
56			int main(int argc, char** argv)
57			{
58
59			vector<vector<string> > inputFiles;
60			inputFiles.push_back(vector<string>());
61
62			//
63			// args
64			//--------------------------------------------------------------------------------------------------------------
65			ControlFlags ctrl;
66			parse_args( argc, argv, ctrl );
67			if( ctrl.inputfile.empty() \|\|
68			ctrl.inputfile.empty() ) {
69			cerr << "usage: applySelection.exe <flags> " << endl;
70			cerr << "\tmandoatory flags : " << endl;
71			cerr << "\t--inputfile \| file containing a list of ntuples to run over" << endl;
72			cerr << "\t--outputfile \| file to store selected evet" << endl;
73			return 1;
74			}
75			ctrl.dump();
76
77	dkralph	1.4	map<string,double> entrymap; // number of unskimmed entries in each file
78
79	khahn	1.1	//
80			// File I/O
81			//--------------------------------------------------------------------------------------------------------------
82			TChain * chain = new TChain("Events");
83			ifstream f(ctrl.inputfile.c_str());
84			string fname;
85			while (f >> fname) {
86	dkralph	1.4	if( !(strncmp( fname.c_str(), "#", 1 ) ) ) continue; // skip commented lines
87	khahn	1.1	cout << "adding inputfile : " << fname.c_str() << endl;
88	dkralph	1.4	entrymap[string(fname.c_str())] = unskimmedEntries(fname.c_str());
89	khahn	1.1	chain->AddFile(fname.c_str());
90			}
91
92	dkralph	1.2	// table of cross section values
93			SimpleTable xstab("./data/xs.dat");
94
95	khahn	1.1	//
96			// Setup
97			//--------------------------------------------------------------------------------------------------------------
98			TH1F * h_evt = new TH1F("hevt", "hevt", 2, 0, 2 );
99	dkralph	1.4	TH1F * h_evtfail = new TH1F("hevtfail", "hevtfail", 1024, 0, 1024 );
100	khahn	1.6	// TNtuple * passtuple = new TNtuple( "passtuple", "passtuple", "run:evt:lumi:channel:mZ1:mZ2:m4l:pt4l:w" );
101			TTree * passtuple = new TTree("passtuple", "passtuple" );
102			unsigned run, evt, lumi, channel;
103	khahn	1.7	float mZ1, mZ2, m4l, pt4l;
104			double w;
105	khahn	1.6	passtuple->Branch("run", &run);
106			passtuple->Branch("evt", &evt);
107			passtuple->Branch("lumi", &lumi);
108			passtuple->Branch("mZ1", &mZ1);
109			passtuple->Branch("mZ2", &mZ2);
110			passtuple->Branch("m4l", &m4l);
111			passtuple->Branch("pt4l", &pt4l);
112			passtuple->Branch("w", &w);
113	khahn	1.1	initCiCSelection();
114	dkralph	1.5	if(ctrl.eleSele=="lik") initLikSelection();
115			if(ctrl.eleSele=="bdt") initBDTSelection();
116	khahn	1.1	initRunLumiRangeMap();
117
118			//
119			// Loop
120			//--------------------------------------------------------------------------------------------------------------
121
122			//
123			// Access samples and fill histograms
124			TFile *inputFile=0;
125			TTree *eventTree=0;
126
127			// Data structures to store info from TTrees
128			mithep::TEventInfo *info = new mithep::TEventInfo();
129			TClonesArray *electronArr = new TClonesArray("mithep::TElectron");
130			TClonesArray *muonArr = new TClonesArray("mithep::TMuon");
131
132
133			TBranch *infoBr;
134			TBranch *electronBr;
135			TBranch *muonBr;
136
137			// chain->SetBranchStatus("*",0); //disable all branches
138			// chain->SetBranchStatus("Info", 1);
139			// chain->SetBranchStatus("Muon", 1);
140			// chain->SetBranchStatus("Electron", 1);
141
142			chain->SetBranchAddress("Info", &info);
143			chain->SetBranchAddress("Electron", &electronArr);
144			chain->SetBranchAddress("Muon", &muonArr);
145
146
147
148			cout << "nEntries: " << chain->GetEntries() << endl;
149	dkralph	1.5	for(UInt_t ientry=0; ientry<chain->GetEntries(); ientry++) {
150	dkralph	1.4
151	khahn	1.1	chain->GetEntry(ientry);
152	dkralph	1.4	if(!(ientry%100000)) cout << "entry: " << ientry << endl;
153	khahn	1.1
154	dkralph	1.2	// get event weight for this file
155	dkralph	1.4	double eventweight=1;
156			if(ctrl.mc) eventweight = getWeight(xstab,entrymap,chain);
157	dkralph	1.2
158	khahn	1.1	#ifdef THEIR_EVENTS
159			if( !( info->evtNum == 504867308 \|\|
160			info->evtNum == 368148849 \|\|
161			info->evtNum == 129514273 \|\|
162			info->evtNum == 286336207 \|\|
163			info->evtNum == 344708580 \|\|
164			info->evtNum == 30998576 \|\|
165			info->evtNum == 155679852 \|\|
166			info->evtNum == 394010457 \|\|
167			info->evtNum == 917379387 \|\|
168			info->evtNum == 78213037 \|\|
169			info->evtNum == 337493970 \|\|
170			info->evtNum == 1491724484 \|\|
171			info->evtNum == 480301165 \|\|
172			info->evtNum == 1038911933 \|\|
173			info->evtNum == 876658967 \|\|
174			info->evtNum == 966824024 \|\|
175			info->evtNum == 141954801 \|\|
176			info->evtNum == 160966858 \|\|
177			info->evtNum == 191231387 \|\|
178			info->evtNum == 66033190 \|\|
179			info->evtNum == 10347106 \|\|
180			info->evtNum == 107360878 ) ) continue;
181			#endif
182
183	dkralph	1.4	unsigned evtfail = fails_HZZ4L_selection(ctrl, info, electronArr, muonArr, eventweight, passtuple);
184	khahn	1.1	h_evtfail->Fill( evtfail, eventweight );
185	dkralph	1.4	if( ctrl.debug ) cout << endl << endl;
186	khahn	1.1
187			} //end loop over data
188
189
190
191			delete info;
192			delete electronArr;
193			delete muonArr;
194
195
196
197			//--------------------------------------------------------------------------------------------------------------
198			// Save Histograms;
199			//==============================================================================================================
200			const char * ofname;
201			if( strcmp( ctrl.outputfile.c_str(), "") ) {
202			ofname = ctrl.outputfile.c_str();
203			} else {
204			ofname = "tmp.root";
205			}
206
207			TFile *file = new TFile(ofname, "RECREATE");
208			h_evt->Write();
209			h_evtfail->Write();
210			passtuple->Write();
211			file->Close();
212			delete file;
213
214	dkralph	1.5	// map<TString,TMVA::Reader*>::iterator it;
215			// for(it=readers.begin(); it!=readers.end(); it++) delete (*it).second;
216
217	khahn	1.1	cerr << "done!" << endl;
218			}
219
220
221