ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/MuonIDMod.cc
Revision: 1.51
Committed: Mon May 16 18:21:42 2011 UTC (13 years, 11 months ago) by ceballos
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_021
Changes since 1.50: +13 -5 lines
Log Message: 
new id

File Contents

# User Rev Content
1 ceballos 1.51 // $Id: MuonIDMod.cc,v 1.50 2011/05/13 08:27:15 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.36 fMuonIDType("WWMuId"),
26 ceballos 1.49 fMuonIsoType("PFIso"),
27 loizides 1.6 fMuonClassType("Global"),
28 ceballos 1.2 fTrackIsolationCut(3.0),
29     fCaloIsolationCut(3.0),
30 ceballos 1.49 fCombIsolationCut(-1.0),
31 ceballos 1.3 fMuonPtMin(10),
32 loizides 1.21 fApplyD0Cut(kTRUE),
33 ceballos 1.42 fApplyDZCut(kTRUE),
34 ceballos 1.30 fD0Cut(0.020),
35 ceballos 1.50 fDZCut(0.10),
36 ceballos 1.42 fWhichVertex(-1),
37 ceballos 1.30 fEtaCut(2.4),
38 loizides 1.15 fMuIDType(kIdUndef),
39     fMuIsoType(kIsoUndef),
40     fMuClassType(kClassUndef),
41 ceballos 1.14 fMuons(0),
42 loizides 1.15 fVertices(0),
43 ceballos 1.42 fBeamSpot(0),
44 ceballos 1.38 fTracks(0),
45     fPFCandidates(0),
46 ceballos 1.39 fNonIsolatedMuons(0),
47 ceballos 1.41 fNonIsolatedElectrons(0),
48     fPileupEnergyDensityName(Names::gkPileupEnergyDensityBrn),
49     fPileupEnergyDensity(0)
50 loizides 1.1 {
51     // Constructor.
52     }
53    
54     //--------------------------------------------------------------------------------------------------
55     void MuonIDMod::Process()
56     {
57     // Process entries of the tree.
58    
59 ceballos 1.38 if(fMuIsoType != kPFIsoNoL) {
60     LoadEventObject(fMuonBranchName, fMuons);
61     }
62     else {
63     fMuons = GetObjThisEvt<MuonOArr>(fMuonBranchName);
64     }
65 ceballos 1.42 LoadEventObject(fBeamSpotName, fBeamSpot);
66 ceballos 1.38 LoadEventObject(fTrackName, fTracks);
67     LoadEventObject(fPFCandidatesName, fPFCandidates);
68 ceballos 1.41 if(fMuIsoType == kTrackCaloSliding) {
69     LoadEventObject(fPileupEnergyDensityName, fPileupEnergyDensity);
70     }
71 loizides 1.7 MuonOArr *CleanMuons = new MuonOArr;
72     CleanMuons->SetName(fCleanMuonsName);
73 loizides 1.1
74 ceballos 1.38 fVertices = GetObjThisEvt<VertexOArr>(fVertexName);
75    
76 loizides 1.1 for (UInt_t i=0; i<fMuons->GetEntries(); ++i) {
77 loizides 1.8 const Muon *mu = fMuons->At(i);
78 loizides 1.6
79     Bool_t pass = kFALSE;
80 ceballos 1.30 Double_t pt = 0; // make sure pt is taken from the correct track!
81     Double_t eta = 0; // make sure eta is taken from the correct track!
82 loizides 1.6 switch (fMuClassType) {
83     case kAll:
84     pass = kTRUE;
85 ceballos 1.30 if (mu->HasTrk()) {
86     pt = mu->Pt();
87     eta = TMath::Abs(mu->Eta());
88     }
89 loizides 1.6 break;
90     case kGlobal:
91 ceballos 1.31 pass = mu->HasGlobalTrk() && mu->IsTrackerMuon();
92 ceballos 1.32 if (pass && mu->TrackerTrk()) {
93 ceballos 1.30 pt = mu->TrackerTrk()->Pt();
94     eta = TMath::Abs(mu->TrackerTrk()->Eta());
95     }
96 ceballos 1.32 else {
97     pt = mu->Pt();
98     eta = TMath::Abs(mu->Eta());
99     }
100 ceballos 1.30 break;
101 loizides 1.6 case kSta:
102 loizides 1.19 pass = mu->HasStandaloneTrk();
103 ceballos 1.30 if (pass) {
104     pt = mu->StandaloneTrk()->Pt();
105     eta = TMath::Abs(mu->StandaloneTrk()->Eta());
106     }
107 loizides 1.6 break;
108 loizides 1.24 case kTrackerMuon:
109 ceballos 1.28 pass = mu->HasTrackerTrk() && mu->IsTrackerMuon() &&
110     mu->Quality().Quality(MuonQuality::TrackerMuonArbitrated);
111 ceballos 1.30 if (pass) {
112     pt = mu->TrackerTrk()->Pt();
113     eta = TMath::Abs(mu->TrackerTrk()->Eta());
114     }
115 loizides 1.24 break;
116     case kCaloMuon:
117     pass = mu->HasTrackerTrk() && mu->IsCaloMuon();
118 ceballos 1.30 if (pass) {
119     pt = mu->TrackerTrk()->Pt();
120     eta = TMath::Abs(mu->TrackerTrk()->Eta());
121     }
122 loizides 1.24 break;
123     case kTrackerBased:
124 loizides 1.19 pass = mu->HasTrackerTrk();
125 ceballos 1.30 if (pass) {
126     pt = mu->TrackerTrk()->Pt();
127     eta = TMath::Abs(mu->TrackerTrk()->Eta());
128     }
129 loizides 1.6 break;
130     default:
131     break;
132 ceballos 1.5 }
133 loizides 1.6
134     if (!pass)
135     continue;
136    
137     if (pt <= fMuonPtMin)
138     continue;
139    
140 ceballos 1.30 if (eta >= fEtaCut)
141     continue;
142    
143 ceballos 1.34 Double_t RChi2 = 0.0;
144     if (mu->HasGlobalTrk()) {
145     RChi2 = mu->GlobalTrk()->Chi2()/mu->GlobalTrk()->Ndof();
146     }
147     else if(mu->BestTrk() != 0){
148     RChi2 = mu->BestTrk()->Chi2()/mu->BestTrk()->Ndof();
149     }
150 loizides 1.6 Bool_t idpass = kFALSE;
151     switch (fMuIDType) {
152 ceballos 1.34 case kWMuId:
153     idpass = mu->BestTrk() != 0 &&
154     mu->BestTrk()->NHits() > 10 &&
155     RChi2 < 10.0 &&
156 ceballos 1.36 (mu->NSegments() > 1 || mu->NMatches() > 1) &&
157 ceballos 1.34 mu->BestTrk()->NPixelHits() > 0 &&
158     mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
159     break;
160     case kZMuId:
161     idpass = mu->BestTrk() != 0 &&
162     mu->BestTrk()->NHits() > 10 &&
163 ceballos 1.36 (mu->NSegments() > 1 || mu->NMatches() > 1) &&
164 ceballos 1.34 mu->BestTrk()->NPixelHits() > 0 &&
165     mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
166     break;
167 loizides 1.6 case kLoose:
168 ceballos 1.32 idpass = mu->BestTrk() != 0 &&
169     mu->Quality().Quality(MuonQuality::TMOneStationLoose) &&
170 ceballos 1.28 mu->Quality().Quality(MuonQuality::TM2DCompatibilityLoose) &&
171 ceballos 1.29 mu->BestTrk()->NHits() > 10 &&
172 ceballos 1.34 RChi2 < 10.0 &&
173 ceballos 1.30 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
174 loizides 1.6 break;
175     case kTight:
176 ceballos 1.32 idpass = mu->BestTrk() != 0 &&
177     mu->Quality().Quality(MuonQuality::TMOneStationTight) &&
178 ceballos 1.28 mu->Quality().Quality(MuonQuality::TM2DCompatibilityTight) &&
179 ceballos 1.29 mu->BestTrk()->NHits() > 10 &&
180 ceballos 1.34 RChi2 < 10.0 &&
181 ceballos 1.30 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
182 ceballos 1.28 break;
183 ceballos 1.36 case kWWMuId:
184 ceballos 1.32 idpass = mu->BestTrk() != 0 &&
185     mu->BestTrk()->NHits() > 10 &&
186 ceballos 1.34 RChi2 < 10.0 &&
187 ceballos 1.36 (mu->NSegments() > 1 || mu->NMatches() > 1) &&
188     mu->BestTrk()->NPixelHits() > 0 &&
189     mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight) &&
190     mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1;
191 loizides 1.6 break;
192 loizides 1.18 case kNoId:
193     idpass = kTRUE;
194     break;
195 loizides 1.6 default:
196     break;
197 ceballos 1.4 }
198 loizides 1.6
199     if (!idpass)
200     continue;
201    
202 ceballos 1.22 Bool_t isocut = kFALSE;
203 loizides 1.6 switch (fMuIsoType) {
204     case kTrackCalo:
205 ceballos 1.22 isocut = (mu->IsoR03SumPt() < fTrackIsolationCut) &&
206 loizides 1.6 (mu->IsoR03EmEt() + mu->IsoR03HadEt() < fCaloIsolationCut);
207     break;
208     case kTrackCaloCombined:
209 ceballos 1.38 isocut = (1.0 * mu->IsoR03SumPt() +
210     1.0 * mu->IsoR03EmEt() +
211     1.0 * mu->IsoR03HadEt() < fCombIsolationCut);
212 loizides 1.6 break;
213     case kTrackCaloSliding:
214     {
215 ceballos 1.41 const PileupEnergyDensity *rho = fPileupEnergyDensity->At(0);
216     Double_t totalIso = mu->IsoR03SumPt() + TMath::Max(mu->IsoR03EmEt() + mu->IsoR03HadEt() - rho->Rho() * TMath::Pi() * 0.3 * 0.3, 0.0);
217 ceballos 1.46 // trick to change the signal region cut
218     double theIsoCut = fCombIsolationCut;
219     if(theIsoCut < 0.20){
220     if(mu->Pt() > 20.0) theIsoCut = 0.15;
221     else theIsoCut = 0.10;
222     }
223     if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
224 ceballos 1.22 }
225 loizides 1.6 break;
226 ceballos 1.41 case kTrackCaloSlidingNoCorrection:
227 ceballos 1.40 {
228     Double_t totalIso = 1.0 * mu->IsoR03SumPt() +
229 ceballos 1.41 1.0 * mu->IsoR03EmEt() +
230     1.0 * mu->IsoR03HadEt();
231 ceballos 1.46 // trick to change the signal region cut
232     double theIsoCut = fCombIsolationCut;
233     if(theIsoCut < 0.20){
234     if(mu->Pt() > 20.0) theIsoCut = 0.15;
235     else theIsoCut = 0.10;
236     }
237     if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
238 ceballos 1.40 }
239     break;
240 ceballos 1.38 case kPFIso:
241     {
242 ceballos 1.49 Double_t pfIsoCutValue = 9999;
243     if(fCombIsolationCut > 0){
244     pfIsoCutValue = fCombIsolationCut;
245     } else {
246 ceballos 1.51 if (mu->AbsEta() < 1.479) {
247     if (mu->Pt() > 20) {
248     pfIsoCutValue = 0.22;
249     } else {
250     pfIsoCutValue = 0.11;
251     }
252 ceballos 1.49 } else {
253 ceballos 1.51 if (mu->Pt() > 20) {
254     pfIsoCutValue = 0.20;
255     } else {
256     pfIsoCutValue = 0.11;
257     }
258     }
259 ceballos 1.49 }
260 ceballos 1.50 Double_t totalIso = IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.4, 0.0);
261 ceballos 1.49 if (totalIso < (mu->Pt()*pfIsoCutValue) )
262 ceballos 1.38 isocut = kTRUE;
263     }
264     break;
265     case kPFIsoNoL:
266     {
267 ceballos 1.39 fNonIsolatedMuons = GetObjThisEvt<MuonCol>(fNonIsolatedMuonsName);
268     fNonIsolatedElectrons = GetObjThisEvt<ElectronCol>(fNonIsolatedElectronsName);
269 ceballos 1.38
270     Double_t beta = IsolationTools::BetaM(fTracks, mu, fVertices->At(0), 0.0, 0.2, 0.3, 0.02);
271 ceballos 1.40 if(beta == 0) beta = 1.0;
272 mzanetti 1.47 Double_t totalIso = IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), fNonIsolatedMuons, fNonIsolatedElectrons, 0.2, 1.0, 0.4, 0.0, 3);
273 ceballos 1.38 if (totalIso < (mu->Pt()*fCombIsolationCut) )
274     isocut = kTRUE;
275     }
276     break;
277 loizides 1.13 case kNoIso:
278 ceballos 1.22 isocut = kTRUE;
279 loizides 1.13 break;
280 loizides 1.6 case kCustomIso:
281     default:
282     break;
283 ceballos 1.4 }
284 ceballos 1.3
285 ceballos 1.22 if (isocut == kFALSE)
286 loizides 1.6 continue;
287    
288 ceballos 1.42 // apply d0 cut
289 loizides 1.21 if (fApplyD0Cut) {
290 ceballos 1.42 Bool_t passD0cut = kTRUE;
291 ceballos 1.46 if(fD0Cut < 0.05) { // trick to change the signal region cut
292     if (mu->Pt() > 20.0) fD0Cut = 0.02;
293     else if (mu->Pt() <= 20.0) fD0Cut = 0.01;
294     }
295 ceballos 1.42 if(fWhichVertex >= -1) passD0cut = MuonTools::PassD0Cut(mu, fVertices, fD0Cut, fWhichVertex);
296     else passD0cut = MuonTools::PassD0Cut(mu, fBeamSpot, fD0Cut);
297 ceballos 1.31 if (!passD0cut)
298 loizides 1.21 continue;
299 ceballos 1.14 }
300    
301 ceballos 1.42 // apply dz cut
302     if (fApplyDZCut) {
303 ceballos 1.45 Bool_t passDZcut = MuonTools::PassDZCut(mu, fVertices, fDZCut, fWhichVertex);
304 ceballos 1.42 if (!passDZcut)
305     continue;
306     }
307    
308 loizides 1.6 // add good muon
309     CleanMuons->Add(mu);
310 loizides 1.1 }
311    
312 loizides 1.10 // sort according to pt
313     CleanMuons->Sort();
314    
315 loizides 1.6 // add objects for other modules to use
316 loizides 1.7 AddObjThisEvt(CleanMuons);
317 loizides 1.1 }
318    
319     //--------------------------------------------------------------------------------------------------
320     void MuonIDMod::SlaveBegin()
321     {
322     // Run startup code on the computer (slave) doing the actual analysis. Here,
323 loizides 1.6 // we just request the muon collection branch.
324    
325 ceballos 1.38 // In this case we cannot have a branch
326     if (fMuonIsoType.CompareTo("PFIsoNoL") != 0) {
327     ReqEventObject(fMuonBranchName, fMuons, kTRUE);
328     }
329 ceballos 1.42 ReqEventObject(fBeamSpotName, fBeamSpot, kTRUE);
330 ceballos 1.38 ReqEventObject(fTrackName, fTracks, kTRUE);
331     ReqEventObject(fPFCandidatesName, fPFCandidates, kTRUE);
332 ceballos 1.41 if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0) {
333     ReqEventObject(fPileupEnergyDensityName, fPileupEnergyDensity, kTRUE);
334     }
335 loizides 1.1
336 ceballos 1.34 if (fMuonIDType.CompareTo("WMuId") == 0)
337     fMuIDType = kWMuId;
338     else if (fMuonIDType.CompareTo("ZMuId") == 0)
339     fMuIDType = kZMuId;
340     else if (fMuonIDType.CompareTo("Tight") == 0)
341 loizides 1.6 fMuIDType = kTight;
342     else if (fMuonIDType.CompareTo("Loose") == 0)
343     fMuIDType = kLoose;
344 ceballos 1.36 else if (fMuonIDType.CompareTo("WWMuId") == 0)
345     fMuIDType = kWWMuId;
346 loizides 1.12 else if (fMuonIDType.CompareTo("NoId") == 0)
347     fMuIDType = kNoId;
348 loizides 1.6 else if (fMuonIDType.CompareTo("Custom") == 0) {
349     fMuIDType = kCustomId;
350     SendError(kWarning, "SlaveBegin",
351     "Custom muon identification is not yet implemented.");
352     } else {
353     SendError(kAbortAnalysis, "SlaveBegin",
354     "The specified muon identification %s is not defined.",
355     fMuonIDType.Data());
356     return;
357     }
358 loizides 1.1
359 loizides 1.6 if (fMuonIsoType.CompareTo("TrackCalo") == 0)
360     fMuIsoType = kTrackCalo;
361     else if (fMuonIsoType.CompareTo("TrackCaloCombined") == 0)
362     fMuIsoType = kTrackCaloCombined;
363     else if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0)
364     fMuIsoType = kTrackCaloSliding;
365 ceballos 1.41 else if (fMuonIsoType.CompareTo("TrackCaloSlidingNoCorrection") == 0)
366     fMuIsoType = kTrackCaloSlidingNoCorrection;
367 ceballos 1.38 else if (fMuonIsoType.CompareTo("PFIso") == 0)
368     fMuIsoType = kPFIso;
369     else if (fMuonIsoType.CompareTo("PFIsoNoL") == 0)
370     fMuIsoType = kPFIsoNoL;
371 loizides 1.6 else if (fMuonIsoType.CompareTo("NoIso") == 0)
372     fMuIsoType = kNoIso;
373     else if (fMuonIsoType.CompareTo("Custom") == 0) {
374     fMuIsoType = kCustomIso;
375     SendError(kWarning, "SlaveBegin",
376     "Custom muon isolation is not yet implemented.");
377     } else {
378     SendError(kAbortAnalysis, "SlaveBegin",
379     "The specified muon isolation %s is not defined.",
380     fMuonIsoType.Data());
381     return;
382     }
383 loizides 1.1
384 loizides 1.6 if (fMuonClassType.CompareTo("All") == 0)
385     fMuClassType = kAll;
386     else if (fMuonClassType.CompareTo("Global") == 0)
387     fMuClassType = kGlobal;
388     else if (fMuonClassType.CompareTo("Standalone") == 0)
389     fMuClassType = kSta;
390 loizides 1.24 else if (fMuonClassType.CompareTo("TrackerMuon") == 0)
391     fMuClassType = kTrackerMuon;
392     else if (fMuonClassType.CompareTo("CaloMuon") == 0)
393     fMuClassType = kCaloMuon;
394     else if (fMuonClassType.CompareTo("TrackerBased") == 0)
395     fMuClassType = kTrackerBased;
396 loizides 1.6 else {
397     SendError(kAbortAnalysis, "SlaveBegin",
398     "The specified muon class %s is not defined.",
399     fMuonClassType.Data());
400     return;
401     }
402 loizides 1.1 }