ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/MitHzz4l/Selection/src/applySelection.cc

Comparing UserCode/MitHzz4l/Selection/src/applySelection.cc (file contents):
Revision 1.26 by khahn, Wed May 23 13:33:15 2012 UTC vs.
Revision 1.52 by khahn, Wed Jul 11 18:24:39 2012 UTC

#	Line 33 | Line 33 | using namespace std;
33		#include "Electron.h"
34		#include "Muon.h"
35		#include "Vertex.h"
36	+	#include "PFMet.h"
37		#include "PFCandidate.h"
38		#include "PFCandidateCol.h"
39		#include "PileupInfoCol.h"
#	Line 40 | Line 41 | using namespace std;
41		#include "MCParticle.h"
42		#include "TriggerMask.h"
43		#include "TriggerTable.h"
44	+	#include "TriggerObject.h"
45	+	#include "TriggerObjectRel.h"
46		#include "Names.h"
47		#include "BaseMod.h"
48	+	#include "TriggerObjectsTable.h"
49
50		//
51		// our headers
#	Line 50 | Line 54 | using namespace std;
54		#include "MuonSelection.h"
55		#include "ElectronSelection.h"
56		#include "IsolationSelection.h"
57	<
54	<	//#include "Selection.h"
55	<	//#include "ReferenceSelection.h"
56	<	#include "ReferenceSelection2.h"
57	>	#include "ReferenceSelection.h"
58
59		#include "TriggerUtils.h"
60		#include "PassHLT.h"
61		#include "Angles.h"
62		#include "KinematicsStruct.h"
63		#include "InfoStruct.h"
64	<	//#include "GenInfoStruct.h"
64	>	#include "GenInfoStruct.h"
65		#include "WeightStruct.h"
65	–	//#include "GlobalDataAndFuncs.h"
66	–	#include "RunLumiRangeMap.h"
67	–
66		#include "AngleTuple.h"
67		#include "FOTuple.h"
68
69	+	#include "RunLumiRangeMap.h"
70		#include "SampleWeight.h"
71		#include "EfficiencyWeightsInterface.h"
72	+	#include "SelectionFuncs.h"
73
74		#include "SimpleLepton.h"
75
#	Line 82 | Line 82 | using namespace mithep;
82		//
83		// globals
84		//----------------------------------------------------------------------------
85	<	TH1D * hpu;
85	>	TH1D *hpu_2011, *hpu_2012;
86		RunLumiRangeMap rlrm;
87		vector<SimpleLepton> failingLeptons ; // for fake estimation
88		vector<SimpleLepton> passingLeptons;      // for fake estimation
#	Line 97 | Line 97 | vector<bool>   PFnoPUflag;;
97		//
98		// prototypes
99		//----------------------------------------------------------------------------
100	<	void initPUWeights();
101	<	double getPUWeight(unsigned npu);
102	<	void setEffiencyWeights(EventData&, WeightStruct& );
103	<	void initRunLumiRangeMap();
104	<	void initEvtRhoMap(map<unsigned,float> &);
105	<	unsigned makePFnoPUArray(mithep::Array<PFCandidate> * fPFCandidates,
106	<	vector<bool> &pfNoPileUpFlag,
107	<	const mithep::VertexCol     * vtxArr );
100	>	void initRunLumiRangeMap(ControlFlags &ctrl);
101		//----------------------------------------------------------------------------
102
103
#	Line 175 | Line 168 | int main(int argc, char** argv)
168		TChain * hltchain = new TChain("HLT");
169
170		string fname;
171	+	unsigned total_unskimmed=0;
172		if( !ctrl.inputfiles.empty() ) {
173		ifstream f(ctrl.inputfiles.c_str());
174		while (f >> fname) {
#	Line 182 | Line 176 | int main(int argc, char** argv)
176		cout << "adding inputfile : " << fname.c_str() << endl;
177		entrymap[string(fname.c_str())] = unskimmedEntries(fname.c_str());
178		cout << "unskimmed entries: " << entrymap[string(fname.c_str())] << endl;
179	+	total_unskimmed += entrymap[string(fname.c_str())];
180		chain->AddFile(fname.c_str());
181		hltchain->AddFile(fname.c_str());
182		}
183		} else {
184		cout << "adding inputfile : " << ctrl.inputfile.c_str() << endl;
185	<	unsigned tmpent = unskimmedEntries(ctrl.inputfile.c_str());
186	<	cout << "tmpent: " << tmpent << endl;
187	<	entrymap[string(ctrl.inputfile.c_str())] = unskimmedEntries(ctrl.inputfile.c_str());
188	<	cout << "unskimmed entries: " << entrymap[string(ctrl.inputfile.c_str())] << endl;
185	>	unsigned unsk_ents = unskimmedEntries(ctrl.inputfile.c_str());
186	>	entrymap[string(ctrl.inputfile.c_str())] = unsk_ents;
187	>	cout << "unskimmed entries: " << unsk_ents << endl;
188	>	total_unskimmed += unsk_ents;
189		chain->AddFile(ctrl.inputfile.c_str());
190		hltchain->AddFile(ctrl.inputfile.c_str());
191		}
192	+	// // write the total number of unskimmed events that went into making this output file to a text file
193	+	// writeEntries(ctrl,total_unskimmed);
194
195		const char * ofname;
196		if( strcmp( ctrl.outputfile.c_str(), "") ) {
#	Line 207 | Line 204 | int main(int argc, char** argv)
204		SimpleTable xstab(xspath.c_str());
205
206
207	+
208		//
209		// Setup
210		//--------------------------------------------------------------------------------------------------------------
#	Line 214 | Line 212 | int main(int argc, char** argv)
212		KinematicsStruct kinematics;
213		InfoStruct evtinfo;
214		WeightStruct weights;
215	<	//  GenInfoStruct geninfo;
215	>	GenInfoStruct geninfo;
216
217		AngleTuple nt( (const char*)ofname, (const char*)"zznt");
218		nt.makeAngleBranch(angles);
219		nt.makeKinematicsBranch(kinematics);
220		nt.makeInfoBranch(evtinfo);
221
222	+
223		FOTuple foTree( nt.getFile(), (const char*)"FOtree");
224	+	foTree.makeInfoBranch(evtinfo);
225		foTree.makeFailedLeptonBranch(failingLeptons);
226	<	foTree.makeZLeptonBranch(passingLeptons);
226	>	foTree.makePassedLeptonBranch(passingLeptons);
227
228	<	TH1F * h_w_hpt;
228	>	TH1F * h_wf11_hpt;
229		if(ctrl.mc) {
230		nt.makeWeightBranch(weights);
231	<	//    nt.makeGenInfoBranch(geninfo);
231	>	if(ctrl.fillGen ) nt.makeGenInfoBranch(geninfo);
232		initEfficiencyWeights();
233		initPUWeights();
234	<	/*
234	>
235		// Higgs only, pt reweighting
236		if( ctrl.mH <= 140 ) {
237		char wptname[256];
#	Line 239 | Line 239 | int main(int argc, char** argv)
239		TFile * f_hpt = new TFile(wptname);
240		f_hpt->Print();
241		sprintf(wptname, "weight_hqt_fehipro_fit_%i", ctrl.mH);
242	<	h_w_hpt  = (TH1F*)(f_hpt->FindObjectAny(wptname));
243	<	h_w_hpt->Print();
242	>	h_wf11_hpt  = (TH1F*)(f_hpt->FindObjectAny(wptname));
243		}
244	<	*/
244	>
245		} else {
246	<	initRunLumiRangeMap();
246	>	initRunLumiRangeMap(ctrl);
247		}
248
249		//  initMuonIDMVA();
250	<	initMuonIsoMVA();
252	<	initElectronIDMVA();
253	<	initElectronIsoMVA();
250	>	initElectronIDMVAV1();
251		initTrigger();
252
256	–	// hack
257	–	initEvtRhoMap(evtrhoMap);
258	–
259	–
253		//
254		// Setup tree I/O
255		//--------------------------------------------------------------------------------------------------------------
#	Line 264 | Line 257 | int main(int argc, char** argv)
257		TTree *eventTree=0;
258		string currentFile("");
259
260	+	UInt_t fNMaxTriggers = 1024;
261		mithep::EventHeader *info    = new mithep::EventHeader();
262	<	//  mithep::TGenInfo    *ginfo  = new mithep::TGenInfo();
262	>	mithep::Array<mithep::PFMet>                *metArr        = new mithep::Array<mithep::PFMet>();
263		mithep::Array<mithep::Electron>             *electronArr   = new mithep::Array<mithep::Electron>();
264		mithep::Array<mithep::Muon>                 *muonArr       = new mithep::Array<mithep::Muon>();
265	<	const mithep::VertexCol              *vtxArr        = 0;//new mithep::Array<mithep::Vertex>();
265	>	mithep::Array<mithep::Vertex>               *vtxArr        = new mithep::Array<mithep::Vertex>();
266		mithep::Array<mithep::PFCandidate>          *pfArr         = new mithep::Array<mithep::PFCandidate>();
267		mithep::Array<mithep::PileupInfo>           *puArr         = new mithep::Array<mithep::PileupInfo>();
268		mithep::Array<mithep::PileupEnergyDensity>  *puDArr        = new mithep::Array<mithep::PileupEnergyDensity>();
269		mithep::Array<mithep::Track>                *trkArr        = new mithep::Array<mithep::Track>();
270		mithep::Array<mithep::MCParticle>           *mcArr         = new mithep::Array<mithep::MCParticle>();
271	+	mithep::MCEventInfo                         *mcEvtInfo     = new mithep::MCEventInfo();
272		mithep::TriggerMask                         *trigMask      = new mithep::TriggerMask();
273	<	mithep::TriggerTable                        *hltTable      = new mithep::TriggerTable();
273	>	mithep::TriggerTable                        *hltTable      = new mithep::TriggerTable(fNMaxTriggers);
274		vector<string>                              *hltTableStrings  = new vector<string>();
275	+	vector<string>                              *hltLabelStrings  = new vector<string>();
276	+	mithep::Array<mithep::TriggerObject>        *hltObjArr     = new mithep::Array<mithep::TriggerObject>();
277	+	mithep::Array<mithep::TriggerObjectRel>     *hltRelsArr    = new mithep::Array<mithep::TriggerObjectRel>();
278	+	std::vector<std::string>                    *hltTab        = new vector<string>();
279	+	std::vector<std::string>                    *hltLab        = new vector<string>();
280	+	TriggerObjectsTable                         *fTrigObjs     = new TriggerObjectsTable(hltTable,fNMaxTriggers);
281	+
282
283		TString fElectronName(Names::gkElectronBrn);
284		TString fMuonName(Names::gkMuonBrn);
285		TString fInfoName(Names::gkEvtHeaderBrn);
286	+	TString fPFMetName("PFMet");
287		TString fPrimVtxName(Names::gkPVBrn);
288		TString fPFCandidateName(Names::gkPFCandidatesBrn);
289		TString fPileupInfoName(Names::gkPileupInfoBrn);
290		TString fPileupEnergyDensityName(Names::gkPileupEnergyDensityBrn);
291		TString fTrackName(Names::gkTrackBrn);
292		TString fMCParticleName(Names::gkMCPartBrn);
293	+	TString fMCEvtInfoName(Names::gkMCEvtInfoBrn);
294		TString fTriggerMaskName(Names::gkHltBitBrn);
295		TString fTriggerTableName(Names::gkHltTableBrn);
296	+	TString fTriggerLabelName(Names::gkHltLabelBrn);
297	+	TString fTriggerObjectName(Names::gkHltObjBrn);
298	+	TString fTriggerObjectRelsName(Form("HLTObjectsRelation"));
299	+
300	+	//   TString fHLTLabName(Names::gkHltLabelBrn);
301	+	//   TString fHLTTabName(Names::gkHltTableBrn);
302	+	//   TString fHLTTabNamePub(Form("%sFwk",fHLTTabName.Data()));
303	+	//   TString fHLTLabNamePub(Form("%sFwk",fHLTLabName.Data()));
304	+	//   TString fObjsNamePub(Form("%sFwk",fTriggerObjectName.Data()));
305	+
306	+
307
308		chain->SetBranchAddress(fInfoName,        &info);
309	+	chain->SetBranchAddress(fPFMetName,       &metArr);
310		chain->SetBranchAddress(fElectronName,    &electronArr);
311		chain->SetBranchAddress(fMuonName,        &muonArr);
312		chain->SetBranchAddress(fPrimVtxName,     &vtxArr);
#	Line 298 | Line 314 | int main(int argc, char** argv)
314		if( ctrl.mc ) {
315		chain->SetBranchAddress(fPileupInfoName,  &puArr);
316		chain->SetBranchAddress(fMCParticleName,  &mcArr);
317	+	chain->SetBranchAddress(fMCEvtInfoName,  &mcEvtInfo);
318		}
319		chain->SetBranchAddress(fPileupEnergyDensityName, &puDArr);
320		chain->SetBranchAddress(fTrackName, &trkArr);
321		chain->SetBranchAddress(fTriggerMaskName, &trigMask);
322	+	chain->SetBranchAddress(fTriggerObjectName,  &hltObjArr);
323	+	chain->SetBranchAddress(fTriggerObjectRelsName,  &hltRelsArr);
324	+	//  chain->SetBranchAddress(fHLTTabNamePub,  &hltTab);
325	+	//  chain->SetBranchAddress(fHLTLabNamePub,  &hltLab);
326	+
327		cout << "hlttable: " << fTriggerTableName << endl;
328		hltchain->SetBranchAddress(fTriggerTableName, &hltTableStrings);
329	+	hltchain->SetBranchAddress(fTriggerLabelName, &hltLabelStrings);
330
331		mithep::Vertex              vtx;          // best primary vertex in the event
332
#	Line 318 | Line 341 | int main(int argc, char** argv)
341
342		// only 1 HLT table / file ???
343		hltchain->GetEntry(0);
321	–	cerr << "here... " << endl;
344
345		int imax = chain->GetEntries();
346		cout << "nEntries: " << imax << endl;
347
326	–
348		//
349		// Loop !!!!!!!!!
350		//--------------------------------------------------------------------------------------------------------------
351		for(UInt_t ientry=0; ientry<imax; ientry++)
352		{
353		chain->GetEntry(ientry);
354	+	if(!(ientry%1000)) cerr << "entry: " << ientry << endl;
355
356	<	/*
357	<	bool found_tau=false;
358	<	if( ctrl.mc) {  // && ctrl.debug
359	<	cout << "nparticles: " << mcArr->GetEntries() << endl;
360	<	cout << "-----------------------------------" << endl;
361	<	for( int i=0; i<mcArr->GetEntries(); i++ ) {
362	<	const mithep::MCParticle *mc = mcArr->At(i);
363	<	cout << "i: " << i << "\t"
342	<	<< "id: " << mc->AbsPdgId()
343	<	<< endl;
344	<	if( mc->Is(MCParticle::kTau) ) {
345	<	found_tau=true;
346	<	}
347	<	}
356	>	if( ctrl.debug ) {
357	>	cout << "-----------------------------------------------------------------" << endl;
358	>	cout << "-----------------------------------------------------------------" << endl;
359	>	cout << "Run: " << info->RunNum()
360	>	<< "\tEvt: " << info->EvtNum()
361	>	<< "\tLumi: " << info->LumiSec()
362	>	<< endl;
363	>	cout << "-----------------------------------------------------------------" << endl;
364		}
349	–	cout << "-----------------------------------" << endl;
350	–	cout << "found tau: " << found_tau << endl;
351	–	*/
352	–
353	–	// events that fail bc ele overlaps w/ reco mu
354	–	/*
355	–	if( info->EvtNum() != 67315 &&
356	–	info->EvtNum() != 288693 ) continue;
357	–	*/
358	–
359	–	// events where they apparently don't veto electrons ...
360	–	/*
361	–	if( info->EvtNum() != 141713 &&
362	–	info->EvtNum() != 178019 &&
363	–	info->EvtNum() != 60823 &&
364	–	info->EvtNum() != 192795 ) continue;
365	–	*/
365
367	–	if(!(ientry%1000)) cerr << "entry: " << ientry << endl;
368	–	//if (!(setPV( ctrl, vtxArr, vtx)) ) continue;
366
367	<	float rho = evtrhoMap[info->EvtNum()];
367	>	string fname = string(chain->GetFile()->GetEndpointUrl()->GetFile());
368	>	if(ctrl.debug) cout << "era is " << ctrl.era << endl;
369	>	if( ctrl.era == 0  ) {
370	>	setEra( fname, ctrl );
371	>	if(ctrl.debug) cout << "era is now " << ctrl.era << endl;
372	>	}
373
374		//
375		// pfNoPU
#	Line 382 | Line 384 | int main(int argc, char** argv)
384		//
385		if(ctrl.mc) {
386		//
387	+	// gen info
388	+	//
389	+	if( ctrl.fillGen )
390	+	fillGenInfo( mcArr, mcEvtInfo, geninfo, ESampleType::kHZZ, ctrl);
391	+
392	+	//
393		// xsec & PU weights
394		//
395	<	weights.w = 1;//getWeight(xstab,entrymap,chain)/ctrl.totalMC;
395	>	weights.w = getWeight(xstab,entrymap,chain)/ctrl.totalMC;
396		cout << "xsec weight: " << weights.w << endl;
397	<	int npu = -1;
398	<	for(int i=0;i<puArr->GetEntries();i++) {
399	<	if(puArr->At(i)->GetBunchCrossing() == 0)
400	<	npu = puArr->At(i)->GetPU_NumInteractions();
397	>
398	>	unsigned npu = getNPU(puArr, 0);
399	>	weights.npu  = npu;
400	>	weights.npuw = getPUWeight(ctrl.era,npu);
401	>	cerr << "weights.npuw: " <<  weights.npuw << endl;
402	>
403	>	//
404	>	// pt reweighting for Higgs < 140, F11
405	>	//
406	>	if( ctrl.fillGen ) {
407	>	if( (fname.find("f11-h115zz4l") != string::npos) ||
408	>	(fname.find("f11-h120zz4l") != string::npos) ||
409	>	(fname.find("f11-h130zz4l") != string::npos)  ) {
410	>	weights.w *= h_wf11_hpt->GetBinContent(h_wf11_hpt->FindBin(geninfo.pt_zz));
411	>	}
412		}
413	<	assert(npu>=0);
395	<	evtinfo.npu;
396	<	weights.npuw = getPUWeight(evtinfo.npu);
397	<	//      cout << "weight: " << weights.w << endl;
413	>
414		//
415		// trigger
416		//
417		if( string(chain->GetFile()->GetEndpointUrl()->GetFile()) != currentFile ) {
418		currentFile = string(chain->GetFile()->GetEndpointUrl()->GetFile());
419		hltchain->SetBranchAddress(fTriggerTableName, &hltTableStrings);
420	+	hltchain->SetBranchAddress(fTriggerLabelName, &hltLabelStrings);
421		hltchain->GetEntry(0);
422		hltTable->Clear();
423	<	fillTriggerBits(hltTable, hltTableStrings );
423	>	fillTriggerNames(hltTable, hltTableStrings );
424	>	hltTable->Print();
425		}
426		if( ctrl.debug ) cout << "file is : " << currentFile  << endl;
427		fillTriggerBits( hltTable, trigMask, triggerBits );
428		passes_HLT_MC = true; // passed to apply as extern global, just for sync
429	<	if( !passHLTMC(triggerBits, info->RunNum(), info->EvtNum() ) ) {
429	>	if( !passHLTMC(triggerBits, info->RunNum(), info->EvtNum(), k2012_MC ) ) {
430		passes_HLT_MC = false;
431		}
432	+
433	+	setHLTObjectRelations( hltObjArr, hltRelsArr, hltTableStrings, hltLabelStrings, fTrigObjs );
434	+
435		} else {
436		//
437		// JSON
438		//
439	<	/*
440	<	RunLumiRangeMap::RunLumiPairType rl(info->RunNum(), info->LumiSec());
441	<	if( !(rlrm.HasRunLumi(rl)) )  {
442	<	if( ctrl.debug ) cout << "\tfails JSON" << endl;
443	<	continue;
439	>	if(!(ctrl.noJSON) ) {
440	>	RunLumiRangeMap::RunLumiPairType rl(info->RunNum(), info->LumiSec());
441	>	if( !(rlrm.HasRunLumi(rl)) )  {
442	>	if( ctrl.debug ) cout << "fails JSON" << endl;
443	>	continue;
444	>	}
445		}
446	<	*/
446	>
447		//
448		// trigger
449		//
450		if( string(chain->GetFile()->GetEndpointUrl()->GetFile()) != currentFile ) {
451		currentFile = string(chain->GetFile()->GetEndpointUrl()->GetFile());
452		hltchain->SetBranchAddress(fTriggerTableName, &hltTableStrings);
453	+	hltchain->SetBranchAddress(fTriggerLabelName, &hltLabelStrings);
454		hltchain->GetEntry(0);
455		hltTable->Clear();
456	<	fillTriggerBits(hltTable, hltTableStrings );
456	>	fillTriggerNames(hltTable, hltTableStrings );
457		}
458	+
459	+	setHLTObjectRelations( hltObjArr, hltRelsArr, hltTableStrings, hltLabelStrings, fTrigObjs );
460	+
461		if( ctrl.debug ) cout << "file is : " << currentFile  << endl;
462		/*
463		fillTriggerBits( hltTable, trigMask, triggerBits );
#	Line 451 | Line 477 | int main(int argc, char** argv)
477
478
479		// reference
480	<	apply_HZZ4L_reference2_selection(ctrl, info,
481	<	vtxArr,
482	<	pfArr,
483	<	#ifdef HACKED_RHOS
484	<	rho,
485	<	#else
486	<	puDArr,
487	<	#endif
488	<	electronArr,
489	<	&electronReferencePreSelection,
490	<	&electronReferenceIDMVASelectionV1,
491	<	&electronReferenceIsoSelection,
492	<	muonArr,
467	<	&muonReferencePreSelection,
468	<	&muonIDPFSelection,
469	<	&muonReferenceIsoSelection);
470	<
471	<	/*
472	<	// evolved
473	<	apply_HZZ4L_reference2_selection(ctrl, info,
474	<	vtx,
475	<	pfArr,
476	<	//rho,
477	<	puDArr,
478	<	electronArr,
479	<	&electronReferencePreSelection,
480	<	//&electronPreSelection,
481	<	&electronReferenceIDMVASelection,
482	<	&electronIsoMVASelection,
483	<	muonArr,
484	<	&muonReferencePreSelection,
485	<	//&muonPreSelection,
486	<	&muonIDPFSelection,
487	<	&muonIsoMVASelection);
488	<	*/
489	<
480	>	apply_HZZ4L_reference_selection(ctrl, info,
481	>	vtxArr,
482	>	pfArr,
483	>	puDArr,
484	>	electronArr,
485	>	&electronReferencePreSelection,
486	>	&electronReferenceIDMVASelectionV1,
487	>	&electronReferenceIsoSelection,
488	>	muonArr,
489	>	&muonReferencePreSelection,
490	>	&muonIDPFSelection,
491	>	&muonReferenceIsoSelection);
492	>
493		if( ctrl.debug ) cout << endl;
494
495
#	Line 500 | Line 503 | int main(int argc, char** argv)
503		if( abs(ret4l.Z2leptons[0].type) == 13  ) Z2type = 13;
504		}
505
503	–	#ifdef SYNC
504	–	cout  << "bits: " << ret4l.status.selectionBits  << "\t"
505	–	<< "Z1: " << Z1type << "\t"
506	–	<< "Z2: " << Z2type << endl;
507	–	cout << endl;
508	–	#endif
509	–
506		if( ret4l.status.pass() ) {
507		fillAngles(ret4l,angles);
508		fillKinematics(ret4l,kinematics);
509	<	fillEventInfo(info,evtinfo);
509	>	fillMassErrors(ret4l,muonArr,electronArr,kinematics);
510	>	TLorentzVector pfmet; pfmet.SetPxPyPzE(metArr->At(0)->Mex(),metArr->At(0)->Mey(),0,0);
511	>	fillEventInfo( info, pfmet, evtinfo, ctrl.mc ? getNPU(puArr) : 0);
512		if( ctrl.mc) {
513	<	// fillGenInfo(ginfo,geninfo);
514	<	setEffiencyWeights(ret4l, weights);
515	<	if(ctrl.debug)
516	<	cout << "w: " << weights.w << "\t"
517	<	<< "won: " << weights.won << "\t"
518	<	<< "woff: " << weights.woff << endl;
513	>	setEffiencyWeights(ctrl.era, ret4l, triggerBits, hltTable, hltObjArr, fTrigObjs, weights);
514	>	if(ctrl.debug)
515	>	cout << "w: " << weights.w << "\t"
516	>	<< "npuw: " << weights.npuw << "\t"
517	>	<< "won: " << weights.won << "\t"
518	>	<< "woff: " << weights.woff << endl;
519		}
520
521	<	/*
524	<	// only for Higgs < 140
525	<	geninfo.weight *= h_w_hpt->GetBinContent(h_w_hpt->FindBin(geninfo.pt_zz));
526	<	angles.bdt = reader->EvaluateMVA("BDTG");
527	<	*/
521	>
522		nt.Fill();
523
524		cerr  << "PASS:: "
525		<< "\trun: " << evtinfo.run
526		<< "\tevt: " << evtinfo.evt
527		<< "\tlumi: " << evtinfo.lumi
528	<	<< "\tcostheta1: " << angles.costheta1
529	<	<< "\tcostheta2: " << angles.costheta2
530	<	<< "\tcostheta*: " << angles.costhetastar
531	<	<< "\tbdt: " << angles.bdt
528	>	<< "\tchannel: " << kinematics.channel
529	>	<< "\tm4l: " << kinematics.m4l
530	>	<< "\tmZ1: " << kinematics.mZ1
531	>	<< "\tmZ2: " << kinematics.mZ2
532		<< endl;
533		pass++;
534		//      if( pass > 3 ) break;
535
536	<	} else if( ret4l.status.selectionBits.test(PASS_GOODZ1) ) { // save for fakes ...
537	<	//      if(ctrl.debug) {
536	>	} else if( ret4l.status.selectionBits.test(PASS_GOODZ1) &&
537	>	ret4l.status.selectionBits.to_ulong() < (1UL<<PASS_4L)) { // save for fakes ...
538	>	if(ctrl.debug) {
539		cout << "failingLeptons : " << failingLeptons.size() << endl;
540		for( int f=0; f<failingLeptons.size(); f++ ) {
541		SimpleLepton  l = failingLeptons[f];
#	Line 549 | Line 544 | int main(int argc, char** argv)
544		<< "pT: " << l.vec.Pt()  << "\t"
545		<< "eta: " << l.vec.Eta()
546		<< endl;
547	<	//    }
547	>	}
548		}
549	+	TLorentzVector pfmet; pfmet.SetPxPyPzE(metArr->At(0)->Mex(),metArr->At(0)->Mey(),0,0);
550	+	fillEventInfo( info, pfmet, evtinfo, ctrl.mc ? getNPU(puArr) : 0);
551		foTree.Fill();
552	+	} else {
553	+	if(ctrl.debug)
554	+	cout << "failing with some other code : " << ret4l.status.selectionBits << endl;
555		}
556		}
557
#	Line 560 | Line 560 | int main(int argc, char** argv)
560		foTree.getTree()->Write();
561		nt.WriteClose();
562
563	<	cout << endl;
564	<	cout << endl;
563	>	if(ctrl.debug)
564	>	cout << endl << endl;
565		for( int i=0; i<cutvec.size(); i++ ) {
566		cout << "cut: " << i << "\t"
567		<< "pass: " << cutvec[i];
#	Line 583 | Line 583 | int main(int argc, char** argv)
583		}
584
585
586	–
586		//----------------------------------------------------------------------------
587	<	void setEffiencyWeights(EventData & evtdat, WeightStruct & weights )
589	<	//----------------------------------------------------------------------------
590	<	{
591	<	vector<SimpleLepton> lepvec = evtdat.Z1leptons;
592	<	lepvec.insert(lepvec.end(), evtdat.Z2leptons.begin(), evtdat.Z2leptons.end());
593	<	double w_offline=-1, werr_offline=0;
594	<	double w_online=-1, werr_online=0;
595	<
596	<	vector< pair <double,double> > wlegs; // pair here is eff & err
597	<	vector< pair <float,float> > mvec;  // pair here is eta & pt
598	<
599	<	//      vector< pair <float,float> > evec;  // pair here is eta & pt
600	<	// now deal with medium vs loose
601	<	vector< pair< bool, pair <float,float> > > evec;  // pair here is eta & pt
602	<
603	<	for( int k=0; k<lepvec.size(); k++ ) {
604	<	if( abs(lepvec[k].type) == 13 ) {
605	<	mvec.push_back( std::pair<float,float> (fabs(lepvec[k].vec.Eta()), lepvec[k].vec.Pt()) );
606	<	wlegs.push_back( muonPerLegOfflineEfficiencyWeight( fabs(lepvec[k].vec.Eta()),
607	<	lepvec[k].vec.Pt() ) );
608	<	} else {
609	<
610	<	// now deal with medium vs loose
611	<	//              evec.push_back( std::pair<float,float> (fabs(lepvec[k].vec.Eta()), lepvec[k].vec.Pt()) );
612	<
613	<	std::pair<float,float> tmppair(fabs(lepvec[k].vec.Eta()), lepvec[k].vec.Pt());
614	<	evec.push_back( std::pair<bool, std::pair<float,float> > (lepvec[k].isTight, tmppair) );
615	<
616	<	//              wlegs.push_back( elePerLegOfflineEfficiencyWeight(  fabs(lepvec[k].vec.Eta()),
617	<	//                                                                 lepvec[k].vec.Pt() ) );
618	<
619	<	wlegs.push_back( elePerLegOfflineEfficiencyWeight(  fabs(lepvec[k].vec.Eta()),
620	<	lepvec[k].vec.Pt(),
621	<	lepvec[k].isTight ) );
622	<
623	<	}
624	<	}
625	<
626	<	pair<double,double> offpair = getOfflineEfficiencyWeight( wlegs );
627	<	weights.woff    = offpair.first;
628	<	weights.werroff = offpair.second;
629	<
630	<	pair<double,double> onpair  = getOnlineEfficiencyWeight( mvec, evec );
631	<	weights.won    = onpair.first;
632	<	weights.werron = onpair.second;
633	<	}
634	<
635	<
636	<	//----------------------------------------------------------------------------
637	<	void initRunLumiRangeMap()
587	>	void initRunLumiRangeMap(ControlFlags &ctrl)
588		//----------------------------------------------------------------------------
589		{
590		/*
#	Line 646 | Line 596 | void initRunLumiRangeMap()
596		// r11b
597		rlrm.AddJSONFile(std::string("./data/Cert_160404-180252_7TeV_PromptReco_Collisions11_JSON.txt"));
598		*/
649	–	};
650	–
651	–
652	–	//----------------------------------------------------------------------------
653	–	void initPUWeights()
654	–	//----------------------------------------------------------------------------
655	–	{
656	–	TFile * puf = new TFile("data/PileupReweighting.Summer11DYmm_To_Full2011.root");
657	–	hpu = (TH1D*)(puf->Get("puWeights"));
658	–	}
659	–
660	–	//----------------------------------------------------------------------------
661	–	double getPUWeight(unsigned npu)
662	–	//----------------------------------------------------------------------------
663	–	{
664	–	return hpu->GetBinContent(hpu->FindBin(npu));
665	–	}
666	–
667	–	//----------------------------------------------------------------------------
668	–	void initEvtRhoMap( map<unsigned, float> & evtrhoMap )
669	–	//----------------------------------------------------------------------------
670	–	{
671	–	unsigned evt;
672	–	float rho;
673	–
674	–	cout << "initialzing EvtRhoMap ";
675	–	//FILE * f = fopen("evtrho.dat", "r");
676	–	//  FILE * f = fopen("allrho.cali.uniq.dat", "r");
677	–	FILE * f = fopen("/data/blue/khahn/CMSSW_4_4_2/src/MitHzz4l/allrhoAA.cali.uniq.dat", "r");
678	–	int lcount=0;
679	–	while( fscanf( f, "%u:%f", &evt, &rho ) ) {
680	–	if( feof(f) ) break;
681	–	evtrhoMap[evt] = rho;
682	–	lcount++;
683	–	if( !(lcount%1000) ) cout << "."; flush(cout);
684	–	}
685	–	cout << " done" << endl;
686	–	}
687	–
688	–	//----------------------------------------------------------------------------
689	–	unsigned makePFnoPUArray(mithep::Array<PFCandidate> * fPFCandidates,
690	–	vector<bool> &pfNoPileUpflag,
691	–	const mithep::VertexCol      * vtxArr )
692	–	//----------------------------------------------------------------------------
693	–	{
694	–	for(UInt_t i = 0; i < fPFCandidates->GetEntries(); i++) {
695	–	const PFCandidate *pf = fPFCandidates->At(i);
696	–	assert(pf);
697	–
698	–	if(pf->PFType() == PFCandidate::eHadron) {
699	–	if(pf->HasTrackerTrk() &&
700	–	vtxArr->At(0)->HasTrack(pf->TrackerTrk()) &&
701	–	vtxArr->At(0)->TrackWeight(pf->TrackerTrk()) > 0) {
702	–
703	–	pfNoPileUpflag.push_back(1);
704	–
705	–	} else {
706	–
707	–	Bool_t vertexFound = kFALSE;
708	–	const Vertex *closestVtx = 0;
709	–	Double_t dzmin = 10000;
710	–
711	–	// loop over vertices
712	–	for(UInt_t j = 0; j < vtxArr->GetEntries(); j++) {
713	–	const Vertex *vtx = vtxArr->At(j);
714	–	assert(vtx);
715	–
716	–	if(pf->HasTrackerTrk() &&
717	–	vtx->HasTrack(pf->TrackerTrk()) &&
718	–	vtx->TrackWeight(pf->TrackerTrk()) > 0) {
719	–	vertexFound = kTRUE;
720	–	closestVtx = vtx;
721	–	break;
722	–	}
723	–	Double_t dz = fabs(pf->SourceVertex().Z() - vtx->Z());
724	–	if(dz < dzmin) {
725	–	closestVtx = vtx;
726	–	dzmin = dz;
727	–	}
728	–	}
729	–
730	–	//      if(fCheckClosestZVertex) {
731	–	if(1) {
732	–	// Fallback: if track is not associated with any vertex,
733	–	// associate it with the vertex closest in z
734	–	if(vertexFound || closestVtx != vtxArr->At(0)) {
735	–	//      pfPileUp->Add(pf);
736	–	pfNoPileUpflag.push_back(0);
737	–	} else {
738	–	pfNoPileUpflag.push_back(1);
739	–	}
740	–	} else {
741	–	if(vertexFound && closestVtx != vtxArr->At(0)) {
742	–	//      pfPileUp->Add(pf);
743	–	pfNoPileUpflag.push_back(0);
744	–	} else {
745	–	//      PFCandidate * pfNoPileUp->AddNew(); // Ridiculous but that's how it is
746	–	pfNoPileUpflag.push_back(1);
747	–	}
748	–	}
749	–	} //hadron & trk stuff
750	–	} else { // hadron
751	–	//      PFCandidate * ptr = pfNoPileUp->AddNew();
752	–	pfNoPileUpflag.push_back(1);
753	–	}
754	–	} // Loop over PF candidates
755	–
756	–	return pfNoPileUpflag.size();
757	–	}
599
600	+	// 2012 only for now ...
601	+	std::string jsonFile;
602	+	if(ctrl.jsonFile=="")
603	+	jsonFile = "./data/Cert_190456-194479_8TeV_PromptReco_Collisions12_JSON.txt";
604	+	else
605	+	jsonFile = ctrl.jsonFile;
606	+	rlrm.AddJSONFile(jsonFile);
607	+	};

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines
>	Changed lines