/* $Id$ */
-//_________________________________________________________________________
-// Algorithm class for the identification of particles detected in EMCAL
-// base class of identified particle
-// Why should I put meaningless comments
-// just to satisfy
-// the code checker
-
-//
-//*-- Author: Yves Schutz (SUBATECH) & Dmitri Peressounko
+// Compute PID weights for all the clusters that are in AliESDs.root file
+// the AliESDs.root have to be in the same directory as the class
+//
+// and do:
+// AliEMCALPID *pid = new AliEMCALPID(kFALSE); // this calls the constructor which avoids the call to recparam
+// pid->SetReconstructor(kFALSE);
+// pid->SetPrintInfo(kTRUE);
+// pid->SetHighFluxParam(); // pid->SetLowFluxParam();
+//
+// then in cluster loop do
+// pid->ComputePID(energy, lambda0);
+//
+// Compute PID Weight for all clusters in AliESDs.root file
+// keep this function for the moment for a simple verification, could be removed
+//
+// pid->GetPIDFinal(idx) gives the probabilities
+//
+// Double_t PIDFinal[AliPID::kSPECIESN] is the standard PID for :
+//
+// kElectron : fPIDFinal[0]
+// kMuon : fPIDFinal[1]
+// kPion : fPIDFinal[2]
+// kKaon : fPIDFinal[3]
+// kProton : fPIDFinal[4]
+// kPhoton : fPIDFinal[5]
+// kPi0 : fPIDFinal[6]
+// kNeutron : fPIDFinal[7]
+// kKaon0 : fPIDFinal[8]
+// kEleCon : fPIDFinal[9]
+// kUnknown : fPIDFinal[10]
+//
+//
+// PID[3] is a simple PID for
+// Electron & Photon PID[0]
+// Pi0 PID[1]
+// Hadron PID[2]
+//
+// Author: Genole Bourdaud 2007 (SUBATECH)
+// Marie Germain 07/2009 (SUBATECH), new parametrization for low and high flux environment
+// Gustavo Conesa 08/2009 (LNF), divide class in AliEMCALPID and AliEMCALPIDUtils, PIDUtils belong to library EMCALUtils
+// --- standard c ---
+// standard C++ includes
+//#include <Riostream.h>
-// --- ROOT system ---
+// ROOT includes
-// --- Standard library ---
-
-// --- AliRoot header files ---
+// STEER includes
+#include "AliESDEvent.h"
#include "AliEMCALPID.h"
-#include "AliEMCALGetter.h"
+#include "AliESDCaloCluster.h"
+#include "AliEMCALReconstructor.h"
+
ClassImp(AliEMCALPID)
-
-//____________________________________________________________________________
- AliEMCALPID::AliEMCALPID():TTask("","")
+
+//______________________________________________
+ AliEMCALPID::AliEMCALPID()
+ : AliEMCALPIDUtils(), fReconstructor(kTRUE)
{
- // ctor
- fEventFolderName = "" ;
- fFirstEvent = 0 ;
- fLastEvent = -1 ;
+ //
+ // Constructor.
+ // Initialize all constant values which have to be used
+ // during PID algorithm execution
+ //
+
+ InitParameters();
+
+
}
+//______________________________________________
+AliEMCALPID::AliEMCALPID(Bool_t reconstructor)
+: AliEMCALPIDUtils(), fReconstructor(reconstructor)
+{
+ //
+ // Constructor.
+ // Initialize all constant values which have to be used
+ // during PID algorithm execution called when used in standalone mode
+ //
+
+ InitParameters();
+
+}
-//____________________________________________________________________________
-AliEMCALPID::AliEMCALPID(const TString alirunFileName, const TString eventFolderName):
- TTask("EMCAL"+AliConfig::Instance()->GetPIDTaskName(), alirunFileName),
- fEventFolderName(eventFolderName)
+//______________________________________________
+void AliEMCALPID::RunPID(AliESDEvent *esd)
{
- // ctor
- fFirstEvent = 0 ;
- fLastEvent = -1 ;
+ //
+ // Make the PID for all the EMCAL clusters containedin the ESDs File
+ // but just gamma/PiO/Hadron
+ //
+ // trivial check against NULL object passed
+
+ if (esd == 0x0) {
+ AliInfo("NULL ESD object passed !!" );
+ return ;
+ }
+
+ Int_t nClusters = esd->GetNumberOfCaloClusters();
+ Int_t firstCluster = 0;
+ Double_t energy=0., lambda0=0.;
+ for (Int_t iCluster = firstCluster; iCluster < (nClusters + firstCluster); iCluster++) {
+
+ AliESDCaloCluster *clust = esd->GetCaloCluster(iCluster);
+ if (!clust->IsEMCAL()) continue ;
+
+ energy = clust->E();
+ lambda0 = clust->GetM02();
+
+ if (lambda0 != 0 && energy < 1000) {
+
+ // reject clusters with lambda0 = 0
+
+
+ ComputePID(energy, lambda0);
+
+
+ if (fPrintInfo) {
+ AliInfo("___________________________________________________");
+ AliInfo(Form( "Particle Energy = %f",energy));
+ AliInfo(Form( "Particle Lambda0 of the particle = %f", lambda0) );
+ AliInfo("PIDWeight of the particle :" );
+ AliInfo(Form( " GAMMA : %f",fPID[0] ));
+ AliInfo(Form( " PiZero : %f",fPID[1] ));
+ AliInfo(Form( " HADRON : %f", fPID[2] ));
+ AliInfo("_________________________________________");
+ AliInfo(Form( " kElectron : %f", fPIDFinal[0]) );
+ AliInfo(Form( " kMuon : %f", fPIDFinal[1] ));
+ AliInfo(Form( " kPion : %f", fPIDFinal[2] ));
+ AliInfo(Form( " kKaon : %f", fPIDFinal[3] ));
+ AliInfo(Form( " kProton : %f", fPIDFinal[4] ));
+ AliInfo(Form( " kPhoton : %f", fPIDFinal[5] ));
+ AliInfo(Form( " kPi0 : %f", fPIDFinal[6] ));
+ AliInfo(Form( " kNeutron : %f", fPIDFinal[7] ));
+ AliInfo(Form( " kKaon0 : %f", fPIDFinal[8] ));
+ AliInfo(Form( " kEleCon : %f", fPIDFinal[9] ));
+ AliInfo(Form( " kUnknown : %f", fPIDFinal[10] ));
+ AliInfo("___________________________________________________");
+ }
+
+ if(fReconstructor){ // In case it is called during reconstruction.
+ // cout << "############# Fill ESDs with PIDWeight ##########" << endl;
+ clust->SetPID(fPIDFinal);}
+ } // end if (lambda0 != 0 && energy < 1000)
+ } // end for (iCluster...)
}
-//____________________________________________________________________________
-AliEMCALPID::~AliEMCALPID()
+
+//_______________________________________________________
+void AliEMCALPID::InitParameters()
{
- // dtor
- //Remove this from the parental task before destroying
- AliEMCALGetter::Instance()->EmcalLoader()->CleanPIDTask();
-
+ // Initialize PID parameters, depending on the use or not of the reconstructor
+ // and the kind of event type if the reconstructor is not used.
+ // fWeightHadronEnergy=0.;
+ // fWeightPiZeroEnergy=0.;
+ // fWeightGammaEnergy=0.;
+
+ fPIDWeight[0] = -1;
+ fPIDWeight[1] = -1;
+ fPIDWeight[2] = -1;
+
+ for(Int_t i=0; i<AliPID::kSPECIESN+1; i++)
+ fPIDFinal[i]= 0;
+
+ const AliEMCALRecParam* recParam = AliEMCALReconstructor::GetRecParam();
+
+ if(fReconstructor){
+
+ if(!recParam) {
+ AliFatal("Reconstruction parameters for EMCAL not set!");
+ }
+ else {
+
+ for(Int_t i=0; i<6; i++){
+ for(Int_t j=0; j<6; j++){
+ fGamma[i][j] = recParam->GetGamma(i,j);
+ fGamma1to10[i][j] = recParam->GetGamma1to10(i,j);
+ fHadron[i][j] = recParam->GetHadron(i,j);
+ fHadron1to10[i][j] = recParam->GetHadron1to10(i,j);
+ fPiZero[i][j] = recParam->GetPiZero(i,j);
+
+
+ // AliDebug(1,Form("PID parameters (%d, %d): fGamma=%.3f, fPi=%.3f, fHadron=%.3f",
+ // i,j, fGamma[i][j],fPiZero[i][j],fHadron[i][j] ));
+ // cout << "PID parameters (" << i << " ,"<<j<<") fGamma= "<< fGamma[i][j]<<" fPi0 ="<< fPiZero[i][j]<< endl;
+
+ } // end loop j
+ fHadronEnergyProb[i] = recParam->GetHadronEnergyProb(i);
+ fPiZeroEnergyProb[i] = recParam->GetPiZeroEnergyProb(i);
+ fGammaEnergyProb[i] = recParam->GetGammaEnergyProb(i);
+ } //end loop i
+
+
+ } // end if !recparam
+
+ }
+
+ else{
+ // init the parameters here instead of from loading from recparam
+ // default parameters are PbPb parameters.
+ SetHighFluxParam();
+
+ }
+
}
+