NonMCBackground/src/select_fakes.cc

#include "select_fakes.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()) {
    cerr << "usage: select_fakes.exe <flags> " << endl;
    cerr << "\tmandoatory flags : " << endl;
    cerr << "\t--inputfile |  file containing a list of ntuples to run over" << 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())] = ctrl.mc ? unskimmedEntries(fname.c_str()): 0;
    assert(chain->AddFile(fname.c_str()));
  }

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

  const Double_t metMax  = 20;

  UInt_t nEvents=0;
  UInt_t nDenom=0;
  UInt_t nNumer=0;

  //
  // Setup
  //--------------------------------------------------------------------------------------------------------------
  const char * ofname;
  if( strcmp( ctrl.outputfile.c_str(), "") ) {
    ofname = ctrl.outputfile.c_str();
  } else {
    ofname = "fo.root";
  }
  // gSystem->mkdir(outputDir,kTRUE);
  // TString outName = outputDir + TString("/") + TString("fo.root");

  TFile *outFile = new TFile(ofname, "RECREATE");
  TTree *outTree = new TTree("FO","FO");
  Float_t fopt, foeta, fophi, fopass, jpt;
  outTree->Branch("pt",   &fopt,  "pt/F");
  outTree->Branch("eta",  &foeta, "eta/F");
  outTree->Branch("phi",  &fophi, "phi/F");
  outTree->Branch("pass", &fopass,"pass/F");
  outTree->Branch("jpt",  &jpt,"jpt/F");

  bdtiface *bdt = new bdtiface("./data","./data/bdtcuts.txt");

  RunLumiRangeMap rlrm;
  rlrm.AddJSONFile(std::string("./data/Cert_136033-149442_7TeV_Apr21ReReco_Collisions10_JSON.txt"));
  //  rlrm.AddJSONFile(std::string("./data/Cert_160404-173244_7TeV_PromptReco_Collisions11_JSON_v2.txt"));
  rlrm.AddJSONFile(std::string("./data/Cert_160404-178078_7TeV_PromptReco_Collisions11_JSON.txt"));
  rlrm.AddJSONFile(std::string("./data/Cert_160404-163869_7TeV_May10ReReco_Collisions11_JSON_v3.txt"));  
  rlrm.AddJSONFile(std::string("./data/Cert_170249-172619_7TeV_ReReco5Aug_Collisions11_JSON.txt"));
  rlrm.AddJSONFile(std::string("./data/Cert_160404-180252_7TeV_PromptReco_Collisions11_JSON.txt"));

  //
  // 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 *jetArr        = new TClonesArray("mithep::TJet");
  
  chain->SetBranchAddress("Info",     &info);
  chain->SetBranchAddress("Electron", &electronArr);
  chain->SetBranchAddress("PFJet",    &jetArr);

  cout << "nEntries: " << chain->GetEntries() << endl;
  UInt_t imax = chain->GetEntries();
  for(UInt_t ientry=0; ientry<imax; ientry++) {
    chain->GetEntry(ientry);
    if(!(ientry%100000)) cout << "entry: " << ientry << "\tfrac: " << setw(15) << float(ientry)/chain->GetEntries() << endl;

    // infoBr->GetEntry(ientry);
    // check for certified runs
    RunLumiRangeMap::RunLumiPairType rl(info->runNum, info->lumiSec);      
    if(!rlrm.HasRunLumi(rl)) continue;  
    
    nEvents++;

    // MET cut
    if(info->pfMET > metMax) continue;
        
    // electronArr->Clear(); electronBr->GetEntry(ientry);

    // count electrons for z veto
    UInt_t nEle=0;
    for(Int_t i=0; i<electronArr->GetEntries(); i++) {
      const mithep::TElectron *ele = (mithep::TElectron*)((*electronArr)[i]);
      if(ele->pt>10) nEle++;
    }

    // jetArr->Clear(); jetBr->GetEntry(ientry);
              
    // count denominator and numerator objects
    for(Int_t i=0; i<electronArr->GetEntriesFast(); i++) {      
      const mithep::TElectron* electron = (mithep::TElectron*)((*electronArr)[i]);

      //find leading jet in the event
      Double_t leadingJetPt = -1;
      for(Int_t j=0; j<jetArr->GetEntries(); j++) {
        const mithep::TJet *jet = (mithep::TJet*)((*jetArr)[j]);        
        if (jet->pt > leadingJetPt &&
            deltaR(jet->eta, jet->phi, electron->eta, electron->phi ) > 1.0) {
          leadingJetPt = jet->pt;          
        }
      }

      if(nEle > 1)                                      continue; //veto events with more than 1 reco electron
      if(leadingJetPt < 25)                             continue; //Jet Threshold Selection
      if(!isLooseEleFO(electron))                       continue;
      if(!(info->triggerBits & kHLT_Ele8))              continue;
      if(!(electron->hltMatchBits & kHLT_Ele8_EleObj))  continue;

      nDenom++;

      bool passID = !failsEleBDTSelection(ctrl,ctrl.eleSeleScheme,electron,bdt);
      bool pass = passID && passPFiso(electron,0);

      if(pass)
        nNumer++;
        
      // fill output TTree
      fopt   = electron->pt;
      foeta  = electron->eta;
      fophi  = electron->phi;
      fopass = pass ? 1 : 0;
      jpt    = leadingJetPt;
      outTree->Fill();
          
    }
  } //end loop over data

  delete info;
  delete electronArr;
  outFile->Write();
  delete outTree;
  outFile->Close();
  delete outFile;

  //--------------------------------------------------------------------------------------------------------------
  // Summary print out
  //==============================================================================================================
  cout << endl;
  cout << "*" << endl;
  cout << "* SUMMARY" << endl;
  cout << "*--------------------------------------------------" << endl;
  cout << endl;

  cout << " >>>    No. of events processed: " << nEvents << endl;
  cout << " >>> No. of denominator objects: " << nDenom << endl;
  cout << " >>>   No. of numerator objects: " << nNumer << endl;     
  cout << " >>>                      ratio: " << float(nNumer)/nDenom << endl;     
  cout << endl;

  TString outputDir("foo");
  ofstream txtfile;
  char txtfname[100];
  sprintf(txtfname,"%s/summary.txt",outputDir.Data());
  txtfile.open(txtfname);
  txtfile << "*" << endl;
  txtfile << "* SUMMARY" << endl;
  txtfile << "*--------------------------------------------------" << endl;
  txtfile << endl;

  txtfile << " >>>    No. of events processed: " << nEvents << endl;
  txtfile << " >>> No. of denominator objects: " << nDenom << endl;
  txtfile << " >>>   No. of numerator objects: " << nNumer << endl;     
  txtfile << " >>>                      ratio: " << float(nNumer)/nDenom << endl;     
  txtfile << endl;  
  txtfile.close();
  
  cout << " <> Output saved in " << outputDir << "/" << endl;      
  cout << endl;
  
  cout << "done!" << endl;

}
Revision:	1.1
Committed:	Sat Dec 17 21:29:32 2011 UTC (13 years, 5 months ago) by dkralph
Content type:	text/plain
Branch:	MAIN
CVS Tags:	synced_FSR_2, synced_FSR, synched2, synched, AN490
Log Message:	* empty log message *
#	Content
1	#include "select_fakes.h"
2
3	// #define THEIR_EVENTS
4
5	//=== MAIN =================================================================================================
6	int main(int argc, char** argv)
7	{
8
9	vector<vector<string> > inputFiles;
10	inputFiles.push_back(vector<string>());
11
12	//
13	// args
14	//--------------------------------------------------------------------------------------------------------------
15	ControlFlags ctrl;
16	parse_args( argc, argv, ctrl );
17	if( ctrl.inputfile.empty()) {
18	cerr << "usage: select_fakes.exe <flags> " << endl;
19	cerr << "\tmandoatory flags : " << endl;
20	cerr << "\t--inputfile \| file containing a list of ntuples to run over" << endl;
21	return 1;
22	}
23	ctrl.dump();
24
25	map<string,double> entrymap; // number of unskimmed entries in each file
26
27	//
28	// File I/O
29	//--------------------------------------------------------------------------------------------------------------
30	TChain * chain = new TChain("Events");
31	ifstream f(ctrl.inputfile.c_str());
32	string fname;
33	while (f >> fname) {
34	if( !(strncmp( fname.c_str(), "#", 1 ) ) ) continue; // skip commented lines
35	cout << "adding inputfile : " << fname.c_str() << endl;
36	entrymap[string(fname.c_str())] = ctrl.mc ? unskimmedEntries(fname.c_str()): 0;
37	assert(chain->AddFile(fname.c_str()));
38	}
39
40	// table of cross section values
41	SimpleTable xstab("./data/xs.dat");
42
43	const Double_t metMax = 20;
44
45	UInt_t nEvents=0;
46	UInt_t nDenom=0;
47	UInt_t nNumer=0;
48
49	//
50	// Setup
51	//--------------------------------------------------------------------------------------------------------------
52	const char * ofname;
53	if( strcmp( ctrl.outputfile.c_str(), "") ) {
54	ofname = ctrl.outputfile.c_str();
55	} else {
56	ofname = "fo.root";
57	}
58	// gSystem->mkdir(outputDir,kTRUE);
59	// TString outName = outputDir + TString("/") + TString("fo.root");
60
61	TFile *outFile = new TFile(ofname, "RECREATE");
62	TTree *outTree = new TTree("FO","FO");
63	Float_t fopt, foeta, fophi, fopass, jpt;
64	outTree->Branch("pt", &fopt, "pt/F");
65	outTree->Branch("eta", &foeta, "eta/F");
66	outTree->Branch("phi", &fophi, "phi/F");
67	outTree->Branch("pass", &fopass,"pass/F");
68	outTree->Branch("jpt", &jpt,"jpt/F");
69
70	bdtiface *bdt = new bdtiface("./data","./data/bdtcuts.txt");
71
72	RunLumiRangeMap rlrm;
73	rlrm.AddJSONFile(std::string("./data/Cert_136033-149442_7TeV_Apr21ReReco_Collisions10_JSON.txt"));
74	// rlrm.AddJSONFile(std::string("./data/Cert_160404-173244_7TeV_PromptReco_Collisions11_JSON_v2.txt"));
75	rlrm.AddJSONFile(std::string("./data/Cert_160404-178078_7TeV_PromptReco_Collisions11_JSON.txt"));
76	rlrm.AddJSONFile(std::string("./data/Cert_160404-163869_7TeV_May10ReReco_Collisions11_JSON_v3.txt"));
77	rlrm.AddJSONFile(std::string("./data/Cert_170249-172619_7TeV_ReReco5Aug_Collisions11_JSON.txt"));
78	rlrm.AddJSONFile(std::string("./data/Cert_160404-180252_7TeV_PromptReco_Collisions11_JSON.txt"));
79
80	//
81	// Access samples and fill histograms
82	TFile *inputFile=0;
83	TTree *eventTree=0;
84
85	// Data structures to store info from TTrees
86	mithep::TEventInfo *info = new mithep::TEventInfo();
87	TClonesArray *electronArr = new TClonesArray("mithep::TElectron");
88	TClonesArray *jetArr = new TClonesArray("mithep::TJet");
89
90	chain->SetBranchAddress("Info", &info);
91	chain->SetBranchAddress("Electron", &electronArr);
92	chain->SetBranchAddress("PFJet", &jetArr);
93
94	cout << "nEntries: " << chain->GetEntries() << endl;
95	UInt_t imax = chain->GetEntries();
96	for(UInt_t ientry=0; ientry<imax; ientry++) {
97	chain->GetEntry(ientry);
98	if(!(ientry%100000)) cout << "entry: " << ientry << "\tfrac: " << setw(15) << float(ientry)/chain->GetEntries() << endl;
99
100	// infoBr->GetEntry(ientry);
101	// check for certified runs
102	RunLumiRangeMap::RunLumiPairType rl(info->runNum, info->lumiSec);
103	if(!rlrm.HasRunLumi(rl)) continue;
104
105	nEvents++;
106
107	// MET cut
108	if(info->pfMET > metMax) continue;
109
110	// electronArr->Clear(); electronBr->GetEntry(ientry);
111
112	// count electrons for z veto
113	UInt_t nEle=0;
114	for(Int_t i=0; i<electronArr->GetEntries(); i++) {
115	const mithep::TElectron ele = (mithep::TElectron)((*electronArr)[i]);
116	if(ele->pt>10) nEle++;
117	}
118
119	// jetArr->Clear(); jetBr->GetEntry(ientry);
120
121	// count denominator and numerator objects
122	for(Int_t i=0; i<electronArr->GetEntriesFast(); i++) {
123	const mithep::TElectron* electron = (mithep::TElectron)((electronArr)[i]);
124
125	//find leading jet in the event
126	Double_t leadingJetPt = -1;
127	for(Int_t j=0; j<jetArr->GetEntries(); j++) {
128	const mithep::TJet jet = (mithep::TJet)((*jetArr)[j]);
129	if (jet->pt > leadingJetPt &&
130	deltaR(jet->eta, jet->phi, electron->eta, electron->phi ) > 1.0) {
131	leadingJetPt = jet->pt;
132	}
133	}
134
135	if(nEle > 1) continue; //veto events with more than 1 reco electron
136	if(leadingJetPt < 25) continue; //Jet Threshold Selection
137	if(!isLooseEleFO(electron)) continue;
138	if(!(info->triggerBits & kHLT_Ele8)) continue;
139	if(!(electron->hltMatchBits & kHLT_Ele8_EleObj)) continue;
140
141	nDenom++;
142
143	bool passID = !failsEleBDTSelection(ctrl,ctrl.eleSeleScheme,electron,bdt);
144	bool pass = passID && passPFiso(electron,0);
145
146	if(pass)
147	nNumer++;
148
149	// fill output TTree
150	fopt = electron->pt;
151	foeta = electron->eta;
152	fophi = electron->phi;
153	fopass = pass ? 1 : 0;
154	jpt = leadingJetPt;
155	outTree->Fill();
156
157	}
158	} //end loop over data
159
160	delete info;
161	delete electronArr;
162	outFile->Write();
163	delete outTree;
164	outFile->Close();
165	delete outFile;
166
167	//--------------------------------------------------------------------------------------------------------------
168	// Summary print out
169	//==============================================================================================================
170	cout << endl;
171	cout << "*" << endl;
172	cout << "* SUMMARY" << endl;
173	cout << "*--------------------------------------------------" << endl;
174	cout << endl;
175
176	cout << " >>> No. of events processed: " << nEvents << endl;
177	cout << " >>> No. of denominator objects: " << nDenom << endl;
178	cout << " >>> No. of numerator objects: " << nNumer << endl;
179	cout << " >>> ratio: " << float(nNumer)/nDenom << endl;
180	cout << endl;
181
182	TString outputDir("foo");
183	ofstream txtfile;
184	char txtfname[100];
185	sprintf(txtfname,"%s/summary.txt",outputDir.Data());
186	txtfile.open(txtfname);
187	txtfile << "*" << endl;
188	txtfile << "* SUMMARY" << endl;
189	txtfile << "*--------------------------------------------------" << endl;
190	txtfile << endl;
191
192	txtfile << " >>> No. of events processed: " << nEvents << endl;
193	txtfile << " >>> No. of denominator objects: " << nDenom << endl;
194	txtfile << " >>> No. of numerator objects: " << nNumer << endl;
195	txtfile << " >>> ratio: " << float(nNumer)/nDenom << endl;
196	txtfile << endl;
197	txtfile.close();
198
199	cout << " <> Output saved in " << outputDir << "/" << endl;
200	cout << endl;
201
202	cout << "done!" << endl;
203
204	}