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 *
#	User	Rev	Content
1	dkralph	1.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			}