ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/MitPhysics/Mods/src/MuonIDMod.cc

Comparing UserCode/MitPhysics/Mods/src/MuonIDMod.cc (file contents):
Revision 1.41 by ceballos, Fri Mar 11 15:13:09 2011 UTC vs.
Revision 1.86 by ceballos, Wed Dec 5 19:06:07 2012 UTC

#	Line 14 | Line 14 | ClassImp(mithep::MuonIDMod)
14		//--------------------------------------------------------------------------------------------------
15		MuonIDMod::MuonIDMod(const char *name, const char *title) :
16		BaseMod(name,title),
17	+	fPrintMVADebugInfo(kFALSE),
18		fMuonBranchName(Names::gkMuonBrn),
19		fCleanMuonsName(ModNames::gkCleanMuonsName),
20		fNonIsolatedMuonsName("random"),
21		fNonIsolatedElectronsName("random"),
22		fVertexName(ModNames::gkGoodVertexesName),
23	+	fBeamSpotName(Names::gkBeamSpotBrn),
24		fTrackName(Names::gkTrackBrn),
25		fPFCandidatesName(Names::gkPFCandidatesBrn),
26	<	fMuonIDType("WWMuId"),
27	<	fMuonIsoType("TrackCaloSliding"),
26	>	fPFNoPileUpName("PFNoPileUp"),
27	>	fPFPileUpName("PFPileUp"),
28	>	fMuonIDType("WWMuIdV3"),
29	>	fMuonIsoType("PFIso"),
30		fMuonClassType("Global"),
31		fTrackIsolationCut(3.0),
32		fCaloIsolationCut(3.0),
33		fCombIsolationCut(0.15),
34	+	fCombRelativeIsolationCut(0.15),
35	+	fPFIsolationCut(-999.0),
36		fMuonPtMin(10),
37		fApplyD0Cut(kTRUE),
38	+	fApplyDZCut(kTRUE),
39		fD0Cut(0.020),
40	+	fDZCut(0.10),
41	+	fWhichVertex(-1),
42		fEtaCut(2.4),
34	–	fReverseIsoCut(kFALSE),
35	–	fReverseD0Cut(kFALSE),
43		fMuIDType(kIdUndef),
44		fMuIsoType(kIsoUndef),
45		fMuClassType(kClassUndef),
46		fMuons(0),
47		fVertices(0),
48	+	fBeamSpot(0),
49		fTracks(0),
50		fPFCandidates(0),
51	+	fPFNoPileUpCands(0),
52	+	fPFPileUpCands(0),
53	+	fIntRadius(0.0),
54		fNonIsolatedMuons(0),
55		fNonIsolatedElectrons(0),
56		fPileupEnergyDensityName(Names::gkPileupEnergyDensityBrn),
57	<	fPileupEnergyDensity(0)
57	>	fPileupEnergyDensity(0),
58	>	fMuonTools(0),
59	>	fMuonIDMVA(0),
60	>	fTheRhoType(RhoUtilities::DEFAULT)
61		{
62		// Constructor.
63		}
#	Line 59 | Line 73 | void MuonIDMod::Process()
73		else {
74		fMuons = GetObjThisEvt<MuonOArr>(fMuonBranchName);
75		}
76	+	LoadEventObject(fBeamSpotName, fBeamSpot);
77		LoadEventObject(fTrackName, fTracks);
78		LoadEventObject(fPFCandidatesName, fPFCandidates);
79	<	if(fMuIsoType == kTrackCaloSliding) {
79	>	if(fMuIsoType == kTrackCaloSliding ||
80	>	fMuIsoType == kCombinedRelativeConeAreaCorrected ||
81	>	fMuIsoType == kPFIsoEffectiveAreaCorrected ||
82	>	fMuIsoType == kMVAIso_BDTG_IDIso ||
83	>	fMuIsoType == kIsoRingsV0_BDTG_Iso ||
84	>	fMuIsoType == kIsoDeltaR
85	>	) {
86		LoadEventObject(fPileupEnergyDensityName, fPileupEnergyDensity);
87		}
88	+	if(fMuIsoType == kPFRadialIso || fMuIsoType == kIsoDeltaR || fMuIsoType == kPFIsoBetaPUCorrected){
89	+	// Name is hardcoded, can be changed if someone feels to do it *** did it--Heng
90	+	fPFNoPileUpCands = GetObjThisEvt<PFCandidateCol>(fPFNoPileUpName);
91	+	fPFPileUpCands = GetObjThisEvt<PFCandidateCol>(fPFPileUpName);
92	+	}
93	+
94		MuonOArr *CleanMuons = new MuonOArr;
95		CleanMuons->SetName(fCleanMuonsName);
96
97		fVertices = GetObjThisEvt<VertexOArr>(fVertexName);
98
99	<	for (UInt_t i=0; i<fMuons->GetEntries(); ++i) {
99	>	for (UInt_t i=0; i<fMuons->GetEntries() && fVertices->GetEntries() > 0 ; ++i) {
100		const Muon *mu = fMuons->At(i);
101
102		Bool_t pass = kFALSE;
#	Line 94 | Line 121 | void MuonIDMod::Process()
121		eta = TMath::Abs(mu->Eta());
122		}
123		break;
124	+	case kGlobalTracker:
125	+	pass = (mu->HasGlobalTrk() && mu->GlobalTrk()->Chi2()/mu->GlobalTrk()->Ndof() < 10 &&
126	+	(mu->NSegments() > 1 || mu->NMatches() > 1) && mu->NValidHits() > 0) ||
127	+	(mu->IsTrackerMuon() &&
128	+	mu->Quality().Quality(MuonQuality::TMLastStationTight));
129	+	if (pass) {
130	+	pt  = mu->TrackerTrk()->Pt();
131	+	eta = TMath::Abs(mu->TrackerTrk()->Eta());
132	+	}
133	+	else {
134	+	pt  = mu->Pt();
135	+	eta = TMath::Abs(mu->Eta());
136	+	}
137	+	break;
138		case kSta:
139		pass = mu->HasStandaloneTrk();
140		if (pass) {
#	Line 123 | Line 164 | void MuonIDMod::Process()
164		eta = TMath::Abs(mu->TrackerTrk()->Eta());
165		}
166		break;
167	+	case kGlobalOnly:
168	+	pass = mu->HasGlobalTrk();
169	+	if (pass && mu->TrackerTrk()) {
170	+	pt  = mu->TrackerTrk()->Pt();
171	+	eta = TMath::Abs(mu->TrackerTrk()->Eta());
172	+	}
173	+	else {
174	+	pt  = mu->Pt();
175	+	eta = TMath::Abs(mu->Eta());
176	+	}
177	+	break;
178		default:
179		break;
180		}
#	Line 169 | Line 221 | void MuonIDMod::Process()
221		mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
222		break;
223		case kTight:
224	<	idpass = mu->BestTrk() !=  0 &&
225	<	mu->Quality().Quality(MuonQuality::TMOneStationTight) &&
226	<	mu->Quality().Quality(MuonQuality::TM2DCompatibilityTight) &&
227	<	mu->BestTrk()->NHits() > 10 &&
228	<	RChi2 < 10.0 &&
229	<	mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
224	>	idpass = mu->BestTrk() != 0 &&
225	>	mu->NTrkLayersHit() > 5 &&
226	>	mu->IsPFMuon() == kTRUE &&
227	>	mu->BestTrk()->NPixelHits() > 0 &&
228	>	RChi2 < 10.0;
229	>	break;
230	>	// 2012 WW analysis for 42x (there is no PFMuon link)
231	>	case kWWMuIdV1:
232	>	idpass = mu->BestTrk() != 0 &&
233	>	mu->NTrkLayersHit() > 5 &&
234	>	mu->BestTrk()->NPixelHits() > 0 &&
235	>	mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1 &&
236	>	mu->TrkKink() < 20.0;
237		break;
238	<	case kWWMuId:
238	>	// 2010 WW analysis
239	>	case kWWMuIdV2:
240		idpass = mu->BestTrk() != 0 &&
241		mu->BestTrk()->NHits() > 10 &&
182	–	RChi2 < 10.0 &&
183	–	(mu->NSegments() > 1 || mu->NMatches() > 1) &&
242		mu->BestTrk()->NPixelHits() > 0 &&
185	–	mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight) &&
243		mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1;
244		break;
245	+	// 2011 WW analysis
246	+	case kWWMuIdV3:
247	+	idpass = mu->BestTrk() != 0 &&
248	+	mu->BestTrk()->NHits() > 10 &&
249	+	mu->BestTrk()->NPixelHits() > 0 &&
250	+	mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1 &&
251	+	mu->TrkKink() < 20.0;
252	+	break;
253	+	// 2012 WW analysis
254	+	case kWWMuIdV4:
255	+	idpass = mu->BestTrk() != 0 &&
256	+	mu->NTrkLayersHit() > 5 &&
257	+	mu->IsPFMuon() == kTRUE &&
258	+	mu->BestTrk()->NPixelHits() > 0 &&
259	+	mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1 &&
260	+	mu->TrkKink() < 20.0;
261	+	break;
262	+	case kMVAID_BDTG_IDIso:
263	+	{
264	+	Bool_t passDenominatorM2 = (mu->BestTrk() != 0 &&
265	+	mu->BestTrk()->NHits() > 10 &&
266	+	mu->BestTrk()->NPixelHits() > 0 &&
267	+	mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1 &&
268	+	MuonTools::PassD0Cut(mu, fVertices, 0.20, 0) &&
269	+	MuonTools::PassDZCut(mu, fVertices, 0.10, 0) &&
270	+	mu->TrkKink() < 20.0
271	+	);
272	+	idpass =  passDenominatorM2;
273	+	//only evaluate MVA if muon passes M2 denominator to save time
274	+	if (idpass) idpass = PassMuonMVA_BDTG_IdIso(mu, fVertices->At(0), fPileupEnergyDensity);
275	+	}
276	+	break;
277		case kNoId:
278		idpass = kTRUE;
279		break;
#	Line 195 | Line 284 | void MuonIDMod::Process()
284		if (!idpass)
285		continue;
286
287	+	Double_t Rho = 0.;
288	+	if(fPileupEnergyDensity){
289	+	const PileupEnergyDensity *rho =  fPileupEnergyDensity->At(0);
290	+
291	+	switch (fTheRhoType) {
292	+	case RhoUtilities::MIT_RHO_VORONOI_LOW_ETA:
293	+	Rho = rho->RhoLowEta();
294	+	break;
295	+	case RhoUtilities::MIT_RHO_VORONOI_HIGH_ETA:
296	+	Rho = rho->Rho();
297	+	break;
298	+	case RhoUtilities::MIT_RHO_RANDOM_LOW_ETA:
299	+	Rho = rho->RhoRandomLowEta();
300	+	break;
301	+	case RhoUtilities::MIT_RHO_RANDOM_HIGH_ETA:
302	+	Rho = rho->RhoRandom();
303	+	break;
304	+	case RhoUtilities::CMS_RHO_RHOKT6PFJETS:
305	+	Rho = rho->RhoKt6PFJets();
306	+	break;
307	+	default:
308	+	Rho = rho->Rho();
309	+	}
310	+
311	+	if ((TMath::IsNaN(fPileupEnergyDensity->At(0)->Rho()) || std::isinf(fPileupEnergyDensity->At(0)->Rho()))) Rho = 0.;
312	+	}
313	+
314		Bool_t isocut = kFALSE;
315		switch (fMuIsoType) {
316		case kTrackCalo:
#	Line 208 | Line 324 | void MuonIDMod::Process()
324		break;
325		case kTrackCaloSliding:
326		{
327	<	//Double_t beta = IsolationTools::BetaM(fTracks, mu, fVertices->At(0), 0.0, 0.2, 0.3, 0.02);
328	<	//if(beta == 0) beta = 1.0;
329	<	const PileupEnergyDensity *rho =  fPileupEnergyDensity->At(0);
330	<	Double_t totalIso =  mu->IsoR03SumPt() + TMath::Max(mu->IsoR03EmEt() + mu->IsoR03HadEt() - rho->Rho() * TMath::Pi() * 0.3 * 0.3, 0.0);
331	<	if (totalIso < (mu->Pt()*fCombIsolationCut) )
332	<	isocut = kTRUE;
333	<
334	<	if     (fReverseIsoCut == kTRUE &&
219	<	isocut == kFALSE && totalIso < 10)
220	<	isocut = kTRUE;
221	<	else if(fReverseIsoCut == kTRUE)
222	<	isocut = kFALSE;
327	>	Double_t totalIso =  mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt() - Rho * TMath::Pi() * 0.3 * 0.3 ;
328	>	// trick to change the signal region cut
329	>	double theIsoCut = fCombIsolationCut;
330	>	if(theIsoCut < 0.20){
331	>	if(mu->Pt() >  20.0) theIsoCut = 0.15;
332	>	else                 theIsoCut = 0.10;
333	>	}
334	>	if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
335		}
336		break;
337	<	case kTrackCaloSlidingNoCorrection:
337	>	case kTrackCaloSlidingNoCorrection:
338		{
339		Double_t totalIso =  1.0 * mu->IsoR03SumPt() +
340		1.0 * mu->IsoR03EmEt()  +
341		1.0 * mu->IsoR03HadEt();
342	<	if (totalIso < (mu->Pt()*fCombIsolationCut) )
343	<	isocut = kTRUE;
342	>	// trick to change the signal region cut
343	>	double theIsoCut = fCombIsolationCut;
344	>	if(theIsoCut < 0.20){
345	>	if(mu->Pt() >  20.0) theIsoCut = 0.15;
346	>	else                 theIsoCut = 0.10;
347	>	}
348	>	if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
349		}
350		break;
351	<	case kPFIso:
352	<	{
353	<	Double_t beta = IsolationTools::BetaM(fTracks, mu, fVertices->At(0), 0.0, 0.2, 0.3, 0.02);
354	<	if(beta == 0) beta = 1.0;
355	<	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.2, 0.5, 0.3, 0.02, 0, beta, fNonIsolatedMuons, fNonIsolatedElectrons);
356	<	if (totalIso < (mu->Pt()*fCombIsolationCut) )
351	>	case kCombinedRelativeConeAreaCorrected:
352	>	{
353	>	//const PileupEnergyDensity *rho =  fPileupEnergyDensity->At(0); // Fabian: made Rho customable
354	>	Double_t totalIso =  mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt() - Rho * TMath::Pi() * 0.3 * 0.3 ;
355	>	double theIsoCut = fCombRelativeIsolationCut;
356	>	if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
357	>	}
358	>	break;
359	>	case kCombinedRelativeEffectiveAreaCorrected:
360	>	{
361	>	Double_t tmpRho = Rho;   // Fabian: made the Rho type customable.
362	>	//if (!(TMath::IsNaN(fPileupEnergyDensity->At(0)->Rho()) || std::isinf(fPileupEnergyDensity->At(0)->Rho())))
363	>	//tmpRho = fPileupEnergyDensity->At(0)->Rho();
364	>
365	>	isocut = ( mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt()
366	>	-  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuEMIso03, mu->Eta())
367	>	-  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuHadIso03, mu->Eta())
368	>	) < (mu->Pt()* 0.40);
369	>	}
370	>	break;
371	>	case kPFIso:
372	>	{
373	>	Double_t pfIsoCutValue = 9999;
374	>	if(fPFIsolationCut > 0){
375	>	pfIsoCutValue = fPFIsolationCut;
376	>	} else {
377	>	if (mu->AbsEta() < 1.479) {
378	>	if (mu->Pt() > 20) {
379	>	pfIsoCutValue = 0.13;
380	>	} else {
381	>	pfIsoCutValue = 0.06;
382	>	}
383	>	} else {
384	>	if (mu->Pt() > 20) {
385	>	pfIsoCutValue = 0.09;
386	>	} else {
387	>	pfIsoCutValue = 0.05;
388	>	}
389	>	}
390	>	}
391	>	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius);
392	>	if (totalIso < (mu->Pt()*pfIsoCutValue) )
393	>	isocut = kTRUE;
394	>	}
395	>	break;
396	>	case kPFRadialIso:
397	>	{
398	>	Double_t pfIsoCutValue = 9999;
399	>	if(fPFIsolationCut > 0){
400	>	pfIsoCutValue = fPFIsolationCut;
401	>	} else {
402	>	if (mu->Pt() > 20) {
403	>	pfIsoCutValue = 0.10;
404	>	} else {
405	>	pfIsoCutValue = 0.05;
406	>	}
407	>	}
408	>	Double_t totalIso =  IsolationTools::PFRadialMuonIsolation(mu, fPFNoPileUpCands, 1.0, 0.3);
409	>	if (totalIso < (mu->Pt()*pfIsoCutValue) )
410		isocut = kTRUE;
411	+	}
412	+	break;
413	+	case kPFIsoBetaPUCorrected:
414	+	{
415	+	Double_t pfIsoCutValue = 9999;
416	+	if(fPFIsolationCut > 0){
417	+	pfIsoCutValue = fPFIsolationCut;
418	+	} else {
419	+	if (mu->Pt() > 20) {
420	+	pfIsoCutValue = 0.2;
421	+	} else {
422	+	pfIsoCutValue = 0.2;
423	+	}
424		}
425	<	break;
426	<	case kPFIsoNoL:
425	>	Double_t totalIso =  IsolationTools::BetaMwithPUCorrection(fPFNoPileUpCands, fPFPileUpCands, mu, 0.4);
426	>	if (totalIso < (mu->Pt()*pfIsoCutValue) )
427	>	isocut = kTRUE;
428	>	}
429	>	break;
430	>	case kPFIsoEffectiveAreaCorrected:
431	>	{
432	>	Double_t pfIsoCutValue = 9999;
433	>	if(fPFIsolationCut > 0){
434	>	pfIsoCutValue = fPFIsolationCut;
435	>	} else {
436	>	pfIsoCutValue = fPFIsolationCut; //leave it like this for now
437	>	}
438	>	Double_t EffectiveAreaCorrectedPFIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius)
439	>	- Rho * MuonTools::MuonEffectiveArea(MuonTools::kMuNeutralIso03, mu->Eta());
440	>	//- fPileupEnergyDensity->At(0)->Rho() * MuonTools::MuonEffectiveArea(MuonTools::kMuNeutralIso03, mu->Eta());  // Fabian: made Rho-type customable
441	>	isocut = EffectiveAreaCorrectedPFIso < (mu->Pt() * pfIsoCutValue);
442	>	break;
443	>	}
444	>
445	>
446	>	case kPFIsoNoL:
447		{
448		fNonIsolatedMuons     = GetObjThisEvt<MuonCol>(fNonIsolatedMuonsName);
449		fNonIsolatedElectrons = GetObjThisEvt<ElectronCol>(fNonIsolatedElectronsName);
450
451	<	Double_t beta = IsolationTools::BetaM(fTracks, mu, fVertices->At(0), 0.0, 0.2, 0.3, 0.02);
452	<	if(beta == 0) beta = 1.0;
453	<	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.2, 0.5, 0.3, 0.02, 3, beta, fNonIsolatedMuons, fNonIsolatedElectrons);
454	<	if (totalIso < (mu->Pt()*fCombIsolationCut) )
451	>	Double_t pfIsoCutValue = 9999;
452	>	if(fPFIsolationCut > 0){
453	>	pfIsoCutValue = fPFIsolationCut;
454	>	} else {
455	>	if (mu->AbsEta() < 1.479) {
456	>	if (mu->Pt() > 20) {
457	>	pfIsoCutValue = 0.13;
458	>	} else {
459	>	pfIsoCutValue = 0.06;
460	>	}
461	>	} else {
462	>	if (mu->Pt() > 20) {
463	>	pfIsoCutValue = 0.09;
464	>	} else {
465	>	pfIsoCutValue = 0.05;
466	>	}
467	>	}
468	>	}
469	>	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fNonIsolatedMuons, fNonIsolatedElectrons, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius);
470	>	if (totalIso < (mu->Pt()*pfIsoCutValue) )
471		isocut = kTRUE;
472		}
473		break;
474	+	case kMVAIso_BDTG_IDIso:
475	+	{
476	+
477	+	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius);
478	+	isocut = (totalIso < (mu->Pt()*0.4));
479	+
480	+	}
481	+	break;
482	+	case kIsoRingsV0_BDTG_Iso:
483	+	{
484	+
485	+	isocut = PassMuonIsoRingsV0_BDTG_Iso(mu, fVertices->At(0), fPileupEnergyDensity);
486	+
487	+	}
488	+	break;
489	+	case kIsoDeltaR:
490	+	{
491	+
492	+	isocut = PassMuonIsoDeltaR(mu, fVertices->At(0), fPileupEnergyDensity);
493	+
494	+	}
495	+	break;
496		case kNoIso:
497		isocut = kTRUE;
498		break;
#	Line 263 | Line 504 | void MuonIDMod::Process()
504		if (isocut == kFALSE)
505		continue;
506
507	+	// apply d0 cut
508		if (fApplyD0Cut) {
509	<	Bool_t passD0cut = MuonTools::PassD0Cut(mu, fVertices, fD0Cut);
509	>	Bool_t passD0cut = kTRUE;
510	>	if(fD0Cut < 0.05) { // trick to change the signal region cut
511	>	if      (mu->Pt() >  20.0) fD0Cut = 0.02;
512	>	else if (mu->Pt() <= 20.0) fD0Cut = 0.01;
513	>	}
514	>	if(fWhichVertex >= -1) passD0cut = MuonTools::PassD0Cut(mu, fVertices, fD0Cut, fWhichVertex);
515	>	else                   passD0cut = MuonTools::PassD0Cut(mu, fBeamSpot, fD0Cut);
516		if (!passD0cut)
517		continue;
518		}
519
520	+	// apply dz cut
521	+	if (fApplyDZCut) {
522	+	Bool_t passDZcut = MuonTools::PassDZCut(mu, fVertices, fDZCut, fWhichVertex);
523	+	if (!passDZcut)
524	+	continue;
525	+	}
526	+
527		// add good muon
528		CleanMuons->Add(mu);
529		}
#	Line 290 | Line 545 | void MuonIDMod::SlaveBegin()
545		if (fMuonIsoType.CompareTo("PFIsoNoL") != 0) {
546		ReqEventObject(fMuonBranchName, fMuons, kTRUE);
547		}
548	+	ReqEventObject(fBeamSpotName, fBeamSpot, kTRUE);
549		ReqEventObject(fTrackName, fTracks, kTRUE);
550		ReqEventObject(fPFCandidatesName, fPFCandidates, kTRUE);
551	<	if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0) {
551	>	if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0
552	>	|| fMuonIsoType.CompareTo("CombinedRelativeConeAreaCorrected") == 0
553	>	|| fMuonIsoType.CompareTo("CombinedRelativeEffectiveAreaCorrected") == 0
554	>	|| fMuonIsoType.CompareTo("PFIsoEffectiveAreaCorrected") == 0
555	>	|| fMuonIsoType.CompareTo("MVA_BDTG_IDIso") == 0
556	>	|| fMuonIsoType.CompareTo("IsoRingsV0_BDTG_Iso") == 0
557	>	|| fMuonIsoType.CompareTo("IsoDeltaR") == 0
558	>	) {
559		ReqEventObject(fPileupEnergyDensityName, fPileupEnergyDensity, kTRUE);
560		}
561
562	+
563		if (fMuonIDType.CompareTo("WMuId") == 0)
564		fMuIDType = kWMuId;
565		else if (fMuonIDType.CompareTo("ZMuId") == 0)
#	Line 304 | Line 568 | void MuonIDMod::SlaveBegin()
568		fMuIDType = kTight;
569		else if (fMuonIDType.CompareTo("Loose") == 0)
570		fMuIDType = kLoose;
571	<	else if (fMuonIDType.CompareTo("WWMuId") == 0)
572	<	fMuIDType = kWWMuId;
571	>	else if (fMuonIDType.CompareTo("WWMuIdV1") == 0)
572	>	fMuIDType = kWWMuIdV1;
573	>	else if (fMuonIDType.CompareTo("WWMuIdV2") == 0)
574	>	fMuIDType = kWWMuIdV2;
575	>	else if (fMuonIDType.CompareTo("WWMuIdV3") == 0)
576	>	fMuIDType = kWWMuIdV3;
577	>	else if (fMuonIDType.CompareTo("WWMuIdV4") == 0)
578	>	fMuIDType = kWWMuIdV4;
579		else if (fMuonIDType.CompareTo("NoId") == 0)
580		fMuIDType = kNoId;
581		else if (fMuonIDType.CompareTo("Custom") == 0) {
582		fMuIDType = kCustomId;
583		SendError(kWarning, "SlaveBegin",
584		"Custom muon identification is not yet implemented.");
585	+	} else if (fMuonIDType.CompareTo("MVA_BDTG_IDIso") == 0) {
586	+	fMuIDType = kMVAID_BDTG_IDIso;
587		} else {
588		SendError(kAbortAnalysis, "SlaveBegin",
589		"The specified muon identification %s is not defined.",
#	Line 327 | Line 599 | void MuonIDMod::SlaveBegin()
599		fMuIsoType = kTrackCaloSliding;
600		else if (fMuonIsoType.CompareTo("TrackCaloSlidingNoCorrection") == 0)
601		fMuIsoType = kTrackCaloSlidingNoCorrection;
602	+	else if (fMuonIsoType.CompareTo("CombinedRelativeConeAreaCorrected") == 0)
603	+	fMuIsoType = kCombinedRelativeConeAreaCorrected;
604	+	else if (fMuonIsoType.CompareTo("CombinedRelativeEffectiveAreaCorrected") == 0)
605	+	fMuIsoType = kCombinedRelativeEffectiveAreaCorrected;
606		else if (fMuonIsoType.CompareTo("PFIso") == 0)
607		fMuIsoType = kPFIso;
608	+	else if (fMuonIsoType.CompareTo("PFRadialIso") == 0)
609	+	fMuIsoType = kPFRadialIso;
610	+	else if (fMuonIsoType.CompareTo("PFIsoBetaPUCorrected") == 0)
611	+	fMuIsoType = kPFIsoBetaPUCorrected;
612	+	else if (fMuonIsoType.CompareTo("PFIsoEffectiveAreaCorrected") == 0)
613	+	fMuIsoType = kPFIsoEffectiveAreaCorrected;
614		else if (fMuonIsoType.CompareTo("PFIsoNoL") == 0)
615		fMuIsoType = kPFIsoNoL;
616		else if (fMuonIsoType.CompareTo("NoIso") == 0)
#	Line 337 | Line 619 | void MuonIDMod::SlaveBegin()
619		fMuIsoType = kCustomIso;
620		SendError(kWarning, "SlaveBegin",
621		"Custom muon isolation is not yet implemented.");
622	+	} else if (fMuonIsoType.CompareTo("MVA_BDTG_IDIso") == 0) {
623	+	fMuIsoType = kMVAIso_BDTG_IDIso;
624	+	} else if (fMuonIsoType.CompareTo("IsoRingsV0_BDTG_Iso") == 0) {
625	+	fMuIsoType = kIsoRingsV0_BDTG_Iso;
626	+	} else if (fMuonIsoType.CompareTo("IsoDeltaR") == 0) {
627	+	fMuIsoType = kIsoDeltaR;
628		} else {
629		SendError(kAbortAnalysis, "SlaveBegin",
630		"The specified muon isolation %s is not defined.",
#	Line 348 | Line 636 | void MuonIDMod::SlaveBegin()
636		fMuClassType = kAll;
637		else if (fMuonClassType.CompareTo("Global") == 0)
638		fMuClassType = kGlobal;
639	+	else if (fMuonClassType.CompareTo("GlobalTracker") == 0)
640	+	fMuClassType = kGlobalTracker;
641		else if (fMuonClassType.CompareTo("Standalone") == 0)
642		fMuClassType = kSta;
643		else if (fMuonClassType.CompareTo("TrackerMuon") == 0)
#	Line 356 | Line 646 | void MuonIDMod::SlaveBegin()
646		fMuClassType = kCaloMuon;
647		else if (fMuonClassType.CompareTo("TrackerBased") == 0)
648		fMuClassType = kTrackerBased;
649	+	else if (fMuonClassType.CompareTo("GlobalOnly") == 0)
650	+	fMuClassType = kGlobalOnly;
651		else {
652		SendError(kAbortAnalysis, "SlaveBegin",
653		"The specified muon class %s is not defined.",
654		fMuonClassType.Data());
655		return;
656		}
657	+
658	+
659	+	//If we use MVA ID, need to load MVA weights
660	+	if     (fMuIsoType == kMVAIso_BDTG_IDIso || fMuIDType == kMVAID_BDTG_IDIso) {
661	+	fMuonTools = new MuonTools();
662	+	fMuonIDMVA = new MuonIDMVA();
663	+	fMuonIDMVA->Initialize("BDTG method",
664	+	string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/BarrelPtBin0_IDIsoCombined_BDTG.weights.xml"))),
665	+	string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/EndcapPtBin0_IDIsoCombined_BDTG.weights.xml"))),
666	+	string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/BarrelPtBin1_IDIsoCombined_BDTG.weights.xml"))),
667	+	string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/EndcapPtBin1_IDIsoCombined_BDTG.weights.xml"))),
668	+	string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/BarrelPtBin2_IDIsoCombined_BDTG.weights.xml"))),
669	+	string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/EndcapPtBin2_IDIsoCombined_BDTG.weights.xml"))),
670	+	MuonIDMVA::kIDIsoCombinedDetIso,
671	+	fTheRhoType);
672	+	}
673	+	else if(fMuIsoType == kIsoRingsV0_BDTG_Iso) {
674	+	std::vector<std::string> muonidiso_weightfiles;
675	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_barrel_lowpt.weights.xml"))));
676	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_barrel_highpt.weights.xml"))));
677	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_endcap_lowpt.weights.xml"))));
678	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_endcap_highpt.weights.xml"))));
679	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_tracker.weights.xml"))));
680	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_global.weights.xml"))));
681	+	fMuonTools = new MuonTools();
682	+	fMuonIDMVA = new MuonIDMVA();
683	+	fMuonIDMVA->Initialize("MuonIso_BDTG_IsoRings",
684	+	MuonIDMVA::kIsoRingsV0,
685	+	kTRUE,
686	+	muonidiso_weightfiles,
687	+	fTheRhoType);
688	+	}
689	+	else if(fMuIsoType == kIsoDeltaR) {
690	+	std::vector<std::string> muonidiso_weightfiles;
691	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_santi-V1_LB_BDT.weights.xml"))));
692	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_santi-V1_LE_BDT.weights.xml"))));
693	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_santi-V1_HB_BDT.weights.xml"))));
694	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_santi-V1_HE_BDT.weights.xml"))));
695	+	fMuonTools = new MuonTools();
696	+	fMuonIDMVA = new MuonIDMVA();
697	+	fMuonIDMVA->Initialize("muonHZZ2012IsoDRMVA",
698	+	MuonIDMVA::kIsoDeltaR,
699	+	kTRUE,
700	+	muonidiso_weightfiles,
701	+	fTheRhoType);
702	+	}
703	+
704	+	}
705	+
706	+
707	+	//--------------------------------------------------------------------------------------------------
708	+	Bool_t MuonIDMod::PassMuonMVA_BDTG_IdIso(const Muon *mu, const Vertex *vertex,
709	+	const PileupEnergyDensityCol *PileupEnergyDensity) const
710	+	{
711	+
712	+	const Track *muTrk=0;
713	+	if(mu->HasTrackerTrk())         { muTrk = mu->TrackerTrk();    }
714	+	else if(mu->HasStandaloneTrk()) { muTrk = mu->StandaloneTrk(); }
715	+
716	+	Double_t MVAValue = fMuonIDMVA->MVAValue(mu,vertex,fMuonTools,fPFCandidates,PileupEnergyDensity);
717	+
718	+	Int_t subdet = 0;
719	+	if (fabs(muTrk->Eta()) < 1.479) subdet = 0;
720	+	else subdet = 1;
721	+	Int_t ptBin = 0;
722	+	if (muTrk->Pt() > 14.5) ptBin = 1;
723	+	if (muTrk->Pt() > 20.0) ptBin = 2;
724	+
725	+	Int_t MVABin = -1;
726	+	if      (subdet == 0 && ptBin == 0) MVABin = 0;
727	+	else if (subdet == 1 && ptBin == 0) MVABin = 1;
728	+	else if (subdet == 0 && ptBin == 1) MVABin = 2;
729	+	else if (subdet == 1 && ptBin == 1) MVABin = 3;
730	+	else if (subdet == 0 && ptBin == 2) MVABin = 4;
731	+	else if (subdet == 1 && ptBin == 2) MVABin = 5;
732	+
733	+	Double_t MVACut = -999;
734	+	if      (MVABin == 0) MVACut = -0.5618;
735	+	else if (MVABin == 1) MVACut = -0.3002;
736	+	else if (MVABin == 2) MVACut = -0.4642;
737	+	else if (MVABin == 3) MVACut = -0.2478;
738	+	else if (MVABin == 4) MVACut =  0.1706;
739	+	else if (MVABin == 5) MVACut =  0.8146;
740	+
741	+	if (MVAValue > MVACut) return kTRUE;
742	+	return kFALSE;
743	+	}
744	+
745	+	//--------------------------------------------------------------------------------------------------
746	+	Bool_t MuonIDMod::PassMuonIsoRingsV0_BDTG_Iso(const Muon *mu, const Vertex *vertex,
747	+	const PileupEnergyDensityCol *PileupEnergyDensity) const
748	+	{
749	+
750	+	Bool_t isDebug = kFALSE;
751	+	const Track *muTrk=0;
752	+	if(mu->HasTrackerTrk())         { muTrk = mu->TrackerTrk();    }
753	+	else if(mu->HasStandaloneTrk()) { muTrk = mu->StandaloneTrk(); }
754	+
755	+	ElectronOArr *tempElectrons = new  ElectronOArr;
756	+	MuonOArr     *tempMuons     = new  MuonOArr;
757	+	Double_t MVAValue = fMuonIDMVA->MVAValue(mu,vertex,fMuonTools,fPFCandidates,
758	+	PileupEnergyDensity,MuonTools::kMuEAFall11MC,tempElectrons,tempMuons,isDebug);
759	+	delete tempElectrons;
760	+	delete tempMuons;
761	+
762	+	Int_t MVABin = fMuonIDMVA->GetMVABin(muTrk->Eta(), muTrk->Pt(), mu->IsGlobalMuon(), mu->IsTrackerMuon());
763	+
764	+	Double_t MVACut = -1.0;
765	+	Double_t eta = mu->AbsEta();
766	+	if     (mu->Pt() <  20 && eta <  1.479) MVACut = 0.86;
767	+	else if(mu->Pt() <  20 && eta >= 1.479) MVACut = 0.82;
768	+	else if(mu->Pt() >= 20 && eta <  1.479) MVACut = 0.82;
769	+	else if(mu->Pt() >= 20 && eta >= 1.479) MVACut = 0.86;
770	+
771	+	if(fPFIsolationCut > -1.0) MVACut = fPFIsolationCut;
772	+
773	+	if(isDebug == kTRUE){
774	+	printf("PassMuonIsoRingsV0_BDTG_IsoDebug: %d, pt, eta = %f, %f, rho = %f(%f) : RingsMVA = %f, bin: %d\n",
775	+	GetEventHeader()->EvtNum(),mu->Pt(), mu->Eta(),
776	+	fPileupEnergyDensity->At(0)->Rho(),fPileupEnergyDensity->At(0)->RhoKt6PFJets(),MVAValue,MVABin);
777	+	}
778	+
779	+	if (MVAValue > MVACut) return kTRUE;
780	+	return kFALSE;
781	+	}
782	+
783	+	//--------------------------------------------------------------------------------------------------
784	+	Bool_t MuonIDMod::PassMuonIsoDeltaR(const Muon *mu, const Vertex *vertex,
785	+	const PileupEnergyDensityCol *PileupEnergyDensity) const
786	+	{
787	+
788	+	const Track *muTrk=0;
789	+	if(mu->HasTrackerTrk())         { muTrk = mu->TrackerTrk();    }
790	+	else if(mu->HasStandaloneTrk()) { muTrk = mu->StandaloneTrk(); }
791	+
792	+	ElectronOArr *tempElectrons = new  ElectronOArr;
793	+	MuonOArr     *tempMuons     = new  MuonOArr;
794	+	Double_t MVAValue = fMuonIDMVA->MVAValue(mu,vertex,fMuonTools,fPFNoPileUpCands,
795	+	PileupEnergyDensity,MuonTools::kMuEAFall11MC,tempElectrons,tempMuons,kFALSE);
796	+	delete tempElectrons;
797	+	delete tempMuons;
798	+
799	+	Int_t MVABin = fMuonIDMVA->GetMVABin(muTrk->Eta(), muTrk->Pt(), mu->IsGlobalMuon(), mu->IsTrackerMuon());
800	+
801	+	Double_t MVACut = -999;
802	+	if      (MVABin == 0) MVACut =  0.000;
803	+	else if (MVABin == 1) MVACut =  0.000;
804	+	else if (MVABin == 2) MVACut =  0.000;
805	+	else if (MVABin == 3) MVACut =  0.000;
806	+
807	+	if (MVAValue > MVACut) return kTRUE;
808	+	return kFALSE;
809	+	}
810	+
811	+	//--------------------------------------------------------------------------------------------------
812	+	void MuonIDMod::Terminate()
813	+	{
814	+	// Run finishing code on the computer (slave) that did the analysis
815	+	delete fMuonIDMVA;
816	+
817	+	delete fMuonTools;
818		}

Diff Legend

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