ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/MuonIDMod.cc
Revision: 1.49
Committed: Thu May 12 23:15:48 2011 UTC (13 years, 11 months ago) by ceballos
Content type: text/plain
Branch: MAIN
Changes since 1.48: +15 -7 lines
Log Message: 
small updates

File Contents

# User Rev Content
1 ceballos 1.49 // $Id: MuonIDMod.cc,v 1.48 2011/05/12 21:31:42 sixie 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.42 fDZCut(0.20),
36     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     if (mu->Pt() > 20) {
247     pfIsoCutValue = 0.17;
248     } else {
249     pfIsoCutValue = 0.13;
250     }
251     }
252     Double_t totalIso = IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0),0.2, 1.0, 0.4, 0.0);
253     if (totalIso < (mu->Pt()*pfIsoCutValue) )
254 ceballos 1.38 isocut = kTRUE;
255     }
256     break;
257     case kPFIsoNoL:
258     {
259 ceballos 1.39 fNonIsolatedMuons = GetObjThisEvt<MuonCol>(fNonIsolatedMuonsName);
260     fNonIsolatedElectrons = GetObjThisEvt<ElectronCol>(fNonIsolatedElectronsName);
261 ceballos 1.38
262     Double_t beta = IsolationTools::BetaM(fTracks, mu, fVertices->At(0), 0.0, 0.2, 0.3, 0.02);
263 ceballos 1.40 if(beta == 0) beta = 1.0;
264 mzanetti 1.47 Double_t totalIso = IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), fNonIsolatedMuons, fNonIsolatedElectrons, 0.2, 1.0, 0.4, 0.0, 3);
265 ceballos 1.38 if (totalIso < (mu->Pt()*fCombIsolationCut) )
266     isocut = kTRUE;
267     }
268     break;
269 loizides 1.13 case kNoIso:
270 ceballos 1.22 isocut = kTRUE;
271 loizides 1.13 break;
272 loizides 1.6 case kCustomIso:
273     default:
274     break;
275 ceballos 1.4 }
276 ceballos 1.3
277 ceballos 1.22 if (isocut == kFALSE)
278 loizides 1.6 continue;
279    
280 ceballos 1.42 // apply d0 cut
281 loizides 1.21 if (fApplyD0Cut) {
282 ceballos 1.42 Bool_t passD0cut = kTRUE;
283 ceballos 1.46 if(fD0Cut < 0.05) { // trick to change the signal region cut
284     if (mu->Pt() > 20.0) fD0Cut = 0.02;
285     else if (mu->Pt() <= 20.0) fD0Cut = 0.01;
286     }
287 ceballos 1.42 if(fWhichVertex >= -1) passD0cut = MuonTools::PassD0Cut(mu, fVertices, fD0Cut, fWhichVertex);
288     else passD0cut = MuonTools::PassD0Cut(mu, fBeamSpot, fD0Cut);
289 ceballos 1.31 if (!passD0cut)
290 loizides 1.21 continue;
291 ceballos 1.14 }
292    
293 ceballos 1.42 // apply dz cut
294     if (fApplyDZCut) {
295 ceballos 1.45 Bool_t passDZcut = MuonTools::PassDZCut(mu, fVertices, fDZCut, fWhichVertex);
296 ceballos 1.42 if (!passDZcut)
297     continue;
298     }
299    
300 loizides 1.6 // add good muon
301     CleanMuons->Add(mu);
302 loizides 1.1 }
303    
304 loizides 1.10 // sort according to pt
305     CleanMuons->Sort();
306    
307 loizides 1.6 // add objects for other modules to use
308 loizides 1.7 AddObjThisEvt(CleanMuons);
309 loizides 1.1 }
310    
311     //--------------------------------------------------------------------------------------------------
312     void MuonIDMod::SlaveBegin()
313     {
314     // Run startup code on the computer (slave) doing the actual analysis. Here,
315 loizides 1.6 // we just request the muon collection branch.
316    
317 ceballos 1.38 // In this case we cannot have a branch
318     if (fMuonIsoType.CompareTo("PFIsoNoL") != 0) {
319     ReqEventObject(fMuonBranchName, fMuons, kTRUE);
320     }
321 ceballos 1.42 ReqEventObject(fBeamSpotName, fBeamSpot, kTRUE);
322 ceballos 1.38 ReqEventObject(fTrackName, fTracks, kTRUE);
323     ReqEventObject(fPFCandidatesName, fPFCandidates, kTRUE);
324 ceballos 1.41 if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0) {
325     ReqEventObject(fPileupEnergyDensityName, fPileupEnergyDensity, kTRUE);
326     }
327 loizides 1.1
328 ceballos 1.34 if (fMuonIDType.CompareTo("WMuId") == 0)
329     fMuIDType = kWMuId;
330     else if (fMuonIDType.CompareTo("ZMuId") == 0)
331     fMuIDType = kZMuId;
332     else if (fMuonIDType.CompareTo("Tight") == 0)
333 loizides 1.6 fMuIDType = kTight;
334     else if (fMuonIDType.CompareTo("Loose") == 0)
335     fMuIDType = kLoose;
336 ceballos 1.36 else if (fMuonIDType.CompareTo("WWMuId") == 0)
337     fMuIDType = kWWMuId;
338 loizides 1.12 else if (fMuonIDType.CompareTo("NoId") == 0)
339     fMuIDType = kNoId;
340 loizides 1.6 else if (fMuonIDType.CompareTo("Custom") == 0) {
341     fMuIDType = kCustomId;
342     SendError(kWarning, "SlaveBegin",
343     "Custom muon identification is not yet implemented.");
344     } else {
345     SendError(kAbortAnalysis, "SlaveBegin",
346     "The specified muon identification %s is not defined.",
347     fMuonIDType.Data());
348     return;
349     }
350 loizides 1.1
351 loizides 1.6 if (fMuonIsoType.CompareTo("TrackCalo") == 0)
352     fMuIsoType = kTrackCalo;
353     else if (fMuonIsoType.CompareTo("TrackCaloCombined") == 0)
354     fMuIsoType = kTrackCaloCombined;
355     else if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0)
356     fMuIsoType = kTrackCaloSliding;
357 ceballos 1.41 else if (fMuonIsoType.CompareTo("TrackCaloSlidingNoCorrection") == 0)
358     fMuIsoType = kTrackCaloSlidingNoCorrection;
359 ceballos 1.38 else if (fMuonIsoType.CompareTo("PFIso") == 0)
360     fMuIsoType = kPFIso;
361     else if (fMuonIsoType.CompareTo("PFIsoNoL") == 0)
362     fMuIsoType = kPFIsoNoL;
363 loizides 1.6 else if (fMuonIsoType.CompareTo("NoIso") == 0)
364     fMuIsoType = kNoIso;
365     else if (fMuonIsoType.CompareTo("Custom") == 0) {
366     fMuIsoType = kCustomIso;
367     SendError(kWarning, "SlaveBegin",
368     "Custom muon isolation is not yet implemented.");
369     } else {
370     SendError(kAbortAnalysis, "SlaveBegin",
371     "The specified muon isolation %s is not defined.",
372     fMuonIsoType.Data());
373     return;
374     }
375 loizides 1.1
376 loizides 1.6 if (fMuonClassType.CompareTo("All") == 0)
377     fMuClassType = kAll;
378     else if (fMuonClassType.CompareTo("Global") == 0)
379     fMuClassType = kGlobal;
380     else if (fMuonClassType.CompareTo("Standalone") == 0)
381     fMuClassType = kSta;
382 loizides 1.24 else if (fMuonClassType.CompareTo("TrackerMuon") == 0)
383     fMuClassType = kTrackerMuon;
384     else if (fMuonClassType.CompareTo("CaloMuon") == 0)
385     fMuClassType = kCaloMuon;
386     else if (fMuonClassType.CompareTo("TrackerBased") == 0)
387     fMuClassType = kTrackerBased;
388 loizides 1.6 else {
389     SendError(kAbortAnalysis, "SlaveBegin",
390     "The specified muon class %s is not defined.",
391     fMuonClassType.Data());
392     return;
393     }
394 loizides 1.1 }