Mods/src/MuonIDMod.cc

// $Id: MuonIDMod.cc,v 1.25 2009/08/30 10:36:43 sixie 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 < (pt-10.0)*5.0/15.0 && pt <= 25) ||
              (totalIso < 5.0 && mu->Pt() > 25) ||
               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.26
Committed:	Mon Aug 31 07:39:28 2009 UTC (15 years, 8 months ago) by loizides
Content type:	text/plain
Branch:	MAIN
Changes since 1.25:	+2 -2 lines
Log Message:	Cleanup.
#	User	Rev	Content
1	loizides	1.26	// $Id: MuonIDMod.cc,v 1.25 2009/08/30 10:36:43 sixie 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	loizides	1.19	if ((totalIso < (pt-10.0)*5.0/15.0 && pt <= 25) \|\|
135	ceballos	1.9	(totalIso < 5.0 && mu->Pt() > 25) \|\|
136			totalIso <= 0)
137	ceballos	1.22	isocut = kTRUE;
138
139			if (fReverseIsoCut == kTRUE &&
140			isocut == kFALSE && totalIso < 10)
141			isocut = kTRUE;
142			else if(fReverseIsoCut == kTRUE)
143			isocut = kFALSE;
144			}
145	loizides	1.6	break;
146	loizides	1.13	case kNoIso:
147	ceballos	1.22	isocut = kTRUE;
148	loizides	1.13	break;
149	loizides	1.6	case kCustomIso:
150			default:
151			break;
152	ceballos	1.4	}
153	ceballos	1.3
154	ceballos	1.22	if (isocut == kFALSE)
155	loizides	1.6	continue;
156
157	loizides	1.21	if (fApplyD0Cut) {
158	ceballos	1.22	Bool_t d0cut = kFALSE;
159	loizides	1.21	const Track *mt = mu->BestTrk();
160			if (!mt)
161			continue;
162			Double_t d0_real = 1e30;
163			for(UInt_t i0 = 0; i0 < fVertices->GetEntries(); i0++) {
164			Double_t pD0 = mt->D0Corrected(*fVertices->At(i0));
165			if(TMath::Abs(pD0) < TMath::Abs(d0_real))
166			d0_real = TMath::Abs(pD0);
167			}
168	ceballos	1.22	if(d0_real < fD0Cut) d0cut = kTRUE;
169
170			if (fReverseD0Cut == kTRUE &&
171			d0cut == kFALSE && d0_real < 0.05)
172			d0cut = kTRUE;
173			else if(fReverseD0Cut == kTRUE)
174			d0cut = kFALSE;
175
176			if (d0cut == kFALSE)
177	loizides	1.21	continue;
178	ceballos	1.14	}
179
180	loizides	1.6	// add good muon
181			CleanMuons->Add(mu);
182	loizides	1.1	}
183
184	loizides	1.10	// sort according to pt
185			CleanMuons->Sort();
186
187	loizides	1.6	// add objects for other modules to use
188	loizides	1.7	AddObjThisEvt(CleanMuons);
189	loizides	1.1	}
190
191			//--------------------------------------------------------------------------------------------------
192			void MuonIDMod::SlaveBegin()
193			{
194			// Run startup code on the computer (slave) doing the actual analysis. Here,
195	loizides	1.6	// we just request the muon collection branch.
196
197	loizides	1.21	ReqEventObject(fMuonBranchName, fMuons, kTRUE);
198
199			if (fApplyD0Cut)
200			ReqEventObject(fVertexName, fVertices, kTRUE);
201	loizides	1.1
202	loizides	1.6	fMuonTools = new MuonTools;
203	loizides	1.1
204	loizides	1.6	if (fMuonIDType.CompareTo("Tight") == 0)
205			fMuIDType = kTight;
206			else if (fMuonIDType.CompareTo("Loose") == 0)
207			fMuIDType = kLoose;
208	loizides	1.12	else if (fMuonIDType.CompareTo("NoId") == 0)
209			fMuIDType = kNoId;
210	loizides	1.6	else if (fMuonIDType.CompareTo("Custom") == 0) {
211			fMuIDType = kCustomId;
212			SendError(kWarning, "SlaveBegin",
213			"Custom muon identification is not yet implemented.");
214			} else {
215			SendError(kAbortAnalysis, "SlaveBegin",
216			"The specified muon identification %s is not defined.",
217			fMuonIDType.Data());
218			return;
219			}
220	loizides	1.1
221	loizides	1.6	if (fMuonIsoType.CompareTo("TrackCalo") == 0)
222			fMuIsoType = kTrackCalo;
223			else if (fMuonIsoType.CompareTo("TrackCaloCombined") == 0)
224			fMuIsoType = kTrackCaloCombined;
225			else if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0)
226			fMuIsoType = kTrackCaloSliding;
227			else if (fMuonIsoType.CompareTo("NoIso") == 0)
228			fMuIsoType = kNoIso;
229			else if (fMuonIsoType.CompareTo("Custom") == 0) {
230			fMuIsoType = kCustomIso;
231			SendError(kWarning, "SlaveBegin",
232			"Custom muon isolation is not yet implemented.");
233			} else {
234			SendError(kAbortAnalysis, "SlaveBegin",
235			"The specified muon isolation %s is not defined.",
236			fMuonIsoType.Data());
237			return;
238			}
239	loizides	1.1
240	loizides	1.6	if (fMuonClassType.CompareTo("All") == 0)
241			fMuClassType = kAll;
242			else if (fMuonClassType.CompareTo("Global") == 0)
243			fMuClassType = kGlobal;
244			else if (fMuonClassType.CompareTo("Standalone") == 0)
245			fMuClassType = kSta;
246	loizides	1.24	else if (fMuonClassType.CompareTo("TrackerMuon") == 0)
247			fMuClassType = kTrackerMuon;
248			else if (fMuonClassType.CompareTo("CaloMuon") == 0)
249			fMuClassType = kCaloMuon;
250			else if (fMuonClassType.CompareTo("TrackerBased") == 0)
251			fMuClassType = kTrackerBased;
252	loizides	1.6	else {
253			SendError(kAbortAnalysis, "SlaveBegin",
254			"The specified muon class %s is not defined.",
255			fMuonClassType.Data());
256			return;
257			}
258	loizides	1.1	}