[ViewVC] Diff of: cvsroot/UserCode/MitHzz4l/Selection/src/applySelection.cc

Comparing UserCode/MitHzz4l/Selection/src/applySelection.cc (file contents):
Revision 1.16 by khahn, Wed May 9 14:57:20 2012 UTC vs.
Revision 1.22 by khahn, Mon May 14 17:27:34 2012 UTC

#	Line 5 \| Line 5
5		#include <iostream> // standard I/O
6		#include <iomanip> // functions to format standard I/O
7		#include <bitset>
8	+	#include <map>
9		using namespace std;
10
11		//
#	Line 50 \| Line 51 \| using namespace std;
51		#include "IsolationSelection.h"
52
53		//#include "Selection.h"
54	<	#include "ReferenceSelection.h"
54	>	//#include "ReferenceSelection.h"
55	>	#include "ReferenceSelection2.h"
56
57		#include "TriggerUtils.h"
58		#include "PassHLT.h"
#	Line 83 \| Line 85 \| *TH1D hpu;**
85		RunLumiRangeMap rlrm;
86		vector<SimpleLepton> failingLeptons ; // for fake estimation
87		vector<SimpleLepton> passingLeptons; // for fake estimation
88	+	vector<unsigned> cutvec;
89	+	vector<vector<unsigned> > zcutvec;
90	+	vector<vector<unsigned> > zzcutvec;
91	+	map<unsigned,float> evtrhoMap;
92	+	vector<string> cutstrs;
93	+	bool passes_HLT_MC;
94	+	vector<bool> PFnoPUflag;;
95
96		//
97		// prototypes
#	Line 92 \| Line 101 \| 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	+	mithep::Array<Vertex> * vtxArr );
108		//----------------------------------------------------------------------------
109
110
#	Line 101 \| Line 114 \| void initRunLumiRangeMap();
114		int main(int argc, char** argv)
115		//----------------------------------------------------------------------------
116		{
117	+	cutstrs.push_back(string("skim0")); //0
118	+	cutstrs.push_back(string("skim1")); //1
119	+	cutstrs.push_back(string("trigger")); //2
120	+	// -------------------------------------------------
121	+	cutstrs.push_back(string("Z candidate")); //3
122	+	cutstrs.push_back(string("good Z1")); //4
123	+	// -------------------------------------------------
124	+	cutstrs.push_back(string("4l")); //5
125	+	cutstrs.push_back(string("ZZ candidate")); //6
126	+	cutstrs.push_back(string("good Z2")); //7
127	+	cutstrs.push_back(string("ZZ 20/10")); //8
128	+	cutstrs.push_back(string("resonance")); //9
129	+	cutstrs.push_back(string("m4l>70")); //10
130	+	cutstrs.push_back(string("m4l>100")); //11
131	+
132	+
133	+
134	+
135		string cmsswpath(CMSSW_BASE);
136		cmsswpath.append("/src");
137		std::bitset<1024> triggerBits;
138		vector<vector<string> > inputFiles;
139		inputFiles.push_back(vector<string>());
140		map<string,unsigned> entrymap; // number of unskimmed entries in each file
141	+	for( int i=0; i<cutstrs.size(); i++ ) cutvec.push_back(0);
142	+	for( int i=0; i<2; i++ ) {
143	+	zcutvec.push_back(vector<unsigned>());
144	+	for( int j=0; j<cutstrs.size(); j++ ) zcutvec[i].push_back(0);
145	+	}
146	+	for( int i=0; i<3; i++ ) {
147	+	zzcutvec.push_back(vector<unsigned>());
148	+	for( int j=0; j<cutstrs.size(); j++ ) zzcutvec[i].push_back(0);
149	+	}
150	+
151
152		//
153		// args
#	Line 178 \| Line 219 \| int main(int argc, char argv)**
219		nt.makeAngleBranch(angles);
220		nt.makeKinematicsBranch(kinematics);
221		nt.makeInfoBranch(evtinfo);
222	<	FOTuple foTree( nt.getFile(), (const char*)"fo");
222	>
223	>	FOTuple foTree( nt.getFile(), (const char*)"FOtree");
224		foTree.makeFailedLeptonBranch(failingLeptons);
225		foTree.makeZLeptonBranch(passingLeptons);
226
#	Line 209 \| Line 251 \| int main(int argc, char argv)**
251		initElectronIDMVA();
252		initElectronIsoMVA();
253		initTrigger();
254	<
254	>
255	>	// hack
256	>	initEvtRhoMap(evtrhoMap);
257
258
259		//
#	Line 248 \| Line 292 \| int main(int argc, char argv)**
292		chain->SetBranchAddress(fMuonName, &muonArr);
293		chain->SetBranchAddress(fPrimVtxName, &vtxArr);
294		chain->SetBranchAddress(fPFCandidateName, &pfArr);
295	<	chain->SetBranchAddress(fPileupInfoName, &puArr);
295	>	if( ctrl.mc ) chain->SetBranchAddress(fPileupInfoName, &puArr);
296		chain->SetBranchAddress(fPileupEnergyDensityName, &puDArr);
297		chain->SetBranchAddress(fTrackName, &trkArr);
298		chain->SetBranchAddress(fTriggerMaskName, &trigMask);
#	Line 280 \| Line 324 \| int main(int argc, char argv)**
324		for(UInt_t ientry=0; ientry<imax; ientry++)
325		{
326		chain->GetEntry(ientry);
327	+
328	+	// events that fail bc ele overlaps w/ reco mu
329	+	/*
330	+	if( info->EvtNum() != 67315 &&
331	+	info->EvtNum() != 288693 ) continue;
332	+	*/
333	+
334	+	// events where they apparently don't veto electrons ...
335	+	/*
336	+	if( info->EvtNum() != 141713 &&
337	+	info->EvtNum() != 178019 &&
338	+	info->EvtNum() != 60823 &&
339	+	info->EvtNum() != 192795 ) continue;
340	+	*/
341	+
342		if(!(ientry%1000)) cerr << "entry: " << ientry << endl;
284	–	if( ientry>100 )break;
343		setPV( ctrl, vtxArr, vtx);
344	+	//if (!(setPV( ctrl, vtxArr, vtx)) ) continue;
345
346	+	float rho = evtrhoMap[info->EvtNum()];
347	+
348	+	//
349	+	// pfNoPU
350	+	//
351	+	//mithep::Array<PFCandidate> pfNoPileUp;
352	+	PFnoPUflag.clear();
353	+	int pfnopu_size = makePFnoPUArray( pfArr, PFnoPUflag, vtxArr );
354	+	assert(pfnopu_size == pfArr->GetEntries());
355
356		//
357		// data/MC
#	Line 302 \| Line 370 \| int main(int argc, char argv)**
370		evtinfo.npu;
371		weights.npuw = getPUWeight(evtinfo.npu);
372		// cout << "weight: " << weights.w << endl;
373	+	//
374	+	// trigger
375	+	//
376	+	if( string(chain->GetFile()->GetEndpointUrl()->GetFile()) != currentFile ) {
377	+	currentFile = string(chain->GetFile()->GetEndpointUrl()->GetFile());
378	+	hltchain->SetBranchAddress(fTriggerTableName, &hltTableStrings);
379	+	hltchain->GetEntry(0);
380	+	hltTable->Clear();
381	+	fillTriggerBits(hltTable, hltTableStrings );
382	+	}
383	+	if( ctrl.debug ) cout << "file is : " << currentFile << endl;
384	+	fillTriggerBits( hltTable, trigMask, triggerBits );
385	+	passes_HLT_MC = true; // passed to apply as extern global, just for sync
386	+	if( !passHLTMC(triggerBits, info->RunNum(), info->EvtNum() ) ) {
387	+	passes_HLT_MC = false;
388	+	}
389		} else {
390		//
391		// JSON
392		//
393	+	/*
394		RunLumiRangeMap::RunLumiPairType rl(info->RunNum(), info->LumiSec());
395		if( !(rlrm.HasRunLumi(rl)) ) {
396		if( ctrl.debug ) cout << "\tfails JSON" << endl;
397		continue;
398		}
399	+	*/
400		//
401		// trigger
402		//
#	Line 322 \| Line 408 \| int main(int argc, char argv)**
408		fillTriggerBits(hltTable, hltTableStrings );
409		}
410		if( ctrl.debug ) cout << "file is : " << currentFile << endl;
411	+	/*
412		fillTriggerBits( hltTable, trigMask, triggerBits );
413		if( !passHLT(triggerBits, info->RunNum(), info->EvtNum() ) ) {
414		if( ctrl.debug ) cout << "\tfails trigger ... " << endl;
415		continue;
416		}
417	+	*/
418		}
419
420		//
#	Line 335 \| Line 423 \| int main(int argc, char argv)**
423		failingLeptons.clear();
424		passingLeptons.clear();
425		EventData ret4l =
426	<	apply_HZZ4L_reference_selection(ctrl, info,
339	<	vtx,
340	<	pfArr,
341	<	puDArr,
342	<	electronArr,
343	<	&electronReferencePreSelection,
344	<	&electronReferenceIDMVASelection,
345	<	&electronReferenceIsoSelection,
346	<	muonArr,
347	<	&muonReferencePreSelection,
348	<	&muonIDPFSelection,
349	<	&muonReferenceIsoSelection);
426	>
427		/*
428	<	apply_HZZ4L_reference_selection(ctrl, info,
429	<	vtx,
430	<	pfArr,
431	<	puArr,
432	<	electronArr,
433	<	&electronPreSelection,
434	<	&electronIDMVASelection,
435	<	&electronIsoMVASelection,
436	<	muonArr,
437	<	&muonPreSelection,
438	<	&muonIDPFSelection,
439	<	&muonIsoMVASelection);
428	>	apply_HZZ4L_reference2_selection(ctrl, info,
429	>	vtx,
430	>	pfArr,
431	>	//rho,
432	>	puDArr,
433	>	electronArr,
434	>	&electronReferencePreSelection,
435	>	&electronReferenceIDMVASelection,
436	>	&electronReferenceIsoSelection,
437	>	muonArr,
438	>	&muonReferencePreSelection,
439	>	&muonIDPFSelection,
440	>	&muonReferenceIsoSelection);
441		*/
442	+	apply_HZZ4L_reference2_selection(ctrl, info,
443	+	vtx,
444	+	pfArr,
445	+	//rho,
446	+	puDArr,
447	+	electronArr,
448	+	&electronReferencePreSelection,
449	+	&electronReferenceIDMVASelection,
450	+	&electronIsoMVASelection,
451	+	muonArr,
452	+	&muonReferencePreSelection,
453	+	&muonIDPFSelection,
454	+	&muonIsoMVASelection);
455	+
456	+
457	+
458
459		if( ctrl.debug ) cout << endl;
460
461	<	// cout << "bits: " << ret4l.status.selectionBits << endl;
461	>
462	>	int Z1type=0, Z2type=0;
463	>	if( ret4l.status.selectionBits.to_ulong() >= PASS_ZCANDIDATE ) {
464	>	if( abs(ret4l.Z1leptons[0].type) == 11 ) Z1type = 11;
465	>	if( abs(ret4l.Z1leptons[0].type) == 13 ) Z1type = 13;
466	>	}
467	>	if( ret4l.status.selectionBits.to_ulong() >= PASS_ZZCANDIDATE ) {
468	>	if( abs(ret4l.Z2leptons[0].type) == 11 ) Z2type = 11;
469	>	if( abs(ret4l.Z2leptons[0].type) == 13 ) Z2type = 13;
470	>	}
471	>
472	>	#ifdef SYNC
473	>	cout << "bits: " << ret4l.status.selectionBits << "\t"
474	>	<< "Z1: " << Z1type << "\t"
475	>	<< "Z2: " << Z2type << endl;
476	>	cout << endl;
477	>	#endif
478
479		if( ret4l.status.pass() ) {
370	–
480		fillAngles(ret4l,angles);
481		fillKinematics(ret4l,kinematics);
482		fillEventInfo(info,evtinfo);
#	Line 399 \| Line 508 \| int main(int argc, char argv)**
508		pass++;
509		// if( pass > 3 ) break;
510
511	<	} else if( ret4l.status.selectionBits.test(PASS_ZCANDIDATE) ) { // save for fakes ...
512	<	cout << "failingLeptons : " << failingLeptons.size() << endl;
513	<	for( int f=0; f<failingLeptons.size(); f++ ) {
514	<	SimpleLepton l = failingLeptons[f];
515	<	cout << "f: " << f << "\t"
516	<	<< "type: " << l.type << "\t"
517	<	<< "pT: " << l.vec->Pt()
518	<	<< endl;
519	<	}
520	<	/*
521	<	cout << "passingLeptons : " << passingLeptons.GetSize() << endl;
413	<	for( int f=0; f<passingLeptons.GetSize(); f++ ) {
414	<	cout << "f: " << f << "\t"
415	<	<< "type: " << passingLeptons.At(f).type << "\t"
416	<	<< "pT: " << passingLeptons.At(f).vec->Pt()
417	<	<< endl;
511	>	} else if( ret4l.status.selectionBits.test(PASS_GOODZ1) ) { // save for fakes ...
512	>	// if(ctrl.debug) {
513	>	cout << "failingLeptons : " << failingLeptons.size() << endl;
514	>	for( int f=0; f<failingLeptons.size(); f++ ) {
515	>	SimpleLepton l = failingLeptons[f];
516	>	cout << "f: " << f << "\t"
517	>	<< "type: " << l.type << "\t"
518	>	<< "pT: " << l.vec.Pt() << "\t"
519	>	<< "eta: " << l.vec.Eta()
520	>	<< endl;
521	>	// }
522		}
419	–	*/
420	–	cout << endl;
523		foTree.Fill();
524		}
525		}
#	Line 426 \| Line 528 \| int main(int argc, char argv)**
528		foTree.getFile()->cd();
529		foTree.getTree()->Write();
530		nt.WriteClose();
531	+
532	+	cout << endl;
533	+	cout << endl;
534	+	for( int i=0; i<cutvec.size(); i++ ) {
535	+	cout << "cut: " << i << "\t"
536	+	<< "pass: " << cutvec[i];
537	+
538	+	if( i>PASS_TRIGGER&&i<=PASS_GOODZ1 )
539	+	cout << "\t2e: " << zcutvec[0][i]
540	+	<< "\t2m: " << zcutvec[1][i];
541	+
542	+	if( i>PASS_ZCANDIDATE )
543	+	cout << "\t4e: " << zzcutvec[0][i]
544	+	<< "\t4m: " << zzcutvec[1][i]
545	+	<< "\t2e2m: " << zzcutvec[2][i];
546	+
547	+	cout << "\t" << cutstrs[i] << endl;
548	+
549	+	cout << endl;
550	+	}
551	+
552		}
553
554		//----------------------------------------------------------------------------
#	Line 454 \| Line 577 \| bool setPV(ControlFlags ctrl,
577		if( ctrl.debug ) cout << "vertex :: " << i << "\tntrks: " << pv->NTracks() << endl;
578
579		// Select best PV for corrected d0; if no PV passing cuts, the first PV in the collection will be used
580	<	// if(!pv->IsValid()) continue;
580	>	if(!pv->IsValid()) continue;
581		if(pv->NTracksFit() < fMinNTracksFit) continue;
582		if(pv->Ndof() < fMinNdof) continue;
583		if(fabs(pv->Z()) > fMaxAbsZ) continue;
#	Line 475 \| Line 598 \| bool setPV(ControlFlags ctrl,
598		if( ctrl.debug) cout << "new PV set, pt : " << best_sumpt << endl;
599		}
600		}
601	<
601	>	// sync
602		if(!bestPV) bestPV = vtxArr->At(0);
603	+	//if(!bestPV) return false;
604		vtx.SetPosition(bestPV->X(),bestPV->Y(),bestPV->Z());
605		vtx.SetErrors(bestPV->XErr(),bestPV->YErr(),bestPV->ZErr());
606		return true;
#	Line 502 \| Line 626 \| void setEffiencyWeights(EventData & evtd
626
627		for( int k=0; k<lepvec.size(); k++ ) {
628		if( abs(lepvec[k].type) == 13 ) {
629	<	mvec.push_back( std::pair<float,float> (fabs(lepvec[k].vec->Eta()), lepvec[k].vec->Pt()) );
630	<	wlegs.push_back( muonPerLegOfflineEfficiencyWeight( fabs(lepvec[k].vec->Eta()),
631	<	lepvec[k].vec->Pt() ) );
629	>	mvec.push_back( std::pair<float,float> (fabs(lepvec[k].vec.Eta()), lepvec[k].vec.Pt()) );
630	>	wlegs.push_back( muonPerLegOfflineEfficiencyWeight( fabs(lepvec[k].vec.Eta()),
631	>	lepvec[k].vec.Pt() ) );
632		} else {
633
634		// now deal with medium vs loose
635		// evec.push_back( std::pair<float,float> (fabs(lepvec[k].vec.Eta()), lepvec[k].vec.Pt()) );
636
637	<	std::pair<float,float> tmppair(fabs(lepvec[k].vec->Eta()), lepvec[k].vec->Pt());
637	>	std::pair<float,float> tmppair(fabs(lepvec[k].vec.Eta()), lepvec[k].vec.Pt());
638		evec.push_back( std::pair<bool, std::pair<float,float> > (lepvec[k].isTight, tmppair) );
639
640	<	// wlegs.push_back( elePerLegOfflineEfficiencyWeight( fabs(lepvec[k].vec->Eta()),
641	<	// lepvec[k].vec->Pt() ) );
640	>	// wlegs.push_back( elePerLegOfflineEfficiencyWeight( fabs(lepvec[k].vec.Eta()),
641	>	// lepvec[k].vec.Pt() ) );
642
643	<	wlegs.push_back( elePerLegOfflineEfficiencyWeight( fabs(lepvec[k].vec->Eta()),
644	<	lepvec[k].vec->Pt(),
643	>	wlegs.push_back( elePerLegOfflineEfficiencyWeight( fabs(lepvec[k].vec.Eta()),
644	>	lepvec[k].vec.Pt(),
645		lepvec[k].isTight ) );
646
647		}
#	Line 563 \| Line 687 \| double getPUWeight(unsigned npu)
687		{
688		return hpu->GetBinContent(hpu->FindBin(npu));
689		}
690	+
691	+	//----------------------------------------------------------------------------
692	+	void initEvtRhoMap( map<unsigned, float> & evtrhoMap )
693	+	//----------------------------------------------------------------------------
694	+	{
695	+	unsigned evt;
696	+	float rho;
697	+
698	+	cout << "initialzing EvtRhoMap ";
699	+	//FILE * f = fopen("evtrho.dat", "r");
700	+	// FILE * f = fopen("allrho.cali.uniq.dat", "r");
701	+	FILE * f = fopen("allrhoAA.cali.uniq.dat", "r");
702	+	int lcount=0;
703	+	while( fscanf( f, "%u:%f", &evt, &rho ) ) {
704	+	if( feof(f) ) break;
705	+	evtrhoMap[evt] = rho;
706	+	lcount++;
707	+	if( !(lcount%1000) ) cout << "."; flush(cout);
708	+	}
709	+	cout << " done" << endl;
710	+	}
711	+
712	+	//----------------------------------------------------------------------------
713	+	unsigned makePFnoPUArray(mithep::Array<PFCandidate> * fPFCandidates,
714	+	vector<bool> &pfNoPileUpflag,
715	+	mithep::Array<Vertex> * vtxArr )
716	+	//----------------------------------------------------------------------------
717	+	{
718	+	for(UInt_t i = 0; i < fPFCandidates->GetEntries(); i++) {
719	+	const PFCandidate *pf = fPFCandidates->At(i);
720	+	assert(pf);
721	+
722	+	if(pf->PFType() == PFCandidate::eHadron) {
723	+	if(pf->HasTrackerTrk() &&
724	+	vtxArr->At(0)->HasTrack(pf->TrackerTrk()) &&
725	+	vtxArr->At(0)->TrackWeight(pf->TrackerTrk()) > 0) {
726	+
727	+	pfNoPileUpflag.push_back(1);
728	+
729	+	} else {
730	+
731	+	Bool_t vertexFound = kFALSE;
732	+	const Vertex *closestVtx = 0;
733	+	Double_t dzmin = 10000;
734	+
735	+	// loop over vertices
736	+	for(UInt_t j = 0; j < vtxArr->GetEntries(); j++) {
737	+	const Vertex *vtx = vtxArr->At(j);
738	+	assert(vtx);
739	+
740	+	if(pf->HasTrackerTrk() &&
741	+	vtx->HasTrack(pf->TrackerTrk()) &&
742	+	vtx->TrackWeight(pf->TrackerTrk()) > 0) {
743	+	vertexFound = kTRUE;
744	+	closestVtx = vtx;
745	+	break;
746	+	}
747	+	Double_t dz = fabs(pf->SourceVertex().Z() - vtx->Z());
748	+	if(dz < dzmin) {
749	+	closestVtx = vtx;
750	+	dzmin = dz;
751	+	}
752	+	}
753	+
754	+	// if(fCheckClosestZVertex) {
755	+	if(1) {
756	+	// Fallback: if track is not associated with any vertex,
757	+	// associate it with the vertex closest in z
758	+	if(vertexFound \|\| closestVtx != vtxArr->At(0)) {
759	+	// pfPileUp->Add(pf);
760	+	pfNoPileUpflag.push_back(0);
761	+	} else {
762	+	pfNoPileUpflag.push_back(1);
763	+	}
764	+	} else {
765	+	if(vertexFound && closestVtx != vtxArr->At(0)) {
766	+	// pfPileUp->Add(pf);
767	+	pfNoPileUpflag.push_back(0);
768	+	} else {
769	+	// PFCandidate * pfNoPileUp->AddNew(); // Ridiculous but that's how it is
770	+	pfNoPileUpflag.push_back(1);
771	+	}
772	+	}
773	+	} //hadron & trk stuff
774	+	} else { // hadron
775	+	// PFCandidate * ptr = pfNoPileUp->AddNew();
776	+	pfNoPileUpflag.push_back(1);
777	+	}
778	+	} // Loop over PF candidates
779	+
780	+	return pfNoPileUpflag.size();
781	+	}
782	+

Diff Legend

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

Comparing UserCode/MitHzz4l/Selection/src/applySelection.cc (file contents): Revision 1.16 by khahn, Wed May 9 14:57:20 2012 UTC vs. Revision 1.22 by khahn, Mon May 14 17:27:34 2012 UTC

Diff Legend

Comparing UserCode/MitHzz4l/Selection/src/applySelection.cc (file contents):
Revision 1.16 by khahn, Wed May 9 14:57:20 2012 UTC vs.
Revision 1.22 by khahn, Mon May 14 17:27:34 2012 UTC