1 |
+ |
// $Id$ |
2 |
+ |
|
3 |
|
#include "MitPhysics/Utils/interface/MuonTools.h" |
4 |
+ |
#include <TFile.h> |
5 |
+ |
|
6 |
+ |
ClassImp(mithep::MuonTools) |
7 |
+ |
|
8 |
|
using namespace mithep; |
9 |
|
|
10 |
< |
MuonTools::MuonTools() : |
11 |
< |
fpion_em_etaEmi(0), |
12 |
< |
fpion_had_etaEmi(0), |
7 |
< |
fpion_had_etaTmi(0), |
8 |
< |
fpion_em_etaB(0), |
9 |
< |
fpion_had_etaB(0), |
10 |
< |
fpion_ho_etaB(0), |
11 |
< |
fpion_had_etaTpl(0), |
12 |
< |
fpion_em_etaEpl(0), |
13 |
< |
fpion_had_etaEpl(0), |
10 |
> |
//-------------------------------------------------------------------------------------------------- |
11 |
> |
MuonTools::MuonTools(const char *mutemp, const char *pitemp) : |
12 |
> |
fIsInit(kFALSE), |
13 |
|
fmuon_em_etaEmi(0), |
14 |
|
fmuon_had_etaEmi(0), |
15 |
|
fmuon_had_etaTmi(0), |
18 |
|
fmuon_ho_etaB(0), |
19 |
|
fmuon_had_etaTpl(0), |
20 |
|
fmuon_em_etaEpl(0), |
21 |
< |
fmuon_had_etaEpl(0) {} |
22 |
< |
MuonTools::~MuonTools() { |
23 |
< |
delete fpion_em_etaEmi; |
24 |
< |
delete fpion_had_etaEmi; |
25 |
< |
delete fpion_had_etaTmi; |
26 |
< |
delete fpion_em_etaB; |
27 |
< |
delete fpion_had_etaB; |
28 |
< |
delete fpion_ho_etaB; |
29 |
< |
delete fpion_had_etaTpl; |
30 |
< |
delete fpion_em_etaEpl; |
31 |
< |
delete fpion_had_etaEpl; |
32 |
< |
delete fmuon_em_etaEmi; |
33 |
< |
delete fmuon_had_etaEmi; |
34 |
< |
delete fmuon_had_etaTmi; |
35 |
< |
delete fmuon_em_etaB; |
36 |
< |
delete fmuon_had_etaB; |
37 |
< |
delete fmuon_ho_etaB; |
38 |
< |
delete fmuon_had_etaTpl; |
39 |
< |
delete fmuon_em_etaEpl; |
40 |
< |
delete fmuon_had_etaEpl; |
41 |
< |
} |
42 |
< |
|
43 |
< |
double MuonTools::getCaloCompatability(mithep::Muon* iMuon,bool iEMSpecial, bool iCorrectedHCAL) { |
44 |
< |
if(fpion_em_etaEmi == 0) { |
45 |
< |
TFile* fPion_templates = new TFile("$CMSSW_BASE/src/MitPhysics/Init/PionCaloTemplate.root","READ"); |
46 |
< |
TFile* fMuon_templates = new TFile("$CMSSW_BASE/src/MitPhysics/Init/MuonCaloTemplate.root","READ"); |
47 |
< |
fpion_em_etaEmi = (TH2D*) fPion_templates->Get("em_etaEmi"); |
48 |
< |
fpion_had_etaEmi = (TH2D*) fPion_templates->Get("had_etaEmi"); |
49 |
< |
fpion_had_etaTmi = (TH2D*) fPion_templates->Get("had_etaTmi"); |
50 |
< |
fpion_em_etaB = (TH2D*) fPion_templates->Get("em_etaB") ; |
51 |
< |
fpion_had_etaB = (TH2D*) fPion_templates->Get("had_etaB") ; |
52 |
< |
fpion_ho_etaB = (TH2D*) fPion_templates->Get("ho_etaB") ; |
53 |
< |
fpion_had_etaTpl = (TH2D*) fPion_templates->Get("had_etaTpl"); |
54 |
< |
fpion_em_etaEpl = (TH2D*) fPion_templates->Get("em_etaEpl") ; |
55 |
< |
fpion_had_etaEpl = (TH2D*) fPion_templates->Get("had_etaEpl"); |
56 |
< |
fmuon_em_etaEmi = (TH2D*) fMuon_templates->Get("em_etaEmi") ; |
57 |
< |
fmuon_had_etaEmi = (TH2D*) fMuon_templates->Get("had_etaEmi"); |
58 |
< |
fmuon_had_etaTmi = (TH2D*) fMuon_templates->Get("had_etaTmi"); |
59 |
< |
fmuon_em_etaB = (TH2D*) fMuon_templates->Get("em_etaB") ; |
60 |
< |
fmuon_had_etaB = (TH2D*) fMuon_templates->Get("had_etaB") ; |
61 |
< |
fmuon_ho_etaB = (TH2D*) fMuon_templates->Get("ho_etaB") ; |
62 |
< |
fmuon_had_etaTpl = (TH2D*) fMuon_templates->Get("had_etaTpl"); |
63 |
< |
fmuon_em_etaEpl = (TH2D*) fMuon_templates->Get("em_etaEpl"); |
64 |
< |
fmuon_had_etaEpl = (TH2D*) fMuon_templates->Get("had_etaEpl"); |
65 |
< |
|
66 |
< |
fpion_em_etaEmi ->SetDirectory(0); |
67 |
< |
fpion_had_etaEmi->SetDirectory(0); |
68 |
< |
fpion_had_etaTmi->SetDirectory(0); |
69 |
< |
fpion_em_etaB ->SetDirectory(0); |
70 |
< |
fpion_had_etaB ->SetDirectory(0); |
71 |
< |
fpion_ho_etaB ->SetDirectory(0); |
72 |
< |
fpion_had_etaTpl->SetDirectory(0); |
73 |
< |
fpion_em_etaEpl ->SetDirectory(0); |
74 |
< |
fpion_had_etaEpl->SetDirectory(0); |
75 |
< |
fmuon_em_etaEmi ->SetDirectory(0); |
76 |
< |
fmuon_had_etaEmi->SetDirectory(0); |
77 |
< |
fmuon_had_etaTmi->SetDirectory(0); |
78 |
< |
fmuon_em_etaB ->SetDirectory(0); |
79 |
< |
fmuon_had_etaB ->SetDirectory(0); |
80 |
< |
fmuon_ho_etaB ->SetDirectory(0); |
81 |
< |
fmuon_had_etaTpl->SetDirectory(0); |
82 |
< |
fmuon_em_etaEpl ->SetDirectory(0); |
83 |
< |
fmuon_had_etaEpl->SetDirectory(0); |
84 |
< |
fPion_templates->Close(); |
85 |
< |
fMuon_templates->Close(); |
86 |
< |
} |
87 |
< |
double lEta = -1.; double lP = -1; |
88 |
< |
double lEM = -5.; double lHad = 0; double lHO = 0; |
89 |
< |
lEta = iMuon->Eta(); |
90 |
< |
lP = iMuon->P(); |
91 |
< |
if(lP >= 2000.) lP = 1999.9; |
92 |
< |
if(!iEMSpecial || iMuon->EmEnergy() != 0.) lEM = iMuon->EmEnergy(); |
93 |
< |
lHad = iMuon->HadEnergy(); |
94 |
< |
lHO = iMuon->HoEnergy(); |
95 |
< |
if(lP < 0. ) return 0.5; |
96 |
< |
if(fabs(lEta) > 2.5 ) return 0.5; |
97 |
< |
TH2D* lTMuonHad = NULL; |
98 |
< |
TH2D* lTPionHad = NULL; |
99 |
< |
TH2D* lTMuonHo = NULL; |
100 |
< |
TH2D* lTPionHo = NULL; |
101 |
< |
TH2D* lTMuonEm = NULL; |
102 |
< |
TH2D* lTPionEm = NULL; |
103 |
< |
|
104 |
< |
if(fabs(lEta) >= 1.27) { |
105 |
< |
if(iCorrectedHCAL) lHad *= 1.8/2.2; |
106 |
< |
if(lEta > 0) { |
21 |
> |
fmuon_had_etaEpl(0), |
22 |
> |
fpion_em_etaEmi(0), |
23 |
> |
fpion_had_etaEmi(0), |
24 |
> |
fpion_had_etaTmi(0), |
25 |
> |
fpion_em_etaB(0), |
26 |
> |
fpion_had_etaB(0), |
27 |
> |
fpion_ho_etaB(0), |
28 |
> |
fpion_had_etaTpl(0), |
29 |
> |
fpion_em_etaEpl(0), |
30 |
> |
fpion_had_etaEpl(0) |
31 |
> |
{ |
32 |
> |
// Constructor. |
33 |
> |
|
34 |
> |
if (mutemp && pitemp) |
35 |
> |
Init(mutemp, pitemp); |
36 |
> |
} |
37 |
> |
|
38 |
> |
//-------------------------------------------------------------------------------------------------- |
39 |
> |
MuonTools::~MuonTools() |
40 |
> |
{ |
41 |
> |
// Destructor. |
42 |
> |
|
43 |
> |
DeleteHistos(); |
44 |
> |
} |
45 |
> |
|
46 |
> |
//-------------------------------------------------------------------------------------------------- |
47 |
> |
void MuonTools::DeleteHistos() |
48 |
> |
{ |
49 |
> |
// Delete histograms. |
50 |
> |
|
51 |
> |
if (fIsInit) { |
52 |
> |
delete fpion_em_etaEmi; |
53 |
> |
delete fpion_had_etaEmi; |
54 |
> |
delete fpion_had_etaTmi; |
55 |
> |
delete fpion_em_etaB; |
56 |
> |
delete fpion_had_etaB; |
57 |
> |
delete fpion_ho_etaB; |
58 |
> |
delete fpion_had_etaTpl; |
59 |
> |
delete fpion_em_etaEpl; |
60 |
> |
delete fpion_had_etaEpl; |
61 |
> |
delete fmuon_em_etaEmi; |
62 |
> |
delete fmuon_had_etaEmi; |
63 |
> |
delete fmuon_had_etaTmi; |
64 |
> |
delete fmuon_em_etaB; |
65 |
> |
delete fmuon_had_etaB; |
66 |
> |
delete fmuon_ho_etaB; |
67 |
> |
delete fmuon_had_etaTpl; |
68 |
> |
delete fmuon_em_etaEpl; |
69 |
> |
delete fmuon_had_etaEpl; |
70 |
> |
fpion_em_etaEmi = 0; |
71 |
> |
fpion_had_etaEmi = 0; |
72 |
> |
fpion_had_etaTmi = 0; |
73 |
> |
fpion_em_etaB = 0; |
74 |
> |
fpion_had_etaB = 0; |
75 |
> |
fpion_ho_etaB = 0; |
76 |
> |
fpion_had_etaTpl = 0; |
77 |
> |
fpion_em_etaEpl = 0; |
78 |
> |
fpion_had_etaEpl = 0; |
79 |
> |
fmuon_em_etaEmi = 0; |
80 |
> |
fmuon_had_etaEmi = 0; |
81 |
> |
fmuon_had_etaTmi = 0; |
82 |
> |
fmuon_em_etaB = 0; |
83 |
> |
fmuon_had_etaB = 0; |
84 |
> |
fmuon_ho_etaB = 0; |
85 |
> |
fmuon_had_etaTpl = 0; |
86 |
> |
fmuon_em_etaEpl = 0; |
87 |
> |
fmuon_had_etaEpl = 0; |
88 |
> |
fIsInit = kFALSE; |
89 |
> |
} |
90 |
> |
} |
91 |
> |
|
92 |
> |
//-------------------------------------------------------------------------------------------------- |
93 |
> |
Double_t MuonTools::GetCaloCompatability(const Muon *iMuon, |
94 |
> |
Bool_t iEMSpecial, Bool_t iCorrectedHCAL) const |
95 |
> |
{ |
96 |
> |
// Get calo compatibility value for given muon based on calorimeter templates. |
97 |
> |
// If iEMSpecial is true, then a use different arrangement of ECAL for compatibility. |
98 |
> |
|
99 |
> |
Double_t lEta = iMuon->Eta(); |
100 |
> |
Double_t aEta = TMath::Abs(lEta); |
101 |
> |
if (aEta > 2.5) |
102 |
> |
return 0.5; |
103 |
> |
|
104 |
> |
Double_t lP = iMuon->P(); |
105 |
> |
if (lP >= 2000.) |
106 |
> |
lP = 1999.9; |
107 |
> |
if(lP < 0. ) |
108 |
> |
return 0.5; |
109 |
> |
|
110 |
> |
Double_t lEM = -5.; |
111 |
> |
if (!iEMSpecial || iMuon->EmEnergy() != 0.) |
112 |
> |
lEM = iMuon->EmEnergy(); |
113 |
> |
|
114 |
> |
Double_t lHad = iMuon->HadEnergy(); |
115 |
> |
Double_t lHO = iMuon->HoEnergy();; |
116 |
> |
|
117 |
> |
TH2D *lTMuonHad = 0; |
118 |
> |
TH2D *lTPionHad = 0; |
119 |
> |
TH2D *lTMuonHo = 0; |
120 |
> |
TH2D *lTPionHo = 0; |
121 |
> |
TH2D *lTMuonEm = 0; |
122 |
> |
TH2D *lTPionEm = 0; |
123 |
> |
|
124 |
> |
if (aEta >= 1.27) { |
125 |
> |
if (iCorrectedHCAL) |
126 |
> |
lHad *= 1.8/2.2; |
127 |
> |
if (lEta > 0) { |
128 |
|
lTPionHad = fpion_had_etaEpl; |
129 |
|
lTMuonHad = fmuon_had_etaEpl; |
130 |
|
} else { |
132 |
|
lTMuonHad = fmuon_had_etaEmi; |
133 |
|
} |
134 |
|
} |
135 |
< |
if(fabs(lEta) < 1.27 && fabs(lEta) >= 1.1 ) { |
136 |
< |
if(iCorrectedHCAL) lHad *= (1.8/(-2.2*fabs(lEta)+5.5)); |
137 |
< |
if(lEta > 0) { |
135 |
> |
|
136 |
> |
if (aEta < 1.27 && aEta >= 1.1) { |
137 |
> |
if (iCorrectedHCAL) |
138 |
> |
lHad *= (1.8/(-2.2*aEta+5.5)); |
139 |
> |
if (lEta > 0) { |
140 |
|
lTPionHad = fpion_had_etaTpl; |
141 |
|
lTMuonHad = fmuon_had_etaTpl; |
142 |
|
} else { |
144 |
|
lTMuonHad = fmuon_had_etaTmi; |
145 |
|
} |
146 |
|
} |
147 |
< |
if(fabs(lEta) < 1.1) { |
148 |
< |
if(iCorrectedHCAL) lHad *= sin(2*atan(exp(iMuon->Eta()))); |
147 |
> |
|
148 |
> |
if (aEta < 1.1) { |
149 |
> |
if(iCorrectedHCAL) |
150 |
> |
lHad *= TMath::Sin(2*TMath::ATan(TMath::Exp(lEta))); |
151 |
|
lTPionHad = fpion_had_etaB; |
152 |
|
lTMuonHad = fmuon_had_etaB; |
153 |
|
} |
154 |
< |
if(lEta > 1.479 ) { |
155 |
< |
lTPionEm = fpion_em_etaEpl; |
156 |
< |
lTMuonEm = fmuon_em_etaEpl; |
154 |
> |
if (lEta > 1.479) { |
155 |
> |
lTPionEm = fpion_em_etaEpl; |
156 |
> |
lTMuonEm = fmuon_em_etaEpl; |
157 |
|
} |
158 |
< |
if(fabs(lEta) <= 1.479) { |
158 |
> |
if (aEta <= 1.479) { |
159 |
|
lTPionEm = fpion_em_etaB; |
160 |
|
lTMuonEm = fmuon_em_etaB; |
161 |
|
} |
162 |
< |
if(lEta < -1.479 ) { |
162 |
> |
if (lEta < -1.479) { |
163 |
|
lTPionEm = fpion_em_etaEmi; |
164 |
|
lTMuonEm = fmuon_em_etaEmi; |
165 |
|
} |
166 |
< |
if(fabs(lEta) < 1.28) { |
166 |
> |
if (aEta < 1.28) { |
167 |
|
lTPionHo = fpion_ho_etaB; |
168 |
|
lTMuonHo = fmuon_ho_etaB; |
169 |
|
} |
170 |
|
|
171 |
< |
double lPBX = 1.; double lPSX = 1.; |
172 |
< |
double lPBY = 1.; double lPSY = 1.; |
173 |
< |
double lPBZ = 1.; double lPSZ = 1.; |
174 |
< |
if(!overflow(lTPionEm, lP,lEM)) lPBX = lTPionEm ->GetBinContent(lTPionEm ->GetXaxis()->FindBin(lP),lTPionEm ->GetYaxis()->FindBin(lEM) ); |
175 |
< |
if(!overflow(lTPionHad,lP,lHad)) lPBY = lTPionHad->GetBinContent(lTPionHad->GetXaxis()->FindBin(lP),lTPionHad->GetYaxis()->FindBin(lHad)); |
176 |
< |
if(!overflow(lTPionHo, lP,lHO)) lPBZ = lTPionHo ->GetBinContent(lTPionHo ->GetXaxis()->FindBin(lP),lTPionHo ->GetYaxis()->FindBin(lHO) ); |
177 |
< |
if(!overflow(lTMuonEm, lP,lEM )) lPSX = lTMuonEm ->GetBinContent(lTMuonEm ->GetXaxis()->FindBin(lP),lTMuonEm ->GetYaxis()->FindBin(lEM) ); |
178 |
< |
if(!overflow(lTMuonHad,lP,lHad)) lPSY = lTMuonHad->GetBinContent(lTMuonHad->GetXaxis()->FindBin(lP),lTMuonHad->GetYaxis()->FindBin(lHad)); |
179 |
< |
if(!overflow(lTMuonHo ,lP,lHO)) lPSZ = lTMuonHo ->GetBinContent(lTMuonHo ->GetXaxis()->FindBin(lP),lTMuonHo ->GetYaxis()->FindBin(lHO) ); |
180 |
< |
|
181 |
< |
if(lPSX == 0. || lPBX == 0. || (lEM <= 0. && !iEMSpecial)) {lPSX = 1.; lPBX = 1.;} |
182 |
< |
if(lPSY == 0. || lPBY == 0. || lHad == 0.) {lPSY = 1.; lPBY = 1.;} |
183 |
< |
if(lPSZ == 0. || lPBZ == 0. || lHO == 0.) {lPSZ = 1.; lPBZ = 1.;} |
184 |
< |
if((lPSX*lPSY*lPSZ+lPBX*lPBY*lPBZ) > 0.) return lPSX*lPSY*lPSZ/(lPSX*lPSY*lPSZ+lPBX*lPBY*lPBZ); |
171 |
> |
Double_t lPBX = 1.; |
172 |
> |
Double_t lPSX = 1.; |
173 |
> |
Double_t lPBY = 1.; |
174 |
> |
Double_t lPSY = 1.; |
175 |
> |
Double_t lPBZ = 1.; |
176 |
> |
Double_t lPSZ = 1.; |
177 |
> |
if (!Overflow(lTPionEm, lP,lEM)) |
178 |
> |
lPBX = lTPionEm ->GetBinContent(lTPionEm ->GetXaxis()->FindBin(lP), |
179 |
> |
lTPionEm ->GetYaxis()->FindBin(lEM)); |
180 |
> |
if (!Overflow(lTPionHad,lP,lHad)) |
181 |
> |
lPBY = lTPionHad->GetBinContent(lTPionHad->GetXaxis()->FindBin(lP), |
182 |
> |
lTPionHad->GetYaxis()->FindBin(lHad)); |
183 |
> |
if (!Overflow(lTPionHo, lP,lHO)) |
184 |
> |
lPBZ = lTPionHo ->GetBinContent(lTPionHo ->GetXaxis()->FindBin(lP), |
185 |
> |
lTPionHo ->GetYaxis()->FindBin(lHO)); |
186 |
> |
if (!Overflow(lTMuonEm, lP,lEM )) |
187 |
> |
lPSX = lTMuonEm ->GetBinContent(lTMuonEm ->GetXaxis()->FindBin(lP), |
188 |
> |
lTMuonEm ->GetYaxis()->FindBin(lEM)); |
189 |
> |
if (!Overflow(lTMuonHad,lP,lHad)) |
190 |
> |
lPSY = lTMuonHad->GetBinContent(lTMuonHad->GetXaxis()->FindBin(lP), |
191 |
> |
lTMuonHad->GetYaxis()->FindBin(lHad)); |
192 |
> |
if (!Overflow(lTMuonHo ,lP,lHO)) |
193 |
> |
lPSZ = lTMuonHo ->GetBinContent(lTMuonHo ->GetXaxis()->FindBin(lP), |
194 |
> |
lTMuonHo ->GetYaxis()->FindBin(lHO)); |
195 |
> |
|
196 |
> |
if (lPSX == 0. || lPBX == 0. || (lEM <= 0. && !iEMSpecial)) { |
197 |
> |
lPSX = 1.; |
198 |
> |
lPBX = 1.; |
199 |
> |
} |
200 |
> |
if (lPSY == 0. || lPBY == 0. || lHad == 0.) { |
201 |
> |
lPSY = 1.; |
202 |
> |
lPBY = 1.; |
203 |
> |
} |
204 |
> |
if (lPSZ == 0. || lPBZ == 0. || lHO == 0.) { |
205 |
> |
lPSZ = 1.; |
206 |
> |
lPBZ = 1.; |
207 |
> |
} |
208 |
> |
if ((lPSX*lPSY*lPSZ+lPBX*lPBY*lPBZ) > 0.) |
209 |
> |
return lPSX*lPSY*lPSZ/(lPSX*lPSY*lPSZ+lPBX*lPBY*lPBZ); |
210 |
> |
|
211 |
|
return 0.5; |
212 |
|
} |
213 |
|
|
214 |
< |
bool MuonTools::isGood(mithep::Muon *iMuon,selection iSelection) { |
215 |
< |
double lVal = 0; |
216 |
< |
switch(iSelection) { |
217 |
< |
case AllArbitrated: |
218 |
< |
if(iMuon->StandaloneTrk() != 0 || iMuon->GlobalTrk()!= 0) return true; |
219 |
< |
if(iMuon->NSegments() > 0) return true; |
220 |
< |
return false; |
221 |
< |
break; |
222 |
< |
case PromptTight: |
223 |
< |
return promptTight(-1,iMuon); |
224 |
< |
break; |
225 |
< |
case TMOneStationLoose: |
226 |
< |
return TMOneStation(iMuon,99999,999999); |
227 |
< |
break; |
228 |
< |
case TMOneStationTight: |
229 |
< |
return TMOneStation(iMuon); |
230 |
< |
break; |
231 |
< |
case TMLastStationLoose: |
232 |
< |
return TMLastStation(iMuon,999999,999999); |
233 |
< |
break; |
234 |
< |
case TMLastStationTight: |
235 |
< |
return TMLastStation(iMuon); |
236 |
< |
break; |
237 |
< |
case TM2DCompatibilityLoose: |
238 |
< |
lVal = 1.2*getSegmentCompatability(iMuon); |
239 |
< |
if(lVal/1.2 == 0.5) return false; |
240 |
< |
lVal += 0.8*getCaloCompatability(iMuon,true,true); |
241 |
< |
if(lVal > 0.7) return true; |
242 |
< |
return false; |
243 |
< |
break; |
244 |
< |
case TM2DCompatibilityTight: |
245 |
< |
lVal = 1.2*getSegmentCompatability(iMuon); |
196 |
< |
if(lVal/1.2 == 0.5) return false; |
197 |
< |
lVal += 0.8*getCaloCompatability(iMuon,true,true); |
198 |
< |
if(lVal > 1.0) return true; |
199 |
< |
return false; |
200 |
< |
break; |
201 |
< |
default: |
202 |
< |
return false; |
214 |
> |
//-------------------------------------------------------------------------------------------------- |
215 |
> |
Bool_t MuonTools::Init(const char *mutemp, const char *pitemp) |
216 |
> |
{ |
217 |
> |
// Read histograms from given files. |
218 |
> |
|
219 |
> |
if (fIsInit) { |
220 |
> |
DeleteHistos(); |
221 |
> |
} |
222 |
> |
|
223 |
> |
TDirectory::TContext context(0); |
224 |
> |
|
225 |
> |
TFile *muon_templates = TFile::Open(mutemp); |
226 |
> |
if (!muon_templates) { |
227 |
> |
Fatal("Init", "Could not open file %s", mutemp); |
228 |
> |
return kFALSE; |
229 |
> |
} |
230 |
> |
fmuon_em_etaEmi = LoadHisto("em_etaEmi", muon_templates); |
231 |
> |
fmuon_had_etaEmi = LoadHisto("had_etaEmi", muon_templates); |
232 |
> |
fmuon_had_etaTmi = LoadHisto("had_etaTmi", muon_templates); |
233 |
> |
fmuon_em_etaB = LoadHisto("em_etaB", muon_templates); |
234 |
> |
fmuon_had_etaB = LoadHisto("had_etaB", muon_templates); |
235 |
> |
fmuon_ho_etaB = LoadHisto("ho_etaB", muon_templates); |
236 |
> |
fmuon_had_etaTpl = LoadHisto("had_etaTpl", muon_templates); |
237 |
> |
fmuon_em_etaEpl = LoadHisto("em_etaEpl", muon_templates); |
238 |
> |
fmuon_had_etaEpl = LoadHisto("had_etaEpl", muon_templates); |
239 |
> |
muon_templates->Close(); |
240 |
> |
delete muon_templates; |
241 |
> |
|
242 |
> |
TFile *pion_templates = TFile::Open(pitemp); |
243 |
> |
if (!pion_templates) { |
244 |
> |
Fatal("Init", "Could not open file %s", pitemp); |
245 |
> |
return kFALSE; |
246 |
|
} |
247 |
< |
return false; |
247 |
> |
|
248 |
> |
fpion_em_etaEmi = LoadHisto("em_etaEmi", pion_templates); |
249 |
> |
fpion_had_etaEmi = LoadHisto("had_etaEmi", pion_templates); |
250 |
> |
fpion_had_etaTmi = LoadHisto("had_etaTmi", pion_templates); |
251 |
> |
fpion_em_etaB = LoadHisto("em_etaB", pion_templates); |
252 |
> |
fpion_had_etaB = LoadHisto("had_etaB", pion_templates); |
253 |
> |
fpion_ho_etaB = LoadHisto("ho_etaB", pion_templates); |
254 |
> |
fpion_had_etaTpl = LoadHisto("had_etaTpl", pion_templates); |
255 |
> |
fpion_em_etaEpl = LoadHisto("em_etaEpl", pion_templates); |
256 |
> |
fpion_had_etaEpl = LoadHisto("had_etaEpl", pion_templates); |
257 |
> |
pion_templates->Close(); |
258 |
> |
delete pion_templates; |
259 |
> |
|
260 |
> |
fIsInit = kTRUE; |
261 |
> |
return kTRUE; |
262 |
> |
} |
263 |
> |
|
264 |
> |
//-------------------------------------------------------------------------------------------------- |
265 |
> |
Bool_t MuonTools::IsGood(const mithep::Muon *iMuon, ESelType iSel) const |
266 |
> |
{ |
267 |
> |
// Return true if given muon qualifies given selection criterium. |
268 |
> |
|
269 |
> |
Double_t tm2dcut = 0.; |
270 |
> |
|
271 |
> |
switch(iSel) { |
272 |
> |
case kAllArbitrated: |
273 |
> |
if (iMuon->StandaloneTrk() != 0 || iMuon->GlobalTrk()!= 0) |
274 |
> |
return kTRUE; |
275 |
> |
if (iMuon->NSegments() > 0) |
276 |
> |
return kTRUE; |
277 |
> |
break; |
278 |
> |
case kPromptTight: |
279 |
> |
return iMuon->PromptTight(Muon::kAny); |
280 |
> |
break; |
281 |
> |
case kTMOneStationLoose: |
282 |
> |
return iMuon->TMOneStation(999999,999999); |
283 |
> |
break; |
284 |
> |
case kTMOneStationTight: |
285 |
> |
return iMuon->TMOneStation(); |
286 |
> |
break; |
287 |
> |
case kTMLastStationLoose: |
288 |
> |
return iMuon->TMLastStation(999999,999999); |
289 |
> |
break; |
290 |
> |
case kTMLastStationTight: |
291 |
> |
return iMuon->TMLastStation(); |
292 |
> |
break; |
293 |
> |
case kTM2DCompatibilityLoose: |
294 |
> |
tm2dcut = 0.7; |
295 |
> |
break; |
296 |
> |
case kTM2DCompatibilityTight: |
297 |
> |
tm2dcut = 1.0; |
298 |
> |
break; |
299 |
> |
default: |
300 |
> |
return kFALSE; |
301 |
> |
break; |
302 |
> |
} |
303 |
> |
|
304 |
> |
Double_t lVal = GetSegmentCompatability(iMuon); |
305 |
> |
if (lVal == 0.5) // exclude this border case |
306 |
> |
return kFALSE; |
307 |
> |
|
308 |
> |
lVal *= 1.2; |
309 |
> |
lVal += 0.8*GetCaloCompatability(iMuon,kTRUE,kTRUE); |
310 |
> |
if (lVal > tm2dcut) |
311 |
> |
return kTRUE; |
312 |
> |
|
313 |
> |
return kFALSE; |
314 |
> |
} |
315 |
> |
|
316 |
> |
//-------------------------------------------------------------------------------------------------- |
317 |
> |
Double_t MuonTools::GetSegmentCompatability(const mithep::Muon *iMuon) const |
318 |
> |
{ |
319 |
> |
// Get segment compatability for given muon based on likelihood of well defined |
320 |
> |
// track through chambers. |
321 |
> |
|
322 |
> |
Int_t lNStationsCrossed = 0; |
323 |
> |
Int_t lNStationsSegment = 0; |
324 |
> |
|
325 |
> |
Int_t lStSegmentmatch[8]; |
326 |
> |
Int_t lStCrossed[8]; |
327 |
> |
Double_t lStBoundary[8]; |
328 |
> |
|
329 |
> |
Double_t lWeight = 0.; |
330 |
> |
Bool_t lAdjust = kTRUE; |
331 |
> |
for (Int_t i0 = 0; i0 < 8; ++i0) { |
332 |
> |
lStBoundary[i0] = 0.; |
333 |
> |
if(iMuon->GetTrackDist(i0) < 999999. ) { |
334 |
> |
lNStationsCrossed++; |
335 |
> |
lStCrossed[i0] = 1; |
336 |
> |
if (iMuon->GetTrackDist(i0) > -10. ) |
337 |
> |
lStBoundary[i0] = iMuon->GetTrackDist(i0); |
338 |
> |
} else |
339 |
> |
lStCrossed[i0] = 0; |
340 |
> |
|
341 |
> |
if(iMuon->GetDX(i0) < 999999.) { |
342 |
> |
lNStationsSegment++; |
343 |
> |
lStSegmentmatch[i0] = 1; |
344 |
> |
} else |
345 |
> |
lStSegmentmatch[i0] = 0; |
346 |
> |
|
347 |
> |
if(iMuon->GetDY(i0) < 999999.) |
348 |
> |
lAdjust = kFALSE; |
349 |
> |
} |
350 |
> |
|
351 |
> |
if (lNStationsCrossed == 0) |
352 |
> |
return 0.5; |
353 |
> |
|
354 |
> |
Double_t lStWeight[8]; |
355 |
> |
Int_t lPCross = -1; |
356 |
> |
const Double_t lAtWeight = 0.5; |
357 |
> |
for (Int_t i0 = 0; i0< 8; ++i0) { |
358 |
> |
lStWeight[i0] = 0; |
359 |
> |
if (lStCrossed[i0] > 0) { |
360 |
> |
lPCross++; |
361 |
> |
|
362 |
> |
switch (lNStationsCrossed) { |
363 |
> |
case 1 : |
364 |
> |
lStWeight[i0] = 1.; |
365 |
> |
break; |
366 |
> |
case 2 : |
367 |
> |
if (lPCross == 0 ) |
368 |
> |
lStWeight[i0] = 0.33; |
369 |
> |
else |
370 |
> |
lStWeight[i0] = 0.67; |
371 |
> |
break; |
372 |
> |
case 3 : |
373 |
> |
if (lPCross == 0) |
374 |
> |
lStWeight[i0] = 0.23; |
375 |
> |
else if (lPCross == 1) |
376 |
> |
lStWeight[i0] = 0.33; |
377 |
> |
else |
378 |
> |
lStWeight[i0] = 0.44; |
379 |
> |
break; |
380 |
> |
case 4 : |
381 |
> |
if (lPCross == 0) |
382 |
> |
lStWeight[i0] = 0.10; |
383 |
> |
else if (lPCross == 1) |
384 |
> |
lStWeight[i0] = 0.20; |
385 |
> |
else if (lPCross == 2) |
386 |
> |
lStWeight[i0] = 0.30; |
387 |
> |
else |
388 |
> |
lStWeight[i0] = 0.40; |
389 |
> |
break; |
390 |
> |
default : |
391 |
> |
lStWeight[i0] = 1./lNStationsCrossed; |
392 |
> |
} |
393 |
> |
|
394 |
> |
if (lStSegmentmatch[i0] <= 0 && lStBoundary[i0] != 0.) |
395 |
> |
lStWeight[i0] *= lAtWeight*0.5*(TMath::Erf(lStBoundary[i0]/6.)+1.); |
396 |
> |
else if (lStSegmentmatch[i0] <= 0 && lStBoundary[i0] == 0) |
397 |
> |
lStWeight[i0] = 0.; |
398 |
> |
|
399 |
> |
if (lStSegmentmatch[i0] > 0) { |
400 |
> |
Double_t lP2X = TMath::Power(iMuon->GetPullX(i0),2.); |
401 |
> |
Double_t lP2Y = TMath::Power(iMuon->GetPullY(i0),2.); |
402 |
> |
Double_t lD2X = TMath::Power(iMuon->GetDX(i0),2.); |
403 |
> |
Double_t lD2Y = TMath::Power(iMuon->GetDY(i0),2.); |
404 |
> |
if (iMuon->GetDY(i0) < 999999 && iMuon->GetDX(i0) < 999999) |
405 |
> |
lStWeight[i0] *= SigWeight(TMath::Sqrt(lD2X+lD2Y),TMath::Sqrt(lP2X+lP2Y)); |
406 |
> |
else if (iMuon->GetDY(i0) >= 999999 && i0 < 4) |
407 |
> |
lStWeight[i0] *= SigWeight(iMuon->GetDX(i0),iMuon->GetPullX(i0)); |
408 |
> |
else if(i0 < 4) |
409 |
> |
lStWeight[i0] *= SigWeight(iMuon->GetDY(i0),iMuon->GetPullY(i0)); |
410 |
> |
} |
411 |
> |
} |
412 |
> |
lWeight += lStWeight[i0]; |
413 |
> |
} |
414 |
> |
|
415 |
> |
return lWeight; |
416 |
> |
} |
417 |
> |
|
418 |
> |
//-------------------------------------------------------------------------------------------------- |
419 |
> |
TH2D *MuonTools::LoadHisto(const char *name, TFile *file) const |
420 |
> |
{ |
421 |
> |
// Load histogram with given name from given file and return it. |
422 |
> |
|
423 |
> |
TH2D *ret = dynamic_cast<TH2D*>(file->Get(name)); |
424 |
> |
if (!ret) { |
425 |
> |
Fatal("LoadHisto", "Could not load histogram %s from file %s", name, file->GetName()); |
426 |
> |
return 0; |
427 |
> |
} |
428 |
> |
ret->SetDirectory(0); |
429 |
> |
return ret; |
430 |
> |
} |
431 |
> |
//-------------------------------------------------------------------------------------------------- |
432 |
> |
Bool_t MuonTools::PassD0Cut(const Muon *mu, const VertexCol *vertices, Double_t fD0Cut) |
433 |
> |
{ |
434 |
> |
Bool_t d0cut = kFALSE; |
435 |
> |
const Track *mt = mu->BestTrk(); |
436 |
> |
if (!mt) return kFALSE; |
437 |
> |
|
438 |
> |
Double_t d0_real = 1e30; |
439 |
> |
for(UInt_t i0 = 0; i0 < vertices->GetEntries(); i0++) { |
440 |
> |
if(vertices->At(i0)->NTracks() > 0){ |
441 |
> |
Double_t pD0 = mt->D0Corrected(*vertices->At(i0)); |
442 |
> |
d0_real = TMath::Abs(pD0); |
443 |
> |
break; |
444 |
> |
} |
445 |
> |
} |
446 |
> |
if(d0_real < fD0Cut) d0cut = kTRUE; |
447 |
> |
|
448 |
> |
return d0cut; |
449 |
> |
} |
450 |
> |
|
451 |
> |
//-------------------------------------------------------------------------------------------------- |
452 |
> |
Bool_t MuonTools::PassD0Cut(const Muon *mu, const BeamSpotCol *beamspots, Double_t fD0Cut) |
453 |
> |
{ |
454 |
> |
Bool_t d0cut = kFALSE; |
455 |
> |
const Track *mt = mu->BestTrk(); |
456 |
> |
if (!mt) return kFALSE; |
457 |
> |
|
458 |
> |
// d0 cut |
459 |
> |
Double_t d0_real = 99999; |
460 |
> |
for(UInt_t i0 = 0; i0 < beamspots->GetEntries(); i0++) { |
461 |
> |
Double_t pD0 = mt->D0Corrected(*beamspots->At(i0)); |
462 |
> |
if(TMath::Abs(pD0) < TMath::Abs(d0_real)) d0_real = TMath::Abs(pD0); |
463 |
> |
} |
464 |
> |
if(d0_real < fD0Cut) d0cut = kTRUE; |
465 |
> |
|
466 |
> |
return d0cut; |
467 |
> |
} |
468 |
> |
|
469 |
> |
//-------------------------------------------------------------------------------------------------- |
470 |
> |
Bool_t MuonTools::PassSoftMuonCut(const Muon *mu, const VertexCol *vertices) |
471 |
> |
{ |
472 |
> |
if(mu->Pt() <= 3.0) return kFALSE; |
473 |
> |
|
474 |
> |
if(!mu->IsTrackerMuon()) return kFALSE; |
475 |
> |
|
476 |
> |
if(!mu->Quality().Quality(MuonQuality::TMLastStationAngTight)) return kFALSE; |
477 |
> |
|
478 |
> |
if(mu->BestTrk()->NHits() <= 10) return kFALSE; |
479 |
> |
|
480 |
> |
if(!PassD0Cut(mu, vertices, 0.2)) return kFALSE; |
481 |
> |
|
482 |
> |
Double_t totalIso = 1.0 * mu->IsoR03SumPt() + |
483 |
> |
1.0 * mu->IsoR03EmEt() + |
484 |
> |
1.0 * mu->IsoR03HadEt(); |
485 |
> |
if (totalIso < (mu->Pt()*0.10) && mu->Pt() > 20.0) return kFALSE; |
486 |
> |
|
487 |
> |
return kTRUE; |
488 |
|
} |