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.3 by ceballos, Wed Nov 5 14:06:09 2008 UTC vs.
Revision 1.35 by ceballos, Wed Oct 20 02:44:16 2010 UTC

#	Line 1 | Line 1
1	<	// $Id$
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 "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
#	Line 13 | Line 13 | ClassImp(mithep::MuonIDMod)
13		//--------------------------------------------------------------------------------------------------
14		MuonIDMod::MuonIDMod(const char *name, const char *title) :
15		BaseMod(name,title),
16	<	fPrintDebug(false),
17	<	fMuonName(Names::gkMuonBrn),
18	<	fCleanMuonsName(Names::gkCleanMuonsName),
19	<	fMuonIDType("Tight"),
20	<	fMuonIsoType("TrackCalo"),
21	<	fMuons(0),
16	>	fMuonBranchName(Names::gkMuonBrn),
17	>	fCleanMuonsName(ModNames::gkCleanMuonsName),
18	>	fVertexName(ModNames::gkGoodVertexesName),
19	>	fMuonIDType("Minimal"),
20	>	fMuonIsoType("TrackCaloSliding"),
21	>	fMuonClassType("Global"),
22		fTrackIsolationCut(3.0),
23		fCaloIsolationCut(3.0),
24	+	fCombIsolationCut(5.0),
25		fMuonPtMin(10),
26	<	fNEventsProcessed(0)
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		//--------------------------------------------------------------------------------------------------
31	–	void MuonIDMod::Begin()
32	–	{
33	–	// Run startup code on the client machine. For this module, we dont do
34	–	// anything here.
35	–	}
36	–
37	–	//--------------------------------------------------------------------------------------------------
42		void MuonIDMod::Process()
43		{
44		// Process entries of the tree.
45
46	<	fNEventsProcessed++;
47	<
48	<	if (fNEventsProcessed % 1000 == 0 || fPrintDebug) {
49	<	time_t systime;
50	<	systime = time(NULL);
47	<
48	<	cerr << endl << "MuonIDMod : Process Event " << fNEventsProcessed << "  Time: " << ctime(&systime) << endl;
49	<	}
50	<
51	<	//Get Muons
52	<	LoadBranch(fMuonName);
53	<	ObjArray<Muon> *CleanMuons = new ObjArray<Muon>;
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	<	Muon *mu = fMuons->At(i);
53	<
54	<	Double_t MuonClass = -1;
55	<	if (mu->GlobalTrk())
56	<	MuonClass = 0;
57	<	else if (mu->StandaloneTrk())
58	<	MuonClass = 1;
59	<	else if (mu->TrackerTrk())
60	<	MuonClass = 2;
61	<
62	<	bool allCuts = false;
63	<
64	<	if(MuonClass == 0) allCuts = true;
65	<
66	<	if(mu->IsoR03SumPt() >= fTrackIsolationCut) allCuts = false;
67	<
68	<	if(mu->IsoR03EmEt() +
69	<	mu->IsoR03HadEt() >= fCaloIsolationCut) allCuts = false;
70	<
71	<	if(mu->Pt() <= fMuonPtMin) allCuts = false;
72	<
73	<	if(allCuts) {
74	<	CleanMuons->Add(mu);
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 &&
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 &&
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 kMinimal:
159	+	idpass = mu->BestTrk() != 0 &&
160	+	mu->BestTrk()->NHits() > 10 &&
161	+	RChi2 < 10.0 &&
162	+	mu->Quality().Quality(MuonQuality::GlobalMuonPromptTight);
163	+	break;
164	+	case kNoId:
165	+	idpass = kTRUE;
166	+	break;
167	+	default:
168	+	break;
169	+	}
170	+
171	+	if (!idpass)
172	+	continue;
173	+
174	+	Bool_t isocut = kFALSE;
175	+	switch (fMuIsoType) {
176	+	case kTrackCalo:
177	+	isocut = (mu->IsoR03SumPt() < fTrackIsolationCut) &&
178	+	(mu->IsoR03EmEt() + mu->IsoR03HadEt() < fCaloIsolationCut);
179	+	break;
180	+	case kTrackCaloCombined:
181	+	isocut = (1.0 * mu->IsoR03SumPt() + 1.0 * mu->IsoR03EmEt() +
182	+	1.0 * mu->IsoR03HadEt() < fCombIsolationCut);
183	+	break;
184	+	case kTrackCaloSliding:
185	+	{
186	+	Double_t totalIso = 1.0 * mu->IsoR03SumPt() +
187	+	1.0 * mu->IsoR03EmEt() +
188	+	1.0 * mu->IsoR03HadEt();
189	+	if (totalIso < (mu->Pt()*0.15) )
190	+	isocut = kTRUE;
191	+
192	+	if     (fReverseIsoCut == kTRUE &&
193	+	isocut == kFALSE && totalIso < 10)
194	+	isocut = kTRUE;
195	+	else if(fReverseIsoCut == kTRUE)
196	+	isocut = kFALSE;
197	+	}
198	+	break;
199	+	case kNoIso:
200	+	isocut = kTRUE;
201	+	break;
202	+	case kCustomIso:
203	+	default:
204	+	break;
205	+	}
206	+
207	+	if (isocut == kFALSE)
208	+	continue;
209	+
210	+	if (fApplyD0Cut) {
211	+	fVertices = GetObjThisEvt<VertexOArr>(fVertexName);
212	+	Bool_t passD0cut = MuonTools::PassD0Cut(mu, fVertices, fD0Cut, fReverseD0Cut);
213	+	if (!passD0cut)
214	+	continue;
215	+	}
216	+
217	+	// add good muon
218	+	CleanMuons->Add(mu);
219		}
220
221	<	//Final Summary Debug Output
222	<	if ( fPrintDebug ) {
83	<	cerr << "Event Dump: " << fNEventsProcessed << endl;
84	<	cerr << "Muons" << endl;
85	<	for (UInt_t i = 0; i < CleanMuons->GetEntries(); i++) {
86	<	cerr << i << " " << CleanMuons->At(i)->Pt() << " " << CleanMuons->At(i)->Eta()
87	<	<< " " << CleanMuons->At(i)->Phi() << endl;
88	<	}
89	<	}
90	<
91	<	//Save Objects for Other Modules to use
92	<	AddObjThisEvt(CleanMuons, fCleanMuonsName.Data());
93	<	}
221	>	// sort according to pt
222	>	CleanMuons->Sort();
223
224	+	// add objects for other modules to use
225	+	AddObjThisEvt(CleanMuons);
226	+	}
227
228		//--------------------------------------------------------------------------------------------------
229		void MuonIDMod::SlaveBegin()
230		{
231		// Run startup code on the computer (slave) doing the actual analysis. Here,
232	<	// we typically initialize histograms and other analysis objects and request
101	<	// branches. For this module, we request a branch of the MitTree.
232	>	// we just request the muon collection branch.
233
234	<	ReqBranch(fMuonName,              fMuons);
104	<	}
234	>	ReqEventObject(fMuonBranchName, fMuons, kTRUE);
235
236	<	//--------------------------------------------------------------------------------------------------
107	<	void MuonIDMod::SlaveTerminate()
108	<	{
109	<	// Run finishing code on the computer (slave) that did the analysis. For this
110	<	// module, we dont do anything here.
236	>	fMuonTools = new MuonTools;
237
238	<	}
238	>	if (fMuonIDType.CompareTo("WMuId") == 0)
239	>	fMuIDType = kWMuId;
240	>	else if (fMuonIDType.CompareTo("ZMuId") == 0)
241	>	fMuIDType = kZMuId;
242	>	else if (fMuonIDType.CompareTo("Tight") == 0)
243	>	fMuIDType = kTight;
244	>	else if (fMuonIDType.CompareTo("Loose") == 0)
245	>	fMuIDType = kLoose;
246	>	else if (fMuonIDType.CompareTo("Minimal") == 0)
247	>	fMuIDType = kMinimal;
248	>	else if (fMuonIDType.CompareTo("NoId") == 0)
249	>	fMuIDType = kNoId;
250	>	else if (fMuonIDType.CompareTo("Custom") == 0) {
251	>	fMuIDType = kCustomId;
252	>	SendError(kWarning, "SlaveBegin",
253	>	"Custom muon identification is not yet implemented.");
254	>	} else {
255	>	SendError(kAbortAnalysis, "SlaveBegin",
256	>	"The specified muon identification %s is not defined.",
257	>	fMuonIDType.Data());
258	>	return;
259	>	}
260
261	<	//--------------------------------------------------------------------------------------------------
262	<	void MuonIDMod::Terminate()
263	<	{
264	<	// Run finishing code on the client computer. For this module, we dont do
265	<	// anything here.
261	>	if (fMuonIsoType.CompareTo("TrackCalo") == 0)
262	>	fMuIsoType = kTrackCalo;
263	>	else if (fMuonIsoType.CompareTo("TrackCaloCombined") == 0)
264	>	fMuIsoType = kTrackCaloCombined;
265	>	else if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0)
266	>	fMuIsoType = kTrackCaloSliding;
267	>	else if (fMuonIsoType.CompareTo("NoIso") == 0)
268	>	fMuIsoType = kNoIso;
269	>	else if (fMuonIsoType.CompareTo("Custom") == 0) {
270	>	fMuIsoType = kCustomIso;
271	>	SendError(kWarning, "SlaveBegin",
272	>	"Custom muon isolation is not yet implemented.");
273	>	} else {
274	>	SendError(kAbortAnalysis, "SlaveBegin",
275	>	"The specified muon isolation %s is not defined.",
276	>	fMuonIsoType.Data());
277	>	return;
278	>	}
279	>
280	>	if (fMuonClassType.CompareTo("All") == 0)
281	>	fMuClassType = kAll;
282	>	else if (fMuonClassType.CompareTo("Global") == 0)
283	>	fMuClassType = kGlobal;
284	>	else if (fMuonClassType.CompareTo("Standalone") == 0)
285	>	fMuClassType = kSta;
286	>	else if (fMuonClassType.CompareTo("TrackerMuon") == 0)
287	>	fMuClassType = kTrackerMuon;
288	>	else if (fMuonClassType.CompareTo("CaloMuon") == 0)
289	>	fMuClassType = kCaloMuon;
290	>	else if (fMuonClassType.CompareTo("TrackerBased") == 0)
291	>	fMuClassType = kTrackerBased;
292	>	else {
293	>	SendError(kAbortAnalysis, "SlaveBegin",
294	>	"The specified muon class %s is not defined.",
295	>	fMuonClassType.Data());
296	>	return;
297	>	}
298		}

Diff Legend

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