ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/LJMet/Utils/src/CutFlow.cc

Comparing UserCode/LJMet/Utils/src/CutFlow.cc (file contents):
Revision 1.1 by kukartse, Fri Oct 9 19:56:45 2009 UTC vs.
Revision 1.4 by kukartse, Tue Oct 27 22:39:44 2009 UTC

#	Line 26 | Line 26
26		#include "DataFormats/PatCandidates/interface/MET.h"
27		#include "DataFormats/PatCandidates/interface/TriggerPath.h"
28		#include "DataFormats/BeamSpot/interface/BeamSpot.h"
29	<
29	>	#include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
30	>	#include "DataFormats/Math/interface/deltaR.h"
31
32		using namespace edm;
33
#	Line 42 | Line 43 | CutFlow::CutFlow(const edm::ParameterSet
43		c_3->setPrintCount(false);
44		c_4->setPrintCount(false);
45		c_5->setPrintCount(false);
46	+
47	+	muon_index = -1;
48	+	electron_index = -1;
49		}
50
51
#	Line 60 | Line 64 | void
64		CutFlow::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
65		{
66		c_total->count();
67	+	//
68	+	//_____ mu+jets cut flow
69	+	//
70	+	/*
71		while(1){
72	<	if ( is_trigger("HLT_Mu9", iEvent) ) c_1->count();
72	>	//
73	>	//_____ stage 1
74	>	//
75	>	if ( is_trigger("HLT_Mu9", iEvent) ){
76	>	}
77	>	else break;
78	>	if ( has_global_muon("selectedLayer1Muons", iEvent) ){
79	>	c_1->count();
80	>	}
81		else break;
82	+
83		if ( is_muon("selectedLayer1Muons",
84		20.0,
85		2.1,
#	Line 76 | Line 93 | CutFlow::analyze(const edm::Event& iEven
93		c_2->count();
94		}
95		else break;
96	<	cout << "Cut Flow: passed muon" << endl;
96	>	if ( is_jets("selectedLayer1Jets",
97	>	30.0,
98	>	2.4,
99	>	iEvent) ){
100	>	c_3->count();
101	>	}
102	>	else break;
103	>	if ( no_loose_muon("selectedLayer1Muons",
104	>	2.5,
105	>	10.0,
106	>	0.2,
107	>	iEvent) ){
108	>	c_4->count();
109	>	}
110	>	else break;
111	>	if ( no_loose_electron("selectedLayer1Electrons",
112	>	2.5,
113	>	15.0,
114	>	0.2,
115	>	iEvent) ){
116	>	c_5->count();
117	>	}
118	>	else break;
119	>	cout << "Cut Flow: passed stage 5 " << endl;
120	>	break;
121	>	}
122	>	*/
123	>	//
124	>	//_____ e+jets cut flow
125	>	//
126	>	while(1){
127	>	//
128	>	//_____ stage 1
129	>	//
130	>	if ( is_trigger("HLT_Ele15_LW_L1R", iEvent) ){
131	>	}
132	>	else break;
133	>	if ( has_electron("selectedLayer1Electrons", iEvent) ){
134	>	c_1->count();
135	>	}
136	>	else break;
137	>	//
138	>	//_____ stage 2
139	>	//
140	>	cout << "CutFlow: passed stage 1: ";
141	>	if ( is_electron("selectedLayer1Electrons",
142	>	30.0,
143	>	2.4,
144	>	0.02,
145	>	0.1,
146	>	iEvent) ){
147	>	c_2->count();
148	>	cout << endl;
149	>	}
150	>	else{
151	>	cout << endl;
152	>	break;
153	>	}
154	>	//
155	>	//_____ stage 3
156	>	//
157	>	cout << "CutFlow: passed stage 2: ";
158	>	if ( !is_muon("selectedLayer1Muons",
159	>	30.0,
160	>	2.1,
161	>	11,
162	>	0.02,
163	>	10,
164	>	4000000.0,
165	>	6000000.0,
166	>	0.05,
167	>	iEvent) ){
168	>	c_3->count();
169	>	cout << endl;
170	>	}
171	>	else{
172	>	cout << endl;
173	>	break;
174	>	}
175	>	//
176	>	//_____ stage 4
177	>	//
178	>	cout << "CutFlow: passed stage 3: ";
179	>	int n_jets = count_jets("selectedLayer1Jets",
180	>	30.0,
181	>	2.4,
182	>	"selectedLayer1Electrons",
183	>	30.0,
184	>	2.4,
185	>	0.02,
186	>	0.3,
187	>	iEvent);
188	>	if ( n_jets >= 0){
189	>	c_4->count();
190	>	cout << endl;
191	>	}
192	>	else{
193	>	cout << endl;
194	>	break;
195	>	}
196	>	//
197	>	//_____ stage 5
198	>	//
199	>	cout << "CutFlow: passed stage 4: ";
200	>	if ( n_jets >= 4 ){
201	>	c_5->count();
202	>	cout << endl;
203	>	}
204	>	else{
205	>	cout << endl;
206	>	break;
207	>	}
208	>	cout << "Cut Flow: passed stage 5 " << endl;
209		break;
210		}
211		}
#	Line 92 | Line 221 | void
221		CutFlow::endJob() {
222		cout << "========= Cut Flow Report ==========================================" << endl;
223		cout << "Total events: " << c_total->getCount() << endl;
224	<	cout << "Pass stage 1 (trigger+global): " << c_1->getCount() << endl;
225	<	cout << "Pass stage 2 (muon kinematics): " << c_2->getCount() << endl;
226	<	cout << "Pass stage 3 (jets): " << c_3->getCount() << endl;
227	<	cout << "Pass stage 4 (no loose muon): " << c_4->getCount() << endl;
228	<	cout << "Pass stage 5 (no loose electron): " << c_5->getCount() << endl;
224	>	cout << "Pass stage 1: " << c_1->getCount() << endl;
225	>	cout << "Pass stage 2: " << c_2->getCount() << endl;
226	>	cout << "Pass stage 3: " << c_3->getCount() << endl;
227	>	cout << "Pass stage 4: " << c_4->getCount() << endl;
228	>	cout << "Pass stage 5: " << c_5->getCount() << endl;
229		cout << "================================================================" << endl;
230		}
231
232	<	bool CutFlow::is_trigger(std::string mLabel, const edm::Event& iEvent){
232	>	bool CutFlow::is_trigger(std::string mTriggerName, const edm::Event& iEvent){
233		bool result = false;
234		Handle<std::vector<pat::TriggerPath> > trigs;
235		iEvent.getByLabel("patTrigger", trigs);
#	Line 108 | Line 237 | bool CutFlow::is_trigger(std::string mLa
237		trig!=trigs->end();
238		++trig){
239		//cout << "Trigger name and number: " << trig->name() << ",  " << trig->index() << endl;
240	<	if( trig->name().find(mLabel)!=std::string::npos ){
240	>	if( trig->name().find(mTriggerName)!=std::string::npos ){
241		result = trig->wasAccept();
242		}
243		}
244		return result;
245		}
246
247	+	bool CutFlow::has_global_muon(std::string mLabel,
248	+	const edm::Event& iEvent){
249	+	bool result = false;
250	+	//cout << "Cut Flow: passed trigger ";
251	+	Handle< vector< pat::Muon > >     muons;
252	+	iEvent . getByLabel( mLabel, muons );
253	+	//
254	+	//_____ loop over muons _______________________________________________
255	+	//
256	+	MeanCounter passed_muons("");
257	+	for ( vector<pat::Muon>::const_iterator mu = muons -> begin(); mu != muons -> end(); mu++){
258	+	if ( mu->isGlobalMuon() ){
259	+	passed_muons.count();
260	+	}
261	+	}
262	+	result = passed_muons.getCount()>0;
263	+	//cout << endl;
264	+	return result;
265	+	}
266	+
267	+
268		bool CutFlow::is_muon(std::string mLabel,
269		double mPt,
270		double mEta,
#	Line 126 | Line 276 | bool CutFlow::is_muon(std::string mLabel
276		double mRelIso,
277		const edm::Event& iEvent){
278		bool result = false;
279	+	//cout << "Cut Flow: passed stage 1 (mu+jets) or stage 2 (e+jets): ";
280		Handle< vector< pat::Muon > >     muons;
281		iEvent . getByLabel( "selectedLayer1Muons", muons );
282		Handle< reco::BeamSpot >          beamSpotHandle;
#	Line 155 | Line 306 | bool CutFlow::is_muon(std::string mLabel
306		double _trackiso = 10000.0;
307		double _caloiso = 10000.0;
308		double _reliso = 10000.0;
309	+	MeanCounter passed_muons("");
310		for ( vector<pat::Muon>::const_iterator mu = muons -> begin(); mu != muons -> end(); mu++){
311	+	//_____ check that the muon is global
312	+	if( !(mu->isGlobalMuon()) ) continue;
313		_pt = mu->pt();
314		_abseta = fabs(mu->eta());
315	<	_d0 = mu->dB();
315	>	//_d0 = mu->dB(); // global track
316	>	_d0 = -(mu->innerTrack()->dxy(_bs));
317		_chi2ndof = mu->normChi2();
318		_nhits = mu->numberOfValidHits();
319		reco::TrackRef _track = mu->track();
#	Line 178 | Line 333 | bool CutFlow::is_muon(std::string mLabel
333		if (_pt > mPt &&
334		_abseta < mEta &&
335		_nhits >= mNHits &&
336	<	_d0 < mD0 &&
336	>	fabs(_d0) < mD0 &&
337		_chi2ndof < mChi2Ndof &&
338	<	_emveto < mEmVeto &&
339	<	_hadveto < mHadVeto &&
340	<	_reliso < mRelIso
338	>	fabs(_emveto) < mEmVeto &&
339	>	fabs(_hadveto) < mHadVeto &&
340	>	fabs(_reliso) < mRelIso
341		){
342	<	result = true;
343	<	break;
342	>	passed_muons.count();
343	>	muon_index = (int)(mu - muons->begin());
344	>	if (passed_muons.getCount()==1){
345	>	cout << _pt << ", ";
346	>	cout << _abseta << ", ";
347	>	cout << _nhits << ", ";
348	>	cout << _d0 << ", ";
349	>	cout << _chi2ndof << ", ";
350	>	cout << _emveto << ", ";
351	>	cout << _hadveto << ", ";
352	>	cout << _trackiso << ", ";
353	>	cout << _caloiso << ", ";
354	>	cout << _reliso;
355	>	//cout << endl;
356	>	}
357	>	}
358	>	}
359	>	cout << "Number of good muons: " << passed_muons.getCount();
360	>	result = passed_muons.getCount()==1;
361	>	return result;
362	>	}
363	>
364	>
365	>	bool CutFlow::is_jets(std::string mLabel,
366	>	double mPt,
367	>	double mEta,
368	>	const edm::Event& iEvent){
369	>	bool result = false;
370	>
371	>	//cout << "Cut Flow: passed stage 2: ";
372	>
373	>	Handle< vector< pat::Jet > >      jets;
374	>	iEvent . getByLabel( mLabel, jets );
375	>	MeanCounter c_jets("");
376	>	MeanCounter passed_jets("");
377	>	for ( vector<pat::Jet>::const_iterator jet = jets -> begin(); jet != jets -> end(); jet++ ){
378	>	c_jets.count();
379	>	double _pt = jet->pt();
380	>	double _eta = jet->eta();
381	>	if (
382	>	( fabs( _eta ) < mEta ) &&
383	>	( _pt > mPt )
384	>	)
385	>	{
386	>	passed_jets.count();
387	>	}
388	>	if (c_jets.getCount()<=4){
389	>	cout << _pt << ", ";
390	>	cout << _eta;
391	>	}
392	>	}
393	>	//cout << endl;
394	>	result = passed_jets.getCount() >=4;
395	>	return result;
396	>	}
397	>
398	>
399	>	bool CutFlow::no_loose_muon(std::string mLabel,
400	>	double mEta,
401	>	double mPt,
402	>	double mRelIso,
403	>	const edm::Event& iEvent){
404	>	bool result = false;
405	>	//cout << "Cut Flow: passed stage 3: ";
406	>	Handle< vector< pat::Muon > >     muons;
407	>	iEvent . getByLabel( "selectedLayer1Muons", muons );
408	>	//
409	>	//_____ loop over muons _______________________________________________
410	>	//
411	>	double _pt = 0.0;
412	>	double _abseta = 0.0;
413	>	double _trackiso = 10000.0;
414	>	double _caloiso = 10000.0;
415	>	double _reliso = 10000.0;
416	>	MeanCounter passed_muons("");
417	>	for ( vector<pat::Muon>::const_iterator mu = muons -> begin(); mu != muons -> end(); mu++){
418	>	//_____ check that the muon is global
419	>	if( !(mu->isGlobalMuon()) ) continue;
420	>	_pt = mu->pt();
421	>	_abseta = fabs(mu->eta());
422	>	//_____ check that isolation is valid
423	>	if ( !(mu->isIsolationValid()) ) continue;
424	>	_trackiso = mu->isolationR03().sumPt;
425	>	_caloiso = mu->isolationR03().emEt + mu->isolationR03().hadEt;
426	>	_reliso = (_trackiso + _caloiso)/_pt;
427	>	int _index = (int)(mu - muons->begin());
428	>	if (_pt > mPt &&
429	>	_abseta < mEta &&
430	>	_reliso < mRelIso &&
431	>	_index != muon_index
432	>	){
433	>	passed_muons.count();
434	>	}
435	>	}
436	>	result = passed_muons.getCount()==0;
437	>	//cout << endl;
438	>	return result;
439	>	}
440	>
441	>
442	>	bool CutFlow::no_loose_electron(std::string mLabel,
443	>	double mEta,
444	>	double mPt,
445	>	double mRelIso,
446	>	const edm::Event& iEvent){
447	>	bool result = false;
448	>	//cout << "Cut Flow: passed stage 4: ";
449	>	Handle< vector< pat::Electron > >     electrons;
450	>	iEvent . getByLabel( "selectedLayer1Electrons", electrons );
451	>	//
452	>	//_____ loop over electrons _______________________________________________
453	>	//
454	>	double _pt = 0.0;
455	>	double _abseta = 0.0;
456	>	double _trackiso = 10000.0;
457	>	double _caloiso = 10000.0;
458	>	double _reliso = 10000.0;
459	>	MeanCounter passed_electrons("");
460	>	for ( vector<pat::Electron>::const_iterator e = electrons -> begin(); e != electrons -> end(); e++){
461	>	_pt = e->pt();
462	>	_abseta = fabs(e->eta());
463	>	_trackiso = e->dr03IsolationVariables().tkSumPt;
464	>	_caloiso = e->dr03IsolationVariables().ecalRecHitSumEt + e->dr03IsolationVariables().hcalDepth1TowerSumEt + e->dr03IsolationVariables().hcalDepth2TowerSumEt;
465	>	_reliso = (_trackiso + _caloiso)/_pt;
466	>	int _index = (int)(e - electrons->begin());
467	>	if (_pt > mPt &&
468	>	_abseta < mEta &&
469	>	_reliso < mRelIso &&
470	>	_index != electron_index
471	>	){
472	>	passed_electrons.count();
473	>	}
474	>	}
475	>	result = passed_electrons.getCount()==0;
476	>	//cout << endl;
477	>	return result;
478	>	}
479	>
480	>
481	>	bool CutFlow::has_electron(std::string mLabel,
482	>	const edm::Event& iEvent){
483	>	bool result = false;
484	>	Handle< vector< pat::Electron > >     electrons;
485	>	iEvent . getByLabel( mLabel, electrons );
486	>	result = electrons->size()>0;
487	>	return result;
488	>	}
489	>
490	>
491	>	bool CutFlow::is_electron(std::string mLabel,
492	>	double mPt,
493	>	double mEta,
494	>	double mD0,
495	>	double mRelIso,
496	>	const edm::Event& iEvent){
497	>	bool result = false;
498	>	//cout << "Cut Flow: passed stage 1: ";
499	>	Handle< vector< pat::Electron > >     electrons;
500	>	iEvent . getByLabel( "selectedLayer1Electrons", electrons );
501	>	Handle< reco::BeamSpot >          beamSpotHandle;
502	>	iEvent . getByLabel( "offlineBeamSpot", beamSpotHandle);
503	>	//
504	>	//_____ Beam spot _____________________________________________________
505	>	//
506	>	reco::BeamSpot beamSpot;
507	>	if ( beamSpotHandle.isValid() ){
508	>	beamSpot = *beamSpotHandle;
509	>	}
510	>	else{
511	>	edm::LogInfo("TtLJetsAnalyzer")
512	>	<< "No beam spot available from EventSetup \n";
513	>	}
514	>	const reco::BeamSpot::Point _bs(beamSpot.x0(), beamSpot.y0(),beamSpot.z0());
515	>	//
516	>	//_____ loop over electrons _______________________________________________
517	>	//
518	>	double _pt = 0.0;
519	>	double _abseta = 0.0;
520	>	double _d0 = 0.0;
521	>	double _d0_pat = 0.0;
522	>	double _trackiso = 10000.0;
523	>	double _ecaloiso = 10000.0;
524	>	double _ecaloiso2 = 10000.0;
525	>	double _pat_ecaloiso = 10000.0;
526	>	double _hcaloiso1 = 10000.0;
527	>	double _hcaloiso2 = 10000.0;
528	>	double _caloiso = 10000.0;
529	>	double _reliso = 10000.0;
530	>	MeanCounter passed_electrons("");
531	>	for ( vector<pat::Electron>::const_iterator e = electrons -> begin(); e != electrons -> end(); e++){
532	>	//_____ check that the electron is robust-tight
533	>	if( !(e->electronID("eidRobustTight")) ) continue;
534	>	_pt = e->pt();
535	>	_abseta = fabs(e->eta());
536	>	_d0_pat = e->dB();
537	>	_d0 = -(e->gsfTrack()->dxy(_bs));
538	>	_trackiso = e->dr03IsolationVariables().tkSumPt;
539	>	_ecaloiso = e->dr04IsolationVariables().ecalRecHitSumEt;
540	>	_ecaloiso2 = e->dr03IsolationVariables().ecalRecHitSumEt;
541	>	_pat_ecaloiso = e->isolation(pat::ECalIso);
542	>	_hcaloiso1 = e->dr04IsolationVariables().hcalDepth1TowerSumEt;
543	>	_hcaloiso2 = e->dr04IsolationVariables().hcalDepth2TowerSumEt;
544	>	_caloiso = _ecaloiso + _hcaloiso1 + _hcaloiso2;
545	>	//_caloiso = _pat_ecaloiso +_hcaloiso1 + _hcaloiso2;
546	>	_reliso = (_trackiso + _caloiso)/_pt;
547	>	if (_pt > mPt &&
548	>	_abseta < mEta &&
549	>	fabs(_d0) < mD0 &&
550	>	fabs(_reliso) < mRelIso
551	>	){
552	>	passed_electrons.count();
553	>	electron_index = (int)(e - electrons->begin());
554	>	if (passed_electrons.getCount()==1){
555	>	cout << "pt=" << _pt << ", abseta=";
556	>	cout << _abseta << ", d0=";
557	>	cout << _d0 << ", trackiso=";
558	>	cout << _trackiso << ", ecaloiso04=";
559	>	cout << _ecaloiso << ", ecaloiso03=";
560	>	cout << _ecaloiso2 << ", pat_ecaloiso=";
561	>	cout << _pat_ecaloiso << ", hcaloiso1=";
562	>	cout << _hcaloiso1 << ", hcaloiso2=";
563	>	cout << _hcaloiso2 << ", caloiso04=";
564	>	cout << _caloiso << ", reliso=";
565	>	cout << _reliso;
566	>	//cout << endl;
567	>	}
568		}
569		}
570	<	cout << "Cut Flow: passed trigger: " << _pt << ", ";
571	<	cout << _abseta << ", ";
193	<	cout << _nhits << ", ";
194	<	cout << _d0 << ", ";
195	<	cout << _chi2ndof << ", ";
196	<	cout << _emveto << ", ";
197	<	cout << _hadveto << ", ";
198	<	cout << _trackiso << ", ";
199	<	cout << _caloiso << ", ";
200	<	cout << _reliso << ", ";
201	<	cout << endl;
570	>	//cout << "Number of good electrons: " << passed_electrons.getCount();
571	>	result = passed_electrons.getCount()==1;
572		return result;
573		}
574
575
576	+	int CutFlow::count_jets(std::string jLabel,
577	+	double jPt,
578	+	double jEta,
579	+	std::string lLabel,
580	+	double lPt,
581	+	double lEta,
582	+	double lD0,
583	+	double lDR,
584	+	const edm::Event & iEvent){
585	+	MeanCounter _njets("");
586	+	Handle< vector< pat::Jet > >     jets;
587	+	iEvent . getByLabel( jLabel, jets );
588	+	Handle< vector< pat::Electron > >     electrons;
589	+	iEvent . getByLabel( lLabel, electrons );
590	+	Handle< reco::BeamSpot >          beamSpotHandle;
591	+	iEvent . getByLabel( "offlineBeamSpot", beamSpotHandle);
592	+	//
593	+	//_____ Beam spot _____________________________________________________
594	+	//
595	+	reco::BeamSpot beamSpot;
596	+	if ( beamSpotHandle.isValid() ){
597	+	beamSpot = *beamSpotHandle;
598	+	}
599	+	else{
600	+	edm::LogInfo("TtLJetsAnalyzer")
601	+	<< "No beam spot available from EventSetup \n";
602	+	}
603	+	const reco::BeamSpot::Point _bs(beamSpot.x0(), beamSpot.y0(),beamSpot.z0());
604	+	//
605	+	//_____ loop over jets _______________________________________________
606	+	//
607	+	double j_pt = 0.0;
608	+	double j_eta = 0.0;
609	+	for ( vector<pat::Jet>::const_iterator jet = jets -> begin(); jet != jets -> end(); ++jet){
610	+	j_pt = jet->pt();
611	+	j_eta = jet->eta();
612	+	//
613	+	//_____ loop over electrons _______________________________________________
614	+	//
615	+	double e_pt = 0.0;
616	+	double e_abseta = 0.0;
617	+	double e_d0 = 0.0;
618	+	double e_dr = 0.0;
619	+	bool jet_removed = false;
620	+	for ( vector<pat::Electron>::const_iterator e = electrons -> begin(); e != electrons -> end(); e++){
621	+	//_____ check that the electron is robust-tight
622	+	if( !(e->electronID("eidRobustTight")) ) continue;
623	+	e_pt = e->pt();
624	+	e_abseta = fabs(e->eta());
625	+	e_d0 = -(e->gsfTrack()->dxy(_bs));
626	+	e_dr = reco::deltaR(*jet, *e);
627	+	if ( e_pt > lPt && e_abseta < lEta && fabs(e_d0)<lD0 && e_dr<lDR){
628	+	jet_removed = true;
629	+	}
630	+	}
631	+	if (j_pt>jPt && fabs(j_eta)<jEta && !jet_removed) _njets.count();
632	+	}
633	+	cout << "Number of jets after cleaning: " << _njets.getCount();
634	+	return _njets.getCount();
635	+	}
636	+
637		//define this as a plug-in
638		DEFINE_FWK_MODULE(CutFlow);

Diff Legend

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