Mods/src/MuonIDMod.cc

// $Id: MuonIDMod.cc,v 1.29 2010/05/21 05:59:55 ceballos Exp $

#include "MitPhysics/Mods/interface/MuonIDMod.h"
#include "MitCommon/MathTools/interface/MathUtils.h"
#include "MitAna/DataTree/interface/MuonCol.h"
#include "MitAna/DataTree/interface/VertexCol.h"
#include "MitPhysics/Init/interface/ModNames.h"

using namespace mithep;

ClassImp(mithep::MuonIDMod)

//--------------------------------------------------------------------------------------------------
  MuonIDMod::MuonIDMod(const char *name, const char *title) : 
  BaseMod(name,title),
  fMuonBranchName(Names::gkMuonBrn),
  fCleanMuonsName(ModNames::gkCleanMuonsName),  
  fVertexName("PrimaryVertexes"),
  fMuonIDType("Minimal"),
  fMuonIsoType("TrackCaloSliding"),  
  fMuonClassType("Global"),  
  fTrackIsolationCut(3.0),
  fCaloIsolationCut(3.0),
  fCombIsolationCut(5.0),
  fMuonPtMin(10),
  fApplyD0Cut(kTRUE),
  fD0Cut(0.020),
  fEtaCut(2.4),
  fReverseIsoCut(kFALSE),
  fReverseD0Cut(kFALSE),
  fMuIDType(kIdUndef),
  fMuIsoType(kIsoUndef),
  fMuClassType(kClassUndef),
  fMuons(0),
  fVertices(0),
  fMuonTools(0)
{
  // Constructor.
}

