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.13 by loizides, Thu Dec 11 15:53:03 2008 UTC vs.
Revision 1.68 by paus, Fri Mar 30 01:08:39 2012 UTC

#	Line 2 | Line 2
2
3		#include "MitPhysics/Mods/interface/MuonIDMod.h"
4		#include "MitCommon/MathTools/interface/MathUtils.h"
5	+	#include "MitAna/DataTree/interface/MuonFwd.h"
6	+	#include "MitAna/DataTree/interface/ElectronFwd.h"
7	+	#include "MitAna/DataTree/interface/VertexCol.h"
8		#include "MitPhysics/Init/interface/ModNames.h"
9
10		using namespace mithep;
#	Line 11 | 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	<	fMuonIDType("Loose"),
21	<	fMuonIsoType("TrackCaloSliding"),
20	>	fNonIsolatedMuonsName("random"),
21	>	fNonIsolatedElectronsName("random"),
22	>	fVertexName(ModNames::gkGoodVertexesName),
23	>	fBeamSpotName(Names::gkBeamSpotBrn),
24	>	fTrackName(Names::gkTrackBrn),
25	>	fPFCandidatesName(Names::gkPFCandidatesBrn),
26	>	fMuonIDType("WWMuIdV3"),
27	>	fMuonIsoType("PFIso"),
28		fMuonClassType("Global"),
29		fTrackIsolationCut(3.0),
30		fCaloIsolationCut(3.0),
31	<	fCombIsolationCut(5.0),
31	>	fCombIsolationCut(0.15),
32	>	fCombRelativeIsolationCut(0.15),
33	>	fPFIsolationCut(-1.0),
34		fMuonPtMin(10),
35	<	fMuons(0)
35	>	fApplyD0Cut(kTRUE),
36	>	fApplyDZCut(kTRUE),
37	>	fD0Cut(0.020),
38	>	fDZCut(0.10),
39	>	fWhichVertex(-1),
40	>	fEtaCut(2.4),
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	>	fIntRadius(0.0),
50	>	fNonIsolatedMuons(0),
51	>	fNonIsolatedElectrons(0),
52	>	fPileupEnergyDensityName(Names::gkPileupEnergyDensityBrn),
53	>	fPileupEnergyDensity(0),
54	>	fMuonTools(0),
55	>	fMuonIDMVA(0),
56	>	fMuonMVAWeights_Subdet0Pt10To14p5(""),
57	>	fMuonMVAWeights_Subdet1Pt10To14p5(""),
58	>	fMuonMVAWeights_Subdet0Pt14p5To20(""),
59	>	fMuonMVAWeights_Subdet1Pt14p5To20(""),
60	>	fMuonMVAWeights_Subdet0Pt20ToInf(""),
61	>	fMuonMVAWeights_Subdet1Pt20ToInf("")
62		{
63		// Constructor.
64		}
#	Line 30 | Line 68 | void MuonIDMod::Process()
68		{
69		// Process entries of the tree.
70
71	<	LoadBranch(fMuonBranchName);
71	>	if(fMuIsoType != kPFIsoNoL) {
72	>	LoadEventObject(fMuonBranchName, fMuons);
73	>	}
74	>	else {
75	>	fMuons = GetObjThisEvt<MuonOArr>(fMuonBranchName);
76	>	}
77	>	LoadEventObject(fBeamSpotName, fBeamSpot);
78	>	LoadEventObject(fTrackName, fTracks);
79	>	LoadEventObject(fPFCandidatesName, fPFCandidates);
80	>	if(fMuIsoType == kTrackCaloSliding ||
81	>	fMuIsoType == kCombinedRelativeConeAreaCorrected ||
82	>	fMuIsoType == kPFIsoEffectiveAreaCorrected ||
83	>	fMuIsoType == kMVAIso_BDTG_IDIso
84	>	) {
85	>	LoadEventObject(fPileupEnergyDensityName, fPileupEnergyDensity);
86	>	}
87
88		MuonOArr *CleanMuons = new MuonOArr;
89		CleanMuons->SetName(fCleanMuonsName);
90
91	+	fVertices = GetObjThisEvt<VertexOArr>(fVertexName);
92	+
93		for (UInt_t i=0; i<fMuons->GetEntries(); ++i) {
94		const Muon *mu = fMuons->At(i);
95
96		Bool_t pass = kFALSE;
97	<	Double_t pt = -1; // make sure pt is taken from the correct track!
97	>	Double_t pt = 0;  // make sure pt is taken from the correct track!
98	>	Double_t eta = 0; // make sure eta is taken from the correct track!
99		switch (fMuClassType) {
100		case kAll:
101		pass = kTRUE;
102	<	pt = mu->Pt();
102	>	if (mu->HasTrk()) {
103	>	pt  = mu->Pt();
104	>	eta = TMath::Abs(mu->Eta());
105	>	}
106		break;
107		case kGlobal:
108	<	pass = (mu->GlobalTrk() != 0);
109	<	if (pass)
110	<	pt = mu->GlobalTrk()->Pt();
111	<	break;
108	>	pass = mu->HasGlobalTrk() && mu->IsTrackerMuon();
109	>	if (pass && mu->TrackerTrk()) {
110	>	pt  = mu->TrackerTrk()->Pt();
111	>	eta = TMath::Abs(mu->TrackerTrk()->Eta());
112	>	}
113	>	else {
114	>	pt  = mu->Pt();
115	>	eta = TMath::Abs(mu->Eta());
116	>	}
117	>	break;
118	>	case kGlobalTracker:
119	>	pass = (mu->HasGlobalTrk() && mu->GlobalTrk()->Chi2()/mu->GlobalTrk()->Ndof() < 10 &&
120	>	(mu->NSegments() > 1 || mu->NMatches() > 1) && mu->NValidHits() > 0) ||
121	>	(mu->IsTrackerMuon() &&
122	>	mu->Quality().Quality(MuonQuality::TMLastStationTight));
123	>	if (pass) {
124	>	pt  = mu->TrackerTrk()->Pt();
125	>	eta = TMath::Abs(mu->TrackerTrk()->Eta());
126	>	}
127	>	else {
128	>	pt  = mu->Pt();
129	>	eta = TMath::Abs(mu->Eta());
130	>	}
131	>	break;
132		case kSta:
133	<	pass = (mu->StandaloneTrk() != 0);
134	<	if (pass)
135	<	pt = mu->StandaloneTrk()->Pt();
136	<	break;
137	<	case kTrackerOnly:
138	<	pass = (mu->TrackerTrk() != 0);
139	<	if (pass)
140	<	pt = mu->TrackerTrk()->Pt();
133	>	pass = mu->HasStandaloneTrk();
134	>	if (pass) {
135	>	pt  = mu->StandaloneTrk()->Pt();
136	>	eta = TMath::Abs(mu->StandaloneTrk()->Eta());
137	>	}
138	>	break;
139	>	case kTrackerMuon:
140	>	pass = mu->HasTrackerTrk() && mu->IsTrackerMuon() &&
141	>	mu->Quality().Quality(MuonQuality::TrackerMuonArbitrated);
142	>	if (pass) {
143	>	pt  = mu->TrackerTrk()->Pt();
144	>	eta = TMath::Abs(mu->TrackerTrk()->Eta());
145	>	}
146	>	break;
147	>	case kCaloMuon:
148	>	pass = mu->HasTrackerTrk() && mu->IsCaloMuon();
149	>	if (pass) {
150	>	pt  = mu->TrackerTrk()->Pt();
151	>	eta = TMath::Abs(mu->TrackerTrk()->Eta());
152	>	}
153	>	break;
154	>	case kTrackerBased:
155	>	pass = mu->HasTrackerTrk();
156	>	if (pass) {
157	>	pt  = mu->TrackerTrk()->Pt();
158	>	eta = TMath::Abs(mu->TrackerTrk()->Eta());
159	>	}
160		break;
161		default:
162		break;
#	Line 70 | Line 168 | void MuonIDMod::Process()
168		if (pt <= fMuonPtMin)
169		continue;
170
171	+	if (eta >= fEtaCut)
172	+	continue;
173	+
174	+
175	+	//***********************************************************************************************
176	+	//Debug Info For Lepton MVA
177	+	//***********************************************************************************************
178	+	if( fPrintMVADebugInfo &&
179	+	(fMuIsoType == kMVAIso_BDTG_IDIso || fMuIDType == kMVAID_BDTG_IDIso)
180	+	) {
181	+	cout << "Event: " << GetEventHeader()->RunNum() << " " << GetEventHeader()->LumiSec() << " "
182	+	<< GetEventHeader()->EvtNum() << " : Rho = " << fPileupEnergyDensity->At(0)->Rho()
183	+	<< " : Muon " << i << " "
184	+	<< endl;
185	+	fMuonIDMVA->MVAValue(mu,fVertices->At(0),fMuonTools,fPFCandidates,fPileupEnergyDensity,kTRUE);
186	+	}
187	+	//***********************************************************************************************
188	+
189	+
190	+	Double_t RChi2 = 0.0;
191	+	if     (mu->HasGlobalTrk()) {
192	+	RChi2 = mu->GlobalTrk()->Chi2()/mu->GlobalTrk()->Ndof();
193	+	}
194	+	else if(mu->BestTrk() != 0){
195	+	RChi2 = mu->BestTrk()->Chi2()/mu->BestTrk()->Ndof();
196	+	}
197		Bool_t idpass = kFALSE;
198		switch (fMuIDType) {
199	+	case kWMuId:
200	+	idpass = mu->BestTrk() != 0 &&
201	+	mu->BestTrk()->NHits() > 10 &&
202	+	RChi2 < 10.0 &&
203	+	(mu->NSegments() > 1 || mu->NMatches() > 1) &&
204	+	mu->BestTrk()->NPixelHits() > 0 &&
205	+	mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
206	+	break;
207	+	case kZMuId:
208	+	idpass = mu->BestTrk() != 0 &&
209	+	mu->BestTrk()->NHits() > 10 &&
210	+	(mu->NSegments() > 1 || mu->NMatches() > 1) &&
211	+	mu->BestTrk()->NPixelHits() > 0 &&
212	+	mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
213	+	break;
214		case kLoose:
215	<	idpass = fMuonTools->IsGood(mu, MuonTools::kTMOneStationLoose) &&
216	<	fMuonTools->IsGood(mu, MuonTools::kTM2DCompatibilityLoose);
215	>	idpass = mu->BestTrk() != 0 &&
216	>	mu->Quality().Quality(MuonQuality::TMOneStationLoose) &&
217	>	mu->Quality().Quality(MuonQuality::TM2DCompatibilityLoose) &&
218	>	mu->BestTrk()->NHits() > 10 &&
219	>	RChi2 < 10.0 &&
220	>	mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
221		break;
222		case kTight:
223	<	idpass = fMuonTools->IsGood(mu, MuonTools::kTMOneStationTight) &&
224	<	fMuonTools->IsGood(mu, MuonTools::kTM2DCompatibilityTight);
223	>	idpass = mu->BestTrk() !=  0 &&
224	>	mu->Quality().Quality(MuonQuality::TMOneStationTight) &&
225	>	mu->Quality().Quality(MuonQuality::TM2DCompatibilityTight) &&
226	>	mu->BestTrk()->NHits() > 10 &&
227	>	RChi2 < 10.0 &&
228	>	mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
229	>	break;
230	>	case kWWMuIdV1:
231	>	idpass = mu->BestTrk() != 0 &&
232	>	mu->BestTrk()->NHits() > 10 &&
233	>	mu->BestTrk()->NPixelHits() > 0 &&
234	>	mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1 &&
235	>	RChi2 < 10.0 &&
236	>	(mu->NSegments() > 1 || mu->NMatches() > 1) &&
237	>	mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
238	>	break;
239	>	case kWWMuIdV2:
240	>	idpass = mu->BestTrk() != 0 &&
241	>	mu->BestTrk()->NHits() > 10 &&
242	>	mu->BestTrk()->NPixelHits() > 0 &&
243	>	mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1;
244	>	break;
245	>	case kWWMuIdV3:
246	>	idpass = mu->BestTrk() != 0 &&
247	>	mu->BestTrk()->NHits() > 10 &&
248	>	mu->BestTrk()->NPixelHits() > 0 &&
249	>	mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1 &&
250	>	mu->TrkKink() < 20.0;
251	>	break;
252	>	case kMVAID_BDTG_IDIso:
253	>	{
254	>	Bool_t passDenominatorM2 = (mu->BestTrk() != 0 &&
255	>	mu->BestTrk()->NHits() > 10 &&
256	>	mu->BestTrk()->NPixelHits() > 0 &&
257	>	mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1 &&
258	>	MuonTools::PassD0Cut(mu, fVertices, 0.20, 0) &&
259	>	MuonTools::PassDZCut(mu, fVertices, 0.10, 0) &&
260	>	mu->TrkKink() < 20.0
261	>	);
262	>	idpass =  passDenominatorM2;
263	>	//only evaluate MVA if muon passes M2 denominator to save time
264	>	if (idpass) idpass = PassMuonMVA_BDTG_IdIso(mu, fVertices->At(0), fPileupEnergyDensity);
265	>	}
266	>	break;
267	>	case kNoId:
268	>	idpass = kTRUE;
269		break;
270		default:
271		break;
#	Line 87 | Line 274 | void MuonIDMod::Process()
274		if (!idpass)
275		continue;
276
277	<	Bool_t isopass = kFALSE;
277	>	Bool_t isocut = kFALSE;
278		switch (fMuIsoType) {
279		case kTrackCalo:
280	<	isopass = (mu->IsoR03SumPt() < fTrackIsolationCut) &&
280	>	isocut = (mu->IsoR03SumPt() < fTrackIsolationCut) &&
281		(mu->IsoR03EmEt() + mu->IsoR03HadEt() < fCaloIsolationCut);
282		break;
283		case kTrackCaloCombined:
284	<	isopass = (1.0 * mu->IsoR03SumPt() + 1.0 * mu->IsoR03EmEt() +
285	<	1.0 * mu->IsoR03HadEt() < fCombIsolationCut);
284	>	isocut = (1.0 * mu->IsoR03SumPt() +
285	>	1.0 * mu->IsoR03EmEt()  +
286	>	1.0 * mu->IsoR03HadEt() < fCombIsolationCut);
287		break;
288		case kTrackCaloSliding:
289		{
290	<	Double_t totalIso = 1.0 * mu->IsoR03SumPt() +
291	<	1.0 * mu->IsoR03EmEt() +
292	<	1.0 * mu->IsoR03HadEt();
293	<	if ((totalIso < (mu->Pt()-10.0)*5.0/15.0 && mu->Pt() <= 25) ||
294	<	(totalIso < 5.0 && mu->Pt() > 25) ||
295	<	totalIso <= 0)
296	<	isopass = kTRUE;
290	>	const PileupEnergyDensity *rho =  fPileupEnergyDensity->At(0);
291	>	Double_t totalIso =  mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt() - rho->Rho() * TMath::Pi() * 0.3 * 0.3 ;
292	>	// trick to change the signal region cut
293	>	double theIsoCut = fCombIsolationCut;
294	>	if(theIsoCut < 0.20){
295	>	if(mu->Pt() >  20.0) theIsoCut = 0.15;
296	>	else                 theIsoCut = 0.10;
297	>	}
298	>	if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
299	>	}
300	>	break;
301	>	case kTrackCaloSlidingNoCorrection:
302	>	{
303	>	Double_t totalIso =  1.0 * mu->IsoR03SumPt() +
304	>	1.0 * mu->IsoR03EmEt()  +
305	>	1.0 * mu->IsoR03HadEt();
306	>	// trick to change the signal region cut
307	>	double theIsoCut = fCombIsolationCut;
308	>	if(theIsoCut < 0.20){
309	>	if(mu->Pt() >  20.0) theIsoCut = 0.15;
310	>	else                 theIsoCut = 0.10;
311	>	}
312	>	if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
313	>	}
314	>	break;
315	>	case kCombinedRelativeConeAreaCorrected:
316	>	{
317	>	const PileupEnergyDensity *rho =  fPileupEnergyDensity->At(0);
318	>	Double_t totalIso =  mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt() - rho->Rho() * TMath::Pi() * 0.3 * 0.3 ;
319	>	double theIsoCut = fCombRelativeIsolationCut;
320	>	if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
321	>	}
322	>	break;
323	>	case kCombinedRelativeEffectiveAreaCorrected:
324	>	{
325	>	Double_t tmpRho = 0;
326	>	if (!(TMath::IsNaN(fPileupEnergyDensity->At(0)->Rho()) || std::isinf(fPileupEnergyDensity->At(0)->Rho())))
327	>	tmpRho = fPileupEnergyDensity->At(0)->Rho();
328	>
329	>	isocut = ( mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt()
330	>	-  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuEMIso03, mu->Eta())
331	>	-  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuHadIso03, mu->Eta())
332	>	) < (mu->Pt()* 0.40);
333		}
334	+	break;
335	+	case kPFIso:
336	+	{
337	+	Double_t pfIsoCutValue = 9999;
338	+	if(fPFIsolationCut > 0){
339	+	pfIsoCutValue = fPFIsolationCut;
340	+	} else {
341	+	if (mu->AbsEta() < 1.479) {
342	+	if (mu->Pt() > 20) {
343	+	pfIsoCutValue = 0.13;
344	+	} else {
345	+	pfIsoCutValue = 0.06;
346	+	}
347	+	} else {
348	+	if (mu->Pt() > 20) {
349	+	pfIsoCutValue = 0.09;
350	+	} else {
351	+	pfIsoCutValue = 0.05;
352	+	}
353	+	}
354	+	}
355	+	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius);
356	+	if (totalIso < (mu->Pt()*pfIsoCutValue) )
357	+	isocut = kTRUE;
358	+	}
359	+	break;
360	+	case kPFIsoEffectiveAreaCorrected:
361	+	{
362	+	Double_t pfIsoCutValue = 9999;
363	+	if(fPFIsolationCut > 0){
364	+	pfIsoCutValue = fPFIsolationCut;
365	+	} else {
366	+	pfIsoCutValue = fPFIsolationCut; //leave it like this for now
367	+	}
368	+	Double_t EffectiveAreaCorrectedPFIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius)
369	+	- fPileupEnergyDensity->At(0)->Rho() * MuonTools::MuonEffectiveArea(MuonTools::kMuNeutralIso03, mu->Eta());
370	+	isocut = EffectiveAreaCorrectedPFIso < (mu->Pt() * pfIsoCutValue);
371	+	break;
372	+	}
373	+	case kPFIsoNoL:
374	+	{
375	+	fNonIsolatedMuons     = GetObjThisEvt<MuonCol>(fNonIsolatedMuonsName);
376	+	fNonIsolatedElectrons = GetObjThisEvt<ElectronCol>(fNonIsolatedElectronsName);
377	+
378	+	Double_t pfIsoCutValue = 9999;
379	+	if(fPFIsolationCut > 0){
380	+	pfIsoCutValue = fPFIsolationCut;
381	+	} else {
382	+	if (mu->AbsEta() < 1.479) {
383	+	if (mu->Pt() > 20) {
384	+	pfIsoCutValue = 0.13;
385	+	} else {
386	+	pfIsoCutValue = 0.06;
387	+	}
388	+	} else {
389	+	if (mu->Pt() > 20) {
390	+	pfIsoCutValue = 0.09;
391	+	} else {
392	+	pfIsoCutValue = 0.05;
393	+	}
394	+	}
395	+	}
396	+	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fNonIsolatedMuons, fNonIsolatedElectrons, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius);
397	+	if (totalIso < (mu->Pt()*pfIsoCutValue) )
398	+	isocut = kTRUE;
399	+	}
400	+	break;
401	+	case kMVAIso_BDTG_IDIso:
402	+	{
403	+
404	+	// **************************************************************************
405	+	// Don't use effective area correction denominator. Instead use the old one.
406	+	// **************************************************************************
407	+
408	+	//         Double_t tmpRho = 0;
409	+	//         if (!(TMath::IsNaN(fPileupEnergyDensity->At(0)->Rho()) || isinf(fPileupEnergyDensity->At(0)->Rho())))
410	+	//           tmpRho = fPileupEnergyDensity->At(0)->Rho();
411	+
412	+	//         isocut = ( mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt()
413	+	//                    -  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuEMIso03, mu->Eta())
414	+	//                    -  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuHadIso03, mu->Eta())
415	+	//           ) < (mu->Pt()* 0.40);
416	+
417	+	Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius);
418	+	isocut = (totalIso < (mu->Pt()*0.4));
419	+
420	+	}
421		break;
422		case kNoIso:
423	<	isopass = kTRUE;
423	>	isocut = kTRUE;
424		break;
425		case kCustomIso:
426		default:
427		break;
428		}
429
430	<	if (!isopass)
430	>	if (isocut == kFALSE)
431		continue;
432
433	+	// apply d0 cut
434	+	if (fApplyD0Cut) {
435	+	Bool_t passD0cut = kTRUE;
436	+	if(fD0Cut < 0.05) { // trick to change the signal region cut
437	+	if      (mu->Pt() >  20.0) fD0Cut = 0.02;
438	+	else if (mu->Pt() <= 20.0) fD0Cut = 0.01;
439	+	}
440	+	if(fWhichVertex >= -1) passD0cut = MuonTools::PassD0Cut(mu, fVertices, fD0Cut, fWhichVertex);
441	+	else                   passD0cut = MuonTools::PassD0Cut(mu, fBeamSpot, fD0Cut);
442	+	if (!passD0cut)
443	+	continue;
444	+	}
445	+
446	+	// apply dz cut
447	+	if (fApplyDZCut) {
448	+	Bool_t passDZcut = MuonTools::PassDZCut(mu, fVertices, fDZCut, fWhichVertex);
449	+	if (!passDZcut)
450	+	continue;
451	+	}
452	+
453		// add good muon
454		CleanMuons->Add(mu);
455		}
#	Line 136 | Line 467 | void MuonIDMod::SlaveBegin()
467		// Run startup code on the computer (slave) doing the actual analysis. Here,
468		// we just request the muon collection branch.
469
470	<	ReqBranch(fMuonBranchName, fMuons);
470	>	// In this case we cannot have a branch
471	>	if (fMuonIsoType.CompareTo("PFIsoNoL") != 0) {
472	>	ReqEventObject(fMuonBranchName, fMuons, kTRUE);
473	>	}
474	>	ReqEventObject(fBeamSpotName, fBeamSpot, kTRUE);
475	>	ReqEventObject(fTrackName, fTracks, kTRUE);
476	>	ReqEventObject(fPFCandidatesName, fPFCandidates, kTRUE);
477	>	if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0
478	>	|| fMuonIsoType.CompareTo("CombinedRelativeConeAreaCorrected") == 0
479	>	|| fMuonIsoType.CompareTo("CombinedRelativeEffectiveAreaCorrected") == 0
480	>	|| fMuonIsoType.CompareTo("PFIsoEffectiveAreaCorrected") == 0
481	>	|| fMuonIsoType.CompareTo("MVA_BDTG_IDIso") == 0
482	>	) {
483	>	ReqEventObject(fPileupEnergyDensityName, fPileupEnergyDensity, kTRUE);
484	>	}
485
141	–	fMuonTools = new MuonTools;
486
487	<	if (fMuonIDType.CompareTo("Tight") == 0)
487	>	if (fMuonIDType.CompareTo("WMuId") == 0)
488	>	fMuIDType = kWMuId;
489	>	else if (fMuonIDType.CompareTo("ZMuId") == 0)
490	>	fMuIDType = kZMuId;
491	>	else if (fMuonIDType.CompareTo("Tight") == 0)
492		fMuIDType = kTight;
493		else if (fMuonIDType.CompareTo("Loose") == 0)
494		fMuIDType = kLoose;
495	+	else if (fMuonIDType.CompareTo("WWMuIdV1") == 0)
496	+	fMuIDType = kWWMuIdV1;
497	+	else if (fMuonIDType.CompareTo("WWMuIdV2") == 0)
498	+	fMuIDType = kWWMuIdV2;
499	+	else if (fMuonIDType.CompareTo("WWMuIdV3") == 0)
500	+	fMuIDType = kWWMuIdV3;
501		else if (fMuonIDType.CompareTo("NoId") == 0)
502		fMuIDType = kNoId;
503		else if (fMuonIDType.CompareTo("Custom") == 0) {
504		fMuIDType = kCustomId;
505		SendError(kWarning, "SlaveBegin",
506		"Custom muon identification is not yet implemented.");
507	+	} else if (fMuonIDType.CompareTo("MVA_BDTG_IDIso") == 0) {
508	+	fMuIDType = kMVAID_BDTG_IDIso;
509		} else {
510		SendError(kAbortAnalysis, "SlaveBegin",
511		"The specified muon identification %s is not defined.",
#	Line 163 | Line 519 | void MuonIDMod::SlaveBegin()
519		fMuIsoType = kTrackCaloCombined;
520		else if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0)
521		fMuIsoType = kTrackCaloSliding;
522	+	else if (fMuonIsoType.CompareTo("TrackCaloSlidingNoCorrection") == 0)
523	+	fMuIsoType = kTrackCaloSlidingNoCorrection;
524	+	else if (fMuonIsoType.CompareTo("CombinedRelativeConeAreaCorrected") == 0)
525	+	fMuIsoType = kCombinedRelativeConeAreaCorrected;
526	+	else if (fMuonIsoType.CompareTo("CombinedRelativeEffectiveAreaCorrected") == 0)
527	+	fMuIsoType = kCombinedRelativeEffectiveAreaCorrected;
528	+	else if (fMuonIsoType.CompareTo("PFIso") == 0)
529	+	fMuIsoType = kPFIso;
530	+	else if (fMuonIsoType.CompareTo("PFIsoEffectiveAreaCorrected") == 0)
531	+	fMuIsoType = kPFIsoEffectiveAreaCorrected;
532	+	else if (fMuonIsoType.CompareTo("PFIsoNoL") == 0)
533	+	fMuIsoType = kPFIsoNoL;
534		else if (fMuonIsoType.CompareTo("NoIso") == 0)
535		fMuIsoType = kNoIso;
536		else if (fMuonIsoType.CompareTo("Custom") == 0) {
537		fMuIsoType = kCustomIso;
538		SendError(kWarning, "SlaveBegin",
539		"Custom muon isolation is not yet implemented.");
540	+	} else if (fMuonIDType.CompareTo("MVA_BDTG_IDIso") == 0) {
541	+	fMuIsoType = kMVAIso_BDTG_IDIso;
542		} else {
543		SendError(kAbortAnalysis, "SlaveBegin",
544		"The specified muon isolation %s is not defined.",
#	Line 180 | Line 550 | void MuonIDMod::SlaveBegin()
550		fMuClassType = kAll;
551		else if (fMuonClassType.CompareTo("Global") == 0)
552		fMuClassType = kGlobal;
553	+	else if (fMuonClassType.CompareTo("GlobalTracker") == 0)
554	+	fMuClassType = kGlobalTracker;
555		else if (fMuonClassType.CompareTo("Standalone") == 0)
556		fMuClassType = kSta;
557	<	else if (fMuonClassType.CompareTo("TrackerOnly") == 0)
558	<	fMuClassType = kTrackerOnly;
557	>	else if (fMuonClassType.CompareTo("TrackerMuon") == 0)
558	>	fMuClassType = kTrackerMuon;
559	>	else if (fMuonClassType.CompareTo("CaloMuon") == 0)
560	>	fMuClassType = kCaloMuon;
561	>	else if (fMuonClassType.CompareTo("TrackerBased") == 0)
562	>	fMuClassType = kTrackerBased;
563		else {
564		SendError(kAbortAnalysis, "SlaveBegin",
565		"The specified muon class %s is not defined.",
566		fMuonClassType.Data());
567		return;
568		}
569	+
570	+
571	+	//If we use MVA ID, need to load MVA weights
572	+	if(fMuIsoType == kMVAIso_BDTG_IDIso || fMuIDType == kMVAID_BDTG_IDIso) {
573	+	fMuonTools = new MuonTools();
574	+	fMuonIDMVA = new MuonIDMVA();
575	+	fMuonIDMVA->Initialize("BDTG method",
576	+	fMuonMVAWeights_Subdet0Pt10To14p5,
577	+	fMuonMVAWeights_Subdet1Pt10To14p5,
578	+	fMuonMVAWeights_Subdet0Pt14p5To20,
579	+	fMuonMVAWeights_Subdet1Pt14p5To20,
580	+	fMuonMVAWeights_Subdet0Pt20ToInf,
581	+	fMuonMVAWeights_Subdet1Pt20ToInf,
582	+	MuonIDMVA::kIDIsoCombinedDetIso);
583	+	}
584	+
585	+	}
586	+
587	+
588	+	//--------------------------------------------------------------------------------------------------
589	+	Bool_t MuonIDMod::PassMuonMVA_BDTG_IdIso(const Muon *mu, const Vertex *vertex,
590	+	const PileupEnergyDensityCol *PileupEnergyDensity) const
591	+	{
592	+
593	+	const Track *muTrk=0;
594	+	if(mu->HasTrackerTrk())         { muTrk = mu->TrackerTrk();    }
595	+	else if(mu->HasStandaloneTrk()) { muTrk = mu->StandaloneTrk(); }
596	+
597	+	Double_t MVAValue = fMuonIDMVA->MVAValue(mu,vertex,fMuonTools,fPFCandidates,PileupEnergyDensity);
598	+
599	+	Int_t subdet = 0;
600	+	if (fabs(muTrk->Eta()) < 1.479) subdet = 0;
601	+	else subdet = 1;
602	+	Int_t ptBin = 0;
603	+	if (muTrk->Pt() > 14.5) ptBin = 1;
604	+	if (muTrk->Pt() > 20.0) ptBin = 2;
605	+
606	+	Int_t MVABin = -1;
607	+	if      (subdet == 0 && ptBin == 0) MVABin = 0;
608	+	else if (subdet == 1 && ptBin == 0) MVABin = 1;
609	+	else if (subdet == 0 && ptBin == 1) MVABin = 2;
610	+	else if (subdet == 1 && ptBin == 1) MVABin = 3;
611	+	else if (subdet == 0 && ptBin == 2) MVABin = 4;
612	+	else if (subdet == 1 && ptBin == 2) MVABin = 5;
613	+
614	+	Double_t MVACut = -999;
615	+	if      (MVABin == 0) MVACut = -0.5618;
616	+	else if (MVABin == 1) MVACut = -0.3002;
617	+	else if (MVABin == 2) MVACut = -0.4642;
618	+	else if (MVABin == 3) MVACut = -0.2478;
619	+	else if (MVABin == 4) MVACut =  0.1706;
620	+	else if (MVABin == 5) MVACut =  0.8146;
621	+
622	+	if (MVAValue > MVACut) return kTRUE;
623	+	return kFALSE;
624		}

Diff Legend

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