ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/Selection/src/applySelection.cc
(Generate patch)

Comparing UserCode/MitHzz4l/Selection/src/applySelection.cc (file contents):
Revision 1.3 by khahn, Wed Sep 14 21:43:21 2011 UTC vs.
Revision 1.60 by dkralph, Tue Oct 23 11:39:21 2012 UTC

# Line 4 | Line 4
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>
7 > #include <bitset>  
8 > #include <map>  
9   using namespace std;
10  
11   //
12 < // ROOT headers
12 > // root headers
13   //
14 < #include <TROOT.h>                  // access to gROOT, entry point to ROOT system
15 < #include <TSystem.h>                // interface to OS
16 < #include <TFile.h>                  // file handle class
17 < #include <TNtuple.h>
18 < #include <TTree.h>                  // class to access ntuples
19 < #include <TChain.h>                 //
20 < #include <TBranch.h>                // class to access branches in TTree
21 < #include <TClonesArray.h>           // ROOT array class
22 < #include <TCanvas.h>                // class for drawing
23 < #include <TH1F.h>                   // 1D histograms
24 < #include <TBenchmark.h>             // class to track macro running statistics
25 < #include <TLorentzVector.h>         // 4-vector class
26 < #include <TVector3.h>               // 3D vector class
14 > #include <TFile.h>
15 > #include <TRegexp.h>
16 > #include <TTree.h>                  
17 > #include <TChain.h>                
18 > #include <TBranch.h>                
19 > #include <TClonesArray.h>          
20 >
21 >
22 > //
23 > // TMVA
24 > //
25 > #include "TMVA/Reader.h"
26 > #include "TMVA/Tools.h"
27 > #include "TMVA/Config.h"
28 > #include "TMVA/MethodBDT.h"
29 >
30 > #include "FactorizedJetCorrector.h"
31 > #include "JetCorrectorParameters.h"
32  
33   //
34   // ntuple format headers
35   //
36 < #include "HiggsAnaDefs.hh"
37 < #include "TEventInfo.hh"
38 < #include "TElectron.hh"
39 < #include "TMuon.hh"
40 < #include "TJet.hh"
41 < #include "RunLumiRangeMap.h"
36 > #include "EventHeader.h"
37 > #include "Electron.h"
38 > #include "Muon.h"
39 > #include "Vertex.h"
40 > #include "PFMet.h"
41 > #include "PFJetCol.h"
42 > #include "PFCandidate.h"
43 > #include "PFCandidateCol.h"
44 > #include "PileupInfoCol.h"
45 > #include "PileupEnergyDensity.h"
46 > #include "MCParticle.h"
47 > #include "TriggerMask.h"
48 > #include "TriggerTable.h"
49 > #include "TriggerObject.h"
50 > #include "TriggerObjectRel.h"
51 > #include "Names.h"
52 > #include "BaseMod.h"
53 > #include "TriggerObjectsTable.h"
54 > #include "JetIDMVA.h"
55  
56   //
57 < // utility headers
57 > // our headers
58   //
59   #include "ParseArgs.h"
60 + #include "MuonSelection.h"
61 + #include "ElectronSelection.h"
62 + #include "ElectronMomentumCorrection.h"
63 + #include "JetSelection.h"
64 + #include "IsolationSelection.h"
65 + #include "ReferenceSelection.h"
66 + #include "MuonMomentumCorrection.h"
67 +
68 + #include "TriggerUtils.h"
69 + #include "PassHLT.h"
70 + #include "Angles.h"
71 + #include "KinematicsStruct.h"
72 + #include "InfoStruct.h"
73 + #include "GenInfoStruct.h"
74 + #include "WeightStruct.h"
75 + #include "AngleTuple.h"
76 + #include "FOTuple.h"
77 +
78 + #include "RunLumiRangeMap.h"
79   #include "SampleWeight.h"
80 < #include "Selection.h"
81 < #include "HZZCiCElectronSelection.h"
82 < #include "SampleWeight.h"
80 > #include "EfficiencyWeightsInterface.h"
81 > #include "SelectionFuncs.h"
82 >
83 > #include "SimpleLepton.h"
84 >
85 > #ifndef CMSSW_BASE
86 > #define CMSSW_BASE "../"
87 > #endif
88 >
89 > using namespace mithep;
90 >
91 > //
92 > // globals
93 > //----------------------------------------------------------------------------
94 > TH1D *hpu_2011, *hpu_2012;
95 > RunLumiRangeMap rlrm;
96 > vector<SimpleLepton> failingLeptons ; // for fake estimation
97 > vector<SimpleLepton> passingLeptons;      // for fake estimation
98 > vector<unsigned> cutvec;
99 > vector<vector<unsigned> > zcutvec;
100 > vector<vector<unsigned> > zzcutvec;
101 > map<unsigned,float>       evtrhoMap;
102 > vector<string> cutstrs;
103 > bool passes_HLT;
104 > vector<bool>   PFnoPUflag;;
105 > map<TString, map<TString,int>* > counts;
106 > ElectronMomentumCorrection electron_momentum_correction;
107 > MuCorr *muCorr;
108 > //----------------------------------------------------------------------------
109 > void initRunLumiRangeMap(ControlFlags &ctrl)
110 > //----------------------------------------------------------------------------
111 > {
112 >  if(ctrl.jsonFile==""){
113 >    cout << "WARNING: running with empty jsonFile, so setting noJSON to true" << endl;
114 >    ctrl.noJSON = true;
115 >  }
116 >  rlrm.AddJSONFile(string(ctrl.jsonFile));  
117 > }
118  
119 < //#define THEIR_EVENTS
119 > //
120 > // prototypes
121 > //----------------------------------------------------------------------------
122  
123 < //=== MAIN =================================================================================================
123 > //
124 > // MAIN
125 > //----------------------------------------------------------------------------
126   int main(int argc, char** argv)
127 + //----------------------------------------------------------------------------
128   {
129 +  vector<TString> cuts;
130 +  cuts.push_back("start"); // NOTE: is incremented in ReferenceSelection, so if you turn on JSONS be careful
131 +  cuts.push_back("trigger");
132 +  cuts.push_back("sfOsHiPt");
133 +  cuts.push_back("4<mZ2<120");
134 +  cuts.push_back("pt>20,10");
135 +  cuts.push_back("mll>4");
136 +  cuts.push_back("m4l>70");
137 +  cuts.push_back("mZ2>12");
138 +  cuts.push_back("m4l>100");
139 +
140 +  for(unsigned icut=0; icut<cuts.size(); icut++)
141 +    counts[cuts[icut]] = new map<TString,int>;
142 +
143 +  map<TString,map<TString,int>* >::iterator cit;
144 +  for(cit=counts.begin(); cit!=counts.end(); cit++) {
145 +    map<TString,int> *cts = (*cit).second;
146 +    (*cts)["all"]  = int(0);
147 +    (*cts)["4e"]   = int(0);
148 +    (*cts)["4m"]   = int(0);
149 +    (*cts)["2e2m"] = int(0);
150 +  }
151    
152 +  string cmsswpath(CMSSW_BASE);
153 +  cmsswpath.append("/src");
154 +  bitset<1024> triggerBits;
155    vector<vector<string> > inputFiles;
156    inputFiles.push_back(vector<string>());
157 +  map<string,unsigned> entrymap; // number of unskimmed entries in each file
158  
159    //
160    // args
161    //--------------------------------------------------------------------------------------------------------------
162    ControlFlags ctrl;
163    parse_args( argc, argv, ctrl );
164 <  if( ctrl.inputfile.empty() ||
165 <      ctrl.inputfile.empty() ) {
166 <    cerr << "usage: applySelection.exe <flags> " << endl;
167 <    cerr << "\tmandoatory flags : " << endl;
168 <    cerr << "\t--inputfile |  file containing a list of ntuples to run over" << endl;
169 <    cerr << "\t--outputfile | file to store  selected evet" << endl;
170 <    return 1;
171 <  }
164 >  if( ctrl.inputfiles.empty() &&ctrl.inputfile.empty() )
165 >    {
166 >      cerr << "usage: applySelection.exe <flags> " << endl;
167 >      cerr << "\tmandoatory flags : " << endl;
168 >      cerr << "\t--inputfiles |  file containing a list of ntuples to run over" << endl;
169 >      cerr << "\t--inputfile |  a file to run over" << endl;
170 >      cerr << "\t--outputfile | file to store  selected evet" << endl;
171 >      return 1;
172 >    }
173    ctrl.dump();
174 +  cout << "\n\nNOTE: now setting era with arguments... double check this doesnt mess something else up\n\n" << endl;
175 +  bool doJets = true;//TString(ctrl.mcfmfname).Contains("doJets");
176 +  assert(ctrl.doFSR);
177  
178    //
179    // File I/O
180    //--------------------------------------------------------------------------------------------------------------
181    TChain * chain = new TChain("Events");
182 <  ifstream f(ctrl.inputfile.c_str());
182 >  TChain * hltchain = new TChain("HLT");
183 >  TChain * runchain = new TChain(Names::gkRunTreeName);
184 >
185    string fname;
186 <  while (f >> fname) {
187 <    if( !(strncmp( fname.c_str(), "#", 1 ) ) ) continue; // skip commented lines  
188 <    cout << "adding inputfile : " << fname.c_str() << endl;
189 <    chain->AddFile(fname.c_str());
186 >  unsigned total_unskimmed=0;
187 >  if( !ctrl.inputfiles.empty() ) {
188 >    ifstream f(ctrl.inputfiles.c_str());
189 >    while (f >> fname) {
190 >      if( !(strncmp( fname.c_str(), "#", 1 ) ) ) continue; // skip commented lines
191 >      if(fname.find("/store/") != string::npos) fname = "root://eoscms//"+fname;
192 >      cout << "adding inputfile : " << fname.c_str() << endl;
193 >      entrymap[string(fname.c_str())] = unskimmedEntries(fname.c_str());
194 >      cout << "unskimmed entries: " << entrymap[string(fname.c_str())] << endl;
195 >      total_unskimmed += entrymap[string(fname.c_str())];
196 >      chain->AddFile(fname.c_str());
197 >      hltchain->AddFile(fname.c_str());
198 >      runchain->AddFile(fname.c_str());
199 >    }
200 >  } else {
201 >    if(ctrl.inputfile.find("/store/") != string::npos) ctrl.inputfile = "root://eoscms//"+ctrl.inputfile;
202 >    cout << "adding inputfile : " << ctrl.inputfile.c_str() << endl;
203 >    unsigned unsk_ents = unskimmedEntries(ctrl.inputfile.c_str());
204 >    entrymap[string(ctrl.inputfile.c_str())] = unsk_ents;
205 >    cout << "unskimmed entries: " << unsk_ents << endl;
206 >    total_unskimmed += unsk_ents;
207 >    chain->AddFile(ctrl.inputfile.c_str());
208 >    hltchain->AddFile(ctrl.inputfile.c_str());
209 >    runchain->AddFile(ctrl.inputfile.c_str());
210    }
211 +  // write the total number of unskimmed events that went into making this output file to a text file
212 +  writeEntries(ctrl,total_unskimmed);
213 +  ctrl.totalMC = total_unskimmed;
214  
215 +  const char * ofname;
216 +  if( strcmp( ctrl.outputfile.c_str(), "") ) {
217 +    ofname = ctrl.outputfile.c_str();
218 +  } else {
219 +    ofname = "tmp.root";
220 +  }
221    // table of cross section values  
222 <  SimpleTable xstab("./data/xs.dat");
222 >  string xspath = (cmsswpath+string("/MitPhysics/data/xs.dat"));
223 >  cout << "xspath: " << xspath.c_str() << endl;
224 >  SimpleTable xstab(xspath.c_str());
225  
226    //
227    // Setup
228    //--------------------------------------------------------------------------------------------------------------
229 <  TH1F * h_evt = new TH1F("hevt", "hevt", 2, 0, 2 );
230 <  TH1F * h_evtfail = new TH1F("hevtfail", "hevtfail", 1024, 0, 1024 );
231 <  TNtuple * passtuple = new TNtuple( "passtuple", "passtuple", "run:evt:lumi:channel:mZ1:mZ2:m4l:pt4l:w" );
232 <  initCiCSelection();
233 <  initRunLumiRangeMap();
229 >  Angles angles;
230 >  KinematicsStruct kine;
231 >  InfoStruct evtinfo;
232 >  WeightStruct weights;
233 >  GenInfoStruct geninfo;
234 >  JetInfoStruct ji;
235 >  vector<SimpleLepton> simpleLeptonJets;
236 >
237 >  AngleTuple nt( (const char*)ofname, (const char*)"zznt");
238 >  nt.makeAngleBranch(angles);
239 >  nt.makeKinematicsBranch(kine);
240 >  nt.makeInfoBranch(evtinfo);
241 >  if(doJets)
242 >    nt.makeJetInfoBranch(ji);
243 >
244 >  int genSampleType(0);
245 >  if(ctrl.mc) {
246 >    nt.makeWeightBranch(weights);
247 >    if(ctrl.fillGen) {
248 >      genSampleType = getGenSampleType(ofname);
249 >      if(genSampleType == -1) {
250 >        cout << "\ngenSampleType returned as -1, setting fillGen to false" << endl;
251 >        ctrl.fillGen = false;
252 >      } else
253 >        nt.makeGenInfoBranch(geninfo);
254 >    }
255 >    initEfficiencyWeights();
256 >    // initPUWeights();
257 >  } else {
258 >    if(!(ctrl.noJSON) )
259 >      initRunLumiRangeMap(ctrl);
260 >  }
261 >
262 >  initElectronIDMVAV1();
263 >  initTrigger();
264 >  TrigInfo ti;
265 >  initAnalysisTriggers(ti);
266 >  if(ctrl.do_escale) {
267 >    electron_momentum_correction.initialize_regression("EGamma/EGammaAnalysisTools/data/eleEnergyRegWeights_V1.root");
268 >    if(ctrl.mc) {
269 >      if(ctrl.era == 2011) assert(0); // don't know which file to use yet
270 >      else                 ctrl.smearing_scaling_dataset = "Summer12_DR53X_HCP2012";
271 >    } else {
272 >      if(ctrl.era == 2011) assert(0); // don't know which file to use yet
273 >      else                 ctrl.smearing_scaling_dataset = "2012Jul13ReReco";
274 >    }
275 >  }
276 >  if(ctrl.correct_muon_momentum) {
277 >    muCorr = new MuCorr;
278 >    muCorr->corr2011 = new rochcor;
279 >    muCorr->corr2012 = new rochcor2012;
280 >  }
281  
282 +  vector<TString> fJetCorrParsv;
283 +  vector<JetCorrectorParameters> correctionParameters;
284 +  FactorizedJetCorrector *fJetCorrector;
285 +  JetIDMVA *fJetIDMVA;
286 +  if(doJets) {
287 +    initJets(ctrl, cmsswpath, fJetCorrParsv, correctionParameters, fJetCorrector, fJetIDMVA);
288 +  }
289  
290    //
291 <  // Loop
291 >  // Setup tree I/O
292    //--------------------------------------------------------------------------------------------------------------
103  
104  //
105  // Access samples and fill histograms
293    TFile *inputFile=0;
294    TTree *eventTree=0;  
295 <  
109 <  // Data structures to store info from TTrees
110 <  mithep::TEventInfo *info    = new mithep::TEventInfo();
111 <  TClonesArray *electronArr   = new TClonesArray("mithep::TElectron");
112 <  TClonesArray *muonArr       = new TClonesArray("mithep::TMuon");
113 <  
114 <
115 <  TBranch *infoBr;
116 <  TBranch *electronBr;
117 <  TBranch *muonBr;
118 <
119 < //   chain->SetBranchStatus("*",0); //disable all branches
120 < //   chain->SetBranchStatus("Info", 1);
121 < //   chain->SetBranchStatus("Muon", 1);
122 < //   chain->SetBranchStatus("Electron", 1);
123 <
124 <  chain->SetBranchAddress("Info",       &info);
125 <  chain->SetBranchAddress("Electron", &electronArr);
126 <  chain->SetBranchAddress("Muon", &muonArr);
127 <
128 <
129 <
130 <  cout << "nEntries: " << chain->GetEntries() << endl;
131 <  for(UInt_t ientry=0; ientry<chain->GetEntries(); ientry++) {          
132 <    
133 <    chain->GetEntry(ientry);
134 <
135 <    // get event weight for this file                                                    
136 <    double xs=1,eventweight=1;                                                            
137 <    if(ctrl.mc) eventweight = xstab.Get(getname(chain).c_str()) / chain->GetEntries();    
138 <
139 < #ifdef THEIR_EVENTS
140 <    if( !( info->evtNum == 504867308  ||
141 <           info->evtNum == 368148849  ||
142 <           info->evtNum == 129514273  ||
143 <           info->evtNum == 286336207  ||
144 <           info->evtNum == 344708580  ||
145 <           info->evtNum == 30998576   ||
146 <           info->evtNum == 155679852  ||
147 <           info->evtNum == 394010457  ||
148 <           info->evtNum == 917379387  ||
149 <           info->evtNum == 78213037   ||
150 <           info->evtNum == 337493970  ||
151 <           info->evtNum == 1491724484 ||
152 <           info->evtNum == 480301165  ||
153 <           info->evtNum == 1038911933 ||
154 <           info->evtNum == 876658967  ||
155 <           info->evtNum == 966824024  ||
156 <           info->evtNum == 141954801  ||
157 <           info->evtNum == 160966858  ||
158 <           info->evtNum == 191231387  ||
159 <           info->evtNum == 66033190   ||
160 <           info->evtNum == 10347106   ||
161 <           info->evtNum == 107360878  ) ) continue;
162 < #endif
163 <
164 <    unsigned evtfail = fails_HZZ4L_selection(ctrl, info, electronArr, muonArr, eventweight, passtuple );    
165 <    h_evtfail->Fill( evtfail, eventweight );
166 <    cout << endl << endl;  
167 <          
168 <  } //end loop over data    
169 <
170 <
171 <
172 <  delete info;
173 <  delete electronArr;
174 <  delete muonArr;
295 >  string currentFile("");
296  
297 <
298 <
299 <  //--------------------------------------------------------------------------------------------------------------
300 <  // Save Histograms;
301 <  //==============================================================================================================
302 <  const char * ofname;
303 <  if( strcmp( ctrl.outputfile.c_str(), "") ) {
304 <    ofname = ctrl.outputfile.c_str();
305 <  } else {
306 <    ofname = "tmp.root";
297 >  UInt_t fNMaxTriggers = TRIGGER_BIG_NUMBER;
298 >  EventHeader                   *info                    = new EventHeader();
299 >  Array<PFMet>                  *metArr                  = new Array<PFMet>();
300 >  Array<Electron>               *electronArr             = new Array<Electron>();
301 >  Array<Muon>                   *muonArr                 = new Array<Muon>();
302 >  Array<Vertex>                 *vtxArr                  = new Array<Vertex>();
303 >  Array<PFCandidate>            *pfArr                   = new Array<PFCandidate>();
304 >  Array<PFJet>                  *jetArr                  = new Array<PFJet>();
305 >  Array<PileupInfo>             *puArr                   = new Array<PileupInfo>();
306 >  Array<PileupEnergyDensity>    *puDArr                  = new Array<PileupEnergyDensity>();
307 >  Array<Track>                  *trkArr                  = new Array<Track>();
308 >  Array<MCParticle>             *mcArr                   = new Array<MCParticle>();
309 >  MCEventInfo                   *mcEvtInfo               = new MCEventInfo();
310 >  TriggerMask                   *trigMask                = new TriggerMask();
311 >  TriggerTable                  *hltTable                = new TriggerTable(fNMaxTriggers);
312 >  vector<string>                *hltTableStrings         = new vector<string>();
313 >  vector<string>                *hltLabelStrings         = new vector<string>();
314 >  Array<TriggerObject>          *hltObjArr               = new Array<TriggerObject>();
315 >  Array<TriggerObjectRel>       *hltRelsArr              = new Array<TriggerObjectRel>();
316 >  TriggerObjectsTable           *fTrigObjs               = new TriggerObjectsTable(hltTable,fNMaxTriggers);
317 >  RunInfo                       *runInfo                 = new RunInfo();
318 >  
319 >  chain->SetBranchAddress(Names::gkEvtHeaderBrn,                &info);
320 >  chain->SetBranchAddress("PFMet",                              &metArr);
321 >  chain->SetBranchAddress(Names::gkElectronBrn,                 &electronArr);
322 >  chain->SetBranchAddress(Names::gkMuonBrn,                     &muonArr);
323 >  chain->SetBranchAddress(Names::gkPVBrn,                       &vtxArr);
324 >  chain->SetBranchAddress(Names::gkPFCandidatesBrn,             &pfArr);
325 >  chain->SetBranchAddress(Names::gkPFJetBrn,                    &jetArr);
326 >  if( ctrl.mc ) {
327 >    chain->SetBranchAddress(Names::gkPileupInfoBrn,             &puArr);
328 >    chain->SetBranchAddress(Names::gkMCPartBrn,                 &mcArr);
329 >    chain->SetBranchAddress(Names::gkMCEvtInfoBrn,              &mcEvtInfo);
330 >  }
331 >  chain->SetBranchAddress(Names::gkPileupEnergyDensityBrn,      &puDArr);
332 >  chain->SetBranchAddress(Names::gkTrackBrn,                    &trkArr);
333 >  chain->SetBranchAddress(Names::gkHltBitBrn,                   &trigMask);
334 >  chain->SetBranchAddress(Names::gkHltObjBrn,                   &hltObjArr);
335 >  chain->SetBranchAddress("HLTObjectsRelation",                 &hltRelsArr);
336 >  hltchain->SetBranchAddress(Names::gkHltTableBrn,              &hltTableStrings);
337 >  hltchain->SetBranchAddress(Names::gkHltLabelBrn,              &hltLabelStrings);
338 >
339 >  runchain->SetBranchAddress(Names::gkRunInfoBrn, &runInfo);
340 >
341 >  Vertex              vtx;          // best primary vertex in the event
342 >
343 >  hltchain->GetEntry(0);
344 >
345 >  int imax = (ctrl.n_events_to_process < 0) ? chain->GetEntries() : ctrl.n_events_to_process;
346 >  cerr << "processing: " << imax << " ( / " << chain->GetEntries() << " )" << endl;
347 >
348 >  if(ctrl.mc && ctrl.debugOffMc) // necessary for batch mode -- we need debug output for, data but size is too large to store MC output
349 >    ctrl.debug = false;
350 >
351 >  int lastHltEntry=-2;
352 >  for(UInt_t ientry=0; ientry<imax; ientry++)
353 >    {
354 >      chain->GetEntry(ientry);
355 >      if(!(ientry%1000)) cerr << "entry: " << ientry << " ( / " << imax << " )" << endl;
356 >
357 >      if( ctrl.debug ) {
358 >        cout << "-----------------------------------------------------------------" << endl;
359 >        cout << "-----------------------------------------------------------------" << endl;
360 >        cout << "Run: " << info->RunNum()
361 >             << "\tEvt: " << info->EvtNum()
362 >             << "\tLumi: " << info->LumiSec()
363 >             << endl;
364 >        cerr << "Run: " << info->RunNum()
365 >             << "\tEvt: " << info->EvtNum()
366 >             << "\tLumi: " << info->LumiSec()
367 >             << endl;
368 >        cout << "-----------------------------------------------------------------" << endl;
369 >      }
370 >
371 >      runchain->GetEvent(info->RunEntry());
372 >      hltchain->GetEntry(runInfo->HltEntry());
373 >      hltTable->Clear();
374 >      fillTriggerNames(hltTable, hltTableStrings );
375 >      if(ctrl.debug && (lastHltEntry != runInfo->HltEntry())) hltTable->Print();
376 >      lastHltEntry = runInfo->HltEntry();
377 >      fillTriggerBits( hltTable, trigMask, triggerBits );
378 >      setHLTObjectRelations( hltObjArr, hltRelsArr, hltTableStrings, hltLabelStrings, fTrigObjs );
379 >      // printTriggerObjs( hltTable, fTrigObjs);
380 >      
381 >      string fname = string(chain->GetFile()->GetEndpointUrl()->GetFile());
382 >      if(ctrl.debug) cout << "era is " << ctrl.era << endl;
383 >      if( ctrl.era == 0  ) {
384 >        setEra( fname, ctrl );
385 >        if(ctrl.debug) cout << "era is now " << ctrl.era << endl;
386 >      }
387 >
388 >      // pfNoPU
389 >      PFnoPUflag.clear();
390 >      int pfnopu_size = makePFnoPUArray( pfArr, PFnoPUflag, vtxArr );
391 >      assert(pfnopu_size == pfArr->GetEntries());
392 >
393 >      // trigger
394 >      if( string(chain->GetFile()->GetEndpointUrl()->GetFile()) != currentFile ) {
395 >        currentFile = string(chain->GetFile()->GetEndpointUrl()->GetFile());
396 >      }
397 >      if( ctrl.debug ) cout << "file is : " << currentFile  << endl;
398 >      passes_HLT = passHLT(ctrl, triggerBits, ti, info->RunNum()); // NOTE: run value is not used if it's MC
399 >
400 >      // data/MC
401 >      if(ctrl.mc) {
402 >        if(ctrl.fillGen) fillGenInfo( mcArr, mcEvtInfo, geninfo, genSampleType, ctrl);
403 >        // NOTE: xs weight and npuw are reset in merge.cc
404 >        weights.w = getWeight(xstab,entrymap,chain)/ctrl.totalMC;
405 >        if(ctrl.debug) cout ;
406 >        weights.npu = getNPU(puArr, 0);
407 >        weights.npuw = getPUWeight(ctrl.era, fname, weights.npu);
408 >      } else {
409 >        // JSON
410 >        if(!(ctrl.noJSON) ) {
411 >          RunLumiRangeMap::RunLumiPairType rl(info->RunNum(), info->LumiSec());      
412 >          if( !(rlrm.HasRunLumi(rl)) )  {
413 >            if( ctrl.debug ) cout << "fails JSON" << endl;
414 >            continue;
415 >          }
416 >        }
417 >      }
418 >
419 >      //
420 >      // lepton/kinematic selection ...
421 >      //
422 >      failingLeptons.clear();
423 >      passingLeptons.clear();
424 >      EventData ret4l = apply_HZZ4L_reference_selection(ctrl, info,
425 >                                                        vtxArr,
426 >                                                        pfArr,
427 >                                                        puDArr,
428 >                                                        electronArr,
429 >                                                        &electronReferencePreSelection,
430 >                                                        &electronReferenceIDMVASelectionV1, // not used anymore! (Blame the damn regression implementation, not me!)
431 >                                                        &electronReferenceIsoSelection,
432 >                                                        muonArr,
433 >                                                        &muonReferencePreSelection,
434 >                                                        &muonIDPFSelection,
435 >                                                        &muonReferenceIsoSelection);
436 >      
437 >      if( ctrl.debug ) cout << endl;
438 >
439 >      if( ret4l.status.pass() ) {
440 >        fillAngles(ret4l,angles);
441 >
442 >        fillKinematics(ret4l,kine);
443 >        // fillMassErrors(ret4l,muonArr,electronArr,kine);
444 >
445 >        if(doJets) {
446 >          fillJets(ctrl, jetArr, simpleLeptonJets, fJetCorrector, puDArr, fJetIDMVA, vtxArr, ret4l);
447 >          fillJetInfo(simpleLeptonJets,ji,ctrl);
448 >        }
449 >
450 >        TLorentzVector pfmet; pfmet.SetPxPyPzE(metArr->At(0)->Mex(),metArr->At(0)->Mey(),0,0);
451 >        fillEventInfo( info, pfmet, evtinfo, ctrl.mc ? getNPU(puArr) : 0, vtxArr->GetEntries());
452 >        if( ctrl.mc) {
453 >          setEfficiencyWeights(ctrl, ctrl.era, ret4l, triggerBits, hltTable, hltObjArr, fTrigObjs, weights);
454 >          if(ctrl.debug)
455 >            cout << "w: " << weights.w << "\t" << "npuw: " << weights.npuw << "\t" << "won: " << weights.won << "\t" << "woff: " << weights.woff << endl;
456 >        }
457 >        
458 >        nt.Fill();
459 >        
460 >        cerr  << "PASS:: " << "\trun: " << evtinfo.run << "\tevt: " << evtinfo.evt << "\tlumi: " << evtinfo.lumi << "\tchannel: "
461 >              << kine.channel << "\tm4l: " << kine.m4l << "\tmZ1: " << kine.mZ1 << "\tmZ2: " << kine.mZ2 << endl;
462 >        cout  << "PASS:: " << "\trun: " << evtinfo.run << "\tevt: " << evtinfo.evt << "\tlumi: " << evtinfo.lumi << "\tchannel: "
463 >              << kine.channel << "\tm4l: " << kine.m4l << "\tmZ1: " << kine.mZ1 << "\tmZ2: " << kine.mZ2 << endl;
464 >
465 >      } else {
466 >        if(ctrl.debug)
467 >          cout << "failing with code : " << ret4l.status.selectionBits << endl;
468 >      }
469 >    }
470 >
471 >  nt.WriteClose();
472 >
473 >  for(unsigned icut=0; icut<cuts.size(); icut++) {
474 >    map<TString,int> thisCount(*(counts[cuts[icut]]));
475 >    cout << setw(20) << cuts[icut] << "\t" << setw(9) << thisCount["all"] << "  ("
476 >         << setw(9) << thisCount["4e"] << setw(9) << thisCount["4m"] << setw(9) << thisCount["2e2m"] << ")" << endl;
477    }
187
188  TFile *file = new TFile(ofname, "RECREATE");
189  h_evt->Write();
190  h_evtfail->Write();
191  passtuple->Write();
192  file->Close();
193  delete file;
194
195  cerr << "done!" << endl;
478   }
197
198
199

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines