ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/Selection/src/applySelection.cc
Revision: 1.55
Committed: Wed Sep 12 12:13:27 2012 UTC (12 years, 8 months ago) by anlevin
Content type: text/plain
Branch: MAIN
Changes since 1.54: +2 -2 lines
Log Message:
changed stuff for the electron regression

File Contents

# 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 khahn 1.13 #include <bitset>
8 khahn 1.19 #include <map>
9 khahn 1.1 using namespace std;
10    
11     //
12 khahn 1.13 // root headers
13 khahn 1.1 //
14 khahn 1.13 #include <TFile.h>
15     #include <TTree.h>
16     #include <TChain.h>
17     #include <TBranch.h>
18     #include <TClonesArray.h>
19 khahn 1.1
20 khahn 1.11
21 khahn 1.1 //
22 khahn 1.13 // TMVA
23     //
24     #include "TMVA/Reader.h"
25     #include "TMVA/Tools.h"
26     #include "TMVA/Config.h"
27     #include "TMVA/MethodBDT.h"
28    
29     //
30 khahn 1.1 // ntuple format headers
31     //
32 khahn 1.13 #include "EventHeader.h"
33     #include "Electron.h"
34     #include "Muon.h"
35     #include "Vertex.h"
36 dkralph 1.48 #include "PFMet.h"
37 khahn 1.13 #include "PFCandidate.h"
38     #include "PFCandidateCol.h"
39 khahn 1.16 #include "PileupInfoCol.h"
40     #include "PileupEnergyDensity.h"
41 khahn 1.23 #include "MCParticle.h"
42 khahn 1.13 #include "TriggerMask.h"
43     #include "TriggerTable.h"
44 khahn 1.41 #include "TriggerObject.h"
45     #include "TriggerObjectRel.h"
46 khahn 1.13 #include "Names.h"
47     #include "BaseMod.h"
48 dkralph 1.44 #include "TriggerObjectsTable.h"
49 khahn 1.43
50 khahn 1.1 //
51 khahn 1.13 // our headers
52 khahn 1.1 //
53     #include "ParseArgs.h"
54 khahn 1.13 #include "MuonSelection.h"
55     #include "ElectronSelection.h"
56     #include "IsolationSelection.h"
57 khahn 1.29 #include "ReferenceSelection.h"
58 khahn 1.14
59 khahn 1.13 #include "TriggerUtils.h"
60     #include "PassHLT.h"
61     #include "Angles.h"
62     #include "KinematicsStruct.h"
63     #include "InfoStruct.h"
64 khahn 1.28 #include "GenInfoStruct.h"
65 khahn 1.13 #include "WeightStruct.h"
66     #include "AngleTuple.h"
67 khahn 1.15 #include "FOTuple.h"
68 khahn 1.13
69 khahn 1.33 #include "RunLumiRangeMap.h"
70 dkralph 1.2 #include "SampleWeight.h"
71 khahn 1.13 #include "EfficiencyWeightsInterface.h"
72 dkralph 1.38 #include "SelectionFuncs.h"
73 khahn 1.1
74 khahn 1.16 #include "SimpleLepton.h"
75    
76 khahn 1.13 #ifndef CMSSW_BASE
77     #define CMSSW_BASE "../"
78     #endif
79 khahn 1.1
80 khahn 1.13 using namespace mithep;
81 khahn 1.11
82 khahn 1.16 //
83     // globals
84     //----------------------------------------------------------------------------
85 khahn 1.47 TH1D *hpu_2011, *hpu_2012;
86 khahn 1.16 RunLumiRangeMap rlrm;
87     vector<SimpleLepton> failingLeptons ; // for fake estimation
88     vector<SimpleLepton> passingLeptons; // for fake estimation
89 khahn 1.19 vector<unsigned> cutvec;
90     vector<vector<unsigned> > zcutvec;
91     vector<vector<unsigned> > zzcutvec;
92     map<unsigned,float> evtrhoMap;
93     vector<string> cutstrs;
94 dkralph 1.53 bool passes_HLT_MC,passes_HLT;
95 khahn 1.22 vector<bool> PFnoPUflag;;
96 khahn 1.11
97 khahn 1.16 //
98     // prototypes
99 khahn 1.13 //----------------------------------------------------------------------------
100 dkralph 1.48 void initRunLumiRangeMap(ControlFlags &ctrl);
101 khahn 1.13 //----------------------------------------------------------------------------
102 khahn 1.11
103 khahn 1.16
104     //
105     // MAIN
106     //----------------------------------------------------------------------------
107 khahn 1.1 int main(int argc, char** argv)
108 khahn 1.16 //----------------------------------------------------------------------------
109 khahn 1.1 {
110 khahn 1.19 cutstrs.push_back(string("skim0")); //0
111     cutstrs.push_back(string("skim1")); //1
112 anlevin 1.25 cutstrs.push_back(string("skim2")); //1
113 khahn 1.21 cutstrs.push_back(string("trigger")); //2
114 khahn 1.19 // -------------------------------------------------
115 khahn 1.21 cutstrs.push_back(string("Z candidate")); //3
116     cutstrs.push_back(string("good Z1")); //4
117 khahn 1.19 // -------------------------------------------------
118 khahn 1.21 cutstrs.push_back(string("4l")); //5
119     cutstrs.push_back(string("ZZ candidate")); //6
120     cutstrs.push_back(string("good Z2")); //7
121     cutstrs.push_back(string("ZZ 20/10")); //8
122     cutstrs.push_back(string("resonance")); //9
123     cutstrs.push_back(string("m4l>70")); //10
124 anlevin 1.25 cutstrs.push_back(string("m4l>100, mZ2 > 12")); //11
125 khahn 1.19
126    
127    
128    
129 khahn 1.13 string cmsswpath(CMSSW_BASE);
130     cmsswpath.append("/src");
131     std::bitset<1024> triggerBits;
132 khahn 1.1 vector<vector<string> > inputFiles;
133     inputFiles.push_back(vector<string>());
134 khahn 1.16 map<string,unsigned> entrymap; // number of unskimmed entries in each file
135 khahn 1.21 for( int i=0; i<cutstrs.size(); i++ ) cutvec.push_back(0);
136 khahn 1.19 for( int i=0; i<2; i++ ) {
137     zcutvec.push_back(vector<unsigned>());
138 khahn 1.21 for( int j=0; j<cutstrs.size(); j++ ) zcutvec[i].push_back(0);
139 khahn 1.19 }
140     for( int i=0; i<3; i++ ) {
141     zzcutvec.push_back(vector<unsigned>());
142 khahn 1.21 for( int j=0; j<cutstrs.size(); j++ ) zzcutvec[i].push_back(0);
143 khahn 1.19 }
144    
145 khahn 1.1
146     //
147     // args
148     //--------------------------------------------------------------------------------------------------------------
149     ControlFlags ctrl;
150     parse_args( argc, argv, ctrl );
151 khahn 1.13 if( ctrl.inputfiles.empty() &&ctrl.inputfile.empty() )
152     {
153     cerr << "usage: applySelection.exe <flags> " << endl;
154     cerr << "\tmandoatory flags : " << endl;
155     cerr << "\t--inputfiles | file containing a list of ntuples to run over" << endl;
156     cerr << "\t--inputfile | a file to run over" << endl;
157     cerr << "\t--outputfile | file to store selected evet" << endl;
158     return 1;
159     }
160 khahn 1.1 ctrl.dump();
161 khahn 1.13
162 khahn 1.1
163 dkralph 1.4
164 khahn 1.1 //
165     // File I/O
166     //--------------------------------------------------------------------------------------------------------------
167     TChain * chain = new TChain("Events");
168 khahn 1.13 TChain * hltchain = new TChain("HLT");
169 dkralph 1.53 TChain * runchain = new TChain(Names::gkRunTreeName);
170 khahn 1.13
171 khahn 1.1 string fname;
172 dkralph 1.44 unsigned total_unskimmed=0;
173 khahn 1.13 if( !ctrl.inputfiles.empty() ) {
174     ifstream f(ctrl.inputfiles.c_str());
175     while (f >> fname) {
176     if( !(strncmp( fname.c_str(), "#", 1 ) ) ) continue; // skip commented lines
177     cout << "adding inputfile : " << fname.c_str() << endl;
178     entrymap[string(fname.c_str())] = unskimmedEntries(fname.c_str());
179 khahn 1.16 cout << "unskimmed entries: " << entrymap[string(fname.c_str())] << endl;
180 dkralph 1.44 total_unskimmed += entrymap[string(fname.c_str())];
181 khahn 1.13 chain->AddFile(fname.c_str());
182     hltchain->AddFile(fname.c_str());
183 dkralph 1.53 runchain->AddFile(fname.c_str());
184 khahn 1.13 }
185     } else {
186     cout << "adding inputfile : " << ctrl.inputfile.c_str() << endl;
187 dkralph 1.44 unsigned unsk_ents = unskimmedEntries(ctrl.inputfile.c_str());
188     entrymap[string(ctrl.inputfile.c_str())] = unsk_ents;
189     cout << "unskimmed entries: " << unsk_ents << endl;
190     total_unskimmed += unsk_ents;
191 khahn 1.13 chain->AddFile(ctrl.inputfile.c_str());
192     hltchain->AddFile(ctrl.inputfile.c_str());
193 dkralph 1.53 runchain->AddFile(ctrl.inputfile.c_str());
194 khahn 1.1 }
195 dkralph 1.44 // // write the total number of unskimmed events that went into making this output file to a text file
196 dkralph 1.53 writeEntries(ctrl,total_unskimmed);
197 khahn 1.1
198 dkralph 1.10 const char * ofname;
199     if( strcmp( ctrl.outputfile.c_str(), "") ) {
200     ofname = ctrl.outputfile.c_str();
201     } else {
202     ofname = "tmp.root";
203     }
204 khahn 1.13 // table of cross section values
205     string xspath = (cmsswpath+string("/MitPhysics/data/xs.dat"));
206     cout << "xspath: " << xspath.c_str() << endl;
207     SimpleTable xstab(xspath.c_str());
208 khahn 1.11
209 khahn 1.1
210 dkralph 1.44
211 khahn 1.1 //
212 khahn 1.13 // Setup
213 khahn 1.1 //--------------------------------------------------------------------------------------------------------------
214 khahn 1.13 Angles angles;
215     KinematicsStruct kinematics;
216     InfoStruct evtinfo;
217     WeightStruct weights;
218 khahn 1.28 GenInfoStruct geninfo;
219 khahn 1.13
220     AngleTuple nt( (const char*)ofname, (const char*)"zznt");
221     nt.makeAngleBranch(angles);
222     nt.makeKinematicsBranch(kinematics);
223     nt.makeInfoBranch(evtinfo);
224 khahn 1.17
225 khahn 1.33
226 khahn 1.17 FOTuple foTree( nt.getFile(), (const char*)"FOtree");
227 khahn 1.52 foTree.makeInfoBranch(evtinfo);
228 khahn 1.16 foTree.makeFailedLeptonBranch(failingLeptons);
229 khahn 1.33 foTree.makePassedLeptonBranch(passingLeptons);
230 khahn 1.13
231 khahn 1.28 TH1F * h_wf11_hpt;
232 khahn 1.13 if(ctrl.mc) {
233     nt.makeWeightBranch(weights);
234 khahn 1.35 if(ctrl.fillGen ) nt.makeGenInfoBranch(geninfo);
235 anlevin 1.55 //initEfficiencyWeights();
236 khahn 1.16 initPUWeights();
237 khahn 1.28
238 khahn 1.13 // Higgs only, pt reweighting
239     if( ctrl.mH <= 140 ) {
240     char wptname[256];
241     sprintf( wptname, "/data/blue/pharris/Flat/ntupler/root/weight_ptH_%i.root", ctrl.mH );
242     TFile * f_hpt = new TFile(wptname);
243     f_hpt->Print();
244     sprintf(wptname, "weight_hqt_fehipro_fit_%i", ctrl.mH);
245 khahn 1.28 h_wf11_hpt = (TH1F*)(f_hpt->FindObjectAny(wptname));
246 khahn 1.13 }
247 khahn 1.28
248 khahn 1.13 } else {
249 khahn 1.50 initRunLumiRangeMap(ctrl);
250 khahn 1.13 }
251    
252 khahn 1.14 // initMuonIDMVA();
253 khahn 1.33 initElectronIDMVAV1();
254 khahn 1.13 initTrigger();
255 dkralph 1.53 TrigInfo ti;
256     initAnalysisTriggers(ti);
257 khahn 1.19
258 khahn 1.16 //
259 khahn 1.13 // Setup tree I/O
260 khahn 1.16 //--------------------------------------------------------------------------------------------------------------
261 khahn 1.1 TFile *inputFile=0;
262     TTree *eventTree=0;
263 khahn 1.13 string currentFile("");
264    
265 dkralph 1.44 UInt_t fNMaxTriggers = 1024;
266 dkralph 1.53 mithep::EventHeader *info = new mithep::EventHeader();
267 dkralph 1.48 mithep::Array<mithep::PFMet> *metArr = new mithep::Array<mithep::PFMet>();
268 khahn 1.13 mithep::Array<mithep::Electron> *electronArr = new mithep::Array<mithep::Electron>();
269     mithep::Array<mithep::Muon> *muonArr = new mithep::Array<mithep::Muon>();
270 khahn 1.28 mithep::Array<mithep::Vertex> *vtxArr = new mithep::Array<mithep::Vertex>();
271 khahn 1.13 mithep::Array<mithep::PFCandidate> *pfArr = new mithep::Array<mithep::PFCandidate>();
272 khahn 1.16 mithep::Array<mithep::PileupInfo> *puArr = new mithep::Array<mithep::PileupInfo>();
273     mithep::Array<mithep::PileupEnergyDensity> *puDArr = new mithep::Array<mithep::PileupEnergyDensity>();
274 khahn 1.13 mithep::Array<mithep::Track> *trkArr = new mithep::Array<mithep::Track>();
275 khahn 1.23 mithep::Array<mithep::MCParticle> *mcArr = new mithep::Array<mithep::MCParticle>();
276 khahn 1.28 mithep::MCEventInfo *mcEvtInfo = new mithep::MCEventInfo();
277 khahn 1.13 mithep::TriggerMask *trigMask = new mithep::TriggerMask();
278 dkralph 1.44 mithep::TriggerTable *hltTable = new mithep::TriggerTable(fNMaxTriggers);
279 khahn 1.13 vector<string> *hltTableStrings = new vector<string>();
280 khahn 1.41 vector<string> *hltLabelStrings = new vector<string>();
281     mithep::Array<mithep::TriggerObject> *hltObjArr = new mithep::Array<mithep::TriggerObject>();
282     mithep::Array<mithep::TriggerObjectRel> *hltRelsArr = new mithep::Array<mithep::TriggerObjectRel>();
283     std::vector<std::string> *hltTab = new vector<string>();
284 dkralph 1.44 std::vector<std::string> *hltLab = new vector<string>();
285     TriggerObjectsTable *fTrigObjs = new TriggerObjectsTable(hltTable,fNMaxTriggers);
286 khahn 1.41
287 dkralph 1.53 RunInfo *runInfo = new RunInfo();
288 khahn 1.13
289     TString fElectronName(Names::gkElectronBrn);
290     TString fMuonName(Names::gkMuonBrn);
291     TString fInfoName(Names::gkEvtHeaderBrn);
292 dkralph 1.48 TString fPFMetName("PFMet");
293 khahn 1.13 TString fPrimVtxName(Names::gkPVBrn);
294     TString fPFCandidateName(Names::gkPFCandidatesBrn);
295 khahn 1.16 TString fPileupInfoName(Names::gkPileupInfoBrn);
296 khahn 1.13 TString fPileupEnergyDensityName(Names::gkPileupEnergyDensityBrn);
297     TString fTrackName(Names::gkTrackBrn);
298 khahn 1.23 TString fMCParticleName(Names::gkMCPartBrn);
299 khahn 1.28 TString fMCEvtInfoName(Names::gkMCEvtInfoBrn);
300 khahn 1.13 TString fTriggerMaskName(Names::gkHltBitBrn);
301     TString fTriggerTableName(Names::gkHltTableBrn);
302 khahn 1.41 TString fTriggerLabelName(Names::gkHltLabelBrn);
303     TString fTriggerObjectName(Names::gkHltObjBrn);
304     TString fTriggerObjectRelsName(Form("HLTObjectsRelation"));
305    
306 khahn 1.43 // TString fHLTLabName(Names::gkHltLabelBrn);
307     // TString fHLTTabName(Names::gkHltTableBrn);
308     // TString fHLTTabNamePub(Form("%sFwk",fHLTTabName.Data()));
309     // TString fHLTLabNamePub(Form("%sFwk",fHLTLabName.Data()));
310     // TString fObjsNamePub(Form("%sFwk",fTriggerObjectName.Data()));
311 khahn 1.41
312    
313 dkralph 1.44
314 khahn 1.13 chain->SetBranchAddress(fInfoName, &info);
315 dkralph 1.48 chain->SetBranchAddress(fPFMetName, &metArr);
316 khahn 1.13 chain->SetBranchAddress(fElectronName, &electronArr);
317     chain->SetBranchAddress(fMuonName, &muonArr);
318     chain->SetBranchAddress(fPrimVtxName, &vtxArr);
319     chain->SetBranchAddress(fPFCandidateName, &pfArr);
320 khahn 1.23 if( ctrl.mc ) {
321     chain->SetBranchAddress(fPileupInfoName, &puArr);
322     chain->SetBranchAddress(fMCParticleName, &mcArr);
323 khahn 1.28 chain->SetBranchAddress(fMCEvtInfoName, &mcEvtInfo);
324 khahn 1.23 }
325 khahn 1.16 chain->SetBranchAddress(fPileupEnergyDensityName, &puDArr);
326 khahn 1.13 chain->SetBranchAddress(fTrackName, &trkArr);
327     chain->SetBranchAddress(fTriggerMaskName, &trigMask);
328 khahn 1.41 chain->SetBranchAddress(fTriggerObjectName, &hltObjArr);
329     chain->SetBranchAddress(fTriggerObjectRelsName, &hltRelsArr);
330 khahn 1.43 // chain->SetBranchAddress(fHLTTabNamePub, &hltTab);
331     // chain->SetBranchAddress(fHLTLabNamePub, &hltLab);
332 khahn 1.41
333 khahn 1.13 cout << "hlttable: " << fTriggerTableName << endl;
334     hltchain->SetBranchAddress(fTriggerTableName, &hltTableStrings);
335 khahn 1.41 hltchain->SetBranchAddress(fTriggerLabelName, &hltLabelStrings);
336 khahn 1.1
337 dkralph 1.53 runchain->SetBranchAddress(Names::gkRunInfoBrn, &runInfo);
338    
339 khahn 1.13 mithep::Vertex vtx; // best primary vertex in the event
340 khahn 1.1
341 khahn 1.13 // ginfo = NULL;
342     // if( ctrl.mc ) { chain->SetBranchAddress("Gen", &ginfo);}
343 khahn 1.1
344 khahn 1.9 int count=0, pass=0;
345 khahn 1.13 float passcorr=0., passcorr_errup=0., passcorr_errdown=0.;
346 khahn 1.11 float denom[3]={0.,0.,0.};
347     float numer[3]={0.,0.,0.};
348     float numercorr[3]={0.,0.,0.};
349    
350 khahn 1.13 // only 1 HLT table / file ???
351     hltchain->GetEntry(0);
352    
353 anlevin 1.54 int imax;
354     //-1 means process all of the events
355     if(ctrl.n_events_to_process != -2 && ctrl.n_events_to_process != -1){
356     if(ctrl.n_events_to_process > chain->GetEntries()){
357     cerr << "error: n_events_to_process is larger than chain->GetEntries()" << endl;
358     return 1;
359     }
360     imax = ctrl.n_events_to_process;
361     }
362     else
363     imax = chain->GetEntries();
364    
365     cout << "number of events to process: " << imax << endl;
366 khahn 1.11
367 khahn 1.16 //
368     // Loop !!!!!!!!!
369     //--------------------------------------------------------------------------------------------------------------
370 dkralph 1.53 int lastHltEntry=-2;
371 khahn 1.13 for(UInt_t ientry=0; ientry<imax; ientry++)
372     {
373     chain->GetEntry(ientry);
374     if(!(ientry%1000)) cerr << "entry: " << ientry << endl;
375 dkralph 1.44
376 khahn 1.41 if( ctrl.debug ) {
377     cout << "-----------------------------------------------------------------" << endl;
378     cout << "-----------------------------------------------------------------" << endl;
379     cout << "Run: " << info->RunNum()
380     << "\tEvt: " << info->EvtNum()
381     << "\tLumi: " << info->LumiSec()
382     << endl;
383 dkralph 1.53 cerr << "Run: " << info->RunNum()
384     << "\tEvt: " << info->EvtNum()
385     << "\tLumi: " << info->LumiSec()
386     << endl;
387 khahn 1.41 cout << "-----------------------------------------------------------------" << endl;
388     }
389    
390 dkralph 1.53 runchain->GetEvent(info->RunEntry());
391     hltchain->GetEntry(runInfo->HltEntry());
392     hltTable->Clear();
393     fillTriggerNames(hltTable, hltTableStrings );
394     if(ctrl.debug && (lastHltEntry != runInfo->HltEntry())) hltTable->Print();
395     lastHltEntry = runInfo->HltEntry();
396     fillTriggerBits( hltTable, trigMask, triggerBits );
397     setHLTObjectRelations( hltObjArr, hltRelsArr, hltTableStrings, hltLabelStrings, fTrigObjs );
398     // printTriggerObjs( hltTable, fTrigObjs);
399    
400 khahn 1.28 string fname = string(chain->GetFile()->GetEndpointUrl()->GetFile());
401 khahn 1.49 if(ctrl.debug) cout << "era is " << ctrl.era << endl;
402     if( ctrl.era == 0 ) {
403     setEra( fname, ctrl );
404     if(ctrl.debug) cout << "era is now " << ctrl.era << endl;
405     }
406 khahn 1.41
407 khahn 1.22 //
408     // pfNoPU
409     //
410     //mithep::Array<PFCandidate> pfNoPileUp;
411     PFnoPUflag.clear();
412     int pfnopu_size = makePFnoPUArray( pfArr, PFnoPUflag, vtxArr );
413     assert(pfnopu_size == pfArr->GetEntries());
414 khahn 1.13
415 khahn 1.16 //
416 dkralph 1.53 // trigger
417     //
418     if( string(chain->GetFile()->GetEndpointUrl()->GetFile()) != currentFile ) {
419     currentFile = string(chain->GetFile()->GetEndpointUrl()->GetFile());
420     }
421     if( ctrl.debug ) cout << "file is : " << currentFile << endl;
422    
423     //
424 khahn 1.16 // data/MC
425     //
426 khahn 1.13 if(ctrl.mc) {
427 khahn 1.16 //
428 khahn 1.28 // gen info
429     //
430 khahn 1.31 if( ctrl.fillGen )
431 khahn 1.28 fillGenInfo( mcArr, mcEvtInfo, geninfo, ESampleType::kHZZ, ctrl);
432    
433     //
434 khahn 1.16 // xsec & PU weights
435     //
436 khahn 1.31 weights.w = getWeight(xstab,entrymap,chain)/ctrl.totalMC;
437 khahn 1.49 cout << "xsec weight: " << weights.w << endl;
438    
439     unsigned npu = getNPU(puArr, 0);
440     weights.npu = npu;
441     weights.npuw = getPUWeight(ctrl.era,npu);
442     cerr << "weights.npuw: " << weights.npuw << endl;
443    
444 khahn 1.28 //
445     // pt reweighting for Higgs < 140, F11
446     //
447     if( ctrl.fillGen ) {
448     if( (fname.find("f11-h115zz4l") != string::npos) ||
449     (fname.find("f11-h120zz4l") != string::npos) ||
450     (fname.find("f11-h130zz4l") != string::npos) ) {
451     weights.w *= h_wf11_hpt->GetBinContent(h_wf11_hpt->FindBin(geninfo.pt_zz));
452     }
453     }
454    
455 dkralph 1.53 // trigger
456 khahn 1.21 passes_HLT_MC = true; // passed to apply as extern global, just for sync
457 anlevin 1.27 if( !passHLTMC(triggerBits, info->RunNum(), info->EvtNum(), k2012_MC ) ) {
458 khahn 1.21 passes_HLT_MC = false;
459 khahn 1.20 }
460 khahn 1.41
461 khahn 1.13 } else {
462 khahn 1.16 //
463     // JSON
464     //
465 khahn 1.40 if(!(ctrl.noJSON) ) {
466 khahn 1.36 RunLumiRangeMap::RunLumiPairType rl(info->RunNum(), info->LumiSec());
467     if( !(rlrm.HasRunLumi(rl)) ) {
468 khahn 1.41 if( ctrl.debug ) cout << "fails JSON" << endl;
469 khahn 1.36 continue;
470 khahn 1.40 }
471 khahn 1.16 }
472 dkralph 1.44
473 khahn 1.16 // trigger
474 dkralph 1.53 passes_HLT = true;
475 anlevin 1.54 //if( !passHLT(ctrl, triggerBits,ti,info->RunNum()) ) {
476     // passes_HLT = false;
477     //}
478 khahn 1.41
479 khahn 1.13 }
480 khahn 1.16
481 khahn 1.13 //
482     // lepton/kinematic selection ...
483     //
484 khahn 1.16 failingLeptons.clear();
485     passingLeptons.clear();
486 khahn 1.13 EventData ret4l =
487 khahn 1.21
488 khahn 1.23
489     // reference
490 khahn 1.29 apply_HZZ4L_reference_selection(ctrl, info,
491 khahn 1.33 vtxArr,
492     pfArr,
493     puDArr,
494     electronArr,
495     &electronReferencePreSelection,
496     &electronReferenceIDMVASelectionV1,
497     &electronReferenceIsoSelection,
498     muonArr,
499     &muonReferencePreSelection,
500     &muonIDPFSelection,
501     &muonReferenceIsoSelection);
502    
503 khahn 1.13 if( ctrl.debug ) cout << endl;
504 khahn 1.14
505 khahn 1.19
506     int Z1type=0, Z2type=0;
507     if( ret4l.status.selectionBits.to_ulong() >= PASS_ZCANDIDATE ) {
508     if( abs(ret4l.Z1leptons[0].type) == 11 ) Z1type = 11;
509     if( abs(ret4l.Z1leptons[0].type) == 13 ) Z1type = 13;
510     }
511     if( ret4l.status.selectionBits.to_ulong() >= PASS_ZZCANDIDATE ) {
512     if( abs(ret4l.Z2leptons[0].type) == 11 ) Z2type = 11;
513     if( abs(ret4l.Z2leptons[0].type) == 13 ) Z2type = 13;
514     }
515    
516 khahn 1.13 if( ret4l.status.pass() ) {
517     fillAngles(ret4l,angles);
518     fillKinematics(ret4l,kinematics);
519 khahn 1.31 fillMassErrors(ret4l,muonArr,electronArr,kinematics);
520 dkralph 1.48 TLorentzVector pfmet; pfmet.SetPxPyPzE(metArr->At(0)->Mex(),metArr->At(0)->Mey(),0,0);
521     fillEventInfo( info, pfmet, evtinfo, ctrl.mc ? getNPU(puArr) : 0);
522 khahn 1.16 if( ctrl.mc) {
523 anlevin 1.55 //setEffiencyWeights(ctrl.era, ret4l, triggerBits, hltTable, hltObjArr, fTrigObjs, weights);
524 dkralph 1.44 if(ctrl.debug)
525     cout << "w: " << weights.w << "\t"
526 khahn 1.49 << "npuw: " << weights.npuw << "\t"
527 dkralph 1.44 << "won: " << weights.won << "\t"
528     << "woff: " << weights.woff << endl;
529 khahn 1.16 }
530 khahn 1.13
531 khahn 1.28
532 khahn 1.13 nt.Fill();
533    
534     cerr << "PASS:: "
535     << "\trun: " << evtinfo.run
536     << "\tevt: " << evtinfo.evt
537     << "\tlumi: " << evtinfo.lumi
538 khahn 1.45 << "\tchannel: " << kinematics.channel
539     << "\tm4l: " << kinematics.m4l
540     << "\tmZ1: " << kinematics.mZ1
541     << "\tmZ2: " << kinematics.mZ2
542 khahn 1.13 << endl;
543     pass++;
544     // if( pass > 3 ) break;
545 khahn 1.8
546 khahn 1.33 } else if( ret4l.status.selectionBits.test(PASS_GOODZ1) &&
547 khahn 1.36 ret4l.status.selectionBits.to_ulong() < (1UL<<PASS_4L)) { // save for fakes ...
548     if(ctrl.debug) {
549 khahn 1.18 cout << "failingLeptons : " << failingLeptons.size() << endl;
550     for( int f=0; f<failingLeptons.size(); f++ ) {
551     SimpleLepton l = failingLeptons[f];
552     cout << "f: " << f << "\t"
553     << "type: " << l.type << "\t"
554     << "pT: " << l.vec.Pt() << "\t"
555     << "eta: " << l.vec.Eta()
556     << endl;
557 khahn 1.36 }
558 khahn 1.16 }
559 khahn 1.52 TLorentzVector pfmet; pfmet.SetPxPyPzE(metArr->At(0)->Mex(),metArr->At(0)->Mey(),0,0);
560     fillEventInfo( info, pfmet, evtinfo, ctrl.mc ? getNPU(puArr) : 0);
561 khahn 1.16 foTree.Fill();
562 khahn 1.36 } else {
563 dkralph 1.48 if(ctrl.debug)
564     cout << "failing with some other code : " << ret4l.status.selectionBits << endl;
565 khahn 1.15 }
566 dkralph 1.53 if(ctrl.debug)
567     cout << "passes HLT: " << passes_HLT << endl;
568     }
569 khahn 1.1
570 khahn 1.16
571     foTree.getFile()->cd();
572     foTree.getTree()->Write();
573 khahn 1.13 nt.WriteClose();
574 khahn 1.19
575 dkralph 1.48 if(ctrl.debug)
576     cout << endl << endl;
577 khahn 1.19 for( int i=0; i<cutvec.size(); i++ ) {
578     cout << "cut: " << i << "\t"
579     << "pass: " << cutvec[i];
580    
581 khahn 1.21 if( i>PASS_TRIGGER&&i<=PASS_GOODZ1 )
582 khahn 1.19 cout << "\t2e: " << zcutvec[0][i]
583     << "\t2m: " << zcutvec[1][i];
584    
585 khahn 1.21 if( i>PASS_ZCANDIDATE )
586 khahn 1.19 cout << "\t4e: " << zzcutvec[0][i]
587     << "\t4m: " << zzcutvec[1][i]
588     << "\t2e2m: " << zzcutvec[2][i];
589    
590     cout << "\t" << cutstrs[i] << endl;
591    
592     cout << endl;
593     }
594    
595 khahn 1.13 }
596 dkralph 1.10
597 khahn 1.16
598     //----------------------------------------------------------------------------
599 dkralph 1.48 void initRunLumiRangeMap(ControlFlags &ctrl)
600 khahn 1.16 //----------------------------------------------------------------------------
601     {
602     /*
603     rlrm.AddJSONFile(std::string("./data/Cert_136033-149442_7TeV_Apr21ReReco_Collisions10_JSON.txt"));
604     // rlrm.AddJSONFile(std::string("./data/Cert_160404-173244_7TeV_PromptReco_Collisions11_JSON_v2.txt"));
605     rlrm.AddJSONFile(std::string("./data/Cert_160404-178078_7TeV_PromptReco_Collisions11_JSON.txt"));
606     rlrm.AddJSONFile(std::string("./data/Cert_160404-163869_7TeV_May10ReReco_Collisions11_JSON_v3.txt"));
607     rlrm.AddJSONFile(std::string("./data/Cert_170249-172619_7TeV_ReReco5Aug_Collisions11_JSON.txt"));
608     // r11b
609     rlrm.AddJSONFile(std::string("./data/Cert_160404-180252_7TeV_PromptReco_Collisions11_JSON.txt"));
610     */
611 khahn 1.36
612     // 2012 only for now ...
613 dkralph 1.48 std::string jsonFile;
614     if(ctrl.jsonFile=="")
615     jsonFile = "./data/Cert_190456-194479_8TeV_PromptReco_Collisions12_JSON.txt";
616     else
617     jsonFile = ctrl.jsonFile;
618     rlrm.AddJSONFile(jsonFile);
619 khahn 1.16 };