]> git.uio.no Git - u/mrichter/AliRoot.git/blame - EMCAL/AliEMCALPID.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Stuff from aldaqdqm09
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALPID.cxx
CommitLineData
dc293ae9 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
16/* $Id$ */
1e7c9b89 17
18// Compute PID weights for all the clusters that are in AliESDs.root file
19// the AliESDs.root have to be in the same directory as the class
dc293ae9 20//
1e7c9b89 21// and do:
22// AliEMCALPID *pid = new AliEMCALPID(kFALSE); // this calls the constructor which avoids the call to recparam
23// pid->SetReconstructor(kFALSE);
24// pid->SetPrintInfo(kTRUE);
25// pid->SetHighFluxParam(); // pid->SetLowFluxParam();
26//
27// then in cluster loop do
28// pid->ComputePID(energy, lambda0);
29//
30// Compute PID Weight for all clusters in AliESDs.root file
dc293ae9 31// keep this function for the moment for a simple verification, could be removed
32//
1e7c9b89 33// pid->GetPIDFinal(idx) gives the probabilities
dc293ae9 34//
1e7c9b89 35// Double_t PIDFinal[AliPID::kSPECIESN] is the standard PID for :
dc293ae9 36//
37//
38//
39// kElectron : fPIDFinal[0]
40// kMuon : fPIDFinal[1]
41// kPion : fPIDFinal[2]
42// kKaon : fPIDFinal[3]
43// kProton : fPIDFinal[4]
44// kPhoton : fPIDFinal[5]
45// kPi0 : fPIDFinal[6]
46// kNeutron : fPIDFinal[7]
47// kKaon0 : fPIDFinal[8]
48// kEleCon : fPIDFinal[9]
49// kUnknown : fPIDFinal[10]
50//
51//
52// PID[3] is a simple PID for
53// Electron & Photon PID[0]
54// Pi0 PID[1]
55// Hadron PID[2]
56//
1e7c9b89 57// --- standard c ---
dc293ae9 58
59// standard C++ includes
1e7c9b89 60//#include <Riostream.h>
dc293ae9 61
62// ROOT includes
1e7c9b89 63//#include "TTree.h"
64//#include "TVector3.h"
65//#include "TBranch.h"
66//#include "TClonesArray.h"
67//#include "TLorentzVector.h"
dc293ae9 68#include "TMath.h"
1e7c9b89 69//#include "TRefArray.h"
70#include "TArrayD.h"
dc293ae9 71
46363b15 72// STEER includes
1e7c9b89 73#include "AliESDEvent.h"
74//#include "AliLog.h"
dc293ae9 75#include "AliEMCALPID.h"
89ffc0b0 76#include "AliESDCaloCluster.h"
1e7c9b89 77//#include "AliEMCALRecParam.h"
8ba062b1 78#include "AliEMCALReconstructor.h"
79
dc293ae9 80
81ClassImp(AliEMCALPID)
1e7c9b89 82
4cfd8ef9 83//______________________________________________
84 AliEMCALPID::AliEMCALPID():
1e7c9b89 85 fPrintInfo(kFALSE), fProbGamma(0.),fProbPiZero(0.),fProbHadron(0.), fWeightHadronEnergy(1.), fWeightGammaEnergy(1.),fWeightPiZeroEnergy(1.),fReconstructor(kTRUE)
dc293ae9 86{
4cfd8ef9 87 //
88 // Constructor.
89 // Initialize all constant values which have to be used
90 // during PID algorithm execution
91 //
1e7c9b89 92
93 InitParameters();
94
95
96}
46363b15 97
1e7c9b89 98//______________________________________________
99AliEMCALPID::AliEMCALPID(Bool_t reconstructor):
100 fPrintInfo(kFALSE), fProbGamma(0.),fProbPiZero(0.),fProbHadron(0.), fWeightHadronEnergy(1.), fWeightGammaEnergy(1.),fWeightPiZeroEnergy(1.),fReconstructor(reconstructor)
101{
102 //
103 // Constructor.
104 // Initialize all constant values which have to be used
105 // during PID algorithm execution called when used in standalone mode
106 //
107
108 InitParameters();
46363b15 109
dc293ae9 110}
4cfd8ef9 111
112//______________________________________________
af885e0f 113void AliEMCALPID::RunPID(AliESDEvent *esd)
dc293ae9 114{
1e7c9b89 115 //
116 // Make the PID for all the EMCAL clusters containedin the ESDs File
117 // but just gamma/PiO/Hadron
118 //
119 // trivial check against NULL object passed
120
dc293ae9 121 if (esd == 0x0) {
4ce73a00 122 AliInfo("NULL ESD object passed !!" );
dc293ae9 123 return ;
124 }
1e7c9b89 125
3a2a23e1 126 Int_t nClusters = esd->GetNumberOfCaloClusters();
127 Int_t firstCluster = 0;
dc293ae9 128 Double_t energy, lambda0;
129 for (Int_t iCluster = firstCluster; iCluster < (nClusters + firstCluster); iCluster++) {
1e7c9b89 130
dc293ae9 131 AliESDCaloCluster *clust = esd->GetCaloCluster(iCluster);
1e7c9b89 132 if (!clust->IsEMCAL()) continue ;
133
7592dfc4 134 energy = clust->E();
dc293ae9 135 lambda0 = clust->GetM02();
136 // verify cluster type
137 Int_t clusterType= clust->GetClusterType();
8ada0ffe 138 if (clusterType == AliESDCaloCluster::kEMCALClusterv1 && lambda0 != 0 && energy < 1000) {
1e7c9b89 139
140 // if (lambda0 != 0 && energy < 1000) {
141
dc293ae9 142 // reject clusters with lambda0 = 0
1e7c9b89 143
144
dc293ae9 145 ComputePID(energy, lambda0);
1e7c9b89 146
147
dc293ae9 148 if (fPrintInfo) {
149 AliInfo("___________________________________________________");
150 AliInfo(Form( "Particle Energy = %f",energy));
151 AliInfo(Form( "Particle Lambda0 of the particle = %f", lambda0) );
152 AliInfo("PIDWeight of the particle :" );
153 AliInfo(Form( " GAMMA : %f",fPID[0] ));
154 AliInfo(Form( " PiZero : %f",fPID[1] ));
155 AliInfo(Form( " HADRON : %f", fPID[2] ));
156 AliInfo("_________________________________________");
157 AliInfo(Form( " kElectron : %f", fPIDFinal[0]) );
158 AliInfo(Form( " kMuon : %f", fPIDFinal[1] ));
159 AliInfo(Form( " kPion : %f", fPIDFinal[2] ));
160 AliInfo(Form( " kKaon : %f", fPIDFinal[3] ));
161 AliInfo(Form( " kProton : %f", fPIDFinal[4] ));
162 AliInfo(Form( " kPhoton : %f", fPIDFinal[5] ));
163 AliInfo(Form( " kPi0 : %f", fPIDFinal[6] ));
164 AliInfo(Form( " kNeutron : %f", fPIDFinal[7] ));
165 AliInfo(Form( " kKaon0 : %f", fPIDFinal[8] ));
166 AliInfo(Form( " kEleCon : %f", fPIDFinal[9] ));
167 AliInfo(Form( " kUnknown : %f", fPIDFinal[10] ));
168 AliInfo("___________________________________________________");
169 }
1e7c9b89 170
171 if(fReconstructor){ // In case it is called during reconstruction.
172 // cout << "#############On remplit l esd avec les PIDWeight##########" << endl;
173 clust->SetPid(fPIDFinal);}
dc293ae9 174 } // end if (clusterType...)
175 } // end for (iCluster...)
176}
4cfd8ef9 177
178//__________________________________________________________
dc293ae9 179void AliEMCALPID::ComputePID(Double_t energy, Double_t lambda0)
180{
181//
182// This is the main command, which uses the distributions computed and parametrised,
183// and gives the PID by the bayesian method.
184//
1e7c9b89 185// cout << "ENERGY " <<energy <<" lambda0 "<< lambda0<< endl;
186
187 Double_t weightGammaEnergy = DistEnergy(energy, 1);
188 Double_t weightPiZeroEnergy = DistEnergy(energy, 2);
189 Double_t weightHadronEnergy = DistEnergy(energy, 3);
190
191 //Double_t weightHadronEnergy = 1.;
192
193 Double_t energyhadron=energy;
194 if(energyhadron<1.)energyhadron=1.; // no energy dependance of parametrisation for hadrons below 1 GeV
195 if (energy<2){energy =2;} // no energy dependance of parametrisation for gamma and pi0 below 2 GeV
196
197 if (energy>55){
198 energy =55.;
199 energyhadron=55.;
200 } // same parametrisation for gamma and hadrons above 55 GeV
201 // for the pi0 above 55GeV the 2 gammas supperposed no way to distinguish from real gamma PIDWeight[1]=0
202
4cfd8ef9 203 TArrayD paramDistribGamma = DistLambda0(energy, 1);
204 TArrayD paramDistribPiZero = DistLambda0(energy, 2);
1e7c9b89 205 TArrayD paramDistribHadron = DistLambda0(energyhadron, 3);
4cfd8ef9 206
207 Bool_t norm = kFALSE;
208
1e7c9b89 209
4cfd8ef9 210 fProbGamma = TMath::Gaus(lambda0, paramDistribGamma[1], paramDistribGamma[2], norm) * paramDistribGamma[0];
1e7c9b89 211 fProbGamma += TMath::Landau(((1-paramDistribGamma[4])-lambda0),paramDistribGamma[4],paramDistribGamma[5],norm)* paramDistribGamma[3];
212 if(fProbGamma<0.)fProbGamma=0.;
213
214 fProbGamma = fProbGamma*weightGammaEnergy;
215
216 if(energy>10. || energy < 55.){
217 fProbPiZero = TMath::Gaus(lambda0, paramDistribPiZero[1], paramDistribPiZero[2], norm) * paramDistribPiZero[0];
218 fProbPiZero += TMath::Landau(lambda0, paramDistribPiZero[4], paramDistribPiZero[5], norm) * paramDistribPiZero[3];
219 if(fProbPiZero<0. || energy<5.)fProbPiZero=0.;
220 fProbPiZero = fProbPiZero*weightPiZeroEnergy;
221 }
222 else {
223 fProbPiZero = 0.;
224 }
225
4cfd8ef9 226 fProbHadron = TMath::Gaus(lambda0, paramDistribHadron[1], paramDistribHadron[2], norm) * paramDistribHadron[0];
227 fProbHadron += TMath::Landau(lambda0, paramDistribHadron[4], paramDistribHadron[5], norm) * paramDistribHadron[3];
1e7c9b89 228 if(fProbHadron<0.)fProbHadron=0.;
229 fProbHadron = fProbHadron*weightHadronEnergy; // to take into account the probability for a hadron to have a given reconstructed energy
4cfd8ef9 230
231 // compute PID Weight
1e7c9b89 232 if( (fProbGamma + fProbPiZero + fProbHadron)>0.){
233 fPIDWeight[0] = fProbGamma / (fProbGamma + fProbPiZero + fProbHadron);
234 fPIDWeight[1] = fProbPiZero / (fProbGamma+fProbPiZero+fProbHadron);
235 fPIDWeight[2] = fProbHadron / (fProbGamma+fProbPiZero+fProbHadron);
236 }
237 else{
238// cases where energy and lambda0 large, probably du to 2 clusters folded the clusters PID not assigned to hadron nor Pi0 nor gammas
239 fPIDWeight[0] = 0.;
240 fPIDWeight[1] = 0.;
241 fPIDWeight[2] = 0.;
242 }
243
244
245 // cout << " PID[0] "<< fPIDWeight[0] << " PID[1] "<< fPIDWeight[1] << " PID[2] "<< fPIDWeight[2] << endl;
4cfd8ef9 246
247 SetPID(fPIDWeight[0], 0);
248 SetPID(fPIDWeight[1], 1);
249 SetPID(fPIDWeight[2], 2);
250
1e7c9b89 251 // print pid Weight only for control (= in english ???)
4cfd8ef9 252 if (fPrintInfo) {
253 AliInfo(Form( "Energy in loop = %f", energy) );
254 AliInfo(Form( "Lambda0 in loop = %f", lambda0) );
255 AliInfo(Form( "fProbGamma in loop = %f", fProbGamma) );
4cfd8ef9 256 AliInfo(Form( "fProbaPiZero = %f", fProbPiZero ));
257 AliInfo(Form( "fProbaHadron = %f", fProbHadron) );
258 AliInfo(Form( "PIDWeight in loop = %f ||| %f ||| %f", fPIDWeight[0] , fPIDWeight[1], fPIDWeight[2]) );
4cfd8ef9 259 AliInfo("********************************************************" );
260 }
261
1e7c9b89 262 fPIDFinal[0] = fPIDWeight[0]/2; // photon
4cfd8ef9 263 fPIDFinal[1] = fPIDWeight[2]/8;
264 fPIDFinal[2] = fPIDWeight[2]/8;
265 fPIDFinal[3] = fPIDWeight[2]/8;
266 fPIDFinal[4] = fPIDWeight[2]/8;
1e7c9b89 267 fPIDFinal[5] = fPIDWeight[0]/2; // electron
268 fPIDFinal[6] = fPIDWeight[1] ; // Pi0
4cfd8ef9 269 fPIDFinal[7] = fPIDWeight[2]/8;
270 fPIDFinal[8] = fPIDWeight[2]/8;
271 fPIDFinal[9] = fPIDWeight[2]/8;
272 fPIDFinal[10] = fPIDWeight[2]/8;
1e7c9b89 273
dc293ae9 274}
4cfd8ef9 275
1e7c9b89 276
277
278
4cfd8ef9 279//________________________________________________________
1e7c9b89 280TArrayD AliEMCALPID::DistLambda0(const Double_t energy, const Int_t type)
dc293ae9 281{
4cfd8ef9 282 //
283 // Compute the values of the parametrised distributions using the data initialised before.
284 //
285 Double_t constGauss = 0., meanGauss = 0., sigmaGauss = 0.;
286 Double_t constLandau=0., mpvLandau=0., sigmaLandau=0.;
287 TArrayD distributionParam(6);
288
289 switch (type) {
1e7c9b89 290
4cfd8ef9 291 case 1:
46363b15 292
1e7c9b89 293 constGauss = PolynomialMixed2(energy, fGamma[0]);
294 meanGauss = PolynomialMixed2(energy, fGamma[1]);
295 sigmaGauss = PolynomialMixed2(energy, fGamma[2]);
296 constLandau = PolynomialMixed2(energy, fGamma[3]);
297 mpvLandau = PolynomialMixed2(energy, fGamma[4]);
298 sigmaLandau = PolynomialMixed2(energy, fGamma[5]);
299 break;
46363b15 300
4cfd8ef9 301 case 2:
46363b15 302
1e7c9b89 303 constGauss = PolynomialMixed2(energy, fPiZero[0]);
304 meanGauss = PolynomialMixed2(energy, fPiZero[1]);
305 sigmaGauss = PolynomialMixed2(energy, fPiZero[2]);
306 constLandau = PolynomialMixed2(energy, fPiZero[3]);
307 mpvLandau = PolynomialMixed2(energy, fPiZero[4]);
308 sigmaLandau = PolynomialMixed2(energy, fPiZero[5]);
309
4cfd8ef9 310 break;
311 case 3:
1e7c9b89 312
313 constGauss = PolynomialMixed2(energy, fHadron[0]);
314 meanGauss = PolynomialMixed2(energy, fHadron[1]);
315 sigmaGauss = PolynomialMixed2(energy, fHadron[2]);
316 constLandau = PolynomialMixed2(energy, fHadron[3]);
317 mpvLandau = PolynomialMixed2(energy, fHadron[4]);
318 sigmaLandau = PolynomialMixed2(energy, fHadron[5]);
46363b15 319
4cfd8ef9 320 break;
321 }
322
323 distributionParam[0] = constGauss;
324 distributionParam[1] = meanGauss;
325 distributionParam[2] = sigmaGauss;
326 distributionParam[3] = constLandau;
327 distributionParam[4] = mpvLandau;
328 distributionParam[5] = sigmaLandau;
329
330 return distributionParam;
dc293ae9 331}
4cfd8ef9 332
1e7c9b89 333//________________________________________________________
334Double_t AliEMCALPID::DistEnergy(const Double_t energy, const Int_t type)
335{
336 //
337 // Compute the values of the weigh for a given energy the parametrised distribution using the data initialised before.
338 //
339 Double_t constante = 0.;
340 Double_t energyParam;
341
342 switch (type) {
343
344 case 1:
345 constante = 1.;
346 break;
347 case 2:
348 constante = 1.;
349 break;
350 case 3:
351 constante = PowerExp(energy, fHadronEnergyProb);
352 break;
353 }
354
355 energyParam = constante;
356
357 // // cout << "Weight " << constante << " for energy "<< energy<< " GeV "<< endl;
358
359 return energyParam;
360}
361
362
4cfd8ef9 363//_______________________________________________________
1e7c9b89 364Double_t AliEMCALPID::Polynomial(const Double_t x, const Double_t *params) const
dc293ae9 365{
4cfd8ef9 366 //
367 // Compute a polynomial for a given value of 'x'
368 // with the array of parameters passed as the second arg
369 //
46363b15 370
4cfd8ef9 371 Double_t y;
372 y = params[0];
373 y += params[1] * x;
374 y += params[2] * x * x;
375 y += params[3] * x * x * x;
376 y += params[4] * x * x * x * x;
377 y += params[5] * x * x * x * x * x;
378
379 return y;
dc293ae9 380}
1e7c9b89 381//_______________________________________________________
382Double_t AliEMCALPID::Polynomial0(const Double_t *params) const
383{
384 //
385 // Compute a polynomial for a given value of 'x'
386 // with the array of parameters passed as the second arg
387 //
388
389 Double_t y;
390 y = params[0];
391 return y;
392}
393
394//_______________________________________________________
395Double_t AliEMCALPID::Polynomialinv(const Double_t x, const Double_t *params) const
396{
397 //
398 // Compute a polynomial for a given value of 'x'
399 // with the array of parameters passed as the second arg
400 //
401
402 Double_t y;
403 if(x>0){
404 y = params[0];
405 y += params[1] / x;
406 y += params[2] / (x * x);
407 y += params[3] / (x * x * x);
408 y += params[4] / (x * x * x * x);
409 y += params[5] / (x * x * x * x * x);
410 }
411 else
412 y=0.;
413 return y;
414
415}
416//_______________________________________________________
417Double_t AliEMCALPID::PolynomialMixed1(const Double_t x, const Double_t *params) const
418{
419 //
420 // Compute a polynomial for a given value of 'x'
421 // with the array of parameters passed as the second arg
422 //
423
424 Double_t y;
425 if(x>0){
426 y = params[0] / x;
427 y += params[1] ;
428 y += params[2] * x ;
429 // y += params[3] * 0.;
430 // y += params[4] * 0.;
431 // y += params[5] * 0.;
432 }
433 else
434 y=0.;
435
436 return y;
437
438}
439
440//_______________________________________________________
441Double_t AliEMCALPID::PolynomialMixed2(const Double_t x, const Double_t *params) const
442{
443 //
444 // Compute a polynomial for a given value of 'x'
445 // with the array of parameters passed as the second arg
446 //
447
448 Double_t y;
449 if(x>0){
450 y = params[0] / ( x * x);
451 y += params[1] / x;
452 y += params[2] ;
453 y += params[3] * x ;
454 y += params[4] * x * x ;
455 // y += params[5] * 0.;
456 }
457 else
458 y=0.;
459 // cout << "y = " << y << endl;
460 return y;
461
462}
463
464//_______________________________________________________
465Double_t AliEMCALPID::PowerExp(const Double_t x, const Double_t *params) const
466{
467 //
468 // Compute a polynomial for a given value of 'x'
469 // with the array of parameters passed as the second arg
470 // par[0]*TMath::Power(x[0],par[1])
471 // par[0]*TMath::Exp((x[0]-par[1])*par[2]);
472
473 Double_t y;
474
475 y = params[0] *TMath::Power( x,params[1]);
476 y += params[2] *TMath::Exp((x-params[3])*params[4]);
477
478 return y;
479
480}
481
482
483//_______________________________________________________
484void AliEMCALPID::InitParameters()
485{
486 // Initialize PID parameters, depending on the use or not of the reconstructor
487 // and the kind of event type if the reconstructor is not used.
488 // fWeightHadronEnergy=0.;
489 // fWeightPiZeroEnergy=0.;
490 // fWeightGammaEnergy=0.;
491
492 fPIDWeight[0] = -1;
493 fPIDWeight[1] = -1;
494 fPIDWeight[2] = -1;
495
496 for(Int_t i=0; i<AliPID::kSPECIESN+1; i++)
497 fPIDFinal[i]= 0;
498
499 const AliEMCALRecParam* recParam = AliEMCALReconstructor::GetRecParam();
500
501 if(fReconstructor){
502
503 if(!recParam) {
504 AliFatal("Reconstruction parameters for EMCAL not set!");
505 }
506 else {
507
508 for(Int_t i=0; i<6; i++){
509 for(Int_t j=0; j<6; j++){
510 fGamma[i][j] = recParam->GetGamma(i,j);
511 fGamma1to10[i][j] = recParam->GetGamma1to10(i,j);
512 fHadron[i][j] = recParam->GetHadron(i,j);
513 fHadron1to10[i][j] = recParam->GetHadron1to10(i,j);
514 fPiZero[i][j] = recParam->GetPiZero(i,j);
515
516
517 // AliDebug(1,Form("PID parameters (%d, %d): fGamma=%.3f, fPi=%.3f, fHadron=%.3f",
518 // i,j, fGamma[i][j],fPiZero[i][j],fHadron[i][j] ));
519 // cout << "PID parameters (" << i << " ,"<<j<<") fGamma= "<< fGamma[i][j]<<" fPi0 ="<< fPiZero[i][j]<< endl;
520
521 } // end loop j
522 fHadronEnergyProb[i] = recParam->GetHadronEnergyProb(i);
523 fPiZeroEnergyProb[i] = recParam->GetPiZeroEnergyProb(i);
524 fGammaEnergyProb[i] = recParam->GetGammaEnergyProb(i);
525 } //end loop i
526
527
528 } // end if !recparam
529
530 }
531
532 else{
533 // init the parameters here instead of from loading from recparam
534 // default parameters are PbPb parameters.
535 SetHighFluxParam();
536
537 }
538
539}
540
541
542//_______________________________________________________
543void AliEMCALPID::SetLowFluxParam()
544{
545
546 // as a first step, all array elements are initialized to 0.0
547 Int_t i, j;
548
549 for (i = 0; i < 6; i++) {
550 for (j = 0; j < 6; j++) {
551 fGamma[i][j] = fHadron[i][j] = fPiZero[i][j] = 0.;
552 fGamma1to10[i][j] = fHadron1to10[i][j] = 0.;
553 }
554 fGammaEnergyProb[i] = fGammaEnergyProb[i];
555 fPiZeroEnergyProb[i] = fPiZeroEnergyProb[i];
556 fHadronEnergyProb[i] = fHadronEnergyProb[i];
557 }
558
559 // New parametrisation for lambda0^2 (=x): f(x) = normLandau*TMath::Landau(x,mpvLandau,widthLandau)+normgaus*TMath::Gaus(x,meangaus,sigmagaus)
560 // See AliEMCALPid (index j) refers to the polynomial parameters of the fit of each parameter vs energy
561 // pp
562
563 // paramtype[0][j] = norm gauss
564 // paramtype[1][j] = mean gaus
565 // paramtype[2][j] = sigma gaus
566 // paramtype[3][j] = norm landau
567 // paramtype[4][j] = mpv landau
568 // paramtype[5][j] = sigma landau
569
570 fGamma[0][0] = -7.656908e-01;
571 fGamma[0][1] = 2.352536e-01;
572 fGamma[0][2] = 1.555996e-02;
573 fGamma[0][3] = 2.243525e-04;
574 fGamma[0][4] = -2.560087e-06;
575
576 fGamma[1][0] = 6.500216e+00;
577 fGamma[1][1] = -2.564958e-01;
578 fGamma[1][2] = 1.967894e-01;
579 fGamma[1][3] = -3.982273e-04;
580 fGamma[1][4] = 2.797737e-06;
581
582 fGamma[2][0] = 2.416489e+00;
583 fGamma[2][1] = -1.601258e-01;
584 fGamma[2][2] = 3.126839e-02;
585 fGamma[2][3] = 3.387532e-04;
586 fGamma[2][4] = -4.089145e-06;
587
588 fGamma[3][0] = 0.;
589 fGamma[3][1] = -2.696008e+00;
590 fGamma[3][2] = 6.920305e-01;
591 fGamma[3][3] = -2.281122e-03;
592 fGamma[3][4] = 0.;
593
594 fGamma[4][0] = 2.281564e-01;
595 fGamma[4][1] = -7.575040e-02;
596 fGamma[4][2] = 3.813423e-01;
597 fGamma[4][3] = -1.243854e-04;
598 fGamma[4][4] = 1.232045e-06;
599
600 fGamma[5][0] = -3.290107e-01;
601 fGamma[5][1] = 3.707545e-02;
602 fGamma[5][2] = 2.917397e-03;
603 fGamma[5][3] = 4.695306e-05;
604 fGamma[5][4] = -3.572981e-07;
605
606 fHadron[0][0] = 9.482243e-01;
607 fHadron[0][1] = -2.780896e-01;
608 fHadron[0][2] = 2.223507e-02;
609 fHadron[0][3] = 7.294263e-04;
610 fHadron[0][4] = -5.665872e-06;
611
612 fHadron[1][0] = 0.;
613 fHadron[1][1] = 0.;
614 fHadron[1][2] = 2.483298e-01;
615 fHadron[1][3] = 0.;
616 fHadron[1][4] = 0.;
617
618 fHadron[2][0] = -5.601199e+00;
619 fHadron[2][1] = 2.097382e+00;
620 fHadron[2][2] = -2.307965e-01;
621 fHadron[2][3] = 9.206871e-03;
622 fHadron[2][4] = -8.887548e-05;
623
624 fHadron[3][0] = 6.543101e+00;
625 fHadron[3][1] = -2.305203e+00;
626 fHadron[3][2] = 2.761673e-01;
627 fHadron[3][3] = -5.465855e-03;
628 fHadron[3][4] = 2.784329e-05;
629
630 fHadron[4][0] = -2.443530e+01;
631 fHadron[4][1] = 8.902578e+00 ;
632 fHadron[4][2] = -5.265901e-01;
633 fHadron[4][3] = 2.549111e-02;
634 fHadron[4][4] = -2.196801e-04;
635
636 fHadron[5][0] = 2.102007e-01;
637 fHadron[5][1] = -3.844418e-02;
638 fHadron[5][2] = 1.234682e-01;
639 fHadron[5][3] = -3.866733e-03;
640 fHadron[5][4] = 3.362719e-05 ;
641
642 fPiZero[0][0] = 5.072157e-01;
643 fPiZero[0][1] = -5.352747e-01;
644 fPiZero[0][2] = 8.499259e-02;
645 fPiZero[0][3] = -3.687401e-03;
646 fPiZero[0][4] = 5.482280e-05;
647
648 fPiZero[1][0] = 4.590137e+02;
649 fPiZero[1][1] = -7.079341e+01;
650 fPiZero[1][2] = 4.990735e+00;
651 fPiZero[1][3] = -1.241302e-01;
652 fPiZero[1][4] = 1.065772e-03;
653
654 fPiZero[2][0] = 1.376415e+02;
655 fPiZero[2][1] = -3.031577e+01;
656 fPiZero[2][2] = 2.474338e+00;
657 fPiZero[2][3] = -6.903410e-02;
658 fPiZero[2][4] = 6.244089e-04;
659
660 fPiZero[3][0] = 0.;
661 fPiZero[3][1] = 1.145983e+00;
662 fPiZero[3][2] = -2.476052e-01;
663 fPiZero[3][3] = 1.367373e-02;
664 fPiZero[3][4] = 0.;
665
666 fPiZero[4][0] = -2.097586e+02;
667 fPiZero[4][1] = 6.300800e+01;
668 fPiZero[4][2] = -4.038906e+00;
669 fPiZero[4][3] = 1.088543e-01;
670 fPiZero[4][4] = -9.362485e-04;
671
672 fPiZero[5][0] = -1.671477e+01;
673 fPiZero[5][1] = 2.995415e+00;
674 fPiZero[5][2] = -6.040360e-02;
675 fPiZero[5][3] = -6.137459e-04;
676 fPiZero[5][4] = 1.847328e-05;
677
678 fHadronEnergyProb[0] = 4.767543e-02;
679 fHadronEnergyProb[1] = -1.537523e+00;
680 fHadronEnergyProb[2] = 2.956727e-01;
681 fHadronEnergyProb[3] = -3.051022e+01;
682 fHadronEnergyProb[4] =-6.036931e-02;
683
684 Int_t ii= 0;
685 Int_t jj= 3;
686 AliDebug(1,Form("PID parameters (%d, %d): fGamma=%.3f, fPi=%.3f, fHadron=%.3f",
687 ii,jj, fGamma[ii][jj],fPiZero[ii][jj],fHadron[ii][jj] ));
688 //cout << " LowFlux Parameters fGamma [2][2] = " << fGamma[2][2] << endl;
689 //cout << " LowFlux Parameters fHadron [2][2] = " << fHadron[2][2] << endl;
690
691 // end for proton-proton
692
693}
694
695//_______________________________________________________
696void AliEMCALPID::SetHighFluxParam()
697{
698
699 // as a first step, all array elements are initialized to 0.0
700 Int_t i, j;
701 for (i = 0; i < 6; i++) {
702 for (j = 0; j < 6; j++) {
703 fGamma[i][j] = fHadron[i][j] = fPiZero[i][j] = 0.;
704 fGamma1to10[i][j] = fHadron1to10[i][j] = 0.;
705 }
706 fGammaEnergyProb[i] = 0.;
707 fPiZeroEnergyProb[i] = 0.;
708 fHadronEnergyProb[i] = 0.;
709 }
710
711 // Pb Pb this goes with inverted landau + gaussian for gammas, landau+gaussian for Pi0 and hadrons
712
713 fGamma[0][0] = -7.656908e-01;
714 fGamma[0][1] = 2.352536e-01;
715 fGamma[0][2] = 1.555996e-02;
716 fGamma[0][3] = 2.243525e-04;
717 fGamma[0][4] = -2.560087e-06;
718
719 fGamma[1][0] = 6.500216e+00;
720 fGamma[1][1] = -2.564958e-01;
721 fGamma[1][2] = 1.967894e-01;
722 fGamma[1][3] = -3.982273e-04;
723 fGamma[1][4] = 2.797737e-06;
724
725 fGamma[2][0] = 2.416489e+00;
726 fGamma[2][1] = -1.601258e-01;
727 fGamma[2][2] = 3.126839e-02;
728 fGamma[2][3] = 3.387532e-04;
729 fGamma[2][4] = -4.089145e-06;
730
731 fGamma[3][0] = 0.;
732 fGamma[3][1] = -2.696008e+00;
733 fGamma[3][2] = 6.920305e-01;
734 fGamma[3][3] = -2.281122e-03;
735 fGamma[3][4] = 0.;
736
737 fGamma[4][0] = 2.281564e-01;
738 fGamma[4][1] = -7.575040e-02;
739 fGamma[4][2] = 3.813423e-01;
740 fGamma[4][3] = -1.243854e-04;
741 fGamma[4][4] = 1.232045e-06;
742
743 fGamma[5][0] = -3.290107e-01;
744 fGamma[5][1] = 3.707545e-02;
745 fGamma[5][2] = 2.917397e-03;
746 fGamma[5][3] = 4.695306e-05;
747 fGamma[5][4] = -3.572981e-07;
748
749 fHadron[0][0] = 1.519112e-01;
750 fHadron[0][1] = -8.267603e-02;
751 fHadron[0][2] = 1.914574e-02;
752 fHadron[0][3] = -2.677921e-04;
753 fHadron[0][4] = 5.447939e-06;
754
755 fHadron[1][0] = 0.;
756 fHadron[1][1] = -7.549870e-02;
757 fHadron[1][2] = 3.930087e-01;
758 fHadron[1][3] = -2.368500e-03;
759 fHadron[1][4] = 0.;
760
761 fHadron[2][0] = 0.;
762 fHadron[2][1] = -2.463152e-02;
763 fHadron[2][2] = 1.349257e-01;
764 fHadron[2][3] = -1.089440e-03;
765 fHadron[2][4] = 0.;
766
767 fHadron[3][0] = 0.;
768 fHadron[3][1] = 5.101560e-01;
769 fHadron[3][2] = 1.458679e-01;
770 fHadron[3][3] = 4.903068e-04;
771 fHadron[3][4] = 0.;
772
773 fHadron[4][0] = 0.;
774 fHadron[4][1] = -6.693943e-03;
775 fHadron[4][2] = 2.444753e-01;
776 fHadron[4][3] = -5.553749e-05;
777 fHadron[4][4] = 0.;
778
779 fHadron[5][0] = -4.414030e-01;
780 fHadron[5][1] = 2.292277e-01;
781 fHadron[5][2] = -2.433737e-02;
782 fHadron[5][3] = 1.758422e-03;
783 fHadron[5][4] = -3.001493e-05;
784
785 fPiZero[0][0] = 5.072157e-01;
786 fPiZero[0][1] = -5.352747e-01;
787 fPiZero[0][2] = 8.499259e-02;
788 fPiZero[0][3] = -3.687401e-03;
789 fPiZero[0][4] = 5.482280e-05;
790
791 fPiZero[1][0] = 4.590137e+02;
792 fPiZero[1][1] = -7.079341e+01;
793 fPiZero[1][2] = 4.990735e+00;
794 fPiZero[1][3] = -1.241302e-01;
795 fPiZero[1][4] = 1.065772e-03;
796
797 fPiZero[2][0] = 1.376415e+02;
798 fPiZero[2][1] = -3.031577e+01;
799 fPiZero[2][2] = 2.474338e+00;
800 fPiZero[2][3] = -6.903410e-02;
801 fPiZero[2][4] = 6.244089e-04;
802
803 fPiZero[3][0] = 0.;
804 fPiZero[3][1] = 1.145983e+00;
805 fPiZero[3][2] = -2.476052e-01;
806 fPiZero[3][3] = 1.367373e-02;
807 fPiZero[3][4] = 0.;
808
809 fPiZero[4][0] = -2.097586e+02;
810 fPiZero[4][1] = 6.300800e+01;
811 fPiZero[4][2] = -4.038906e+00;
812 fPiZero[4][3] = 1.088543e-01;
813 fPiZero[4][4] = -9.362485e-04;
814
815 fPiZero[5][0] = -1.671477e+01;
816 fPiZero[5][1] = 2.995415e+00;
817 fPiZero[5][2] = -6.040360e-02;
818 fPiZero[5][3] = -6.137459e-04;
819 fPiZero[5][4] = 1.847328e-05;
820
821 // those are the High Flux PbPb ones
822 fHadronEnergyProb[0] = 0.;
823 fHadronEnergyProb[1] = 0.;
824 fHadronEnergyProb[2] = 6.188452e-02;
825 fHadronEnergyProb[3] = 2.030230e+00;
826 fHadronEnergyProb[4] = -6.402242e-02;
827
828 Int_t ii= 0;
829 Int_t jj= 3;
830 AliDebug(1,Form("PID parameters (%d, %d): fGamma=%.3f, fPi=%.3f, fHadron=%.3f",
831 ii,jj, fGamma[ii][jj],fPiZero[ii][jj],fHadron[ii][jj] ));
832 //cout << " HighFlux Parameters fGamma [2][2] = " << fGamma[2][2] << endl;
833 //cout << " HighFlux Parameters fHadron [2][2] = " << fHadron[2][2] << endl;
834
835}