//--------------------------------------------------------------------------------------------------
void MuonIDMod::Process()
{
  // Process entries of the tree. 

  LoadEventObject(fMuonBranchName, fMuons);
  if (fApplyD0Cut) {
    LoadEventObject(fVertexName,     fVertices);
  }

  MuonOArr *CleanMuons = new MuonOArr;
  CleanMuons->SetName(fCleanMuonsName);

  for (UInt_t i=0; i<fMuons->GetEntries(); ++i) {
    const Muon *mu = fMuons->At(i);

    Bool_t pass = kFALSE;
    Double_t pt = 0;  // make sure pt is taken from the correct track!
    Double_t eta = 0; // make sure eta is taken from the correct track!
    switch (fMuClassType) {
      case kAll:
        pass = kTRUE;
        if (mu->HasTrk()) {
          pt  = mu->Pt();
          eta = TMath::Abs(mu->Eta());
        }
        break;
      case kGlobal:
        pass = mu->HasGlobalTrk() && mu->IsTrackerMuon() &&
               mu->Quality().Quality(MuonQuality::TrackerMuonArbitrated);
        if (pass) {
          pt  = mu->TrackerTrk()->Pt();
          eta = TMath::Abs(mu->TrackerTrk()->Eta());
        }
        break;
      case kSta:
        pass = mu->HasStandaloneTrk();
        if (pass) {
          pt  = mu->StandaloneTrk()->Pt();
          eta = TMath::Abs(mu->StandaloneTrk()->Eta());
        }
        break;
      case kTrackerMuon:
        pass = mu->HasTrackerTrk() && mu->IsTrackerMuon() &&
               mu->Quality().Quality(MuonQuality::TrackerMuonArbitrated);
        if (pass) {
          pt  = mu->TrackerTrk()->Pt();
          eta = TMath::Abs(mu->TrackerTrk()->Eta());
        }
        break;
      case kCaloMuon:
        pass = mu->HasTrackerTrk() && mu->IsCaloMuon();
        if (pass) {
          pt  = mu->TrackerTrk()->Pt();
          eta = TMath::Abs(mu->TrackerTrk()->Eta());
        }
        break;
      case kTrackerBased:
        pass = mu->HasTrackerTrk();
        if (pass) {
          pt  = mu->TrackerTrk()->Pt();
          eta = TMath::Abs(mu->TrackerTrk()->Eta());
        }
        break;
      default:
        break;
    }

    if (!pass)
      continue;

    if (pt <= fMuonPtMin) 
      continue;

    if (eta >= fEtaCut) 
      continue;

    Bool_t idpass = kFALSE;
    switch (fMuIDType) {
      case kLoose:
        idpass = mu->Quality().Quality(MuonQuality::TMOneStationLoose) &&
                 mu->Quality().Quality(MuonQuality::TM2DCompatibilityLoose) &&
                 mu->BestTrk()->NHits() > 10 &&
                 mu->BestTrk()->Chi2()/mu->BestTrk()->Ndof() < 10 &&
                 mu->NSegments() > 0 &&
                 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
        break;
      case kTight:
        idpass = mu->Quality().Quality(MuonQuality::TMOneStationTight) &&
                 mu->Quality().Quality(MuonQuality::TM2DCompatibilityTight) &&
                 mu->BestTrk()->NHits() > 10 &&
                 mu->BestTrk()->Chi2()/mu->BestTrk()->Ndof() < 10 &&
                 mu->NSegments() > 0 &&
                 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
        break;
      case kMinimal:
        idpass = mu->BestTrk()->NHits() > 10 &&
                 mu->BestTrk()->Chi2()/mu->BestTrk()->Ndof() < 10 &&
                 mu->NSegments() > 0 &&
                 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
        break;
      case kNoId:
        idpass = kTRUE;
        break;
      default:
        break;
    }

    if (!idpass)
      continue;

    Bool_t isocut = kFALSE;
    switch (fMuIsoType) {
      case kTrackCalo:
        isocut = (mu->IsoR03SumPt() < fTrackIsolationCut) &&
          (mu->IsoR03EmEt() + mu->IsoR03HadEt() < fCaloIsolationCut);
        break;
      case kTrackCaloCombined:
        isocut = (1.0 * mu->IsoR03SumPt() + 1.0 * mu->IsoR03EmEt() + 
                  1.0 * mu->IsoR03HadEt() < fCombIsolationCut);
        break;
      case kTrackCaloSliding:
        { 
          Double_t totalIso = 1.0 * mu->IsoR03SumPt() + 
                              1.0 * mu->IsoR03EmEt() + 
                              1.0 * mu->IsoR03HadEt();
          if (totalIso < (mu->Pt()*0.15) )
            isocut = kTRUE;

          if     (fReverseIsoCut == kTRUE &&
                  isocut == kFALSE && totalIso < 10)
            isocut = kTRUE;
          else if(fReverseIsoCut == kTRUE)
            isocut = kFALSE;
        }
        break;
      case kNoIso:
        isocut = kTRUE;
        break;
      case kCustomIso:
      default:
        break;
    }

    if (isocut == kFALSE)
      continue;

    if (fApplyD0Cut) {
      Bool_t d0cut = kFALSE;
      const Track *mt = mu->BestTrk();
      if (!mt)
        continue;
      Double_t d0_real = 1e30;
      for(UInt_t i0 = 0; i0 < fVertices->GetEntries(); i0++) {
        Double_t pD0 = mt->D0Corrected(*fVertices->At(i0));
        if(TMath::Abs(pD0) < TMath::Abs(d0_real)) 
          d0_real = TMath::Abs(pD0);
      }
      if(d0_real < fD0Cut) d0cut = kTRUE;

      if     (fReverseD0Cut == kTRUE &&
              d0cut == kFALSE && d0_real < 0.05)
        d0cut = kTRUE;
      else if(fReverseD0Cut == kTRUE)
        d0cut = kFALSE;

      if (d0cut == kFALSE)
        continue;
    }

    // add good muon
    CleanMuons->Add(mu);
  }

  // sort according to pt
  CleanMuons->Sort();

  // add objects for other modules to use
  AddObjThisEvt(CleanMuons);  
}

