ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/MuonIDMod.cc
Revision: 1.53
Committed: Fri Jun 10 10:42:29 2011 UTC (13 years, 11 months ago) by ceballos
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_023, Mit_022a, Mit_022
Changes since 1.52: +31 -6 lines
Log Message: 
fixes in the electron id

File Contents

# User Rev Content
1 ceballos 1.53 // $Id: MuonIDMod.cc,v 1.52 2011/05/21 17:05:58 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.53 fMuonIDType("WWMuIdV1"),
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 ceballos 1.53 case kGlobalTracker:
102     pass = (mu->HasGlobalTrk() && mu->GlobalTrk()->Chi2()/mu->GlobalTrk()->Ndof() < 10 &&
103     (mu->NSegments() > 1 || mu->NMatches() > 1) && mu->NValidHits() > 0) ||
104     (mu->IsTrackerMuon() &&
105     mu->Quality().Quality(MuonQuality::TMLastStationTight));
106     if (pass) {
107     pt = mu->TrackerTrk()->Pt();
108     eta = TMath::Abs(mu->TrackerTrk()->Eta());
109     }
110     else {
111     pt = mu->Pt();
112     eta = TMath::Abs(mu->Eta());
113     }
114     break;
115 loizides 1.6 case kSta:
116 loizides 1.19 pass = mu->HasStandaloneTrk();
117 ceballos 1.30 if (pass) {
118     pt = mu->StandaloneTrk()->Pt();
119     eta = TMath::Abs(mu->StandaloneTrk()->Eta());
120     }
121 loizides 1.6 break;
122 loizides 1.24 case kTrackerMuon:
123 ceballos 1.28 pass = mu->HasTrackerTrk() && mu->IsTrackerMuon() &&
124     mu->Quality().Quality(MuonQuality::TrackerMuonArbitrated);
125 ceballos 1.30 if (pass) {
126     pt = mu->TrackerTrk()->Pt();
127     eta = TMath::Abs(mu->TrackerTrk()->Eta());
128     }
129 loizides 1.24 break;
130     case kCaloMuon:
131     pass = mu->HasTrackerTrk() && mu->IsCaloMuon();
132 ceballos 1.30 if (pass) {
133     pt = mu->TrackerTrk()->Pt();
134     eta = TMath::Abs(mu->TrackerTrk()->Eta());
135     }
136 loizides 1.24 break;
137     case kTrackerBased:
138 loizides 1.19 pass = mu->HasTrackerTrk();
139 ceballos 1.30 if (pass) {
140     pt = mu->TrackerTrk()->Pt();
141     eta = TMath::Abs(mu->TrackerTrk()->Eta());
142     }
143 loizides 1.6 break;
144     default:
145     break;
146 ceballos 1.5 }
147 loizides 1.6
148     if (!pass)
149     continue;
150    
151     if (pt <= fMuonPtMin)
152     continue;
153    
154 ceballos 1.30 if (eta >= fEtaCut)
155     continue;
156    
157 ceballos 1.34 Double_t RChi2 = 0.0;
158     if (mu->HasGlobalTrk()) {
159     RChi2 = mu->GlobalTrk()->Chi2()/mu->GlobalTrk()->Ndof();
160     }
161     else if(mu->BestTrk() != 0){
162     RChi2 = mu->BestTrk()->Chi2()/mu->BestTrk()->Ndof();
163     }
164 loizides 1.6 Bool_t idpass = kFALSE;
165     switch (fMuIDType) {
166 ceballos 1.34 case kWMuId:
167     idpass = mu->BestTrk() != 0 &&
168     mu->BestTrk()->NHits() > 10 &&
169     RChi2 < 10.0 &&
170 ceballos 1.36 (mu->NSegments() > 1 || mu->NMatches() > 1) &&
171 ceballos 1.34 mu->BestTrk()->NPixelHits() > 0 &&
172     mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
173     break;
174     case kZMuId:
175     idpass = mu->BestTrk() != 0 &&
176     mu->BestTrk()->NHits() > 10 &&
177 ceballos 1.36 (mu->NSegments() > 1 || mu->NMatches() > 1) &&
178 ceballos 1.34 mu->BestTrk()->NPixelHits() > 0 &&
179     mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
180     break;
181 loizides 1.6 case kLoose:
182 ceballos 1.32 idpass = mu->BestTrk() != 0 &&
183     mu->Quality().Quality(MuonQuality::TMOneStationLoose) &&
184 ceballos 1.28 mu->Quality().Quality(MuonQuality::TM2DCompatibilityLoose) &&
185 ceballos 1.29 mu->BestTrk()->NHits() > 10 &&
186 ceballos 1.34 RChi2 < 10.0 &&
187 ceballos 1.30 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
188 loizides 1.6 break;
189     case kTight:
190 ceballos 1.32 idpass = mu->BestTrk() != 0 &&
191     mu->Quality().Quality(MuonQuality::TMOneStationTight) &&
192 ceballos 1.28 mu->Quality().Quality(MuonQuality::TM2DCompatibilityTight) &&
193 ceballos 1.29 mu->BestTrk()->NHits() > 10 &&
194 ceballos 1.34 RChi2 < 10.0 &&
195 ceballos 1.30 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
196 ceballos 1.28 break;
197 ceballos 1.53 case kWWMuIdV1:
198 ceballos 1.32 idpass = mu->BestTrk() != 0 &&
199     mu->BestTrk()->NHits() > 10 &&
200 ceballos 1.53 mu->BestTrk()->NPixelHits() > 0 &&
201     mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1 &&
202 ceballos 1.34 RChi2 < 10.0 &&
203 ceballos 1.36 (mu->NSegments() > 1 || mu->NMatches() > 1) &&
204 ceballos 1.53 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
205     break;
206     case kWWMuIdV2:
207     idpass = mu->BestTrk() != 0 &&
208     mu->BestTrk()->NHits() > 10 &&
209 ceballos 1.36 mu->BestTrk()->NPixelHits() > 0 &&
210     mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1;
211 loizides 1.6 break;
212 ceballos 1.53
213 loizides 1.18 case kNoId:
214     idpass = kTRUE;
215     break;
216 loizides 1.6 default:
217     break;
218 ceballos 1.4 }
219 loizides 1.6
220     if (!idpass)
221     continue;
222    
223 ceballos 1.22 Bool_t isocut = kFALSE;
224 loizides 1.6 switch (fMuIsoType) {
225     case kTrackCalo:
226 ceballos 1.22 isocut = (mu->IsoR03SumPt() < fTrackIsolationCut) &&
227 loizides 1.6 (mu->IsoR03EmEt() + mu->IsoR03HadEt() < fCaloIsolationCut);
228     break;
229     case kTrackCaloCombined:
230 ceballos 1.38 isocut = (1.0 * mu->IsoR03SumPt() +
231     1.0 * mu->IsoR03EmEt() +
232     1.0 * mu->IsoR03HadEt() < fCombIsolationCut);
233 loizides 1.6 break;
234     case kTrackCaloSliding:
235     {
236 ceballos 1.41 const PileupEnergyDensity *rho = fPileupEnergyDensity->At(0);
237     Double_t totalIso = mu->IsoR03SumPt() + TMath::Max(mu->IsoR03EmEt() + mu->IsoR03HadEt() - rho->Rho() * TMath::Pi() * 0.3 * 0.3, 0.0);
238 ceballos 1.46 // trick to change the signal region cut
239     double theIsoCut = fCombIsolationCut;
240     if(theIsoCut < 0.20){
241     if(mu->Pt() > 20.0) theIsoCut = 0.15;
242     else theIsoCut = 0.10;
243     }
244     if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
245 ceballos 1.22 }
246 loizides 1.6 break;
247 ceballos 1.41 case kTrackCaloSlidingNoCorrection:
248 ceballos 1.40 {
249     Double_t totalIso = 1.0 * mu->IsoR03SumPt() +
250 ceballos 1.41 1.0 * mu->IsoR03EmEt() +
251     1.0 * mu->IsoR03HadEt();
252 ceballos 1.46 // trick to change the signal region cut
253     double theIsoCut = fCombIsolationCut;
254     if(theIsoCut < 0.20){
255     if(mu->Pt() > 20.0) theIsoCut = 0.15;
256     else theIsoCut = 0.10;
257     }
258     if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
259 ceballos 1.40 }
260     break;
261 ceballos 1.38 case kPFIso:
262     {
263 ceballos 1.49 Double_t pfIsoCutValue = 9999;
264     if(fCombIsolationCut > 0){
265     pfIsoCutValue = fCombIsolationCut;
266     } else {
267 ceballos 1.51 if (mu->AbsEta() < 1.479) {
268     if (mu->Pt() > 20) {
269 ceballos 1.52 pfIsoCutValue = 0.13;
270 ceballos 1.51 } else {
271 ceballos 1.52 pfIsoCutValue = 0.06;
272 ceballos 1.51 }
273 ceballos 1.49 } else {
274 ceballos 1.51 if (mu->Pt() > 20) {
275 ceballos 1.52 pfIsoCutValue = 0.09;
276 ceballos 1.51 } else {
277 ceballos 1.52 pfIsoCutValue = 0.05;
278 ceballos 1.51 }
279     }
280 ceballos 1.49 }
281 ceballos 1.52 Double_t totalIso = IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0);
282 ceballos 1.49 if (totalIso < (mu->Pt()*pfIsoCutValue) )
283 ceballos 1.38 isocut = kTRUE;
284     }
285     break;
286     case kPFIsoNoL:
287     {
288 ceballos 1.39 fNonIsolatedMuons = GetObjThisEvt<MuonCol>(fNonIsolatedMuonsName);
289     fNonIsolatedElectrons = GetObjThisEvt<ElectronCol>(fNonIsolatedElectronsName);
290 ceballos 1.38
291     Double_t beta = IsolationTools::BetaM(fTracks, mu, fVertices->At(0), 0.0, 0.2, 0.3, 0.02);
292 ceballos 1.40 if(beta == 0) beta = 1.0;
293 mzanetti 1.47 Double_t totalIso = IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), fNonIsolatedMuons, fNonIsolatedElectrons, 0.2, 1.0, 0.4, 0.0, 3);
294 ceballos 1.38 if (totalIso < (mu->Pt()*fCombIsolationCut) )
295     isocut = kTRUE;
296     }
297     break;
298 loizides 1.13 case kNoIso:
299 ceballos 1.22 isocut = kTRUE;
300 loizides 1.13 break;
301 loizides 1.6 case kCustomIso:
302     default:
303     break;
304 ceballos 1.4 }
305 ceballos 1.3
306 ceballos 1.22 if (isocut == kFALSE)
307 loizides 1.6 continue;
308    
309 ceballos 1.42 // apply d0 cut
310 loizides 1.21 if (fApplyD0Cut) {
311 ceballos 1.42 Bool_t passD0cut = kTRUE;
312 ceballos 1.46 if(fD0Cut < 0.05) { // trick to change the signal region cut
313     if (mu->Pt() > 20.0) fD0Cut = 0.02;
314     else if (mu->Pt() <= 20.0) fD0Cut = 0.01;
315     }
316 ceballos 1.42 if(fWhichVertex >= -1) passD0cut = MuonTools::PassD0Cut(mu, fVertices, fD0Cut, fWhichVertex);
317     else passD0cut = MuonTools::PassD0Cut(mu, fBeamSpot, fD0Cut);
318 ceballos 1.31 if (!passD0cut)
319 loizides 1.21 continue;
320 ceballos 1.14 }
321    
322 ceballos 1.42 // apply dz cut
323     if (fApplyDZCut) {
324 ceballos 1.45 Bool_t passDZcut = MuonTools::PassDZCut(mu, fVertices, fDZCut, fWhichVertex);
325 ceballos 1.42 if (!passDZcut)
326     continue;
327     }
328    
329 loizides 1.6 // add good muon
330     CleanMuons->Add(mu);
331 loizides 1.1 }
332    
333 loizides 1.10 // sort according to pt
334     CleanMuons->Sort();
335    
336 loizides 1.6 // add objects for other modules to use
337 loizides 1.7 AddObjThisEvt(CleanMuons);
338 loizides 1.1 }
339    
340     //--------------------------------------------------------------------------------------------------
341     void MuonIDMod::SlaveBegin()
342     {
343     // Run startup code on the computer (slave) doing the actual analysis. Here,
344 loizides 1.6 // we just request the muon collection branch.
345    
346 ceballos 1.38 // In this case we cannot have a branch
347     if (fMuonIsoType.CompareTo("PFIsoNoL") != 0) {
348     ReqEventObject(fMuonBranchName, fMuons, kTRUE);
349     }
350 ceballos 1.42 ReqEventObject(fBeamSpotName, fBeamSpot, kTRUE);
351 ceballos 1.38 ReqEventObject(fTrackName, fTracks, kTRUE);
352     ReqEventObject(fPFCandidatesName, fPFCandidates, kTRUE);
353 ceballos 1.41 if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0) {
354     ReqEventObject(fPileupEnergyDensityName, fPileupEnergyDensity, kTRUE);
355     }
356 loizides 1.1
357 ceballos 1.34 if (fMuonIDType.CompareTo("WMuId") == 0)
358     fMuIDType = kWMuId;
359     else if (fMuonIDType.CompareTo("ZMuId") == 0)
360     fMuIDType = kZMuId;
361     else if (fMuonIDType.CompareTo("Tight") == 0)
362 loizides 1.6 fMuIDType = kTight;
363     else if (fMuonIDType.CompareTo("Loose") == 0)
364     fMuIDType = kLoose;
365 ceballos 1.53 else if (fMuonIDType.CompareTo("WWMuIdV1") == 0)
366     fMuIDType = kWWMuIdV1;
367     else if (fMuonIDType.CompareTo("WWMuIdV2") == 0)
368     fMuIDType = kWWMuIdV2;
369 loizides 1.12 else if (fMuonIDType.CompareTo("NoId") == 0)
370     fMuIDType = kNoId;
371 loizides 1.6 else if (fMuonIDType.CompareTo("Custom") == 0) {
372     fMuIDType = kCustomId;
373     SendError(kWarning, "SlaveBegin",
374     "Custom muon identification is not yet implemented.");
375     } else {
376     SendError(kAbortAnalysis, "SlaveBegin",
377     "The specified muon identification %s is not defined.",
378     fMuonIDType.Data());
379     return;
380     }
381 loizides 1.1
382 loizides 1.6 if (fMuonIsoType.CompareTo("TrackCalo") == 0)
383     fMuIsoType = kTrackCalo;
384     else if (fMuonIsoType.CompareTo("TrackCaloCombined") == 0)
385     fMuIsoType = kTrackCaloCombined;
386     else if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0)
387     fMuIsoType = kTrackCaloSliding;
388 ceballos 1.41 else if (fMuonIsoType.CompareTo("TrackCaloSlidingNoCorrection") == 0)
389     fMuIsoType = kTrackCaloSlidingNoCorrection;
390 ceballos 1.38 else if (fMuonIsoType.CompareTo("PFIso") == 0)
391     fMuIsoType = kPFIso;
392     else if (fMuonIsoType.CompareTo("PFIsoNoL") == 0)
393     fMuIsoType = kPFIsoNoL;
394 loizides 1.6 else if (fMuonIsoType.CompareTo("NoIso") == 0)
395     fMuIsoType = kNoIso;
396     else if (fMuonIsoType.CompareTo("Custom") == 0) {
397     fMuIsoType = kCustomIso;
398     SendError(kWarning, "SlaveBegin",
399     "Custom muon isolation is not yet implemented.");
400     } else {
401     SendError(kAbortAnalysis, "SlaveBegin",
402     "The specified muon isolation %s is not defined.",
403     fMuonIsoType.Data());
404     return;
405     }
406 loizides 1.1
407 loizides 1.6 if (fMuonClassType.CompareTo("All") == 0)
408     fMuClassType = kAll;
409     else if (fMuonClassType.CompareTo("Global") == 0)
410     fMuClassType = kGlobal;
411 ceballos 1.53 else if (fMuonClassType.CompareTo("GlobalTracker") == 0)
412     fMuClassType = kGlobalTracker;
413 loizides 1.6 else if (fMuonClassType.CompareTo("Standalone") == 0)
414     fMuClassType = kSta;
415 loizides 1.24 else if (fMuonClassType.CompareTo("TrackerMuon") == 0)
416     fMuClassType = kTrackerMuon;
417     else if (fMuonClassType.CompareTo("CaloMuon") == 0)
418     fMuClassType = kCaloMuon;
419     else if (fMuonClassType.CompareTo("TrackerBased") == 0)
420     fMuClassType = kTrackerBased;
421 loizides 1.6 else {
422     SendError(kAbortAnalysis, "SlaveBegin",
423     "The specified muon class %s is not defined.",
424     fMuonClassType.Data());
425     return;
426     }
427 loizides 1.1 }