ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/MuonIDMod.cc
Revision: 1.64
Committed: Mon Jan 23 20:27:14 2012 UTC (13 years, 3 months ago) by sixie
Content type: text/plain
Branch: MAIN
Changes since 1.63: +38 -20 lines
Log Message:
change MVA muon ID to use detector based isolation

File Contents

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