]>
Commit | Line | Data |
---|---|---|
1c5acb87 | 1 | /************************************************************************** |
2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | * * | |
4 | * Author: The ALICE Off-line Project. * | |
5 | * Contributors are mentioned in the code where appropriate. * | |
6 | * * | |
7 | * Permission to use, copy, modify and distribute this software and its * | |
8 | * documentation strictly for non-commercial purposes is hereby granted * | |
9 | * without fee, provided that the above copyright notice appears in all * | |
10 | * copies and that both the copyright notice and this permission notice * | |
11 | * appear in the supporting documentation. The authors make no claims * | |
12 | * about the suitability of this software for any purpose. It is * | |
13 | * provided "as is" without express or implied warranty. * | |
14 | **************************************************************************/ | |
15 | /* $Id: AliCaloPID.cxx 21839 2007-10-29 13:49:42Z gustavo $ */ | |
16 | ||
17 | //_________________________________________________________________________ | |
bdd2a262 | 18 | // Class for PID selection with calorimeters |
21a4b1c0 | 19 | // The Output of the 2 main methods GetIdentifiedParticleType is a PDG number identifying the cluster, |
bdd2a262 | 20 | // being kPhoton, kElectron, kPi0 ... as defined in the header file |
21a4b1c0 | 21 | // - GetIdentifiedParticleType(const TString calo, const Double_t * pid, const Float_t energy) |
bdd2a262 | 22 | // Reads the PID weights array of the ESDs and depending on its magnitude identifies the particle |
21a4b1c0 | 23 | // - GetIdentifiedParticleType(const TString calo,const TLorentzVector mom, const AliVCluster * cluster) |
bdd2a262 | 24 | // Recalcultes PID, the bayesian or any new one to be implemented in the future |
25 | // Right now only the possibility to recalculate EMCAL with bayesian and simple PID. | |
26 | // In order to recalculate Bayesian, it is necessary to load the EMCALUtils library | |
27 | // and do SwitchOnBayesianRecalculation(). | |
28 | // To change the PID parameters from Low to High like the ones by default, use the constructor | |
29 | // AliCaloPID(flux) | |
30 | // where flux is AliCaloPID::kLow or AliCaloPID::kHigh | |
31 | // If it is necessary to change the parameters use the constructor | |
32 | // AliCaloPID(AliEMCALPIDUtils *utils) and set the parameters before. | |
9a6fa057 | 33 | // - SetPIDBits: Simple PID, depending on the thresholds fL0Cut fTOFCut and even the |
bdd2a262 | 34 | // result of the PID bayesian a different PID bit is set. |
35 | // | |
36 | // All these methods can be called in the analysis you are interested. | |
37 | // | |
1c5acb87 | 38 | //*-- Author: Gustavo Conesa (LNF-INFN) |
39 | ////////////////////////////////////////////////////////////////////////////// | |
40 | ||
41 | ||
42 | // --- ROOT system --- | |
43 | #include <TMath.h> | |
1c5acb87 | 44 | #include <TString.h> |
f21fc003 | 45 | #include <TList.h> |
1c5acb87 | 46 | |
47 | //---- ANALYSIS system ---- | |
1c5acb87 | 48 | #include "AliCaloPID.h" |
0ae57829 | 49 | #include "AliVCluster.h" |
d39cba7e | 50 | #include "AliVTrack.h" |
1c5acb87 | 51 | #include "AliAODPWG4Particle.h" |
f2ccb5b8 | 52 | #include "AliCalorimeterUtils.h" |
53 | ||
1c5acb87 | 54 | ClassImp(AliCaloPID) |
55 | ||
56 | ||
57 | //________________________________________________ | |
58 | AliCaloPID::AliCaloPID() : | |
a5fb4114 | 59 | TObject(), |
60 | fEMCALPhotonWeight(0.), fEMCALPi0Weight(0.), | |
61 | fEMCALElectronWeight(0.), fEMCALChargeWeight(0.), fEMCALNeutralWeight(0.), | |
62 | fPHOSPhotonWeight(0.), fPHOSPi0Weight(0.), | |
63 | fPHOSElectronWeight(0.), fPHOSChargeWeight(0.) , fPHOSNeutralWeight(0.), | |
64 | fPHOSWeightFormula(0), fPHOSPhotonWeightFormula(0), fPHOSPi0WeightFormula(0), | |
65 | fPHOSPhotonWeightFormulaExpression(""), fPHOSPi0WeightFormulaExpression(""), | |
9a6fa057 | 66 | fL0CutMax(100.), fL0CutMin(0), fTOFCut(0.), |
ae182e60 | 67 | fRcutPHOS(0), fDebug(-1), |
a5fb4114 | 68 | fRecalculateBayesian(kFALSE), fParticleFlux(kLow), fEMCALPIDUtils(), |
69 | fHistoNEBins(100), fHistoEMax(100.), fHistoEMin(0.), | |
70 | fHistoNDEtaBins(100), fHistoDEtaMax(0.15), fHistoDEtaMin(-0.15), | |
71 | fHistoNDPhiBins(100), fHistoDPhiMax(0.15), fHistoDPhiMin(-0.15), | |
72 | fhTrackMatchedDEta(0x0), fhTrackMatchedDPhi(0x0), fhTrackMatchedDEtaDPhi(0x0) | |
1c5acb87 | 73 | { |
477d6cee | 74 | //Ctor |
75 | ||
76 | //Initialize parameters | |
77 | InitParameters(); | |
1c5acb87 | 78 | } |
79 | ||
bdd2a262 | 80 | //________________________________________________ |
81 | AliCaloPID::AliCaloPID(const Int_t flux) : | |
a5fb4114 | 82 | TObject(), |
83 | fEMCALPhotonWeight(0.), fEMCALPi0Weight(0.), | |
84 | fEMCALElectronWeight(0.), fEMCALChargeWeight(0.), fEMCALNeutralWeight(0.), | |
85 | fPHOSPhotonWeight(0.), fPHOSPi0Weight(0.), | |
86 | fPHOSElectronWeight(0.), fPHOSChargeWeight(0.) , fPHOSNeutralWeight(0.), | |
87 | fPHOSWeightFormula(0), fPHOSPhotonWeightFormula(0), fPHOSPi0WeightFormula(0), | |
88 | fPHOSPhotonWeightFormulaExpression(""), fPHOSPi0WeightFormulaExpression(""), | |
9a6fa057 | 89 | fL0CutMax(100.), fL0CutMin(0), fTOFCut(0.), |
ae182e60 | 90 | fRcutPHOS(0), fDebug(-1), |
a5fb4114 | 91 | fRecalculateBayesian(kFALSE), fParticleFlux(flux), fEMCALPIDUtils(), |
92 | fHistoNEBins(100), fHistoEMax(100.), fHistoEMin(0.), | |
93 | fHistoNDEtaBins(100), fHistoDEtaMax(0.15), fHistoDEtaMin(-0.15), | |
94 | fHistoNDPhiBins(100), fHistoDPhiMax(0.15), fHistoDPhiMin(-0.15), | |
95 | fhTrackMatchedDEta(0x0), fhTrackMatchedDPhi(0x0), fhTrackMatchedDEtaDPhi(0x0) | |
bdd2a262 | 96 | { |
9a6fa057 | 97 | //Ctor |
bdd2a262 | 98 | |
9a6fa057 | 99 | //Initialize parameters |
100 | InitParameters(); | |
101 | ||
bdd2a262 | 102 | } |
103 | ||
104 | //________________________________________________ | |
f21fc003 | 105 | AliCaloPID::AliCaloPID(const TNamed * emcalpid) : |
a5fb4114 | 106 | TObject(), |
107 | fEMCALPhotonWeight(0.), fEMCALPi0Weight(0.), | |
108 | fEMCALElectronWeight(0.), fEMCALChargeWeight(0.), fEMCALNeutralWeight(0.), | |
109 | fPHOSPhotonWeight(0.), fPHOSPi0Weight(0.), | |
110 | fPHOSElectronWeight(0.), fPHOSChargeWeight(0.) , fPHOSNeutralWeight(0.), | |
111 | fPHOSWeightFormula(0), fPHOSPhotonWeightFormula(0), fPHOSPi0WeightFormula(0), | |
112 | fPHOSPhotonWeightFormulaExpression(""), fPHOSPi0WeightFormulaExpression(""), | |
ae182e60 | 113 | fL0CutMax(100.), fL0CutMin(0), fTOFCut(0.), |
114 | fRcutPHOS(0), fDebug(-1), | |
a5fb4114 | 115 | fRecalculateBayesian(kFALSE), fParticleFlux(kLow), fEMCALPIDUtils((AliEMCALPIDUtils*) emcalpid), |
116 | fHistoNEBins(100), fHistoEMax(100.), fHistoEMin(0.), | |
117 | fHistoNDEtaBins(100), fHistoDEtaMax(0.15), fHistoDEtaMin(-0.15), | |
118 | fHistoNDPhiBins(100), fHistoDPhiMax(0.15), fHistoDPhiMin(-0.15), | |
119 | fhTrackMatchedDEta(0x0), fhTrackMatchedDPhi(0x0), fhTrackMatchedDEtaDPhi(0x0) | |
bdd2a262 | 120 | { |
9a6fa057 | 121 | //Ctor |
bdd2a262 | 122 | |
9a6fa057 | 123 | //Initialize parameters |
124 | InitParameters(); | |
bdd2a262 | 125 | } |
126 | ||
1c5acb87 | 127 | //_________________________________ |
128 | AliCaloPID::~AliCaloPID() { | |
477d6cee | 129 | //Dtor |
130 | ||
a5fb4114 | 131 | delete fPHOSPhotonWeightFormula ; |
132 | delete fPHOSPi0WeightFormula ; | |
133 | delete fEMCALPIDUtils ; | |
1c5acb87 | 134 | |
9a6fa057 | 135 | } |
1c5acb87 | 136 | |
d39cba7e | 137 | //________________________________________________________________________ |
138 | TList * AliCaloPID::GetCreateOutputObjects() | |
139 | { | |
140 | // Create histograms to be saved in output file and | |
141 | // store them in outputContainer of the analysis class that calls this class. | |
142 | ||
143 | TList * outputContainer = new TList() ; | |
144 | outputContainer->SetName("CaloPIDHistos") ; | |
145 | ||
146 | outputContainer->SetOwner(kFALSE); | |
147 | ||
148 | fhTrackMatchedDEta = new TH2F | |
149 | ("TrackMatchedDEta", | |
a5fb4114 | 150 | "d#eta of cluster-track vs cluster energy", |
151 | fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNDEtaBins,fHistoDEtaMin,fHistoDEtaMax); | |
d39cba7e | 152 | fhTrackMatchedDEta->SetYTitle("d#eta"); |
153 | fhTrackMatchedDEta->SetXTitle("E_{cluster} (GeV)"); | |
154 | ||
155 | fhTrackMatchedDPhi = new TH2F | |
156 | ("TrackMatchedDPhi", | |
157 | "d#phi of cluster-track vs cluster energy" | |
158 | ,fHistoNEBins,fHistoEMin,fHistoEMax,fHistoNDPhiBins,fHistoDPhiMin,fHistoDPhiMax); | |
159 | fhTrackMatchedDPhi->SetYTitle("d#phi (rad)"); | |
160 | fhTrackMatchedDPhi->SetXTitle("E_{cluster} (GeV)"); | |
161 | ||
162 | fhTrackMatchedDEtaDPhi = new TH2F | |
163 | ("TrackMatchedDEtaDPhi", | |
164 | "d#eta vs d#phi of cluster-track vs cluster energy" | |
165 | ,fHistoNDEtaBins,fHistoDEtaMin,fHistoDEtaMax,fHistoNDPhiBins,fHistoDPhiMin,fHistoDPhiMax); | |
166 | fhTrackMatchedDEtaDPhi->SetYTitle("d#phi (rad)"); | |
167 | fhTrackMatchedDEtaDPhi->SetXTitle("d#eta"); | |
168 | ||
169 | outputContainer->Add(fhTrackMatchedDEta) ; | |
170 | outputContainer->Add(fhTrackMatchedDPhi) ; | |
171 | outputContainer->Add(fhTrackMatchedDEtaDPhi) ; | |
172 | ||
173 | return outputContainer; | |
174 | } | |
175 | ||
176 | ||
177 | ||
1c5acb87 | 178 | //_______________________________________________________________ |
179 | void AliCaloPID::InitParameters() | |
180 | { | |
477d6cee | 181 | //Initialize the parameters of the PID. |
182 | ||
2007809d | 183 | fEMCALPhotonWeight = 0.6 ; |
184 | fEMCALPi0Weight = 0.6 ; | |
185 | fEMCALElectronWeight = 0.6 ; | |
186 | fEMCALChargeWeight = 0.6 ; | |
187 | fEMCALNeutralWeight = 0.6 ; | |
477d6cee | 188 | |
2007809d | 189 | fPHOSPhotonWeight = 0.6 ; |
190 | fPHOSPi0Weight = 0.6 ; | |
191 | fPHOSElectronWeight = 0.6 ; | |
192 | fPHOSChargeWeight = 0.6 ; | |
193 | fPHOSNeutralWeight = 0.6 ; | |
477d6cee | 194 | |
195 | //Formula to set the PID weight threshold for photon or pi0 | |
a5fb4114 | 196 | fPHOSWeightFormula = kFALSE; |
197 | fPHOSPhotonWeightFormulaExpression = "0.98*(x<40)+ 0.68*(x>=100)+(x>=40 && x<100)*(0.98+x*(6e-3)-x*x*(2e-04)+x*x*x*(1.1e-06))"; | |
198 | fPHOSPi0WeightFormulaExpression = "0.98*(x<65)+ 0.915*(x>=100)+(x>=65 && x-x*(1.95e-3)-x*x*(4.31e-05)+x*x*x*(3.61e-07))" ; | |
199 | ||
9a6fa057 | 200 | fL0CutMax = 0.3; |
201 | fL0CutMin = 0.01; | |
202 | fTOFCut = 5.e-9; | |
ae182e60 | 203 | fRcutPHOS = 10.; |
204 | ||
9a6fa057 | 205 | fDebug =-1; |
bdd2a262 | 206 | |
207 | if(fRecalculateBayesian){ | |
bdd2a262 | 208 | if(fParticleFlux == kLow){ |
209 | printf("AliCaloPID::Init() - SetLOWFluxParam\n"); | |
210 | fEMCALPIDUtils->SetLowFluxParam() ; | |
211 | } | |
212 | else if (fParticleFlux == kHigh){ | |
213 | printf("AliCaloPID::Init() - SetHIGHFluxParam\n"); | |
214 | fEMCALPIDUtils->SetHighFluxParam() ; | |
215 | } | |
bdd2a262 | 216 | } |
1c5acb87 | 217 | } |
218 | ||
219 | //_______________________________________________________________ | |
21a4b1c0 | 220 | Int_t AliCaloPID::GetIdentifiedParticleType(const TString calo, const Double_t * pid, const Float_t energy) { |
477d6cee | 221 | //Return most probable identity of the particle. |
222 | ||
223 | if(!pid){ | |
21a4b1c0 | 224 | printf("AliCaloPID::GetIdentifiedParticleType() - pid pointer not initialized!!!\n"); |
477d6cee | 225 | abort(); |
226 | } | |
227 | ||
15800db4 | 228 | Float_t wPh = fPHOSPhotonWeight ; |
477d6cee | 229 | Float_t wPi0 = fPHOSPi0Weight ; |
15800db4 | 230 | Float_t wE = fPHOSElectronWeight ; |
231 | Float_t wCh = fPHOSChargeWeight ; | |
232 | Float_t wNe = fPHOSNeutralWeight ; | |
a5fb4114 | 233 | |
234 | if(calo == "PHOS" && fPHOSWeightFormula){ | |
235 | wPh = GetPHOSPhotonWeightFormula()->Eval(energy) ; | |
236 | wPi0 = GetPHOSPi0WeightFormula() ->Eval(energy); | |
237 | } | |
477d6cee | 238 | |
239 | if(calo == "EMCAL"){ | |
240 | ||
241 | wPh = fEMCALPhotonWeight ; | |
242 | wPi0 = fEMCALPi0Weight ; | |
243 | wE = fEMCALElectronWeight ; | |
244 | wCh = fEMCALChargeWeight ; | |
245 | wNe = fEMCALNeutralWeight ; | |
246 | ||
247 | } | |
248 | ||
21a4b1c0 | 249 | if(fDebug > 0) printf("AliCaloPID::GetIdentifiedParticleType: calo %s, ph %0.2f, pi0 %0.2f, el %0.2f, conv el %0.2f, hadrons: pion %0.2f, kaon %0.2f, proton %0.2f , neutron %0.2f, kaon %0.2f \n", |
c8fe2783 | 250 | calo.Data(),pid[AliVCluster::kPhoton], pid[AliVCluster::kPi0], |
251 | pid[AliVCluster::kElectron], pid[AliVCluster::kEleCon], | |
252 | pid[AliVCluster::kPion], pid[AliVCluster::kKaon], pid[AliVCluster::kProton], | |
253 | pid[AliVCluster::kNeutron], pid[AliVCluster::kKaon0]); | |
477d6cee | 254 | |
255 | Int_t pdg = kNeutralUnknown ; | |
c8fe2783 | 256 | Float_t chargedHadronWeight = pid[AliVCluster::kProton]+pid[AliVCluster::kKaon]+ |
257 | pid[AliVCluster::kPion]+pid[AliVCluster::kMuon]; | |
258 | Float_t neutralHadronWeight = pid[AliVCluster::kNeutron]+pid[AliVCluster::kKaon0]; | |
259 | Float_t allChargedWeight = pid[AliVCluster::kElectron]+pid[AliVCluster::kEleCon]+ chargedHadronWeight; | |
260 | Float_t allNeutralWeight = pid[AliVCluster::kPhoton]+pid[AliVCluster::kPi0]+ neutralHadronWeight; | |
477d6cee | 261 | |
262 | //Select most probable ID | |
263 | if(calo=="PHOS"){ | |
a5fb4114 | 264 | if(pid[AliVCluster::kPhoton] > wPh) pdg = kPhoton ; |
265 | else if(pid[AliVCluster::kPi0] > wPi0) pdg = kPi0 ; | |
c8fe2783 | 266 | else if(pid[AliVCluster::kElectron] > wE) pdg = kElectron ; |
a5fb4114 | 267 | else if(pid[AliVCluster::kEleCon] > wE) pdg = kEleCon ; |
268 | else if(chargedHadronWeight > wCh) pdg = kChargedHadron ; | |
269 | else if(neutralHadronWeight > wNe) pdg = kNeutralHadron ; | |
477d6cee | 270 | else if(allChargedWeight > allNeutralWeight) |
271 | pdg = kChargedUnknown ; | |
272 | else | |
273 | pdg = kNeutralUnknown ; | |
274 | } | |
275 | else{//EMCAL | |
2007809d | 276 | |
277 | if(pid[AliVCluster::kPhoton] > wPh) pdg = kPhoton ; | |
278 | else if(pid[AliVCluster::kElectron] > wE) pdg = kElectron ; | |
279 | else if(pid[AliVCluster::kPhoton]+pid[AliVCluster::kElectron] > wPh) pdg = kPhoton ; //temporal sollution until track matching for electrons is considered | |
280 | else if(pid[AliVCluster::kPi0] > wPi0) pdg = kPi0 ; | |
477d6cee | 281 | else if(chargedHadronWeight + neutralHadronWeight > wCh) pdg = kChargedHadron ; |
282 | else if(neutralHadronWeight + chargedHadronWeight > wNe) pdg = kNeutralHadron ; | |
2007809d | 283 | else pdg = kNeutralUnknown ; |
477d6cee | 284 | } |
285 | ||
21a4b1c0 | 286 | if(fDebug > 0)printf("AliCaloPID::GetIdentifiedParticleType:Final Pdg: %d, cluster energy %2.2f \n", pdg,energy); |
a5fb4114 | 287 | //printf("PDG1\n"); |
477d6cee | 288 | return pdg ; |
289 | ||
1c5acb87 | 290 | } |
291 | ||
292 | //_______________________________________________________________ | |
21a4b1c0 | 293 | Int_t AliCaloPID::GetIdentifiedParticleType(const TString calo,const TLorentzVector mom, const AliVCluster * cluster) { |
477d6cee | 294 | //Recalculated PID with all parameters |
a5fb4114 | 295 | |
bdd2a262 | 296 | Float_t lambda0 = cluster->GetM02(); |
9a6fa057 | 297 | Float_t lambda1 = cluster->GetM20(); |
a5fb4114 | 298 | Float_t energy = mom.E(); |
299 | ||
21a4b1c0 | 300 | if(fDebug > 0) printf("AliCaloPID::GetIdentifiedParticleType: Calorimeter %s, E %3.2f, l0 %3.2f, l1 %3.2f, disp %3.2f, tof %1.11f, distCPV %3.2f, distToBC %1.1f, NMax %d\n", |
a5fb4114 | 301 | calo.Data(),energy,lambda0,cluster->GetM20(),cluster->GetDispersion(),cluster->GetTOF(), |
302 | cluster->GetEmcCpvDistance(), cluster->GetDistanceToBadChannel(),cluster->GetNExMax()); | |
303 | ||
477d6cee | 304 | if(calo == "EMCAL") { |
bdd2a262 | 305 | //Recalculate Bayesian |
306 | if(fRecalculateBayesian){ | |
bdd2a262 | 307 | if(fDebug > 0) { |
c8fe2783 | 308 | const Double_t *pid0 = cluster->GetPID(); |
21a4b1c0 | 309 | printf("AliCaloPID::GetIdentifiedParticleType: BEFORE calo %s, ph %0.2f, pi0 %0.2f, el %0.2f, conv el %0.2f, hadrons: pion %0.2f, kaon %0.2f, proton %0.2f , neutron %0.2f, kaon %0.2f \n", |
a5fb4114 | 310 | calo.Data(),pid0[AliVCluster::kPhoton], pid0[AliVCluster::kPi0], |
311 | pid0[AliVCluster::kElectron], pid0[AliVCluster::kEleCon], | |
312 | pid0[AliVCluster::kPion], pid0[AliVCluster::kKaon], pid0[AliVCluster::kProton], | |
313 | pid0[AliVCluster::kNeutron], pid0[AliVCluster::kKaon0]); | |
bdd2a262 | 314 | } |
315 | ||
a5fb4114 | 316 | fEMCALPIDUtils->ComputePID(energy, lambda0); |
317 | Double_t pid[AliPID::kSPECIESN]; | |
318 | for(Int_t i = 0; i < AliPID::kSPECIESN; i++) pid[i] = fEMCALPIDUtils->GetPIDFinal(i); | |
21a4b1c0 | 319 | return GetIdentifiedParticleType(calo, pid, energy); |
bdd2a262 | 320 | |
a5fb4114 | 321 | } |
9a6fa057 | 322 | |
323 | // If no use of bayesian, simplest PID | |
324 | if(lambda0 < fL0CutMax && lambda0 > fL0CutMin) return kPhoton ; | |
325 | else return kPi0 ; // Wild guess, it can be hadron but not photon | |
477d6cee | 326 | |
9a6fa057 | 327 | }//EMCAL |
328 | else {//PHOS | |
329 | ||
330 | // Do not use bayesian, cut based on shower ellipse | |
ec6d2a72 | 331 | return IsPHOSPhoton(lambda0,lambda1); |
9a6fa057 | 332 | |
333 | } | |
21a4b1c0 | 334 | |
1c5acb87 | 335 | } |
336 | ||
337 | //__________________________________________________ | |
338 | TString AliCaloPID::GetPIDParametersList() { | |
477d6cee | 339 | //Put data member values in string to keep in output container |
340 | ||
341 | TString parList ; //this will be list of parameters used for this analysis. | |
5ae09196 | 342 | const Int_t buffersize = 255; |
343 | char onePar[buffersize] ; | |
344 | snprintf(onePar,buffersize,"--- AliCaloPID ---\n") ; | |
477d6cee | 345 | parList+=onePar ; |
9a6fa057 | 346 | snprintf(onePar,buffersize,"fL0CutMin =%2.2f, fL0CutMax =%2.2f (Cut on Shower Shape, used in PID evaluation) \n",fL0CutMin, fL0CutMax) ; |
477d6cee | 347 | parList+=onePar ; |
5ae09196 | 348 | snprintf(onePar,buffersize,"fTOFCut =%e (Cut on TOF, used in PID evaluation) \n",fTOFCut) ; |
477d6cee | 349 | parList+=onePar ; |
5ae09196 | 350 | snprintf(onePar,buffersize,"fEMCALPhotonWeight =%2.2f (EMCAL bayesian weight for photons)\n",fEMCALPhotonWeight) ; |
477d6cee | 351 | parList+=onePar ; |
5ae09196 | 352 | snprintf(onePar,buffersize,"fEMCALPi0Weight =%2.2f (EMCAL bayesian weight for pi0)\n",fEMCALPi0Weight) ; |
477d6cee | 353 | parList+=onePar ; |
5ae09196 | 354 | snprintf(onePar,buffersize,"fEMCALElectronWeight =%2.2f(EMCAL bayesian weight for electrons)\n",fEMCALElectronWeight) ; |
477d6cee | 355 | parList+=onePar ; |
5ae09196 | 356 | snprintf(onePar,buffersize,"fEMCALChargeWeight =%2.2f (EMCAL bayesian weight for charged hadrons)\n",fEMCALChargeWeight) ; |
477d6cee | 357 | parList+=onePar ; |
5ae09196 | 358 | snprintf(onePar,buffersize,"fEMCALNeutralWeight =%2.2f (EMCAL bayesian weight for neutral hadrons)\n",fEMCALNeutralWeight) ; |
477d6cee | 359 | parList+=onePar ; |
5ae09196 | 360 | snprintf(onePar,buffersize,"fPHOSPhotonWeight =%2.2f (PHOS bayesian weight for photons)\n",fPHOSPhotonWeight) ; |
477d6cee | 361 | parList+=onePar ; |
5ae09196 | 362 | snprintf(onePar,buffersize,"fPHOSPi0Weight =%2.2f (PHOS bayesian weight for pi0)\n",fPHOSPi0Weight) ; |
477d6cee | 363 | parList+=onePar ; |
5ae09196 | 364 | snprintf(onePar,buffersize,"fPHOSElectronWeight =%2.2f(PHOS bayesian weight for electrons)\n",fPHOSElectronWeight) ; |
477d6cee | 365 | parList+=onePar ; |
5ae09196 | 366 | snprintf(onePar,buffersize,"fPHOSChargeWeight =%2.2f (PHOS bayesian weight for charged hadrons)\n",fPHOSChargeWeight) ; |
477d6cee | 367 | parList+=onePar ; |
5ae09196 | 368 | snprintf(onePar,buffersize,"fPHOSNeutralWeight =%2.2f (PHOS bayesian weight for neutral hadrons)\n",fPHOSNeutralWeight) ; |
477d6cee | 369 | parList+=onePar ; |
370 | ||
a5fb4114 | 371 | if(fPHOSWeightFormula){ |
372 | snprintf(onePar,buffersize,"PHOS Photon Weight Formula: %s\n",fPHOSPhotonWeightFormulaExpression.Data() ) ; | |
373 | parList+=onePar; | |
374 | snprintf(onePar,buffersize,"PHOS Pi0 Weight Formula: %s\n",fPHOSPi0WeightFormulaExpression.Data() ) ; | |
375 | parList+=onePar; | |
376 | } | |
477d6cee | 377 | |
378 | return parList; | |
379 | ||
1c5acb87 | 380 | } |
381 | ||
382 | //________________________________________________________________ | |
383 | void AliCaloPID::Print(const Option_t * opt) const | |
384 | { | |
477d6cee | 385 | |
386 | //Print some relevant parameters set for the analysis | |
387 | if(! opt) | |
388 | return; | |
389 | ||
390 | printf("***** Print: %s %s ******\n", GetName(), GetTitle() ) ; | |
391 | ||
392 | printf("PHOS PID weight , photon %0.2f, pi0 %0.2f, e %0.2f, charge %0.2f, neutral %0.2f \n", | |
393 | fPHOSPhotonWeight, fPHOSPi0Weight, | |
394 | fPHOSElectronWeight, fPHOSChargeWeight, fPHOSNeutralWeight) ; | |
395 | printf("EMCAL PID weight, photon %0.2f, pi0 %0.2f, e %0.2f, charge %0.2f, neutral %0.2f\n", | |
396 | fEMCALPhotonWeight, fEMCALPi0Weight, | |
397 | fEMCALElectronWeight, fEMCALChargeWeight, fEMCALNeutralWeight) ; | |
398 | ||
a5fb4114 | 399 | printf("PHOS Parametrized weight on? = %d\n", fPHOSWeightFormula) ; |
400 | if(fPHOSWeightFormula){ | |
401 | printf("Photon weight formula = %s\n", fPHOSPhotonWeightFormulaExpression.Data()); | |
402 | printf("Pi0 weight formula = %s\n", fPHOSPi0WeightFormulaExpression .Data()); | |
403 | } | |
477d6cee | 404 | |
405 | printf("TOF cut = %e\n",fTOFCut); | |
9a6fa057 | 406 | printf("Lambda0 cut min = %2.2f; max = %2.2f\n",fL0CutMin, fL0CutMax); |
477d6cee | 407 | printf("Debug level = %d\n",fDebug); |
bdd2a262 | 408 | printf("Recalculate Bayesian? = %d\n",fRecalculateBayesian); |
409 | if(fRecalculateBayesian) printf("Particle Flux? = %d\n",fParticleFlux); | |
477d6cee | 410 | printf(" \n"); |
411 | ||
1c5acb87 | 412 | } |
413 | ||
414 | //_______________________________________________________________ | |
f2ccb5b8 | 415 | void AliCaloPID::SetPIDBits(const TString calo, const AliVCluster * cluster, AliAODPWG4Particle * ph, const AliCalorimeterUtils* cu) { |
477d6cee | 416 | //Set Bits for PID selection |
417 | ||
418 | //Dispersion/lambdas | |
5ae09196 | 419 | //Double_t disp= cluster->GetDispersion() ; |
420 | Double_t l1 = cluster->GetM20() ; | |
421 | Double_t l0 = cluster->GetM02() ; | |
422 | Bool_t isDispOK = kTRUE ; | |
9a6fa057 | 423 | if(cluster->IsPHOS()){ |
5ae09196 | 424 | |
9a6fa057 | 425 | isDispOK = IsPHOSPhoton(l0,l1); |
5ae09196 | 426 | |
427 | } | |
428 | else{//EMCAL | |
429 | ||
9a6fa057 | 430 | if(l0 > fL0CutMin && l0 < fL0CutMax) isDispOK = kTRUE; |
5ae09196 | 431 | |
432 | } | |
433 | ||
434 | ph->SetDispBit(isDispOK) ; | |
477d6cee | 435 | |
436 | //TOF | |
437 | Double_t tof=cluster->GetTOF() ; | |
438 | ph->SetTOFBit(TMath::Abs(tof)<fTOFCut) ; | |
439 | ||
5ae09196 | 440 | //Charged veto |
441 | //Bool_t isNeutral = kTRUE ; | |
442 | //if(cluster->IsPHOS()) isNeutral = cluster->GetEmcCpvDistance() > 5. ; | |
443 | //else | |
f2ccb5b8 | 444 | Bool_t isNeutral = IsTrackMatched(cluster,cu); |
5ae09196 | 445 | |
446 | ph->SetChargedBit(isNeutral); | |
477d6cee | 447 | |
448 | //Set PID pdg | |
21a4b1c0 | 449 | ph->SetIdentifiedParticleType(GetIdentifiedParticleType(calo,cluster->GetPID(),ph->E())); |
477d6cee | 450 | |
451 | if(fDebug > 0){ | |
5ae09196 | 452 | printf("AliCaloPID::SetPIDBits: TOF %e, Lambda0 %2.2f, Lambda1 %2.2f\n",tof , l0, l1); |
477d6cee | 453 | printf("AliCaloPID::SetPIDBits: pdg %d, bits: TOF %d, Dispersion %d, Charge %d\n", |
21a4b1c0 | 454 | ph->GetIdentifiedParticleType(), ph->GetTOFBit() , ph->GetDispBit() , ph->GetChargedBit()); |
477d6cee | 455 | } |
1c5acb87 | 456 | } |
457 | ||
5ae09196 | 458 | //__________________________________________________________________________ |
f2ccb5b8 | 459 | Bool_t AliCaloPID::IsTrackMatched(const AliVCluster* cluster, const AliCalorimeterUtils * cu) const { |
5ae09196 | 460 | //Check if there is any track attached to this cluster |
461 | ||
462 | Int_t nMatches = cluster->GetNTracksMatched(); | |
ae182e60 | 463 | |
464 | // if(nMatches>0){ | |
465 | // printf("N matches %d, first match (ESD) %d or (AOD) %d\n",nMatches,cluster->GetTrackMatchedIndex(), cluster->GetTrackMatched(0)); | |
466 | // if (cluster->GetTrackMatched(0)) printf("\t matched track id %d\n",((AliVTrack*)cluster->GetTrackMatched(0)) ->GetID() ) ; | |
467 | // } | |
468 | // else { | |
469 | // printf("Not Matched"); | |
470 | // } | |
471 | ||
d39cba7e | 472 | //If EMCAL track matching needs to be recalculated |
473 | if(cluster->IsEMCAL() && cu && cu->IsRecalculationOfClusterTrackMatchingOn()){ | |
474 | Float_t dR = 999., dZ = 999.; | |
475 | cu->GetEMCALRecoUtils()->GetMatchedResiduals(cluster->GetID(),dR,dZ); | |
ae182e60 | 476 | |
d39cba7e | 477 | if(dR < 999) { |
478 | ||
c76f0089 | 479 | if(fhTrackMatchedDEta){ |
480 | fhTrackMatchedDEta->Fill(cluster->E(),dZ); | |
481 | fhTrackMatchedDPhi->Fill(cluster->E(),dR); | |
482 | if(cluster->E() > 0.5)fhTrackMatchedDEtaDPhi->Fill(dZ,dR); | |
483 | } | |
d39cba7e | 484 | //printf("dR %f, dZ %f \n",dR,dZ); |
485 | return kTRUE; | |
486 | } | |
487 | else | |
488 | return kFALSE; | |
489 | }//EMCAL cluster and recalculation of matching on | |
5ae09196 | 490 | |
c76f0089 | 491 | if(fhTrackMatchedDEta){ |
492 | fhTrackMatchedDEta->Fill(cluster->GetTrackDz(),cluster->E()); | |
256822f1 | 493 | fhTrackMatchedDPhi->Fill(cluster->GetTrackDx(),cluster->E()); |
c76f0089 | 494 | if(cluster->E() > 0.5)fhTrackMatchedDEtaDPhi->Fill(cluster->GetTrackDz(),cluster->GetTrackDx()); |
495 | } | |
f2ccb5b8 | 496 | |
5ae09196 | 497 | if(!strcmp("AliESDCaloCluster",Form("%s",cluster->ClassName()))) |
d39cba7e | 498 | { |
5ae09196 | 499 | if (nMatches > 0) { |
500 | if (nMatches == 1 ) { | |
501 | Int_t iESDtrack = cluster->GetTrackMatchedIndex(); | |
502 | //printf("Track Matched index %d\n",iESDtrack); | |
503 | if(iESDtrack==-1) return kFALSE ;// Default value of array, there is no match | |
ae182e60 | 504 | else { |
505 | if(cluster->IsPHOS()) { | |
506 | Double_t r = TMath::Sqrt(cluster->GetTrackDx()*cluster->GetTrackDx()+cluster->GetTrackDz()*cluster->GetTrackDz()); | |
507 | if( r < fRcutPHOS) return kTRUE ;// Default value of array, there is no match | |
508 | else return kFALSE; | |
509 | }//PHOS | |
510 | else return kTRUE; // EMCAL and not default | |
511 | } | |
5ae09196 | 512 | }//Just one, check |
513 | else return kTRUE ;//More than one, there is a match. | |
514 | }// > 0 | |
515 | else return kFALSE; //It does not happen, but in case | |
516 | ||
517 | }//ESDs | |
518 | else | |
a5fb4114 | 519 | { |
d39cba7e | 520 | if(nMatches > 0) return kTRUE; //There is at least one match. |
521 | else return kFALSE; | |
522 | ||
5ae09196 | 523 | }//AODs or MC (copy into AOD) |
524 | ||
525 | return kFALSE; | |
526 | ||
527 | } | |
528 | ||
9a6fa057 | 529 | //__________________________________________________________________________ |
a50dabec | 530 | Bool_t AliCaloPID::IsPHOSPhoton(const Double_t l0, const Double_t l1) { |
9a6fa057 | 531 | // Check the shape of the PHOS cluster |
532 | // Return true if photon like, from Dmitri P. | |
533 | // TO DO, move parameters to data members | |
534 | ||
535 | Double_t l0Mean = 1.22 ; | |
536 | Double_t l1Mean = 2.0 ; | |
537 | Double_t l0Sigma = 0.42 ; | |
538 | Double_t l1Sigma = 0.71 ; | |
539 | Double_t c =-0.59 ; | |
540 | ||
541 | Double_t R2=(l1-l1Mean)*(l0-l0Mean)/l0Sigma/l0Sigma+(l1-l1Mean)*(l1-l1Mean)/l1Sigma/l1Sigma-c*(l0-l0Mean)*(l1-l1Mean)/l0Sigma/l1Sigma ; | |
542 | ||
543 | return (R2<9.) ; | |
544 | ||
545 | } | |
1c5acb87 | 546 |