1 |
+ |
//-------------------------------------------------------------------------------------------------- |
2 |
+ |
// $Id$ |
3 |
+ |
// |
4 |
+ |
// MuonTools |
5 |
+ |
// |
6 |
+ |
// This class allows you to classify a given muon according to defined criteria. For this purpose |
7 |
+ |
// is loads histograms from two ROOT files, specified in the constructor. The main function then |
8 |
+ |
// is "IsGood(*muon, selection)" which returns true if the given muon fulfills the selection |
9 |
+ |
// criteria. |
10 |
+ |
// |
11 |
+ |
// Logically, the code has been put together by Phil who took most of the ideas from |
12 |
+ |
// http://cmslxr.fnal.gov/lxr/source/RecoMuon/MuonIdentification/ |
13 |
+ |
// http://cmslxr.fnal.gov/lxr/source/DataFormats/MuonReco/ |
14 |
+ |
// |
15 |
+ |
// Authors: P.Harris, C.Loizides |
16 |
+ |
//-------------------------------------------------------------------------------------------------- |
17 |
+ |
|
18 |
|
#ifndef MITPHYSICS_UTIL_MUONTOOLS_H |
19 |
|
#define MITPHYSICS_UTIL_MUONTOOLS_H |
3 |
– |
#include <vector> |
4 |
– |
#include <iostream> |
20 |
|
|
21 |
< |
#include "TH2D.h" |
7 |
< |
#include "TFile.h" |
21 |
> |
#include "MitAna/DataTree/interface/Muon.h" |
22 |
|
#include "MitAna/DataTree/interface/Collections.h" |
23 |
|
#include "MitCommon/MathTools/interface/MathUtils.h" |
24 |
+ |
#include "TH2D.h" |
25 |
|
|
26 |
|
namespace mithep { |
27 |
|
class MuonTools { |
28 |
< |
public: |
29 |
< |
MuonTools(); |
30 |
< |
~MuonTools(); |
31 |
< |
double getCaloCompatability(mithep::Muon* iMuon,bool iEMSpecial, bool iCorrectedHCAL); |
32 |
< |
enum selection { |
33 |
< |
AllArbitrated, |
34 |
< |
PromptTight, |
35 |
< |
TMLastStationLoose, |
36 |
< |
TMLastStationTight, |
37 |
< |
TMOneStationLoose, |
38 |
< |
TMOneStationTight, |
39 |
< |
TM2DCompatibilityLoose, |
40 |
< |
TM2DCompatibilityTight |
41 |
< |
}; |
42 |
< |
static inline double sigWeight(double iVal0,double iVal1) { |
43 |
< |
if(iVal1 < 1.) return 1; |
44 |
< |
if(iVal0 < 3. && iVal1 > 3.) { |
45 |
< |
double lVal = TMath::Max(iVal0,1.); |
46 |
< |
return 1/TMath::Power(lVal,0.25); |
47 |
< |
} |
48 |
< |
double lVal = TMath::Max(iVal1,1.); |
49 |
< |
return 1/TMath::Power(lVal,0.25); |
50 |
< |
} |
51 |
< |
static inline int lastHit(mithep::Muon *iMuon) { |
52 |
< |
int lId = -1; |
53 |
< |
for(int i0 = 0; i0 < 8; i0++) { |
54 |
< |
if(iMuon->GetDX(i0) < 99999 || iMuon->GetDY(i0) < 99999) lId = i0; |
55 |
< |
} |
56 |
< |
return lId; |
57 |
< |
} |
58 |
< |
static inline int maxChamberId(mithep::Muon* iMuon,double iMaxD,double iMaxP) { |
59 |
< |
int lId = -1; |
60 |
< |
for(int i0 = 0; i0 < 8; i0++) { |
61 |
< |
if(iMuon->GetTrackDist(i0) < iMaxD && |
62 |
< |
iMuon->GetTrackDist(i0)/iMuon->GetTrackDistErr(i0) < iMaxP) { |
63 |
< |
if(iMuon->GetDX(i0) >= 999999) {lId = i0;} else {lId = -1;} |
64 |
< |
} |
65 |
< |
} |
66 |
< |
return lId; |
67 |
< |
} |
68 |
< |
static inline int lastStation(mithep::Muon* iMuon,double iMaxD,double iMaxP) { |
69 |
< |
int lId = -1; |
70 |
< |
for(int i0 = 0; i0 < 8; i0++) { |
71 |
< |
if((lId % 4) > (i0 % 4)) continue; |
72 |
< |
if(iMuon->GetTrackDist(i0) < iMaxD && |
73 |
< |
iMuon->GetTrackDist(i0)/iMuon->GetTrackDistErr(i0) < iMaxP) lId = i0; |
74 |
< |
} |
75 |
< |
return lId; |
76 |
< |
} |
77 |
< |
static inline int lastStation(mithep::Muon* iMuon,int iMax=8) { |
78 |
< |
int lId = -1; if(iMax > 8) iMax = 8; |
79 |
< |
for(int i0 = 0; i0 < iMax; i0++) {if(iMuon->StationBit(i0) && ((lId % 4) < (i0 % 4)))lId = i0;} |
65 |
< |
return lId; |
66 |
< |
} |
67 |
< |
static inline bool overflow(TH2D *iHist,double lVal0,double lVal1) { |
68 |
< |
if(iHist == 0) return true; |
69 |
< |
if(iHist ->GetXaxis()->FindBin(lVal0) == 0 || |
70 |
< |
iHist ->GetXaxis()->FindBin(lVal0) > iHist->GetNbinsX() || |
71 |
< |
iHist ->GetYaxis()->FindBin(lVal0) == 0 || |
72 |
< |
iHist ->GetYaxis()->FindBin(lVal0) > iHist->GetNbinsY()) return true; |
73 |
< |
return false; |
74 |
< |
} |
75 |
< |
static inline MCParticle* etaPhiMatch(MCParticleCol* iMCParts,Muon *iMuon) { |
76 |
< |
mithep::MCParticle *lMC = 0; double lDR = 1.; double lPt = -1.; |
77 |
< |
for(unsigned int i0 = 0; i0 < iMCParts->GetEntries(); i0++) { |
78 |
< |
mithep::MCParticle* pMC = iMCParts->At(i0); |
79 |
< |
if(pMC->Status() != 1) continue; |
80 |
< |
double pDR = mithep::MathUtils::DeltaR(pMC->Mom(),iMuon->Mom()); |
81 |
< |
if(pDR > lDR && pDR > 0.01) continue; |
82 |
< |
if(fabs(pMC->Pt()) < lPt) continue; |
83 |
< |
lDR = pDR; |
84 |
< |
lMC = pMC; |
85 |
< |
lPt = pMC->Pt(); |
86 |
< |
} |
87 |
< |
return lMC; |
88 |
< |
} |
89 |
< |
static inline Muon* match(MuonCol* iMuons,Track *iTrack) { |
90 |
< |
mithep::Muon *lMuon = 0; double lDR = 1.; double lPt = -1.; |
91 |
< |
for(unsigned int i0 = 0; i0 < iMuons->GetEntries(); i0++) { |
92 |
< |
mithep::Muon* pMuon = iMuons->At(i0); |
93 |
< |
double pDR = mithep::MathUtils::DeltaR(pMuon->Mom(),iTrack->Mom4(0.109)); |
94 |
< |
if(pDR > lDR && pDR > 0.01) continue; |
95 |
< |
if(fabs(pMuon->Pt()) < lPt) continue; |
96 |
< |
lDR = pDR; |
97 |
< |
lMuon = pMuon; |
98 |
< |
lPt = pMuon->Pt(); |
99 |
< |
} |
100 |
< |
return lMuon; |
101 |
< |
} |
102 |
< |
static inline int NSegments(mithep::Muon *iMuon) { |
103 |
< |
int lNSegs = 0; |
104 |
< |
for(int i0 = 0; i0 < 8; i0++) { |
105 |
< |
if(iMuon->StationBit(i0)) lNSegs++; |
106 |
< |
} |
107 |
< |
return lNSegs; |
108 |
< |
} |
109 |
< |
static inline bool promptTight(int iId,Muon *iMuon) { |
110 |
< |
const mithep::Track *lTrack = 0; |
111 |
< |
if(iId == 0 ) lTrack = iMuon->GlobalTrk(); |
112 |
< |
if(iId == 1 || (iId == -1 && lTrack == 0)) lTrack = iMuon->TrackerTrk(); |
113 |
< |
if(iId == 2 || (iId == -1 && lTrack == 0)) lTrack = iMuon->StandaloneTrk(); |
114 |
< |
if(lTrack == 0) return 0; |
115 |
< |
if(lTrack->NHits() < 11) return false; |
116 |
< |
if(lTrack->Chi2()/lTrack->Ndof() > 10.) return false; |
117 |
< |
if(lTrack->D0() > 0.2) return false; |
118 |
< |
if((lastHit(iMuon)% 4) == 0) return false; |
119 |
< |
return true; |
120 |
< |
} |
121 |
< |
static inline bool TMLastStation(Muon *iMuon,double iDYMin = 3., double iPYMin = 3.,double iDXMin = 3., double iPXMin = 3.,int iN = 2) { |
122 |
< |
if(iMuon->NSegments() < iN) return false; //2Last 1One |
123 |
< |
int lLast = lastStation(iMuon,-3.,-3.); //Last Required Station |
124 |
< |
if(lLast < 0) return false; |
125 |
< |
if(iMuon->GetDX(lLast) > 9999.) return false; |
126 |
< |
lLast = lastStation(iMuon); //No Requirements Imply StationMask (Segment and Track Arbitrated) |
127 |
< |
if(lLast < 0) return false; |
128 |
< |
if(!(fabs(iMuon->GetDX(lLast)) < iDXMin || |
129 |
< |
fabs(iMuon->GetPullX(lLast)) < iPXMin)) return false; |
130 |
< |
if(lLast == 3) lLast = lastStation(iMuon,3); |
131 |
< |
if(lLast < 0) return false; |
132 |
< |
if(!(fabs(iMuon->GetDY(lLast)) < iDYMin || //Old Code it was 9999 |
133 |
< |
fabs(iMuon->GetPullY(lLast)) < iPYMin)) return false; |
134 |
< |
return true; |
135 |
< |
} |
136 |
< |
static inline bool TMOneStation(Muon *iMuon,double iDYMin = 3., double iPYMin = 3.,double iDXMin = 3., double iPXMin = 3.,int iN = 1) { |
137 |
< |
if(iMuon->NSegments() < iN) return false; //2Last 1One |
138 |
< |
bool pGoodX = false; bool pBadY = false; |
139 |
< |
for(int i0 = 0; i0 < 8; i0++) { |
140 |
< |
if((fabs(iMuon->GetDX(i0)) < iDXMin || |
141 |
< |
fabs(iMuon->GetPullX(i0)) < iPXMin)) pGoodX = true; |
142 |
< |
if(pGoodX && |
143 |
< |
(fabs(iMuon->GetDY(i0)) < iDYMin || |
144 |
< |
fabs(iMuon->GetPullY(i0)) < iPYMin)) return true; |
145 |
< |
if(fabs(iMuon->GetDY(i0)) < 999999) pBadY = true; |
146 |
< |
if(i0 == 3 && pGoodX && !pBadY) return true; |
147 |
< |
} |
148 |
< |
return false; |
149 |
< |
} |
150 |
< |
static inline double getCaloCompatabilitySlowly(mithep::Muon* iMuon,bool iEMSpecial, bool iCorrectedHCAL) { |
151 |
< |
std::cout << "Warning Loading Compatability Slowly" << std::endl; |
152 |
< |
TFile* lPion_templates = new TFile("MitCommon/MuonTools/data/PionCaloTemplate.root","READ"); |
153 |
< |
TFile* lMuon_templates = new TFile("MitCommon/MuonTools/data/MuonCaloTemplate.root","READ"); |
154 |
< |
TH2D* lpion_em_etaEmi = (TH2D*) lPion_templates->Get("em_etaEmi"); |
155 |
< |
TH2D* lpion_had_etaEmi = (TH2D*) lPion_templates->Get("had_etaEmi"); |
156 |
< |
TH2D* lpion_had_etaTmi = (TH2D*) lPion_templates->Get("had_etaTmi"); |
157 |
< |
TH2D* lpion_em_etaB = (TH2D*) lPion_templates->Get("em_etaB") ; |
158 |
< |
TH2D* lpion_had_etaB = (TH2D*) lPion_templates->Get("had_etaB") ; |
159 |
< |
TH2D* lpion_ho_etaB = (TH2D*) lPion_templates->Get("ho_etaB") ; |
160 |
< |
TH2D* lpion_had_etaTpl = (TH2D*) lPion_templates->Get("had_etaTpl"); |
161 |
< |
TH2D* lpion_em_etaEpl = (TH2D*) lPion_templates->Get("em_etaEpl") ; |
162 |
< |
TH2D* lpion_had_etaEpl = (TH2D*) lPion_templates->Get("had_etaEpl"); |
163 |
< |
TH2D* lmuon_em_etaEmi = (TH2D*) lMuon_templates->Get("em_etaEmi") ; |
164 |
< |
TH2D* lmuon_had_etaEmi = (TH2D*) lMuon_templates->Get("had_etaEmi"); |
165 |
< |
TH2D* lmuon_had_etaTmi = (TH2D*) lMuon_templates->Get("had_etaTmi"); |
166 |
< |
TH2D* lmuon_em_etaB = (TH2D*) lMuon_templates->Get("em_etaB") ; |
167 |
< |
TH2D* lmuon_had_etaB = (TH2D*) lMuon_templates->Get("had_etaB") ; |
168 |
< |
TH2D* lmuon_ho_etaB = (TH2D*) lMuon_templates->Get("ho_etaB") ; |
169 |
< |
TH2D* lmuon_had_etaTpl = (TH2D*) lMuon_templates->Get("had_etaTpl"); |
170 |
< |
TH2D* lmuon_em_etaEpl = (TH2D*) lMuon_templates->Get("em_etaEpl"); |
171 |
< |
TH2D* lmuon_had_etaEpl = (TH2D*) lMuon_templates->Get("had_etaEpl"); |
172 |
< |
|
173 |
< |
double lEta = -1.; double lP = -1; |
174 |
< |
double lEM = -5.; double lHad = 0; double lHO = 0; |
175 |
< |
lEta = iMuon->Eta(); |
176 |
< |
lP = iMuon->P(); |
177 |
< |
if(lP >= 2000.) lP = 1999.9; |
178 |
< |
if(!iEMSpecial || iMuon->EmEnergy() != 0.) lEM = iMuon->EmEnergy(); |
179 |
< |
lHad = iMuon->HadEnergy(); |
180 |
< |
lHO = iMuon->HoEnergy(); |
181 |
< |
if(lP < 0. ) return 0.5; |
182 |
< |
if(fabs(lEta) > 2.5 ) return 0.5; |
183 |
< |
TH2D* lTMuonHad = NULL; |
184 |
< |
TH2D* lTPionHad = NULL; |
185 |
< |
TH2D* lTMuonHo = NULL; |
186 |
< |
TH2D* lTPionHo = NULL; |
187 |
< |
TH2D* lTMuonEm = NULL; |
188 |
< |
TH2D* lTPionEm = NULL; |
189 |
< |
|
190 |
< |
if(fabs(lEta) >= 1.27) { |
191 |
< |
if(iCorrectedHCAL) lHad *= 1.8/2.2; |
192 |
< |
if(lEta > 0) { |
193 |
< |
lTPionHad = lpion_had_etaEpl; |
194 |
< |
lTMuonHad = lmuon_had_etaEpl; |
195 |
< |
} else { |
196 |
< |
lTPionHad = lpion_had_etaEmi; |
197 |
< |
lTMuonHad = lmuon_had_etaEmi; |
198 |
< |
} |
199 |
< |
} |
200 |
< |
if(fabs(lEta) < 1.27 && fabs(lEta) >= 1.1 ) { |
201 |
< |
if(iCorrectedHCAL) lHad *= (1.8/(-2.2*fabs(lEta)+5.5)); |
202 |
< |
if(lEta > 0) { |
203 |
< |
lTPionHad = lpion_had_etaTpl; |
204 |
< |
lTMuonHad = lmuon_had_etaTpl; |
205 |
< |
} else { |
206 |
< |
lTPionHad = lpion_had_etaTmi; |
207 |
< |
lTMuonHad = lmuon_had_etaTmi; |
208 |
< |
} |
209 |
< |
} |
210 |
< |
if(fabs(lEta) < 1.1) { |
211 |
< |
if(iCorrectedHCAL) lHad *= sin(2*atan(exp(iMuon->Eta()))); |
212 |
< |
lTPionHad = lpion_had_etaB; |
213 |
< |
lTMuonHad = lmuon_had_etaB; |
214 |
< |
} |
215 |
< |
if(lEta > 1.479 ) { |
216 |
< |
lTPionEm = lpion_em_etaEpl; |
217 |
< |
lTMuonEm = lmuon_em_etaEpl; |
218 |
< |
} |
219 |
< |
if(fabs(lEta) <= 1.479) { |
220 |
< |
lTPionEm = lpion_em_etaB; |
221 |
< |
lTMuonEm = lmuon_em_etaB; |
222 |
< |
} |
223 |
< |
if(lEta < -1.479 ) { |
224 |
< |
lTPionEm = lpion_em_etaEmi; |
225 |
< |
lTMuonEm = lmuon_em_etaEmi; |
226 |
< |
} |
227 |
< |
if(fabs(lEta) < 1.28) { |
228 |
< |
lTPionHo = lpion_ho_etaB; |
229 |
< |
lTMuonHo = lmuon_ho_etaB; |
230 |
< |
} |
231 |
< |
|
232 |
< |
double lPBX = 1.; double lPSX = 1.; |
233 |
< |
double lPBY = 1.; double lPSY = 1.; |
234 |
< |
double lPBZ = 1.; double lPSZ = 1.; |
235 |
< |
if(!overflow(lTPionEm, lP,lEM)) lPBX = lTPionEm ->GetBinContent(lTPionEm ->GetXaxis()->FindBin(lP),lTPionEm ->GetYaxis()->FindBin(lEM) ); |
236 |
< |
if(!overflow(lTPionHad,lP,lHad)) lPBY = lTPionHad->GetBinContent(lTPionHad->GetXaxis()->FindBin(lP),lTPionHad->GetYaxis()->FindBin(lHad)); |
237 |
< |
if(!overflow(lTPionHo, lP,lHO)) lPBZ = lTPionHo ->GetBinContent(lTPionHo ->GetXaxis()->FindBin(lP),lTPionHo ->GetYaxis()->FindBin(lHO) ); |
238 |
< |
if(!overflow(lTMuonEm, lP,lEM )) lPSX = lTMuonEm ->GetBinContent(lTMuonEm ->GetXaxis()->FindBin(lP),lTMuonEm ->GetYaxis()->FindBin(lEM) ); |
239 |
< |
if(!overflow(lTMuonHad,lP,lHad)) lPSY = lTMuonHad->GetBinContent(lTMuonHad->GetXaxis()->FindBin(lP),lTMuonHad->GetYaxis()->FindBin(lHad)); |
240 |
< |
if(!overflow(lTMuonHo ,lP,lHO)) lPSZ = lTMuonHo ->GetBinContent(lTMuonHo ->GetXaxis()->FindBin(lP),lTMuonHo ->GetYaxis()->FindBin(lHO) ); |
241 |
< |
|
242 |
< |
if(lPSX == 0. || lPBX == 0. || (lEM <= 0. && !iEMSpecial)) {lPSX = 1.; lPBX = 1.;} |
243 |
< |
if(lPSY == 0. || lPBY == 0. || lHad == 0.) {lPSY = 1.; lPBY = 1.;} |
244 |
< |
if(lPSZ == 0. || lPBZ == 0. || lHO == 0.) {lPSZ = 1.; lPBZ = 1.;} |
245 |
< |
lPion_templates->Close(); |
246 |
< |
lMuon_templates->Close(); |
247 |
< |
if((lPSX*lPSY*lPSZ+lPBX*lPBY*lPBZ) > 0.) return lPSX*lPSY*lPSZ/(lPSX*lPSY*lPSZ+lPBX*lPBY*lPBZ); |
248 |
< |
return 0.5; |
249 |
< |
} |
250 |
< |
static inline float getSegmentCompatability(const mithep::Muon* iMuon) { |
251 |
< |
int lNStationsCrossed = 0; |
252 |
< |
int lNStationsSegment = 0; |
253 |
< |
|
254 |
< |
std::vector<int> lStTrack(8); |
255 |
< |
std::vector<int> lStSegmentmatch(8); |
256 |
< |
std::vector<int> lStCrossed(8); |
257 |
< |
std::vector<float> lStBoundary(8); |
258 |
< |
std::vector<float> lStWeight(8); |
259 |
< |
float lWeight = 0.; |
260 |
< |
bool lAdjust = true; |
261 |
< |
for(int i0 = 0; i0 < 8; i0++) { |
262 |
< |
if(iMuon->GetTrackDist(i0) < 999999. ) { |
263 |
< |
lNStationsCrossed++; |
264 |
< |
lStCrossed[i0] = 1; |
265 |
< |
lStBoundary[i0] = 0.; |
266 |
< |
if(iMuon->GetTrackDist(i0) > -10. ) lStBoundary[i0] = iMuon->GetTrackDist(i0); |
267 |
< |
} |
268 |
< |
if(iMuon->GetDX(i0) < 999999.) { //Use iMuon->GetSegmentX--> CHECK |
269 |
< |
lNStationsSegment++; |
270 |
< |
lStSegmentmatch[i0] = 1; |
271 |
< |
} |
272 |
< |
if(iMuon->GetDY(i0) < 999999.) lAdjust = false; |
273 |
< |
} |
274 |
< |
int lPCross = -1; |
275 |
< |
const float lAtWeight = 0.5; |
276 |
< |
for(int i0 = 0; i0< 8; i0++) { |
277 |
< |
lStWeight[i0] = 0; |
278 |
< |
if(lStCrossed[i0] > 0 ) { lPCross++; |
279 |
< |
switch ( lNStationsCrossed ) { |
280 |
< |
case 1 : |
281 |
< |
lStWeight[i0] = 1.; |
282 |
< |
break; |
283 |
< |
case 2 : |
284 |
< |
if ( lPCross == 0 ) lStWeight[i0] = 0.33; |
285 |
< |
else lStWeight[i0] = 0.67; |
286 |
< |
break; |
287 |
< |
case 3 : |
288 |
< |
if ( lPCross == 0 ) lStWeight[i0] = 0.23; |
289 |
< |
else if( lPCross == 1 ) lStWeight[i0] = 0.33; |
290 |
< |
else lStWeight[i0] = 0.44; |
291 |
< |
break; |
292 |
< |
case 4 : |
293 |
< |
if ( lPCross == 0 ) lStWeight[i0] = 0.10; |
294 |
< |
else if( lPCross == 1 ) lStWeight[i0] = 0.20; |
295 |
< |
else if( lPCross == 2 ) lStWeight[i0] = 0.30; |
296 |
< |
else lStWeight[i0] = 0.40; |
297 |
< |
break; |
298 |
< |
|
299 |
< |
default : |
300 |
< |
lStWeight[i0] = 1./lNStationsCrossed; |
301 |
< |
} |
302 |
< |
|
303 |
< |
if(lStSegmentmatch[i0] <= 0 && lStBoundary[i0] != 0. ) { |
304 |
< |
lStWeight[i0] *= lAtWeight*0.5*(TMath::Erf(lStBoundary[i0]/6.)+1.); |
305 |
< |
} else if(lStSegmentmatch[i0] <= 0 && lStBoundary[i0] == 0) {lStWeight[i0] = 0.;} |
306 |
< |
|
307 |
< |
if( lStSegmentmatch[i0] > 0) { |
308 |
< |
double lP2X = TMath::Power(iMuon->GetPullX(i0),2.); |
309 |
< |
double lP2Y = TMath::Power(iMuon->GetPullY(i0),2.); |
310 |
< |
double lD2X = TMath::Power(iMuon->GetDX(i0),2.); |
311 |
< |
double lD2Y = TMath::Power(iMuon->GetDY(i0),2.); |
312 |
< |
if( iMuon->GetDY(i0) < 999999 && iMuon->GetDX(i0) < 999999) { |
313 |
< |
lStWeight[i0] *= sigWeight(sqrt(lD2X+lD2Y),sqrt(lP2X+lP2Y)); |
314 |
< |
} else if (iMuon->GetDY(i0) >= 999999 && i0 < 4) { |
315 |
< |
lStWeight[i0] *= sigWeight(iMuon->GetDX(i0),iMuon->GetPullX(i0)); |
316 |
< |
} else if(i0 < 4) { |
317 |
< |
lStWeight[i0] *= sigWeight(iMuon->GetDY(i0),iMuon->GetPullY(i0)); |
318 |
< |
} |
319 |
< |
} |
320 |
< |
} |
321 |
< |
lWeight += lStWeight[i0]; |
322 |
< |
} |
323 |
< |
if(lNStationsCrossed == 0) lWeight = 0.5; |
324 |
< |
return lWeight; |
325 |
< |
} |
326 |
< |
static inline bool isGoodMuon(mithep::Muon *iMuon,selection iSelection) { |
327 |
< |
double lVal = 0; |
328 |
< |
switch(iSelection) { |
329 |
< |
case AllArbitrated: |
330 |
< |
if(iMuon->StandaloneTrk() != 0 || iMuon->GlobalTrk()!= 0) return true; |
331 |
< |
if(iMuon->NSegments() > 0) return true; |
332 |
< |
return false; |
333 |
< |
break; |
334 |
< |
case PromptTight: |
335 |
< |
return promptTight(-1,iMuon); |
336 |
< |
break; |
337 |
< |
case TMOneStationLoose: |
338 |
< |
return TMOneStation(iMuon,99999,999999); |
339 |
< |
case TMOneStationTight: |
340 |
< |
return TMOneStation(iMuon); |
341 |
< |
case TMLastStationLoose: |
342 |
< |
return TMLastStation(iMuon,999999,999999); |
343 |
< |
break; |
344 |
< |
case TMLastStationTight: |
345 |
< |
return TMLastStation(iMuon); |
346 |
< |
break; |
347 |
< |
case TM2DCompatibilityLoose: |
348 |
< |
lVal = 1.2*getSegmentCompatability(iMuon); |
349 |
< |
if(lVal/1.2 == 0.5) return false; |
350 |
< |
//if(fClass ) lVal += 0.8*getCaloCompatability(iMuon,true,true); |
351 |
< |
lVal += 0.8*getCaloCompatabilitySlowly(iMuon,true,true); |
352 |
< |
if(lVal > 0.7) return true; |
353 |
< |
return false; |
354 |
< |
break; |
355 |
< |
case TM2DCompatibilityTight: |
356 |
< |
lVal = 1.2*getSegmentCompatability(iMuon); |
357 |
< |
if(lVal/1.2 == 0.5) return false; |
358 |
< |
//if(fClass ) lVal += 0.8*getCaloCompatability(iMuon,true,true); |
359 |
< |
lVal += 0.8*getCaloCompatabilitySlowly(iMuon,true,true); |
360 |
< |
if(lVal > 1.0) return true; |
361 |
< |
return false; |
362 |
< |
break; |
363 |
< |
default: |
364 |
< |
return false; |
365 |
< |
} |
366 |
< |
return false; |
367 |
< |
} |
368 |
< |
bool isGood(mithep::Muon *iMuon,selection iSelection); |
369 |
< |
private: |
370 |
< |
TFile *fMuon_templates; |
371 |
< |
TFile *fPion_templates; |
372 |
< |
TH2D* fpion_em_etaEmi; |
373 |
< |
TH2D* fpion_had_etaEmi; |
374 |
< |
TH2D* fpion_had_etaTmi; |
375 |
< |
TH2D* fpion_em_etaB; |
376 |
< |
TH2D* fpion_had_etaB; |
377 |
< |
TH2D* fpion_ho_etaB; |
378 |
< |
TH2D* fpion_had_etaTpl; |
379 |
< |
TH2D* fpion_em_etaEpl; |
380 |
< |
TH2D* fpion_had_etaEpl; |
381 |
< |
TH2D* fmuon_em_etaEmi; |
382 |
< |
TH2D* fmuon_had_etaEmi; |
383 |
< |
TH2D* fmuon_had_etaTmi; |
384 |
< |
TH2D* fmuon_em_etaB; |
385 |
< |
TH2D* fmuon_had_etaB; |
386 |
< |
TH2D* fmuon_ho_etaB; |
387 |
< |
TH2D* fmuon_had_etaTpl; |
388 |
< |
TH2D* fmuon_em_etaEpl; |
389 |
< |
TH2D* fmuon_had_etaEpl; |
28 |
> |
public: |
29 |
> |
MuonTools(const char *mutemp="$CMSSW_BASE/src/MitPhysics/data/MuonCaloTemplate.root", |
30 |
> |
const char *pitemp="$CMSSW_BASE/src/MitPhysics/data/PionCaloTemplate.root"); |
31 |
> |
~MuonTools(); |
32 |
> |
|
33 |
> |
enum ESelType { |
34 |
> |
kAllArbitrated, //All arbitration (DT/CSC/RPC Hits) put on at least one |
35 |
> |
// segments given a global Muon |
36 |
> |
kPromptTight, //Standard global muon identification |
37 |
> |
kTMLastStationLoose, //Loose matching requirements on lastmost muon station of reco |
38 |
> |
kTMLastStationTight, //Tight matching requirements on lastmost muon station of reco |
39 |
> |
kTMOneStationLoose, //Loose matching requirements on at least one muon station of reco |
40 |
> |
kTMOneStationTight, //Tight matching requirements on at least one muon station of reco |
41 |
> |
kTM2DCompatibilityLoose, //Loose requirement on sum of compatabiliity variables |
42 |
> |
// ===> 1.2 Segment compatability + 0.8 calo compatability > 0.8 |
43 |
> |
kTM2DCompatibilityTight //Tight requirement on sum of compatabiliity variables |
44 |
> |
// ===> 1.2 Segment compatability + 0.8 calo compatability > 1.2 |
45 |
> |
}; |
46 |
> |
|
47 |
> |
Bool_t Init(const char *mutemp, const char *pitemp); |
48 |
> |
Bool_t IsGood(const mithep::Muon *iMuon, ESelType iSel) const; |
49 |
> |
Double_t GetCaloCompatability(const mithep::Muon *iMuon, |
50 |
> |
Bool_t iEMSpecial, Bool_t iCorrectedHCAL) const; |
51 |
> |
Double_t GetSegmentCompatability(const mithep::Muon *iMuon) const; |
52 |
> |
|
53 |
> |
protected: |
54 |
> |
void DeleteHistos(); |
55 |
> |
Bool_t Overflow(const TH2D *iHist, Double_t lVal0, Double_t lVal1) const; |
56 |
> |
Double_t SigWeight(Double_t iVal0, Double_t iVal1) const; |
57 |
> |
|
58 |
> |
private: |
59 |
> |
Bool_t fIsInit; //!=true if histograms are loaded |
60 |
> |
TH2D *fmuon_em_etaEmi; //!Neg Endcap EM Calo Deposit Template for Muons |
61 |
> |
TH2D *fmuon_had_etaEmi; //!Neg Endcap Hadronic Calo Deposit Template for Muons |
62 |
> |
TH2D *fmuon_had_etaTmi; //!Neg Transition Hadronic Calo Deposit Template for Muons |
63 |
> |
TH2D *fmuon_em_etaB; //!Barrel EM Calo Deposit Template for Muons |
64 |
> |
TH2D *fmuon_had_etaB; //!Barrel Hadronic Calo Deposit Template for Muons |
65 |
> |
TH2D *fmuon_ho_etaB; //!Barrel HO Calo Deposit Template for Muons |
66 |
> |
TH2D *fmuon_had_etaTpl; //!Plus Transition Hadronic Calo Deposit Template for Muons |
67 |
> |
TH2D *fmuon_em_etaEpl; //!Plus Endcap EM Calo Deposit Template for Muons |
68 |
> |
TH2D *fmuon_had_etaEpl; //!Plus Endcap Hadronic Calo Deposit Template for Muons |
69 |
> |
TH2D *fpion_em_etaEmi; //!Neg Endcap EM Calo Deposit Template for Pions |
70 |
> |
TH2D *fpion_had_etaEmi; //!Neg Endcap Hadronic Calo Deposit Template for Pions |
71 |
> |
TH2D *fpion_had_etaTmi; //!Neg Transition Hadronic Calo Deposit Template for Pions |
72 |
> |
TH2D *fpion_em_etaB; //!Barrel EM Calo Deposit Template for Pions |
73 |
> |
TH2D *fpion_had_etaB; //!Barrel Hadronic Calo Deposit Template for Pions |
74 |
> |
TH2D *fpion_ho_etaB; //!Barrel HO Calo Deposit Template for Pions |
75 |
> |
TH2D *fpion_had_etaTpl; //!Plus Transition Hadronic Calo Deposit Template for Pions |
76 |
> |
TH2D *fpion_em_etaEpl; //!Plus Endcap EM Calo Deposit Template for Pions |
77 |
> |
TH2D *fpion_had_etaEpl; //!Plus Endcap Hadronic Calo Deposit Template for Pions |
78 |
> |
|
79 |
> |
TH2D *LoadHisto(const char *fname, TFile *file) const; |
80 |
|
}; |
81 |
|
} |
82 |
|
|
83 |
+ |
//-------------------------------------------------------------------------------------------------- |
84 |
+ |
inline Double_t mithep::MuonTools::SigWeight(Double_t iVal0, Double_t iVal1) const |
85 |
+ |
{ |
86 |
+ |
// Returns weighted uncertainty given segment matching uncertainty (iVal0) and |
87 |
+ |
// segment matching pull (iVal1). |
88 |
+ |
|
89 |
+ |
if (iVal1 < 1.) //if pull defined and within range |
90 |
+ |
return 1.; |
91 |
+ |
if (iVal0 < 3. && iVal1 > 3.) { //if pull not well defined and uncertainty defined |
92 |
+ |
Double_t lVal = TMath::Max(iVal0,1.); |
93 |
+ |
return 1./TMath::Power(lVal,0.25); |
94 |
+ |
} |
95 |
+ |
|
96 |
+ |
Double_t lVal = TMath::Max(iVal1,1.); //if pull > 1 and pull < 3 return 1/pull^4 |
97 |
+ |
return 1./TMath::Power(lVal,0.25); |
98 |
+ |
} |
99 |
+ |
//-------------------------------------------------------------------------------------------------- |
100 |
+ |
inline Bool_t mithep::MuonTools::Overflow(const TH2D *iHist, Double_t lVal0, Double_t lVal1) const |
101 |
+ |
{ |
102 |
+ |
// Check if values are in overflow bins of given histogram. |
103 |
+ |
|
104 |
+ |
if(iHist == 0) |
105 |
+ |
return kTRUE; |
106 |
+ |
|
107 |
+ |
if (iHist ->GetXaxis()->FindBin(lVal0) == 0 || |
108 |
+ |
iHist ->GetXaxis()->FindBin(lVal0) > iHist->GetNbinsX() || |
109 |
+ |
iHist ->GetYaxis()->FindBin(lVal0) == 0 || |
110 |
+ |
iHist ->GetYaxis()->FindBin(lVal0) > iHist->GetNbinsY()) { |
111 |
+ |
return kTRUE; |
112 |
+ |
} |
113 |
+ |
return kFALSE; |
114 |
+ |
} |
115 |
|
#endif |