Mods/src/MuonIDMod.cc

        // $Id: MuonIDMod.cc,v 1.35 2010/10/20 02:44:16 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(ModNames::gkGoodVertexesName),
  fMuonIDType("WWMuId"),
  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);

  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();
        if (pass && mu->TrackerTrk()) {
          pt  = mu->TrackerTrk()->Pt();
          eta = TMath::Abs(mu->TrackerTrk()->Eta());
        }
        else {
          pt  = mu->Pt();
          eta = TMath::Abs(mu->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;

    Double_t RChi2 = 0.0;
    if     (mu->HasGlobalTrk()) {
      RChi2 = mu->GlobalTrk()->Chi2()/mu->GlobalTrk()->Ndof();
    }
    else if(mu->BestTrk() != 0){
      RChi2 = mu->BestTrk()->Chi2()/mu->BestTrk()->Ndof();
    }
    Bool_t idpass = kFALSE;
    switch (fMuIDType) {
      case kWMuId:
        idpass = mu->BestTrk() != 0 &&
                 mu->BestTrk()->NHits() > 10 &&
                 RChi2 < 10.0 &&
                (mu->NSegments() > 1 || mu->NMatches() > 1) &&
                 mu->BestTrk()->NPixelHits() > 0 &&
                 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
        break;
      case kZMuId:
        idpass = mu->BestTrk() != 0 &&
                 mu->BestTrk()->NHits() > 10 &&
                (mu->NSegments() > 1 || mu->NMatches() > 1) &&
                 mu->BestTrk()->NPixelHits() > 0 &&
                 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
        break;
      case kLoose:
        idpass = mu->BestTrk() != 0 &&
                 mu->Quality().Quality(MuonQuality::TMOneStationLoose) &&
                 mu->Quality().Quality(MuonQuality::TM2DCompatibilityLoose) &&
                 mu->BestTrk()->NHits() > 10 &&
                 RChi2 < 10.0 &&
                 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
        break;
      case kTight:
        idpass = mu->BestTrk() !=  0 &&
                 mu->Quality().Quality(MuonQuality::TMOneStationTight) &&
                 mu->Quality().Quality(MuonQuality::TM2DCompatibilityTight) &&
                 mu->BestTrk()->NHits() > 10 &&
                 RChi2 < 10.0 &&
                 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
        break;
      case kWWMuId:
        idpass = mu->BestTrk() != 0 &&
                 mu->BestTrk()->NHits() > 10 &&
                 RChi2 < 10.0 &&
                (mu->NSegments() > 1 || mu->NMatches() > 1) &&
                 mu->BestTrk()->NPixelHits() > 0 &&
                 mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight) &&
                 mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1;
        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) {
      fVertices = GetObjThisEvt<VertexOArr>(fVertexName);
      Bool_t passD0cut = MuonTools::PassD0Cut(mu, fVertices, fD0Cut, fReverseD0Cut);
      if (!passD0cut)
        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);

  fMuonTools = new MuonTools;

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