//--------------------------------------------------------------------------------------------------
void MuonIDMod::SlaveBegin()
{
  // Run startup code on the computer (slave) doing the actual analysis. Here,
  // we just request the muon collection branch.

  ReqEventObject(fMuonBranchName, fMuons, kTRUE);

  if (fApplyD0Cut)
    ReqEventObject(fVertexName, fVertices, kTRUE);

  fMuonTools = new MuonTools;

  if (fMuonIDType.CompareTo("Tight") == 0) 
    fMuIDType = kTight;
  else if (fMuonIDType.CompareTo("Loose") == 0) 
    fMuIDType = kLoose;
  else if (fMuonIDType.CompareTo("Minimal") == 0) 
    fMuIDType = kMinimal;
  else if (fMuonIDType.CompareTo("NoId") == 0) 
    fMuIDType = kNoId;
  else if (fMuonIDType.CompareTo("Custom") == 0) {
    fMuIDType = kCustomId;
    SendError(kWarning, "SlaveBegin",
              "Custom muon identification is not yet implemented.");
  } else {
    SendError(kAbortAnalysis, "SlaveBegin",
              "The specified muon identification %s is not defined.",
              fMuonIDType.Data());
    return;
  }

  if (fMuonIsoType.CompareTo("TrackCalo") == 0)
    fMuIsoType = kTrackCalo;
  else if (fMuonIsoType.CompareTo("TrackCaloCombined") == 0)
    fMuIsoType = kTrackCaloCombined;
  else if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0)
    fMuIsoType = kTrackCaloSliding;
  else if (fMuonIsoType.CompareTo("NoIso") == 0)
    fMuIsoType = kNoIso;
  else if (fMuonIsoType.CompareTo("Custom") == 0) {
    fMuIsoType = kCustomIso;
    SendError(kWarning, "SlaveBegin",
              "Custom muon isolation is not yet implemented.");
  } else {
    SendError(kAbortAnalysis, "SlaveBegin",
              "The specified muon isolation %s is not defined.",
              fMuonIsoType.Data());
    return;
  }

  if (fMuonClassType.CompareTo("All") == 0) 
    fMuClassType = kAll;
  else if (fMuonClassType.CompareTo("Global") == 0) 
    fMuClassType = kGlobal;
  else if (fMuonClassType.CompareTo("Standalone") == 0) 
    fMuClassType = kSta;
  else if (fMuonClassType.CompareTo("TrackerMuon") == 0) 
    fMuClassType = kTrackerMuon;
  else if (fMuonClassType.CompareTo("CaloMuon") == 0) 
    fMuClassType = kCaloMuon;
  else if (fMuonClassType.CompareTo("TrackerBased") == 0) 
    fMuClassType = kTrackerBased;
  else {
    SendError(kAbortAnalysis, "SlaveBegin",
              "The specified muon class %s is not defined.",
              fMuonClassType.Data());
    return;
  }
}
Revision:	1.30
Committed:	Thu May 27 07:59:03 2010 UTC (14 years, 11 months ago) by ceballos
Content type:	text/plain
Branch:	MAIN
CVS Tags:	Mit_014a, Mit_014, Mit_014pre3
Changes since 1.29:	+40 -20 lines
Log Message:	updates on electrons and muons
#	User	Rev	Content
1	ceballos	1.30	// $Id: MuonIDMod.cc,v 1.29 2010/05/21 05:59:55 ceballos Exp $
2	loizides	1.1
3			#include "MitPhysics/Mods/interface/MuonIDMod.h"
4	loizides	1.6	#include "MitCommon/MathTools/interface/MathUtils.h"
5	loizides	1.23	#include "MitAna/DataTree/interface/MuonCol.h"
6			#include "MitAna/DataTree/interface/VertexCol.h"
7	loizides	1.6	#include "MitPhysics/Init/interface/ModNames.h"
8	loizides	1.1
9			using namespace mithep;
10
11			ClassImp(mithep::MuonIDMod)
12
13			//--------------------------------------------------------------------------------------------------
14			MuonIDMod::MuonIDMod(const char name, const char title) :
15			BaseMod(name,title),
16	loizides	1.6	fMuonBranchName(Names::gkMuonBrn),
17			fCleanMuonsName(ModNames::gkCleanMuonsName),
18	ceballos	1.28	fVertexName("PrimaryVertexes"),
19	ceballos	1.30	fMuonIDType("Minimal"),
20	loizides	1.8	fMuonIsoType("TrackCaloSliding"),
21	loizides	1.6	fMuonClassType("Global"),
22	ceballos	1.2	fTrackIsolationCut(3.0),
23			fCaloIsolationCut(3.0),
24	loizides	1.8	fCombIsolationCut(5.0),
25	ceballos	1.3	fMuonPtMin(10),
26	loizides	1.21	fApplyD0Cut(kTRUE),
27	ceballos	1.30	fD0Cut(0.020),
28			fEtaCut(2.4),
29	loizides	1.20	fReverseIsoCut(kFALSE),
30	ceballos	1.22	fReverseD0Cut(kFALSE),
31	loizides	1.15	fMuIDType(kIdUndef),
32			fMuIsoType(kIsoUndef),
33			fMuClassType(kClassUndef),
34	ceballos	1.14	fMuons(0),
35	loizides	1.15	fVertices(0),
36	loizides	1.20	fMuonTools(0)
37	loizides	1.1	{
38			// Constructor.
39			}
40
41			//--------------------------------------------------------------------------------------------------
42			void MuonIDMod::Process()
43			{
44			// Process entries of the tree.
45
46	loizides	1.21	LoadEventObject(fMuonBranchName, fMuons);
47	sixie	1.25	if (fApplyD0Cut) {
48			LoadEventObject(fVertexName, fVertices);
49			}
50	loizides	1.1
51	loizides	1.7	MuonOArr *CleanMuons = new MuonOArr;
52			CleanMuons->SetName(fCleanMuonsName);
53	loizides	1.1
54			for (UInt_t i=0; i<fMuons->GetEntries(); ++i) {
55	loizides	1.8	const Muon *mu = fMuons->At(i);
56	loizides	1.6
57			Bool_t pass = kFALSE;
58	ceballos	1.30	Double_t pt = 0; // make sure pt is taken from the correct track!
59			Double_t eta = 0; // make sure eta is taken from the correct track!
60	loizides	1.6	switch (fMuClassType) {
61			case kAll:
62			pass = kTRUE;
63	ceballos	1.30	if (mu->HasTrk()) {
64			pt = mu->Pt();
65			eta = TMath::Abs(mu->Eta());
66			}
67	loizides	1.6	break;
68			case kGlobal:
69	ceballos	1.28	pass = mu->HasGlobalTrk() && mu->IsTrackerMuon() &&
70			mu->Quality().Quality(MuonQuality::TrackerMuonArbitrated);
71	ceballos	1.30	if (pass) {
72			pt = mu->TrackerTrk()->Pt();
73			eta = TMath::Abs(mu->TrackerTrk()->Eta());
74			}
75			break;
76	loizides	1.6	case kSta:
77	loizides	1.19	pass = mu->HasStandaloneTrk();
78	ceballos	1.30	if (pass) {
79			pt = mu->StandaloneTrk()->Pt();
80			eta = TMath::Abs(mu->StandaloneTrk()->Eta());
81			}
82	loizides	1.6	break;
83	loizides	1.24	case kTrackerMuon:
84	ceballos	1.28	pass = mu->HasTrackerTrk() && mu->IsTrackerMuon() &&
85			mu->Quality().Quality(MuonQuality::TrackerMuonArbitrated);
86	ceballos	1.30	if (pass) {
87			pt = mu->TrackerTrk()->Pt();
88			eta = TMath::Abs(mu->TrackerTrk()->Eta());
89			}
90	loizides	1.24	break;
91			case kCaloMuon:
92			pass = mu->HasTrackerTrk() && mu->IsCaloMuon();
93	ceballos	1.30	if (pass) {
94			pt = mu->TrackerTrk()->Pt();
95			eta = TMath::Abs(mu->TrackerTrk()->Eta());
96			}
97	loizides	1.24	break;
98			case kTrackerBased:
99	loizides	1.19	pass = mu->HasTrackerTrk();
100	ceballos	1.30	if (pass) {
101			pt = mu->TrackerTrk()->Pt();
102			eta = TMath::Abs(mu->TrackerTrk()->Eta());
103			}
104	loizides	1.6	break;
105			default:
106			break;
107	ceballos	1.5	}
108	loizides	1.6
109			if (!pass)
110			continue;
111
112			if (pt <= fMuonPtMin)
113			continue;
114
115	ceballos	1.30	if (eta >= fEtaCut)
116			continue;
117
118	loizides	1.6	Bool_t idpass = kFALSE;
119			switch (fMuIDType) {
120			case kLoose:
121	ceballos	1.28	idpass = mu->Quality().Quality(MuonQuality::TMOneStationLoose) &&
122			mu->Quality().Quality(MuonQuality::TM2DCompatibilityLoose) &&
123	ceballos	1.29	mu->BestTrk()->NHits() > 10 &&
124	ceballos	1.28	mu->BestTrk()->Chi2()/mu->BestTrk()->Ndof() < 10 &&
125	ceballos	1.30	mu->NSegments() > 0 &&
126			mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
127	loizides	1.6	break;
128			case kTight:
129	ceballos	1.28	idpass = mu->Quality().Quality(MuonQuality::TMOneStationTight) &&
130			mu->Quality().Quality(MuonQuality::TM2DCompatibilityTight) &&
131	ceballos	1.29	mu->BestTrk()->NHits() > 10 &&
132	ceballos	1.28	mu->BestTrk()->Chi2()/mu->BestTrk()->Ndof() < 10 &&
133	ceballos	1.30	mu->NSegments() > 0 &&
134			mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
135	ceballos	1.28	break;
136			case kMinimal:
137	ceballos	1.29	idpass = mu->BestTrk()->NHits() > 10 &&
138	ceballos	1.28	mu->BestTrk()->Chi2()/mu->BestTrk()->Ndof() < 10 &&
139	ceballos	1.30	mu->NSegments() > 0 &&
140			mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
141	loizides	1.6	break;
142	loizides	1.18	case kNoId:
143			idpass = kTRUE;
144			break;
145	loizides	1.6	default:
146			break;
147	ceballos	1.4	}
148	loizides	1.6
149			if (!idpass)
150			continue;
151
152	ceballos	1.22	Bool_t isocut = kFALSE;
153	loizides	1.6	switch (fMuIsoType) {
154			case kTrackCalo:
155	ceballos	1.22	isocut = (mu->IsoR03SumPt() < fTrackIsolationCut) &&
156	loizides	1.6	(mu->IsoR03EmEt() + mu->IsoR03HadEt() < fCaloIsolationCut);
157			break;
158			case kTrackCaloCombined:
159	ceballos	1.22	isocut = (1.0 * mu->IsoR03SumPt() + 1.0 * mu->IsoR03EmEt() +
160	ceballos	1.28	1.0 * mu->IsoR03HadEt() < fCombIsolationCut);
161	loizides	1.6	break;
162			case kTrackCaloSliding:
163			{
164			Double_t totalIso = 1.0 * mu->IsoR03SumPt() +
165			1.0 * mu->IsoR03EmEt() +
166			1.0 * mu->IsoR03HadEt();
167	ceballos	1.29	if (totalIso < (mu->Pt()*0.15) )
168	ceballos	1.22	isocut = kTRUE;
169
170			if (fReverseIsoCut == kTRUE &&
171			isocut == kFALSE && totalIso < 10)
172			isocut = kTRUE;
173			else if(fReverseIsoCut == kTRUE)
174			isocut = kFALSE;
175			}
176	loizides	1.6	break;
177	loizides	1.13	case kNoIso:
178	ceballos	1.22	isocut = kTRUE;
179	loizides	1.13	break;
180	loizides	1.6	case kCustomIso:
181			default:
182			break;
183	ceballos	1.4	}
184	ceballos	1.3
185	ceballos	1.22	if (isocut == kFALSE)
186	loizides	1.6	continue;
187
188	loizides	1.21	if (fApplyD0Cut) {
189	ceballos	1.22	Bool_t d0cut = kFALSE;
190	loizides	1.21	const Track *mt = mu->BestTrk();
191			if (!mt)
192			continue;
193			Double_t d0_real = 1e30;
194			for(UInt_t i0 = 0; i0 < fVertices->GetEntries(); i0++) {
195			Double_t pD0 = mt->D0Corrected(*fVertices->At(i0));
196			if(TMath::Abs(pD0) < TMath::Abs(d0_real))
197			d0_real = TMath::Abs(pD0);
198			}
199	ceballos	1.22	if(d0_real < fD0Cut) d0cut = kTRUE;
200
201			if (fReverseD0Cut == kTRUE &&
202			d0cut == kFALSE && d0_real < 0.05)
203			d0cut = kTRUE;
204			else if(fReverseD0Cut == kTRUE)
205			d0cut = kFALSE;
206
207			if (d0cut == kFALSE)
208	loizides	1.21	continue;
209	ceballos	1.14	}
210
211	loizides	1.6	// add good muon
212			CleanMuons->Add(mu);
213	loizides	1.1	}
214
215	loizides	1.10	// sort according to pt
216			CleanMuons->Sort();
217
218	loizides	1.6	// add objects for other modules to use
219	loizides	1.7	AddObjThisEvt(CleanMuons);
220	loizides	1.1	}
221
222			//--------------------------------------------------------------------------------------------------
223			void MuonIDMod::SlaveBegin()
224			{
225			// Run startup code on the computer (slave) doing the actual analysis. Here,
226	loizides	1.6	// we just request the muon collection branch.
227
228	loizides	1.21	ReqEventObject(fMuonBranchName, fMuons, kTRUE);
229
230			if (fApplyD0Cut)
231			ReqEventObject(fVertexName, fVertices, kTRUE);
232	loizides	1.1
233	loizides	1.6	fMuonTools = new MuonTools;
234	loizides	1.1
235	loizides	1.6	if (fMuonIDType.CompareTo("Tight") == 0)
236			fMuIDType = kTight;
237			else if (fMuonIDType.CompareTo("Loose") == 0)
238			fMuIDType = kLoose;
239	ceballos	1.28	else if (fMuonIDType.CompareTo("Minimal") == 0)
240			fMuIDType = kMinimal;
241	loizides	1.12	else if (fMuonIDType.CompareTo("NoId") == 0)
242			fMuIDType = kNoId;
243	loizides	1.6	else if (fMuonIDType.CompareTo("Custom") == 0) {
244			fMuIDType = kCustomId;
245			SendError(kWarning, "SlaveBegin",
246			"Custom muon identification is not yet implemented.");
247			} else {
248			SendError(kAbortAnalysis, "SlaveBegin",
249			"The specified muon identification %s is not defined.",
250			fMuonIDType.Data());
251			return;
252			}
253	loizides	1.1
254	loizides	1.6	if (fMuonIsoType.CompareTo("TrackCalo") == 0)
255			fMuIsoType = kTrackCalo;
256			else if (fMuonIsoType.CompareTo("TrackCaloCombined") == 0)
257			fMuIsoType = kTrackCaloCombined;
258			else if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0)
259			fMuIsoType = kTrackCaloSliding;
260			else if (fMuonIsoType.CompareTo("NoIso") == 0)
261			fMuIsoType = kNoIso;
262			else if (fMuonIsoType.CompareTo("Custom") == 0) {
263			fMuIsoType = kCustomIso;
264			SendError(kWarning, "SlaveBegin",
265			"Custom muon isolation is not yet implemented.");
266			} else {
267			SendError(kAbortAnalysis, "SlaveBegin",
268			"The specified muon isolation %s is not defined.",
269			fMuonIsoType.Data());
270			return;
271			}
272	loizides	1.1
273	loizides	1.6	if (fMuonClassType.CompareTo("All") == 0)
274			fMuClassType = kAll;
275			else if (fMuonClassType.CompareTo("Global") == 0)
276			fMuClassType = kGlobal;
277			else if (fMuonClassType.CompareTo("Standalone") == 0)
278			fMuClassType = kSta;
279	loizides	1.24	else if (fMuonClassType.CompareTo("TrackerMuon") == 0)
280			fMuClassType = kTrackerMuon;
281			else if (fMuonClassType.CompareTo("CaloMuon") == 0)
282			fMuClassType = kCaloMuon;
283			else if (fMuonClassType.CompareTo("TrackerBased") == 0)
284			fMuClassType = kTrackerBased;
285	loizides	1.6	else {
286			SendError(kAbortAnalysis, "SlaveBegin",
287			"The specified muon class %s is not defined.",
288			fMuonClassType.Data());
289			return;
290			}
291	loizides	1.1	}