ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/MitPhysics/Mods/src/MuonIDMod.cc

Comparing UserCode/MitPhysics/Mods/src/MuonIDMod.cc (file contents):
Revision 1.1 by loizides, Wed Oct 15 06:05:00 2008 UTC vs.
Revision 1.20 by loizides, Mon May 11 08:02:43 2009 UTC

#	Line 1 | Line 1
1		// $Id$
2
3		#include "MitPhysics/Mods/interface/MuonIDMod.h"
4	–	#include "MitAna/DataTree/interface/Names.h"
5	–	#include "MitAna/DataCont/interface/ObjArray.h"
6	–	#include "MitPhysics/Utils/interface/IsolationTools.h"
4		#include "MitCommon/MathTools/interface/MathUtils.h"
5	+	#include "MitPhysics/Init/interface/ModNames.h"
6
7		using namespace mithep;
8
#	Line 13 | Line 11 | ClassImp(mithep::MuonIDMod)
11		//--------------------------------------------------------------------------------------------------
12		MuonIDMod::MuonIDMod(const char *name, const char *title) :
13		BaseMod(name,title),
14	<	fPrintDebug(false),
15	<	fMuonName(Names::gkMuonBrn),
16	<	fCleanMuonsName(Names::gkCleanMuonsName),
17	<	fMuonIDType("Tight"),
18	<	fMuonIsoType("TrackCalo"),
14	>	fMuonBranchName(Names::gkMuonBrn),
15	>	fCleanMuonsName(ModNames::gkCleanMuonsName),
16	>	fVertexName("PrimaryVertexesBeamSpot"),
17	>	fMuonIDType("Loose"),
18	>	fMuonIsoType("TrackCaloSliding"),
19	>	fMuonClassType("Global"),
20	>	fTrackIsolationCut(3.0),
21	>	fCaloIsolationCut(3.0),
22	>	fCombIsolationCut(5.0),
23	>	fMuonPtMin(10),
24	>	fD0Cut(0.025),
25	>	fReverseIsoCut(kFALSE),
26	>	fMuIDType(kIdUndef),
27	>	fMuIsoType(kIsoUndef),
28	>	fMuClassType(kClassUndef),
29		fMuons(0),
30	<	fTrackIsolationCut(5.0),
31	<	fCaloIsolationCut(5.0),
24	<	fNEventsProcessed(0)
30	>	fVertices(0),
31	>	fMuonTools(0)
32		{
33		// Constructor.
34		}
35
36		//--------------------------------------------------------------------------------------------------
30	–	void MuonIDMod::Begin()
31	–	{
32	–	// Run startup code on the client machine. For this module, we dont do
33	–	// anything here.
34	–	}
35	–
36	–	//--------------------------------------------------------------------------------------------------
37		void MuonIDMod::Process()
38		{
39		// Process entries of the tree.
40
41	<	fNEventsProcessed++;
42	<
43	<	if (fNEventsProcessed % 1000 == 0 || fPrintDebug) {
44	<	time_t systime;
45	<	systime = time(NULL);
46	<
47	<	cerr << endl << "MuonIDMod : Process Event " << fNEventsProcessed << "  Time: " << ctime(&systime) << endl;
48	<	}
49	<
50	<	//Get Muons
51	<	LoadBranch(fMuonName);
52	<	ObjArray<Muon> *CleanMuons = new ObjArray<Muon>;
41	>	LoadBranch(fMuonBranchName);
42	>	LoadBranch(fVertexName);
43	>
44	>	MuonOArr *CleanMuons = new MuonOArr;
45	>	CleanMuons->SetName(fCleanMuonsName);
46	>
47		for (UInt_t i=0; i<fMuons->GetEntries(); ++i) {
48	<	Muon *mu = fMuons->At(i);
55	<
56	<	Double_t MuonClass = -1;
57	<	if (mu->GlobalTrk())
58	<	MuonClass = 0;
59	<	else if (mu->StandaloneTrk())
60	<	MuonClass = 1;
61	<	else if (mu->TrackerTrk())
62	<	MuonClass = 2;
63	<
64	<	//These cuts are from the 1.6.X analysis. I'm waiting for Phil to finalize his Muon ID class
65	<	const int nCuts = 4;
66	<	double cutValue[nCuts] = {0.1, 3.0, 3.0, 1.5 };
67	<	bool passCut[nCuts] = {false, false, false, false};
68	<	double muonD0 = fabs(mu->BestTrk()->D0());
69	<	if(muonD0 < cutValue[0] &&  MuonClass == 0 )
70	<	passCut[0] = true;
71	<	if(mu->IsoR03SumPt() < cutValue[1]) passCut[1] = true;
72	<	if(mu->IsoR03EmEt() +
73	<	mu->IsoR03HadEt() < cutValue[2]) passCut[2] = true;
74	<	if(mu->Pt() > 10)
75	<	passCut[3] = true;
76	<
77	<	// Final decision
78	<	bool allCuts = true;
79	<	for(int c=0; c<nCuts; c++) {
80	<	allCuts = allCuts & passCut[c];
81	<	}
82	<
83	<	if ( allCuts
84	<	&& abs(mu->Eta()) < 2.5
85	<	) {
86	<	CleanMuons->Add(mu);
87	<	}
48	>	const Muon *mu = fMuons->At(i);
49
50	+	Bool_t pass = kFALSE;
51	+	Double_t pt = 0; // make sure pt is taken from the correct track!
52	+	switch (fMuClassType) {
53	+	case kAll:
54	+	pass = kTRUE;
55	+	if (mu->HasTrk())
56	+	pt = mu->Pt();
57	+	break;
58	+	case kGlobal:
59	+	pass = mu->HasGlobalTrk();
60	+	if (pass)
61	+	pt = mu->GlobalTrk()->Pt();
62	+	break;
63	+	case kSta:
64	+	pass = mu->HasStandaloneTrk();
65	+	if (pass)
66	+	pt = mu->StandaloneTrk()->Pt();
67	+	break;
68	+	case kTrackerOnly:
69	+	pass = mu->HasTrackerTrk();
70	+	if (pass)
71	+	pt = mu->TrackerTrk()->Pt();
72	+	break;
73	+	default:
74	+	break;
75	+	}
76	+
77	+	if (!pass)
78	+	continue;
79	+
80	+	if (pt <= fMuonPtMin)
81	+	continue;
82	+
83	+	Bool_t idpass = kFALSE;
84	+	switch (fMuIDType) {
85	+	case kLoose:
86	+	idpass = fMuonTools->IsGood(mu, MuonTools::kTMOneStationLoose) &&
87	+	fMuonTools->IsGood(mu, MuonTools::kTM2DCompatibilityLoose);
88	+	break;
89	+	case kTight:
90	+	idpass = fMuonTools->IsGood(mu, MuonTools::kTMOneStationTight) &&
91	+	fMuonTools->IsGood(mu, MuonTools::kTM2DCompatibilityTight);
92	+	break;
93	+	case kNoId:
94	+	idpass = kTRUE;
95	+	break;
96	+	default:
97	+	break;
98	+	}
99	+
100	+	if (!idpass)
101	+	continue;
102	+
103	+	Bool_t isopass = kFALSE;
104	+	switch (fMuIsoType) {
105	+	case kTrackCalo:
106	+	isopass = (mu->IsoR03SumPt() < fTrackIsolationCut) &&
107	+	(mu->IsoR03EmEt() + mu->IsoR03HadEt() < fCaloIsolationCut);
108	+	break;
109	+	case kTrackCaloCombined:
110	+	isopass = (1.0 * mu->IsoR03SumPt() + 1.0 * mu->IsoR03EmEt() +
111	+	1.0 * mu->IsoR03HadEt() < fCombIsolationCut);
112	+	break;
113	+	case kTrackCaloSliding:
114	+	{
115	+	Double_t totalIso = 1.0 * mu->IsoR03SumPt() +
116	+	1.0 * mu->IsoR03EmEt() +
117	+	1.0 * mu->IsoR03HadEt();
118	+	if ((totalIso < (pt-10.0)*5.0/15.0 && pt <= 25) ||
119	+	(totalIso < 5.0 && mu->Pt() > 25) ||
120	+	totalIso <= 0)
121	+	isopass = kTRUE;
122	+	}
123	+	break;
124	+	case kNoIso:
125	+	isopass = kTRUE;
126	+	break;
127	+	case kCustomIso:
128	+	default:
129	+	break;
130	+	}
131	+
132	+	if ((isopass == kFALSE && fReverseIsoCut == kFALSE) ||
133	+	(isopass == kTRUE  && fReverseIsoCut == kTRUE))
134	+	continue;
135	+
136	+	// d0 cut only for muons that have a track
137	+	const Track *mt = mu->BestTrk();
138	+	if (!mt)
139	+	continue;
140	+	Double_t d0_real = 1e30;
141	+	for(UInt_t i0 = 0; i0 < fVertices->GetEntries(); i0++) {
142	+	Double_t pD0 = mt->D0Corrected(*fVertices->At(i0));
143	+	if(TMath::Abs(pD0) < TMath::Abs(d0_real))
144	+	d0_real = TMath::Abs(pD0);
145	+	}
146	+	if(d0_real >= fD0Cut)
147	+	continue;
148
149	+	// add good muon
150	+	CleanMuons->Add(mu);
151		}
152
153	<	//Final Summary Debug Output
154	<	if ( fPrintDebug ) {
94	<	cerr << "Event Dump: " << fNEventsProcessed << endl;
95	<	cerr << "Muons" << endl;
96	<	for (UInt_t i = 0; i < CleanMuons->GetEntries(); i++) {
97	<	cerr << i << " " << CleanMuons->At(i)->Pt() << " " << CleanMuons->At(i)->Eta()
98	<	<< " " << CleanMuons->At(i)->Phi() << endl;
99	<	}
100	<	}
101	<
102	<	//Save Objects for Other Modules to use
103	<	AddObjThisEvt(CleanMuons, fCleanMuonsName.Data());
104	<	}
153	>	// sort according to pt
154	>	CleanMuons->Sort();
155
156	+	// add objects for other modules to use
157	+	AddObjThisEvt(CleanMuons);
158	+	}
159
160		//--------------------------------------------------------------------------------------------------
161		void MuonIDMod::SlaveBegin()
162		{
163		// Run startup code on the computer (slave) doing the actual analysis. Here,
164	<	// we typically initialize histograms and other analysis objects and request
112	<	// branches. For this module, we request a branch of the MitTree.
164	>	// we just request the muon collection branch.
165
166	<	ReqBranch(fMuonName,              fMuons);
167	<	}
166	>	ReqBranch(fMuonBranchName, fMuons);
167	>	ReqBranch(fVertexName, fVertices);
168
169	<	//--------------------------------------------------------------------------------------------------
118	<	void MuonIDMod::SlaveTerminate()
119	<	{
120	<	// Run finishing code on the computer (slave) that did the analysis. For this
121	<	// module, we dont do anything here.
169	>	fMuonTools = new MuonTools;
170
171	<	}
171	>	if (fMuonIDType.CompareTo("Tight") == 0)
172	>	fMuIDType = kTight;
173	>	else if (fMuonIDType.CompareTo("Loose") == 0)
174	>	fMuIDType = kLoose;
175	>	else if (fMuonIDType.CompareTo("NoId") == 0)
176	>	fMuIDType = kNoId;
177	>	else if (fMuonIDType.CompareTo("Custom") == 0) {
178	>	fMuIDType = kCustomId;
179	>	SendError(kWarning, "SlaveBegin",
180	>	"Custom muon identification is not yet implemented.");
181	>	} else {
182	>	SendError(kAbortAnalysis, "SlaveBegin",
183	>	"The specified muon identification %s is not defined.",
184	>	fMuonIDType.Data());
185	>	return;
186	>	}
187
188	<	//--------------------------------------------------------------------------------------------------
189	<	void MuonIDMod::Terminate()
190	<	{
191	<	// Run finishing code on the client computer. For this module, we dont do
192	<	// anything here.
188	>	if (fMuonIsoType.CompareTo("TrackCalo") == 0)
189	>	fMuIsoType = kTrackCalo;
190	>	else if (fMuonIsoType.CompareTo("TrackCaloCombined") == 0)
191	>	fMuIsoType = kTrackCaloCombined;
192	>	else if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0)
193	>	fMuIsoType = kTrackCaloSliding;
194	>	else if (fMuonIsoType.CompareTo("NoIso") == 0)
195	>	fMuIsoType = kNoIso;
196	>	else if (fMuonIsoType.CompareTo("Custom") == 0) {
197	>	fMuIsoType = kCustomIso;
198	>	SendError(kWarning, "SlaveBegin",
199	>	"Custom muon isolation is not yet implemented.");
200	>	} else {
201	>	SendError(kAbortAnalysis, "SlaveBegin",
202	>	"The specified muon isolation %s is not defined.",
203	>	fMuonIsoType.Data());
204	>	return;
205	>	}
206	>
207	>	if (fMuonClassType.CompareTo("All") == 0)
208	>	fMuClassType = kAll;
209	>	else if (fMuonClassType.CompareTo("Global") == 0)
210	>	fMuClassType = kGlobal;
211	>	else if (fMuonClassType.CompareTo("Standalone") == 0)
212	>	fMuClassType = kSta;
213	>	else if (fMuonClassType.CompareTo("TrackerOnly") == 0)
214	>	fMuClassType = kTrackerOnly;
215	>	else {
216	>	SendError(kAbortAnalysis, "SlaveBegin",
217	>	"The specified muon class %s is not defined.",
218	>	fMuonClassType.Data());
219	>	return;
220	>	}
221		}

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines
>	Changed lines