Mods/src/MuonIDMod.cc

// $Id: MuonIDMod.cc,v 1.26 2009/08/31 07:39:28 loizides 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("PrimaryVertexesBeamSpot"),
  fMuonIDType("Loose"),
  fMuonIsoType("TrackCaloSliding"),  
  fMuonClassType("Global"),  
  fTrackIsolationCut(3.0),
  fCaloIsolationCut(3.0),
  fCombIsolationCut(5.0),
  fMuonPtMin(10),
  fApplyD0Cut(kTRUE),
  fD0Cut(0.025),
  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!
    switch (fMuClassType) {
      case kAll:
        pass = kTRUE;
        if (mu->HasTrk())
          pt = mu->Pt();
        break;
      case kGlobal:
        pass = mu->HasGlobalTrk();
        if (pass) 
          pt = mu->GlobalTrk()->Pt();
        break;
      case kSta:
        pass = mu->HasStandaloneTrk();
        if (pass) 
          pt = mu->StandaloneTrk()->Pt();
        break;
      case kTrackerMuon:
        pass = mu->HasTrackerTrk() && mu->IsTrackerMuon();
        if (pass) 
          pt = mu->TrackerTrk()->Pt();
        break;
      case kCaloMuon:
        pass = mu->HasTrackerTrk() && mu->IsCaloMuon();
        if (pass) 
          pt = mu->TrackerTrk()->Pt();
        break;
      case kTrackerBased:
        pass = mu->HasTrackerTrk();
        if (pass) 
          pt = mu->TrackerTrk()->Pt();
        break;
      default:
        break;
    }

    if (!pass)
      continue;

    if (pt <= fMuonPtMin) 
      continue;

    Bool_t idpass = kFALSE;
    switch (fMuIDType) {
      case kLoose:
        idpass = fMuonTools->IsGood(mu, MuonTools::kTMOneStationLoose) &&
                 fMuonTools->IsGood(mu, MuonTools::kTM2DCompatibilityLoose);
        break;
      case kTight:
        idpass = fMuonTools->IsGood(mu, MuonTools::kTMOneStationTight) &&
                 fMuonTools->IsGood(mu, MuonTools::kTM2DCompatibilityTight);
        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()-10.0)*5.0/20.0 ||
              totalIso <= 0)
            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("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.27
Committed:	Thu Sep 3 13:31:35 2009 UTC (15 years, 8 months ago) by ceballos
Content type:	text/plain
Branch:	MAIN
CVS Tags:	Mit_014pre2, Mit_014pre1, Mit_013d, Mit_013c, Mit_013b, Mit_013a, Mit_013, Mit_013pre1, Mit_012i, Mit_012g, Mit_012f, Mit_012e, Mit_012d, Mit_012c, Mit_012b, Mit_012a, Mit_012, Mit_011a
Changes since 1.26:	+3 -4 lines
Log Message:	tighter muon id
#	User	Rev	Content
1	ceballos	1.27	// $Id: MuonIDMod.cc,v 1.26 2009/08/31 07:39:28 loizides 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	loizides	1.26	fVertexName("PrimaryVertexesBeamSpot"),
19	loizides	1.8	fMuonIDType("Loose"),
20			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	loizides	1.15	fD0Cut(0.025),
28	loizides	1.20	fReverseIsoCut(kFALSE),
29	ceballos	1.22	fReverseD0Cut(kFALSE),
30	loizides	1.15	fMuIDType(kIdUndef),
31			fMuIsoType(kIsoUndef),
32			fMuClassType(kClassUndef),
33	ceballos	1.14	fMuons(0),
34	loizides	1.15	fVertices(0),
35	loizides	1.20	fMuonTools(0)
36	loizides	1.1	{
37			// Constructor.
38			}
39
40			//--------------------------------------------------------------------------------------------------
41			void MuonIDMod::Process()
42			{
43			// Process entries of the tree.
44
45	loizides	1.21	LoadEventObject(fMuonBranchName, fMuons);
46	sixie	1.25	if (fApplyD0Cut) {
47			LoadEventObject(fVertexName, fVertices);
48			}
49	loizides	1.1
50	loizides	1.7	MuonOArr *CleanMuons = new MuonOArr;
51			CleanMuons->SetName(fCleanMuonsName);
52	loizides	1.1
53			for (UInt_t i=0; i<fMuons->GetEntries(); ++i) {
54	loizides	1.8	const Muon *mu = fMuons->At(i);
55	loizides	1.6
56			Bool_t pass = kFALSE;
57	loizides	1.19	Double_t pt = 0; // make sure pt is taken from the correct track!
58	loizides	1.6	switch (fMuClassType) {
59			case kAll:
60			pass = kTRUE;
61	loizides	1.19	if (mu->HasTrk())
62			pt = mu->Pt();
63	loizides	1.6	break;
64			case kGlobal:
65	loizides	1.19	pass = mu->HasGlobalTrk();
66	loizides	1.6	if (pass)
67			pt = mu->GlobalTrk()->Pt();
68			break;
69			case kSta:
70	loizides	1.19	pass = mu->HasStandaloneTrk();
71	loizides	1.6	if (pass)
72			pt = mu->StandaloneTrk()->Pt();
73			break;
74	loizides	1.24	case kTrackerMuon:
75			pass = mu->HasTrackerTrk() && mu->IsTrackerMuon();
76			if (pass)
77			pt = mu->TrackerTrk()->Pt();
78			break;
79			case kCaloMuon:
80			pass = mu->HasTrackerTrk() && mu->IsCaloMuon();
81			if (pass)
82			pt = mu->TrackerTrk()->Pt();
83			break;
84			case kTrackerBased:
85	loizides	1.19	pass = mu->HasTrackerTrk();
86	loizides	1.6	if (pass)
87			pt = mu->TrackerTrk()->Pt();
88			break;
89			default:
90			break;
91	ceballos	1.5	}
92	loizides	1.6
93			if (!pass)
94			continue;
95
96			if (pt <= fMuonPtMin)
97			continue;
98
99			Bool_t idpass = kFALSE;
100			switch (fMuIDType) {
101			case kLoose:
102			idpass = fMuonTools->IsGood(mu, MuonTools::kTMOneStationLoose) &&
103			fMuonTools->IsGood(mu, MuonTools::kTM2DCompatibilityLoose);
104			break;
105			case kTight:
106			idpass = fMuonTools->IsGood(mu, MuonTools::kTMOneStationTight) &&
107			fMuonTools->IsGood(mu, MuonTools::kTM2DCompatibilityTight);
108			break;
109	loizides	1.18	case kNoId:
110			idpass = kTRUE;
111			break;
112	loizides	1.6	default:
113			break;
114	ceballos	1.4	}
115	loizides	1.6
116			if (!idpass)
117			continue;
118
119	ceballos	1.22	Bool_t isocut = kFALSE;
120	loizides	1.6	switch (fMuIsoType) {
121			case kTrackCalo:
122	ceballos	1.22	isocut = (mu->IsoR03SumPt() < fTrackIsolationCut) &&
123	loizides	1.6	(mu->IsoR03EmEt() + mu->IsoR03HadEt() < fCaloIsolationCut);
124			break;
125			case kTrackCaloCombined:
126	ceballos	1.22	isocut = (1.0 * mu->IsoR03SumPt() + 1.0 * mu->IsoR03EmEt() +
127	loizides	1.6	1.0 * mu->IsoR03HadEt() < fCombIsolationCut);
128			break;
129			case kTrackCaloSliding:
130			{
131			Double_t totalIso = 1.0 * mu->IsoR03SumPt() +
132			1.0 * mu->IsoR03EmEt() +
133			1.0 * mu->IsoR03HadEt();
134	ceballos	1.27	if (totalIso < (mu->Pt()-10.0)*5.0/20.0 \|\|
135			totalIso <= 0)
136	ceballos	1.22	isocut = kTRUE;
137
138			if (fReverseIsoCut == kTRUE &&
139			isocut == kFALSE && totalIso < 10)
140			isocut = kTRUE;
141			else if(fReverseIsoCut == kTRUE)
142			isocut = kFALSE;
143			}
144	loizides	1.6	break;
145	loizides	1.13	case kNoIso:
146	ceballos	1.22	isocut = kTRUE;
147	loizides	1.13	break;
148	loizides	1.6	case kCustomIso:
149			default:
150			break;
151	ceballos	1.4	}
152	ceballos	1.3
153	ceballos	1.22	if (isocut == kFALSE)
154	loizides	1.6	continue;
155
156	loizides	1.21	if (fApplyD0Cut) {
157	ceballos	1.22	Bool_t d0cut = kFALSE;
158	loizides	1.21	const Track *mt = mu->BestTrk();
159			if (!mt)
160			continue;
161			Double_t d0_real = 1e30;
162			for(UInt_t i0 = 0; i0 < fVertices->GetEntries(); i0++) {
163			Double_t pD0 = mt->D0Corrected(*fVertices->At(i0));
164			if(TMath::Abs(pD0) < TMath::Abs(d0_real))
165			d0_real = TMath::Abs(pD0);
166			}
167	ceballos	1.22	if(d0_real < fD0Cut) d0cut = kTRUE;
168
169			if (fReverseD0Cut == kTRUE &&
170			d0cut == kFALSE && d0_real < 0.05)
171			d0cut = kTRUE;
172			else if(fReverseD0Cut == kTRUE)
173			d0cut = kFALSE;
174
175			if (d0cut == kFALSE)
176	loizides	1.21	continue;
177	ceballos	1.14	}
178
179	loizides	1.6	// add good muon
180			CleanMuons->Add(mu);
181	loizides	1.1	}
182
183	loizides	1.10	// sort according to pt
184			CleanMuons->Sort();
185
186	loizides	1.6	// add objects for other modules to use
187	loizides	1.7	AddObjThisEvt(CleanMuons);
188	loizides	1.1	}
189
190			//--------------------------------------------------------------------------------------------------
191			void MuonIDMod::SlaveBegin()
192			{
193			// Run startup code on the computer (slave) doing the actual analysis. Here,
194	loizides	1.6	// we just request the muon collection branch.
195
196	loizides	1.21	ReqEventObject(fMuonBranchName, fMuons, kTRUE);
197
198			if (fApplyD0Cut)
199			ReqEventObject(fVertexName, fVertices, kTRUE);
200	loizides	1.1
201	loizides	1.6	fMuonTools = new MuonTools;
202	loizides	1.1
203	loizides	1.6	if (fMuonIDType.CompareTo("Tight") == 0)
204			fMuIDType = kTight;
205			else if (fMuonIDType.CompareTo("Loose") == 0)
206			fMuIDType = kLoose;
207	loizides	1.12	else if (fMuonIDType.CompareTo("NoId") == 0)
208			fMuIDType = kNoId;
209	loizides	1.6	else if (fMuonIDType.CompareTo("Custom") == 0) {
210			fMuIDType = kCustomId;
211			SendError(kWarning, "SlaveBegin",
212			"Custom muon identification is not yet implemented.");
213			} else {
214			SendError(kAbortAnalysis, "SlaveBegin",
215			"The specified muon identification %s is not defined.",
216			fMuonIDType.Data());
217			return;
218			}
219	loizides	1.1
220	loizides	1.6	if (fMuonIsoType.CompareTo("TrackCalo") == 0)
221			fMuIsoType = kTrackCalo;
222			else if (fMuonIsoType.CompareTo("TrackCaloCombined") == 0)
223			fMuIsoType = kTrackCaloCombined;
224			else if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0)
225			fMuIsoType = kTrackCaloSliding;
226			else if (fMuonIsoType.CompareTo("NoIso") == 0)
227			fMuIsoType = kNoIso;
228			else if (fMuonIsoType.CompareTo("Custom") == 0) {
229			fMuIsoType = kCustomIso;
230			SendError(kWarning, "SlaveBegin",
231			"Custom muon isolation is not yet implemented.");
232			} else {
233			SendError(kAbortAnalysis, "SlaveBegin",
234			"The specified muon isolation %s is not defined.",
235			fMuonIsoType.Data());
236			return;
237			}
238	loizides	1.1
239	loizides	1.6	if (fMuonClassType.CompareTo("All") == 0)
240			fMuClassType = kAll;
241			else if (fMuonClassType.CompareTo("Global") == 0)
242			fMuClassType = kGlobal;
243			else if (fMuonClassType.CompareTo("Standalone") == 0)
244			fMuClassType = kSta;
245	loizides	1.24	else if (fMuonClassType.CompareTo("TrackerMuon") == 0)
246			fMuClassType = kTrackerMuon;
247			else if (fMuonClassType.CompareTo("CaloMuon") == 0)
248			fMuClassType = kCaloMuon;
249			else if (fMuonClassType.CompareTo("TrackerBased") == 0)
250			fMuClassType = kTrackerBased;
251	loizides	1.6	else {
252			SendError(kAbortAnalysis, "SlaveBegin",
253			"The specified muon class %s is not defined.",
254			fMuonClassType.Data());
255			return;
256			}
257	loizides	1.1	}