ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/MuonIDMod.cc
Revision: 1.46
Committed: Tue Apr 5 06:37:07 2011 UTC (14 years, 1 month ago) by ceballos
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_020d, TMit_020d, Mit_020c, Mit_021pre1, Mit_020b
Changes since 1.45: +19 -9 lines
Log Message: 
new

File Contents

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