-//
-//
-void AliEMCALPID::ComputePID(Double_t energy, Double_t lambda0)
-{
-//
-// This is the main command, which uses the distributions computed and parametrised,
-// and gives the PID by the bayesian method.
-//
- TArrayD paramDistribGamma = DistLambda0(energy, 1);
- TArrayD paramDistribPiZero = DistLambda0(energy, 2);
- TArrayD paramDistribHadron = DistLambda0(energy, 3);
-
- Bool_t norm = kFALSE;
-
- fProbGamma = TMath::Gaus(lambda0, paramDistribGamma[1], paramDistribGamma[2], norm) * paramDistribGamma[0];
- fProbGamma += TMath::Landau(lambda0, paramDistribGamma[4], paramDistribGamma[5], norm) * paramDistribGamma[3];
- fProbPiZero = TMath::Gaus(lambda0, paramDistribPiZero[1], paramDistribPiZero[2], norm) * paramDistribPiZero[0];
- fProbPiZero += TMath::Landau(lambda0, paramDistribPiZero[4], paramDistribPiZero[5], norm) * paramDistribPiZero[3];
- fProbHadron = TMath::Gaus(lambda0, paramDistribHadron[1], paramDistribHadron[2], norm) * paramDistribHadron[0];
- fProbHadron += TMath::Landau(lambda0, paramDistribHadron[4], paramDistribHadron[5], norm) * paramDistribHadron[3];
-
- // compute PID Weight
- fPIDWeight[0] = fProbGamma / (fProbGamma + fProbPiZero + fProbHadron);
- fPIDWeight[1] = fProbPiZero / (fProbGamma+fProbPiZero+fProbHadron);
- fPIDWeight[2] = fProbHadron / (fProbGamma+fProbPiZero+fProbHadron);
-
- SetPID(fPIDWeight[0], 0);
- SetPID(fPIDWeight[1], 1);
- SetPID(fPIDWeight[2], 2);
-
- // sortie ecran pid Weight only for control (= in english ???)
- if (fPrintInfo) {
- AliInfo(Form( "Energy in loop = %f", energy) );
- AliInfo(Form( "Lambda0 in loop = %f", lambda0) );
- AliInfo(Form( "fProbGamma in loop = %f", fProbGamma) );
- // AliInfo(Form( "fParametresDistribGamma[2] = %f", fParamDistribGamma[2]) );
- AliInfo(Form( "fProbaPiZero = %f", fProbPiZero ));
- AliInfo(Form( "fProbaHadron = %f", fProbHadron) );
- AliInfo(Form( "PIDWeight in loop = %f ||| %f ||| %f", fPIDWeight[0] , fPIDWeight[1], fPIDWeight[2]) );
- AliInfo(Form( "fGamma[2][2] = %f", fGamma[2][2] ));
- AliInfo("********************************************************" );
- }
-
- fPIDFinal[0] = fPIDWeight[0]/2;
- fPIDFinal[1] = fPIDWeight[2]/8;
- fPIDFinal[2] = fPIDWeight[2]/8;
- fPIDFinal[3] = fPIDWeight[2]/8;
- fPIDFinal[4] = fPIDWeight[2]/8;
- fPIDFinal[5] = fPIDWeight[0]/2;
- fPIDFinal[6] = fPIDWeight[1]/2;
- fPIDFinal[7] = fPIDWeight[2]/8;
- fPIDFinal[8] = fPIDWeight[2]/8;
- fPIDFinal[9] = fPIDWeight[2]/8;
- fPIDFinal[10] = fPIDWeight[2]/8;
- fPIDFinal[11] = 0;
-}
-//
-//
-TArrayD AliEMCALPID::DistLambda0(Double_t energy, Int_t type)
-{
-//
-// Compute the values of the parametrised distributions using the data initialised before.
-//
- Double_t constGauss = 0., meanGauss = 0., sigmaGauss = 0.;
- Double_t constLandau=0., mpvLandau=0., sigmaLandau=0.;
- TArrayD distributionParam(6);
-
- switch (type) {
- case 1:
- constGauss = Polynomial(energy, fGamma[0]);
- meanGauss = Polynomial(energy, fGamma[1]);
- sigmaGauss = Polynomial(energy, fGamma[2]);
- constLandau = Polynomial(energy, fGamma[3]);
- mpvLandau = Polynomial(energy, fGamma[4]);
- sigmaLandau = Polynomial(energy, fGamma[5]);
- break;
- case 2:
- if (energy < 10) {
- constGauss = Polynomial(energy, fPiZero5to10[0]);
- meanGauss = Polynomial(energy, fPiZero5to10[1]);
- sigmaGauss = Polynomial(energy, fPiZero5to10[2]);
- constLandau = Polynomial(energy, fPiZero5to10[3]);
- mpvLandau = Polynomial(energy, fPiZero5to10[4]);
- sigmaLandau = Polynomial(energy, fPiZero5to10[5]);
- }
- else {
- constGauss = Polynomial(energy, fPiZero10to60[0]);
- meanGauss = Polynomial(energy, fPiZero10to60[1]);
- sigmaGauss = Polynomial(energy, fPiZero10to60[2]);
- constLandau = Polynomial(energy, fPiZero10to60[3]);
- mpvLandau = Polynomial(energy, fPiZero10to60[4]);
- sigmaLandau = Polynomial(energy, fPiZero10to60[5]);
- }
- break;
- case 3:
- constGauss = Polynomial(energy, fHadron[0]);
- meanGauss = Polynomial(energy, fHadron[1]);
- sigmaGauss = Polynomial(energy, fHadron[2]);
- constLandau = Polynomial(energy, fHadron[3]);
- mpvLandau = Polynomial(energy, fHadron[4]);
- sigmaLandau = Polynomial(energy, fHadron[5]);
- break;
- }
-
- distributionParam[0] = constGauss;
- distributionParam[1] = meanGauss;
- distributionParam[2] = sigmaGauss;
- distributionParam[3] = constLandau;
- distributionParam[4] = mpvLandau;
- distributionParam[5] = sigmaLandau;
-
- return distributionParam;
-}
-//
-//
-Double_t AliEMCALPID::Polynomial(Double_t x, Double_t *params)