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.78 by ceballos, Sun May 6 10:32:30 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	>	fMuonIDType("WWMuIdV3"),
27	>	fMuonIsoType("PFIso"),
28		fMuonClassType("Global"),
29		fTrackIsolationCut(3.0),
30		fCaloIsolationCut(3.0),
31		fCombIsolationCut(0.15),
32	+	fCombRelativeIsolationCut(0.15),
33	+	fPFIsolationCut(-1.0),
34		fMuonPtMin(10),
35		fApplyD0Cut(kTRUE),
36	+	fApplyDZCut(kTRUE),
37		fD0Cut(0.020),
38	+	fDZCut(0.10),
39	+	fWhichVertex(-1),
40		fEtaCut(2.4),
34	–	fReverseIsoCut(kFALSE),
35	–	fReverseD0Cut(kFALSE),
41		fMuIDType(kIdUndef),
42		fMuIsoType(kIsoUndef),
43		fMuClassType(kClassUndef),
44		fMuons(0),
45		fVertices(0),
46	+	fBeamSpot(0),
47		fTracks(0),
48		fPFCandidates(0),
49	+	fPFNoPileUpCands(0),
50	+	fIntRadius(0.0),
51		fNonIsolatedMuons(0),
52		fNonIsolatedElectrons(0),
53		fPileupEnergyDensityName(Names::gkPileupEnergyDensityBrn),
54	<	fPileupEnergyDensity(0)
54	>	fPileupEnergyDensity(0),
55	>	fMuonTools(0),
56	>	fMuonIDMVA(0),
57	>	fTheRhoType(RhoUtilities::DEFAULT)
58		{
59		// Constructor.
60		}
#	Line 59 | Line 70 | void MuonIDMod::Process()
70		else {
71		fMuons = GetObjThisEvt<MuonOArr>(fMuonBranchName);
72		}
73	+	LoadEventObject(fBeamSpotName, fBeamSpot);
74		LoadEventObject(fTrackName, fTracks);
75		LoadEventObject(fPFCandidatesName, fPFCandidates);
76	<	if(fMuIsoType == kTrackCaloSliding) {
76	>	if(fMuIsoType == kTrackCaloSliding ||
77	>	fMuIsoType == kCombinedRelativeConeAreaCorrected ||
78	>	fMuIsoType == kPFIsoEffectiveAreaCorrected ||
79	>	fMuIsoType == kMVAIso_BDTG_IDIso ||
80	>	fMuIsoType == kIsoRingsV0_BDTG_Iso ||
81	>	fMuIsoType == kIsoDeltaR
82	>	) {
83		LoadEventObject(fPileupEnergyDensityName, fPileupEnergyDensity);
84		}
85	+	if(fMuIsoType == kPFRadialIso || fMuIsoType == kIsoDeltaR){
86	+	// Name is hardcoded, can be changed if someone feels to do it
87	+	fPFNoPileUpCands = GetObjThisEvt<PFCandidateCol>("PFNoPileUp");
88	+	}
89	+
90		MuonOArr *CleanMuons = new MuonOArr;
91		CleanMuons->SetName(fCleanMuonsName);
92
93		fVertices = GetObjThisEvt<VertexOArr>(fVertexName);
94
95	<	for (UInt_t i=0; i<fMuons->GetEntries(); ++i) {
95	>	for (UInt_t i=0; i<fMuons->GetEntries() && fVertices->GetEntries() > 0 ; ++i) {
96		const Muon *mu = fMuons->At(i);
97
98		Bool_t pass = kFALSE;
#	Line 94 | Line 117 | void MuonIDMod::Process()
117		eta = TMath::Abs(mu->Eta());
118		}
119		break;
120	+	case kGlobalTracker:
121	+	pass = (mu->HasGlobalTrk() && mu->GlobalTrk()->Chi2()/mu->GlobalTrk()->Ndof() < 10 &&
122	+	(mu->NSegments() > 1 || mu->NMatches() > 1) && mu->NValidHits() > 0) ||
123	+	(mu->IsTrackerMuon() &&
124	+	mu->Quality().Quality(MuonQuality::TMLastStationTight));
125	+	if (pass) {
126	+	pt  = mu->TrackerTrk()->Pt();
127	+	eta = TMath::Abs(mu->TrackerTrk()->Eta());
128	+	}
129	+	else {
130	+	pt  = mu->Pt();
131	+	eta = TMath::Abs(mu->Eta());
132	+	}
133	+	break;
134		case kSta:
135		pass = mu->HasStandaloneTrk();
136		if (pass) {
#	Line 136 | Line 173 | void MuonIDMod::Process()
173		if (eta >= fEtaCut)
174		continue;
175
176	+
177	+	//***********************************************************************************************
178	+	//Debug Info For Lepton MVA
179	+	//***********************************************************************************************
180	+	if( fPrintMVADebugInfo &&
181	+	(fMuIsoType == kMVAIso_BDTG_IDIso || fMuIDType == kMVAID_BDTG_IDIso)
182	+	) {
183	+	cout << "Event: " << GetEventHeader()->RunNum() << " " << GetEventHeader()->LumiSec() << " "
184	+	<< GetEventHeader()->EvtNum() << " : Rho = " << fPileupEnergyDensity->At(0)->Rho()
185	+	<< " : Muon " << i << " "
186	+	<< endl;
187	+	fMuonIDMVA->MVAValue(mu,fVertices->At(0),fMuonTools,fPFCandidates,fPileupEnergyDensity,kTRUE);
188	+	}
189	+	//***********************************************************************************************
190	+
191	+
192		Double_t RChi2 = 0.0;
193		if     (mu->HasGlobalTrk()) {
194		RChi2 = mu->GlobalTrk()->Chi2()/mu->GlobalTrk()->Ndof();
#	Line 176 | Line 229 | void MuonIDMod::Process()
229		RChi2 < 10.0 &&
230		mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
231		break;
232	<	case kWWMuId:
232	>	case kWWMuIdV1:
233		idpass = mu->BestTrk() != 0 &&
234		mu->BestTrk()->NHits() > 10 &&
235	+	mu->BestTrk()->NPixelHits() > 0 &&
236	+	mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1 &&
237		RChi2 < 10.0 &&
238		(mu->NSegments() > 1 || mu->NMatches() > 1) &&
239	+	mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
240	+	break;
241	+	case kWWMuIdV2:
242	+	idpass = mu->BestTrk() != 0 &&
243	+	mu->BestTrk()->NHits() > 10 &&
244		mu->BestTrk()->NPixelHits() > 0 &&
185	–	mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight) &&
245		mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1;
246		break;
247	+	case kWWMuIdV3:
248	+	idpass = mu->BestTrk() != 0 &&
249	+	mu->BestTrk()->NHits() > 10 &&
250	+	mu->BestTrk()->NPixelHits() > 0 &&
251	+	mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1 &&
252	+	mu->TrkKink() < 20.0;
253	+	break;
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	+	default:
305	+	Rho = rho->Rho();
306	+	}
307	+
308	+	if ((TMath::IsNaN(fPileupEnergyDensity->At(0)->Rho()) || std::isinf(fPileupEnergyDensity->At(0)->Rho()))) Rho = 0.;
309	+	}
310	+
311		Bool_t isocut = kFALSE;
312		switch (fMuIsoType) {
313		case kTrackCalo:
#	Line 208 | Line 321 | void MuonIDMod::Process()
321		break;
322		case kTrackCaloSliding:
323		{
324	<	//Double_t beta = IsolationTools::BetaM(fTracks, mu, fVertices->At(0), 0.0, 0.2, 0.3, 0.02);
325	<	//if(beta == 0) beta = 1.0;
326	<	const PileupEnergyDensity *rho =  fPileupEnergyDensity->At(0);
327	<	Double_t totalIso =  mu->IsoR03SumPt() + TMath::Max(mu->IsoR03EmEt() + mu->IsoR03HadEt() - rho->Rho() * TMath::Pi() * 0.3 * 0.3, 0.0);
328	<	if (totalIso < (mu->Pt()*fCombIsolationCut) )
329	<	isocut = kTRUE;
330	<
331	<	if     (fReverseIsoCut == kTRUE &&
219	<	isocut == kFALSE && totalIso < 10)
220	<	isocut = kTRUE;
221	<	else if(fReverseIsoCut == kTRUE)
222	<	isocut = kFALSE;
324	>	Double_t totalIso =  mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt() - Rho * TMath::Pi() * 0.3 * 0.3 ;
325	>	// trick to change the signal region cut
326	>	double theIsoCut = fCombIsolationCut;
327	>	if(theIsoCut < 0.20){
328	>	if(mu->Pt() >  20.0) theIsoCut = 0.15;
329	>	else                 theIsoCut = 0.10;
330	>	}
331	>	if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
332		}
333		break;
334		case kTrackCaloSlidingNoCorrection:
#	Line 227 | Line 336 | void MuonIDMod::Process()
336		Double_t totalIso =  1.0 * mu->IsoR03SumPt() +
337		1.0 * mu->IsoR03EmEt()  +
338		1.0 * mu->IsoR03HadEt();
339	<	if (totalIso < (mu->Pt()*fCombIsolationCut) )
340	<	isocut = kTRUE;
339	>	// trick to change the signal region cut
340	>	double theIsoCut = fCombIsolationCut;
341	>	if(theIsoCut < 0.20){
342	>	if(mu->Pt() >  20.0) theIsoCut = 0.15;
343	>	else                 theIsoCut = 0.10;
344	>	}
345	>	if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
346		}
347		break;
348	+	case kCombinedRelativeConeAreaCorrected:
349	+	{
350	+	//const PileupEnergyDensity *rho =  fPileupEnergyDensity->At(0); // Fabian: made Rho customable
351	+	Double_t totalIso =  mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt() - Rho * TMath::Pi() * 0.3 * 0.3 ;
352	+	double theIsoCut = fCombRelativeIsolationCut;
353	+	if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
354	+	}
355	+	break;
356	+	case kCombinedRelativeEffectiveAreaCorrected:
357	+	{
358	+	Double_t tmpRho = Rho;   // Fabian: made the Rho type customable.
359	+	//if (!(TMath::IsNaN(fPileupEnergyDensity->At(0)->Rho()) || std::isinf(fPileupEnergyDensity->At(0)->Rho())))
360	+	//tmpRho = fPileupEnergyDensity->At(0)->Rho();
361	+
362	+	isocut = ( mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt()
363	+	-  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuEMIso03, mu->Eta())
364	+	-  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuHadIso03, mu->Eta())
365	+	) < (mu->Pt()* 0.40);
366	+	}
367	+	break;
368		case kPFIso:
369		{
370	<	Double_t beta = IsolationTools::BetaM(fTracks, mu, fVertices->At(0), 0.0, 0.2, 0.3, 0.02);
371	<	if(beta == 0) beta = 1.0;
372	<	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.2, 0.5, 0.3, 0.02, 0, beta, fNonIsolatedMuons, fNonIsolatedElectrons);
373	<	if (totalIso < (mu->Pt()*fCombIsolationCut) )
370	>	Double_t pfIsoCutValue = 9999;
371	>	if(fPFIsolationCut > 0){
372	>	pfIsoCutValue = fPFIsolationCut;
373	>	} else {
374	>	if (mu->AbsEta() < 1.479) {
375	>	if (mu->Pt() > 20) {
376	>	pfIsoCutValue = 0.13;
377	>	} else {
378	>	pfIsoCutValue = 0.06;
379	>	}
380	>	} else {
381	>	if (mu->Pt() > 20) {
382	>	pfIsoCutValue = 0.09;
383	>	} else {
384	>	pfIsoCutValue = 0.05;
385	>	}
386	>	}
387	>	}
388	>	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius);
389	>	if (totalIso < (mu->Pt()*pfIsoCutValue) )
390	>	isocut = kTRUE;
391	>	}
392	>	break;
393	>	case kPFRadialIso:
394	>	{
395	>	Double_t pfIsoCutValue = 9999;
396	>	if(fPFIsolationCut > 0){
397	>	pfIsoCutValue = fPFIsolationCut;
398	>	} else {
399	>	if (mu->Pt() > 20) {
400	>	pfIsoCutValue = 0.10;
401	>	} else {
402	>	pfIsoCutValue = 0.05;
403	>	}
404	>	}
405	>	Double_t totalIso =  IsolationTools::PFRadialMuonIsolation(mu, fPFNoPileUpCands, 1.0, 0.3);
406	>	if (totalIso < (mu->Pt()*pfIsoCutValue) )
407		isocut = kTRUE;
408		}
409		break;
410	+	case kPFIsoEffectiveAreaCorrected:
411	+	{
412	+	Double_t pfIsoCutValue = 9999;
413	+	if(fPFIsolationCut > 0){
414	+	pfIsoCutValue = fPFIsolationCut;
415	+	} else {
416	+	pfIsoCutValue = fPFIsolationCut; //leave it like this for now
417	+	}
418	+	Double_t EffectiveAreaCorrectedPFIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius)
419	+	- Rho * MuonTools::MuonEffectiveArea(MuonTools::kMuNeutralIso03, mu->Eta());
420	+	//- fPileupEnergyDensity->At(0)->Rho() * MuonTools::MuonEffectiveArea(MuonTools::kMuNeutralIso03, mu->Eta());  // Fabian: made Rho-type customable
421	+	isocut = EffectiveAreaCorrectedPFIso < (mu->Pt() * pfIsoCutValue);
422	+	break;
423	+	}
424		case kPFIsoNoL:
425		{
426		fNonIsolatedMuons     = GetObjThisEvt<MuonCol>(fNonIsolatedMuonsName);
427		fNonIsolatedElectrons = GetObjThisEvt<ElectronCol>(fNonIsolatedElectronsName);
428
429	<	Double_t beta = IsolationTools::BetaM(fTracks, mu, fVertices->At(0), 0.0, 0.2, 0.3, 0.02);
430	<	if(beta == 0) beta = 1.0;
431	<	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.2, 0.5, 0.3, 0.02, 3, beta, fNonIsolatedMuons, fNonIsolatedElectrons);
432	<	if (totalIso < (mu->Pt()*fCombIsolationCut) )
429	>	Double_t pfIsoCutValue = 9999;
430	>	if(fPFIsolationCut > 0){
431	>	pfIsoCutValue = fPFIsolationCut;
432	>	} else {
433	>	if (mu->AbsEta() < 1.479) {
434	>	if (mu->Pt() > 20) {
435	>	pfIsoCutValue = 0.13;
436	>	} else {
437	>	pfIsoCutValue = 0.06;
438	>	}
439	>	} else {
440	>	if (mu->Pt() > 20) {
441	>	pfIsoCutValue = 0.09;
442	>	} else {
443	>	pfIsoCutValue = 0.05;
444	>	}
445	>	}
446	>	}
447	>	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fNonIsolatedMuons, fNonIsolatedElectrons, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius);
448	>	if (totalIso < (mu->Pt()*pfIsoCutValue) )
449		isocut = kTRUE;
450		}
451		break;
452	+	case kMVAIso_BDTG_IDIso:
453	+	{
454	+
455	+	// **************************************************************************
456	+	// Don't use effective area correction denominator. Instead use the old one.
457	+	// **************************************************************************
458	+
459	+	//         Double_t tmpRho = 0;
460	+	//         if (!(TMath::IsNaN(fPileupEnergyDensity->At(0)->Rho()) || isinf(fPileupEnergyDensity->At(0)->Rho())))
461	+	//           tmpRho = fPileupEnergyDensity->At(0)->Rho();
462	+
463	+	//         isocut = ( mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt()
464	+	//                    -  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuEMIso03, mu->Eta())
465	+	//                    -  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuHadIso03, mu->Eta())
466	+	//           ) < (mu->Pt()* 0.40);
467	+
468	+	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius);
469	+	isocut = (totalIso < (mu->Pt()*0.4));
470	+
471	+	}
472	+	break;
473	+	case kIsoRingsV0_BDTG_Iso:
474	+	{
475	+
476	+	isocut = PassMuonIsoRingsV0_BDTG_Iso(mu, fVertices->At(0), fPileupEnergyDensity);
477	+
478	+	}
479	+	break;
480	+	case kIsoDeltaR:
481	+	{
482	+
483	+	isocut = PassMuonIsoDeltaR(mu, fVertices->At(0), fPileupEnergyDensity);
484	+
485	+	}
486	+	break;
487		case kNoIso:
488		isocut = kTRUE;
489		break;
#	Line 263 | Line 495 | void MuonIDMod::Process()
495		if (isocut == kFALSE)
496		continue;
497
498	+	// apply d0 cut
499		if (fApplyD0Cut) {
500	<	Bool_t passD0cut = MuonTools::PassD0Cut(mu, fVertices, fD0Cut);
500	>	Bool_t passD0cut = kTRUE;
501	>	if(fD0Cut < 0.05) { // trick to change the signal region cut
502	>	if      (mu->Pt() >  20.0) fD0Cut = 0.02;
503	>	else if (mu->Pt() <= 20.0) fD0Cut = 0.01;
504	>	}
505	>	if(fWhichVertex >= -1) passD0cut = MuonTools::PassD0Cut(mu, fVertices, fD0Cut, fWhichVertex);
506	>	else                   passD0cut = MuonTools::PassD0Cut(mu, fBeamSpot, fD0Cut);
507		if (!passD0cut)
508		continue;
509		}
510
511	+	// apply dz cut
512	+	if (fApplyDZCut) {
513	+	Bool_t passDZcut = MuonTools::PassDZCut(mu, fVertices, fDZCut, fWhichVertex);
514	+	if (!passDZcut)
515	+	continue;
516	+	}
517	+
518		// add good muon
519		CleanMuons->Add(mu);
520		}
#	Line 290 | Line 536 | void MuonIDMod::SlaveBegin()
536		if (fMuonIsoType.CompareTo("PFIsoNoL") != 0) {
537		ReqEventObject(fMuonBranchName, fMuons, kTRUE);
538		}
539	+	ReqEventObject(fBeamSpotName, fBeamSpot, kTRUE);
540		ReqEventObject(fTrackName, fTracks, kTRUE);
541		ReqEventObject(fPFCandidatesName, fPFCandidates, kTRUE);
542	<	if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0) {
542	>	if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0
543	>	|| fMuonIsoType.CompareTo("CombinedRelativeConeAreaCorrected") == 0
544	>	|| fMuonIsoType.CompareTo("CombinedRelativeEffectiveAreaCorrected") == 0
545	>	|| fMuonIsoType.CompareTo("PFIsoEffectiveAreaCorrected") == 0
546	>	|| fMuonIsoType.CompareTo("MVA_BDTG_IDIso") == 0
547	>	|| fMuonIsoType.CompareTo("IsoRingsV0_BDTG_Iso") == 0
548	>	|| fMuonIsoType.CompareTo("IsoDeltaR") == 0
549	>	) {
550		ReqEventObject(fPileupEnergyDensityName, fPileupEnergyDensity, kTRUE);
551		}
552
553	+
554		if (fMuonIDType.CompareTo("WMuId") == 0)
555		fMuIDType = kWMuId;
556		else if (fMuonIDType.CompareTo("ZMuId") == 0)
#	Line 304 | Line 559 | void MuonIDMod::SlaveBegin()
559		fMuIDType = kTight;
560		else if (fMuonIDType.CompareTo("Loose") == 0)
561		fMuIDType = kLoose;
562	<	else if (fMuonIDType.CompareTo("WWMuId") == 0)
563	<	fMuIDType = kWWMuId;
562	>	else if (fMuonIDType.CompareTo("WWMuIdV1") == 0)
563	>	fMuIDType = kWWMuIdV1;
564	>	else if (fMuonIDType.CompareTo("WWMuIdV2") == 0)
565	>	fMuIDType = kWWMuIdV2;
566	>	else if (fMuonIDType.CompareTo("WWMuIdV3") == 0)
567	>	fMuIDType = kWWMuIdV3;
568		else if (fMuonIDType.CompareTo("NoId") == 0)
569		fMuIDType = kNoId;
570		else if (fMuonIDType.CompareTo("Custom") == 0) {
571		fMuIDType = kCustomId;
572		SendError(kWarning, "SlaveBegin",
573		"Custom muon identification is not yet implemented.");
574	+	} else if (fMuonIDType.CompareTo("MVA_BDTG_IDIso") == 0) {
575	+	fMuIDType = kMVAID_BDTG_IDIso;
576		} else {
577		SendError(kAbortAnalysis, "SlaveBegin",
578		"The specified muon identification %s is not defined.",
#	Line 327 | Line 588 | void MuonIDMod::SlaveBegin()
588		fMuIsoType = kTrackCaloSliding;
589		else if (fMuonIsoType.CompareTo("TrackCaloSlidingNoCorrection") == 0)
590		fMuIsoType = kTrackCaloSlidingNoCorrection;
591	+	else if (fMuonIsoType.CompareTo("CombinedRelativeConeAreaCorrected") == 0)
592	+	fMuIsoType = kCombinedRelativeConeAreaCorrected;
593	+	else if (fMuonIsoType.CompareTo("CombinedRelativeEffectiveAreaCorrected") == 0)
594	+	fMuIsoType = kCombinedRelativeEffectiveAreaCorrected;
595		else if (fMuonIsoType.CompareTo("PFIso") == 0)
596		fMuIsoType = kPFIso;
597	+	else if (fMuonIsoType.CompareTo("PFRadialIso") == 0)
598	+	fMuIsoType = kPFRadialIso;
599	+	else if (fMuonIsoType.CompareTo("PFIsoEffectiveAreaCorrected") == 0)
600	+	fMuIsoType = kPFIsoEffectiveAreaCorrected;
601		else if (fMuonIsoType.CompareTo("PFIsoNoL") == 0)
602		fMuIsoType = kPFIsoNoL;
603		else if (fMuonIsoType.CompareTo("NoIso") == 0)
#	Line 337 | Line 606 | void MuonIDMod::SlaveBegin()
606		fMuIsoType = kCustomIso;
607		SendError(kWarning, "SlaveBegin",
608		"Custom muon isolation is not yet implemented.");
609	+	} else if (fMuonIsoType.CompareTo("MVA_BDTG_IDIso") == 0) {
610	+	fMuIsoType = kMVAIso_BDTG_IDIso;
611	+	} else if (fMuonIsoType.CompareTo("IsoRingsV0_BDTG_Iso") == 0) {
612	+	fMuIsoType = kIsoRingsV0_BDTG_Iso;
613	+	} else if (fMuonIsoType.CompareTo("IsoDeltaR") == 0) {
614	+	fMuIsoType = kIsoDeltaR;
615		} else {
616		SendError(kAbortAnalysis, "SlaveBegin",
617		"The specified muon isolation %s is not defined.",
#	Line 348 | Line 623 | void MuonIDMod::SlaveBegin()
623		fMuClassType = kAll;
624		else if (fMuonClassType.CompareTo("Global") == 0)
625		fMuClassType = kGlobal;
626	+	else if (fMuonClassType.CompareTo("GlobalTracker") == 0)
627	+	fMuClassType = kGlobalTracker;
628		else if (fMuonClassType.CompareTo("Standalone") == 0)
629		fMuClassType = kSta;
630		else if (fMuonClassType.CompareTo("TrackerMuon") == 0)
#	Line 362 | Line 639 | void MuonIDMod::SlaveBegin()
639		fMuonClassType.Data());
640		return;
641		}
642	+
643	+
644	+	//If we use MVA ID, need to load MVA weights
645	+	if     (fMuIsoType == kMVAIso_BDTG_IDIso || fMuIDType == kMVAID_BDTG_IDIso) {
646	+	fMuonTools = new MuonTools();
647	+	fMuonIDMVA = new MuonIDMVA();
648	+	fMuonIDMVA->Initialize("BDTG method",
649	+	string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/BarrelPtBin0_IDIsoCombined_BDTG.weights.xml"))),
650	+	string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/EndcapPtBin0_IDIsoCombined_BDTG.weights.xml"))),
651	+	string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/BarrelPtBin1_IDIsoCombined_BDTG.weights.xml"))),
652	+	string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/EndcapPtBin1_IDIsoCombined_BDTG.weights.xml"))),
653	+	string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/BarrelPtBin2_IDIsoCombined_BDTG.weights.xml"))),
654	+	string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/EndcapPtBin2_IDIsoCombined_BDTG.weights.xml"))),
655	+	MuonIDMVA::kIDIsoCombinedDetIso);
656	+	}
657	+	else if(fMuIsoType == kIsoRingsV0_BDTG_Iso) {
658	+	std::vector<std::string> muonidiso_weightfiles;
659	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_barrel_lowpt.weights.xml"))));
660	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_barrel_highpt.weights.xml"))));
661	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_endcap_lowpt.weights.xml"))));
662	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_endcap_highpt.weights.xml"))));
663	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_tracker.weights.xml"))));
664	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_global.weights.xml"))));
665	+	fMuonTools = new MuonTools();
666	+	fMuonIDMVA = new MuonIDMVA();
667	+	fMuonIDMVA->Initialize("MuonIso_BDTG_IsoRings",
668	+	MuonIDMVA::kIsoRingsV0,
669	+	kTRUE,
670	+	muonidiso_weightfiles);
671	+	}
672	+	else if(fMuIsoType == kIsoDeltaR) {
673	+	std::vector<std::string> muonidiso_weightfiles;
674	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_santi-V1_LB_BDT.weights.xml"))));
675	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_santi-V1_LE_BDT.weights.xml"))));
676	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_santi-V1_HB_BDT.weights.xml"))));
677	+	muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_santi-V1_HE_BDT.weights.xml"))));
678	+	fMuonTools = new MuonTools();
679	+	fMuonIDMVA = new MuonIDMVA();
680	+	fMuonIDMVA->Initialize("muonHZZ2012IsoDRMVA",
681	+	MuonIDMVA::kIsoDeltaR,
682	+	kTRUE,
683	+	muonidiso_weightfiles);
684	+	}
685	+
686	+	}
687	+
688	+
689	+	//--------------------------------------------------------------------------------------------------
690	+	Bool_t MuonIDMod::PassMuonMVA_BDTG_IdIso(const Muon *mu, const Vertex *vertex,
691	+	const PileupEnergyDensityCol *PileupEnergyDensity) const
692	+	{
693	+
694	+	const Track *muTrk=0;
695	+	if(mu->HasTrackerTrk())         { muTrk = mu->TrackerTrk();    }
696	+	else if(mu->HasStandaloneTrk()) { muTrk = mu->StandaloneTrk(); }
697	+
698	+	Double_t MVAValue = fMuonIDMVA->MVAValue(mu,vertex,fMuonTools,fPFCandidates,PileupEnergyDensity);
699	+
700	+	Int_t subdet = 0;
701	+	if (fabs(muTrk->Eta()) < 1.479) subdet = 0;
702	+	else subdet = 1;
703	+	Int_t ptBin = 0;
704	+	if (muTrk->Pt() > 14.5) ptBin = 1;
705	+	if (muTrk->Pt() > 20.0) ptBin = 2;
706	+
707	+	Int_t MVABin = -1;
708	+	if      (subdet == 0 && ptBin == 0) MVABin = 0;
709	+	else if (subdet == 1 && ptBin == 0) MVABin = 1;
710	+	else if (subdet == 0 && ptBin == 1) MVABin = 2;
711	+	else if (subdet == 1 && ptBin == 1) MVABin = 3;
712	+	else if (subdet == 0 && ptBin == 2) MVABin = 4;
713	+	else if (subdet == 1 && ptBin == 2) MVABin = 5;
714	+
715	+	Double_t MVACut = -999;
716	+	if      (MVABin == 0) MVACut = -0.5618;
717	+	else if (MVABin == 1) MVACut = -0.3002;
718	+	else if (MVABin == 2) MVACut = -0.4642;
719	+	else if (MVABin == 3) MVACut = -0.2478;
720	+	else if (MVABin == 4) MVACut =  0.1706;
721	+	else if (MVABin == 5) MVACut =  0.8146;
722	+
723	+	if (MVAValue > MVACut) return kTRUE;
724	+	return kFALSE;
725	+	}
726	+
727	+	//--------------------------------------------------------------------------------------------------
728	+	Bool_t MuonIDMod::PassMuonIsoRingsV0_BDTG_Iso(const Muon *mu, const Vertex *vertex,
729	+	const PileupEnergyDensityCol *PileupEnergyDensity) const
730	+	{
731	+
732	+	const Track *muTrk=0;
733	+	if(mu->HasTrackerTrk())         { muTrk = mu->TrackerTrk();    }
734	+	else if(mu->HasStandaloneTrk()) { muTrk = mu->StandaloneTrk(); }
735	+
736	+	ElectronOArr *tempElectrons = new  ElectronOArr;
737	+	MuonOArr     *tempMuons     = new  MuonOArr;
738	+	Double_t MVAValue = fMuonIDMVA->MVAValue(mu,vertex,fMuonTools,fPFCandidates,
739	+	PileupEnergyDensity,MuonTools::kMuEAFall11MC,tempElectrons,tempMuons,kFALSE);
740	+	delete tempElectrons;
741	+	delete tempMuons;
742	+
743	+	Int_t MVABin = fMuonIDMVA->GetMVABin(muTrk->Eta(), muTrk->Pt(), mu->IsGlobalMuon(), mu->IsTrackerMuon());
744	+
745	+	Double_t MVACut = -999;
746	+	if      (MVABin == 0) MVACut = -0.593;
747	+	else if (MVABin == 1) MVACut =  0.337;
748	+	else if (MVABin == 2) MVACut = -0.767;
749	+	else if (MVABin == 3) MVACut =  0.410;
750	+	else if (MVABin == 4) MVACut = -0.989;
751	+	else if (MVABin == 5) MVACut = -0.995;
752	+
753	+	if (MVAValue > MVACut) return kTRUE;
754	+	return kFALSE;
755	+	}
756	+
757	+	//--------------------------------------------------------------------------------------------------
758	+	Bool_t MuonIDMod::PassMuonIsoDeltaR(const Muon *mu, const Vertex *vertex,
759	+	const PileupEnergyDensityCol *PileupEnergyDensity) const
760	+	{
761	+
762	+	const Track *muTrk=0;
763	+	if(mu->HasTrackerTrk())         { muTrk = mu->TrackerTrk();    }
764	+	else if(mu->HasStandaloneTrk()) { muTrk = mu->StandaloneTrk(); }
765	+
766	+	ElectronOArr *tempElectrons = new  ElectronOArr;
767	+	MuonOArr     *tempMuons     = new  MuonOArr;
768	+	Double_t MVAValue = fMuonIDMVA->MVAValue(mu,vertex,fMuonTools,fPFNoPileUpCands,
769	+	PileupEnergyDensity,MuonTools::kMuEAFall11MC,tempElectrons,tempMuons,kTRUE);
770	+	delete tempElectrons;
771	+	delete tempMuons;
772	+
773	+	Int_t MVABin = fMuonIDMVA->GetMVABin(muTrk->Eta(), muTrk->Pt(), mu->IsGlobalMuon(), mu->IsTrackerMuon());
774	+
775	+	Double_t MVACut = -999;
776	+	if      (MVABin == 0) MVACut =  0.000;
777	+	else if (MVABin == 1) MVACut =  0.000;
778	+	else if (MVABin == 2) MVACut =  0.000;
779	+	else if (MVABin == 3) MVACut =  0.000;
780	+
781	+	if (MVAValue > MVACut) return kTRUE;
782	+	return kFALSE;
783	+	}
784	+
785	+	//--------------------------------------------------------------------------------------------------
786	+	void MuonIDMod::Terminate()
787	+	{
788	+	// Run finishing code on the computer (slave) that did the analysis
789	+	delete fMuonIDMVA;
790	+
791	+	delete fMuonTools;
792		}

Diff Legend

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