[u/mrichter/AliRoot.git] / PWG4 / AliGammaDataReader.cxx


/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/
/* $Id$ */

/* History of cvs commits:
 *
 * $Log$
 * Revision 1.2  2007/08/17 12:40:04  schutz
 * New analysis classes by Gustavo Conesa
 *
 * Revision 1.1.2.1  2007/07/26 10:32:09  schutz
 * new analysis classes in the the new analysis framework
 *
 *
 */

//_________________________________________________________________________
// Class for reading data (ESDs) in order to do prompt gamma correlations
//  Class created from old AliPHOSGammaJet 
//  (see AliRoot versions previous Release 4-09)
//
//*-- Author: Gustavo Conesa (LNF-INFN) 
//////////////////////////////////////////////////////////////////////////////


// --- ROOT system ---
#include "Riostream.h"
#include <TParticle.h>

//---- ANALYSIS system ----
#include "AliGammaDataReader.h" 
#include "AliLog.h"
#include "AliESDEvent.h"
#include "AliESDVertex.h"
#include "AliESDCaloCluster.h"

ClassImp(AliGammaDataReader)

//____________________________________________________________________________
AliGammaDataReader::AliGammaDataReader() : 
  AliGammaReader()
{
  //Default Ctor
  
  //Initialize parameters
  fDataType=kData;
  
}

//____________________________________________________________________________
AliGammaDataReader::AliGammaDataReader(const AliGammaDataReader & g) :   
  AliGammaReader(g)
{
  // cpy ctor
}

//_________________________________________________________________________
AliGammaDataReader & AliGammaDataReader::operator = (const AliGammaDataReader & source)
{
  // assignment operator

  if(&source == this) return *this;

  return *this;

}

//____________________________________________________________________________
void AliGammaDataReader::CreateParticleList(TObject * data, TObject *,
					    TClonesArray * plCTS, 
					    TClonesArray * plEMCAL,  
					    TClonesArray * plPHOS, 
					    TClonesArray * , 
					    TClonesArray * ,  
					    TClonesArray * ){
  
  //Create a list of particles from the ESD. These particles have been measured 
  //by the Central Tracking system (TPC+ITS), PHOS and EMCAL 
  //Also create particle list with mothers.

  AliESDEvent* esd = (AliESDEvent*) data;

  Int_t npar  = 0 ;
  Double_t *pid = new Double_t[AliPID::kSPECIESN];  
   AliDebug(3,"Fill particle lists");
  
  //Get vertex for momentum calculation  
  Double_t v[3] ; //vertex ;
  esd->GetVertex()->GetXYZ(v) ; 
  
  //########### PHOS ##############
  
  Int_t begphos = 0;  
  Int_t endphos = esd->GetNumberOfCaloClusters() ;  
  Int_t indexPH = plPHOS->GetEntries() ;
  Int_t begem = 0;
  Int_t endem = endphos;

  AliDebug(3,Form("First Calorimeter particle %d, last particle %d", begphos,endphos));
  
  for (npar =  begphos; npar <  endphos; npar++) {//////////////PHOS track loop
    AliESDCaloCluster * clus = esd->GetCaloCluster(npar) ; // retrieve track from esd

    //Check if it is PHOS cluster
    if(!clus->IsEMCAL())
      begem=npar+1;
    else
      continue;

    AliDebug(4,Form("PHOS clusters: E %f, match %d", clus->E(),clus->GetTrackMatched()));
    if(clus->GetTrackMatched()==-1){
      //Create a TParticle to fill the particle list
      TLorentzVector momentum ;
      clus->GetMomentum(momentum, v);      
      Double_t phi = momentum.Phi();
      if(phi<0) phi+=TMath::TwoPi() ;
      if(momentum.Pt() > fNeutralPtCut &&  TMath::Abs(momentum.Eta()) < fPHOSEtaCut &&
	 phi > fPhiPHOSCut[0] && phi < fPhiPHOSCut[1] ) {
      	
	pid=clus->GetPid();	
	Int_t pdg = 22;

	if(IsPHOSPIDOn()){
	  AliDebug(5,Form("E %1.2f; PID: ph %0.2f, pi0 %0.2f, el %0.2f, conv el %0.2f,pi %0.2f, k %0.2f, p %0.2f, k0 %0.2f, n %0.2f, mu %0.2f ",
			  momentum.E(),pid[AliPID::kPhoton],pid[AliPID::kPi0],pid[AliPID::kElectron],pid[AliPID::kEleCon],pid[AliPID::kPion],
			  pid[AliPID::kKaon],pid[AliPID::kProton], pid[AliPID::kKaon0],pid[AliPID::kNeutron], pid[AliPID::kMuon]));
	  
	  Float_t wPhoton =  fPHOSPhotonWeight;
	  Float_t wPi0 =  fPHOSPi0Weight;
	  
	  if(fPHOSWeightFormula){
	    wPhoton = fPHOSPhotonWeightFormula->Eval(momentum.E()) ;
	    wPi0 =    fPHOSPi0WeightFormula->Eval(momentum.E());
	  }
	  
	  if(pid[AliPID::kPhoton] > wPhoton) 
	    pdg = kPhoton ;
	  else if(pid[AliPID::kPi0] > wPi0) 
	    pdg = kPi0 ; 
	  else if(pid[AliPID::kElectron] > fPHOSElectronWeight)  
	    pdg = kElectron ;
	  else if(pid[AliPID::kEleCon] > fPHOSElectronWeight) 
	    pdg = kEleCon ;
	  else if(pid[AliPID::kPion]+pid[AliPID::kKaon]+pid[AliPID::kProton] > fPHOSChargeWeight) 
	    pdg = kChargedHadron ;  
	  else if(pid[AliPID::kKaon0]+pid[AliPID::kNeutron] > fPHOSNeutralWeight) 
	    pdg = kNeutralHadron ; 
	  
	  else if(pid[AliPID::kElectron]+pid[AliPID::kEleCon]+pid[AliPID::kPion]+pid[AliPID::kKaon]+pid[AliPID::kProton]  >  
		  pid[AliPID::kPhoton] + pid[AliPID::kPi0]+pid[AliPID::kKaon0]+pid[AliPID::kNeutron]) 
	    pdg = kChargedUnknown  ; 
	  else 
	    pdg = kNeutralUnknown ; 
	  //neutral cluster, unidentifed.
	}
	
	if(pdg != kElectron && pdg != kEleCon && pdg !=kChargedHadron && pdg !=kChargedUnknown ){//keep only neutral particles in the array
	  TParticle * particle = new TParticle(pdg, 1, -1, -1, -1, -1, 
					       momentum.Px(), momentum.Py(), momentum.Pz(), momentum.E(), v[0], v[1], v[2], 0);
	  
	  AliDebug(4,Form("PHOS added: pdg %d, pt %f, phi %f, eta %f", pdg, particle->Pt(),particle->Phi(),particle->Eta()));
	  
	  new((*plPHOS)[indexPH++])   TParticle(*particle) ;
	}
	else AliDebug(4,Form("PHOS charged cluster NOT added: pdg %d, pt %f, phi %f, eta %f\n", 
			     pdg, momentum.Pt(),momentum.Phi(),momentum.Eta()));	

      }//pt, eta, phi cut
      else 	AliDebug(4,"Particle not added");
    }//not charged
  }//cluster loop
  
  //########### CTS (TPC+ITS) #####################
  Int_t begtpc   = 0 ;  
  Int_t endtpc   = esd->GetNumberOfTracks() ;
  Int_t indexCh  = plCTS->GetEntries() ;
  AliDebug(3,Form("First CTS particle %d, last particle %d", begtpc,endtpc));
  
  for (npar =  begtpc; npar <  endtpc; npar++) {////////////// track loop
    AliESDtrack * track = esd->GetTrack(npar) ; // retrieve track from esd
    
    //We want tracks fitted in the detectors:
    ULong_t status=AliESDtrack::kTPCrefit;
    status|=AliESDtrack::kITSrefit;
    
    //We want tracks whose PID bit is set:
    //     ULong_t status =AliESDtrack::kITSpid;
    //     status|=AliESDtrack::kTPCpid;
    
    if ( (track->GetStatus() & status) == status) {//Check if the bits we want are set
      // Do something with the tracks which were successfully
      // re-fitted 
      Double_t en = 0; //track ->GetTPCsignal() ;
      Double_t mom[3];
      track->GetPxPyPz(mom) ;
      Double_t px = mom[0];
      Double_t py = mom[1];
      Double_t pz = mom[2]; //Check with TPC people if this is correct.
      Int_t pdg = 11; //Give any charged PDG code, in this case electron.
      //I just want to tag the particle as charged
       TParticle * particle = new TParticle(pdg, 1, -1, -1, -1, -1, 
					    px, py, pz, en, v[0], v[1], v[2], 0);
  
      //TParticle * particle = new TParticle() ;
      //particle->SetMomentum(px,py,pz,en) ;
       if(particle->Pt() > fChargedPtCut && TMath::Abs(particle->Eta())<fCTSEtaCut)
	 new((*plCTS)[indexCh++])       TParticle(*particle) ;    
    }
  }
  
  //################ EMCAL ##############
  
  Int_t indexEM  = plEMCAL->GetEntries() ; 
  
  for (npar =  begem; npar <  endem; npar++) {//////////////EMCAL track loop
    AliESDCaloCluster * clus = esd->GetCaloCluster(npar) ; // retrieve track from esd
    Int_t clustertype= clus->GetClusterType();
    AliDebug(4,Form("EMCAL clusters: E %f, match %d", clus->E(),clus->GetTrackMatched()));

    if(clustertype == AliESDCaloCluster::kEMCALClusterv1 && clus->GetTrackMatched()==-1 ){
      
      TLorentzVector momentum ;
      clus->GetMomentum(momentum, v); 
      Double_t phi = momentum.Phi();
      if(phi<0) phi+=TMath::TwoPi() ;
      if(momentum.Pt() > fNeutralPtCut &&  TMath::Abs(momentum.Eta()) < fEMCALEtaCut &&
	phi > fPhiEMCALCut[0] && phi < fPhiEMCALCut[1] ) {
      
	pid=clus->GetPid();	
	Int_t pdg = 22;

	if(IsEMCALPIDOn()){
	  AliDebug(5,Form("E %1.2f; PID: ph %0.2f, pi0 %0.2f, el %0.2f, conv el %0.2f,pi %0.2f, k %0.2f, p %0.2f, k0 %0.2f, n %0.2f, mu %0.2f ",
			  momentum.E(),pid[AliPID::kPhoton],pid[AliPID::kPi0],pid[AliPID::kElectron],pid[AliPID::kEleCon],pid[AliPID::kPion],
			  pid[AliPID::kKaon],pid[AliPID::kProton], pid[AliPID::kKaon0],pid[AliPID::kNeutron], pid[AliPID::kMuon]));
	  
	  if(pid[AliPID::kPhoton] > fEMCALPhotonWeight) 
	    pdg = kPhoton ;
	  else if(pid[AliPID::kPi0] > fEMCALPi0Weight) 
	    pdg = kPi0 ; 
	  else if(pid[AliPID::kElectron] > fEMCALElectronWeight)  
	    pdg = kElectron ;
	  else if(pid[AliPID::kEleCon] > fEMCALElectronWeight) 
	    pdg = kEleCon ;
	  else if(pid[AliPID::kPion]+pid[AliPID::kKaon]+pid[AliPID::kProton] > fEMCALChargeWeight) 
	    pdg = kChargedHadron ;  
	  else if(pid[AliPID::kKaon0]+pid[AliPID::kNeutron] > fEMCALNeutralWeight) 
	    pdg = kNeutralHadron ; 
	  else if(pid[AliPID::kElectron]+pid[AliPID::kEleCon]+pid[AliPID::kPion]+pid[AliPID::kKaon]+pid[AliPID::kProton]  >  
		  pid[AliPID::kPhoton] + pid[AliPID::kPi0]+pid[AliPID::kKaon0]+pid[AliPID::kNeutron]) 
	    pdg = kChargedUnknown ; 
	  else 
	    pdg = kNeutralUnknown ;
	}
	
	if(pdg != kElectron && pdg != kEleCon && pdg !=kChargedHadron && pdg !=kChargedUnknown){//keep only neutral particles in the array
	  
	  TParticle * particle = new TParticle(pdg, 1, -1, -1, -1, -1, 
					       momentum.Px(), momentum.Py(), momentum.Pz(), momentum.E(), v[0], v[1], v[2], 0);
	  AliDebug(4,Form("EMCAL cluster added: pdg %f, pt %f, phi %f, eta %f", pdg, particle->Pt(),particle->Phi(),particle->Eta()));
	  
	  new((*plEMCAL)[indexEM++])   TParticle(*particle) ;
	}
	else AliDebug(4,Form("EMCAL charged cluster NOT added: pdg %d, pt %f, phi %f, eta %f", 
			     pdg, momentum.Pt(),momentum.Phi(),momentum.Eta()));
	
      }//pt, phi, eta cut
      else 	AliDebug(4,"Particle not added");
    }//not charged, not pseudocluster
  }//Cluster loop
  
  AliDebug(3,Form("Particle lists filled, PHOS clusters %d, EMCAL clusters %d, CTS tracks %d", 
		  indexPH,indexEM,indexCh));

}
Commit	Line	Data
bdcfac30	1
	2	/**************************************************************************
	3	* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	4	* *
	5	* Author: The ALICE Off-line Project. *
	6	* Contributors are mentioned in the code where appropriate. *
	7	* *
	8	* Permission to use, copy, modify and distribute this software and its *
	9	* documentation strictly for non-commercial purposes is hereby granted *
	10	* without fee, provided that the above copyright notice appears in all *
	11	* copies and that both the copyright notice and this permission notice *
	12	* appear in the supporting documentation. The authors make no claims *
	13	* about the suitability of this software for any purpose. It is *
	14	* provided "as is" without express or implied warranty. *
	15	**************************************************************************/
	16	/* $Id$ */
	17
	18	/* History of cvs commits:
	19	*
	20	* $Log$
4b707925	21	* Revision 1.2 2007/08/17 12:40:04 schutz
	22	* New analysis classes by Gustavo Conesa
	23	*
bdcfac30	24	* Revision 1.1.2.1 2007/07/26 10:32:09 schutz
	25	* new analysis classes in the the new analysis framework
	26	*
	27	*
	28	*/
	29
	30	//_________________________________________________________________________
	31	// Class for reading data (ESDs) in order to do prompt gamma correlations
	32	// Class created from old AliPHOSGammaJet
	33	// (see AliRoot versions previous Release 4-09)
	34	//
	35	//*-- Author: Gustavo Conesa (LNF-INFN)
	36	//////////////////////////////////////////////////////////////////////////////
	37
	38
	39	// --- ROOT system ---
	40	#include "Riostream.h"
	41	#include <TParticle.h>
	42
	43	//---- ANALYSIS system ----
	44	#include "AliGammaDataReader.h"
	45	#include "AliLog.h"
	46	#include "AliESDEvent.h"
	47	#include "AliESDVertex.h"
	48	#include "AliESDCaloCluster.h"
	49
	50	ClassImp(AliGammaDataReader)
	51
	52	//____________________________________________________________________________
	53	AliGammaDataReader::AliGammaDataReader() :
4b707925	54	AliGammaReader()
bdcfac30	55	{
	56	//Default Ctor
	57
	58	//Initialize parameters
	59	fDataType=kData;
bdcfac30	60
	61	}
	62
	63	//____________________________________________________________________________
	64	AliGammaDataReader::AliGammaDataReader(const AliGammaDataReader & g) :
4b707925	65	AliGammaReader(g)
bdcfac30	66	{
	67	// cpy ctor
	68	}
	69
	70	//_________________________________________________________________________
	71	AliGammaDataReader & AliGammaDataReader::operator = (const AliGammaDataReader & source)
	72	{
	73	// assignment operator
	74
	75	if(&source == this) return *this;
	76
bdcfac30	77	return *this;
	78
	79	}
	80
	81	//____________________________________________________________________________
	82	void AliGammaDataReader::CreateParticleList(TObject * data, TObject *,
	83	TClonesArray * plCTS,
	84	TClonesArray * plEMCAL,
4b707925	85	TClonesArray * plPHOS,
	86	TClonesArray * ,
	87	TClonesArray * ,
	88	TClonesArray * ){
bdcfac30	89
	90	//Create a list of particles from the ESD. These particles have been measured
	91	//by the Central Tracking system (TPC+ITS), PHOS and EMCAL
4b707925	92	//Also create particle list with mothers.
bdcfac30	93
	94	AliESDEvent* esd = (AliESDEvent*) data;
	95
	96	Int_t npar = 0 ;
4b707925	97	Double_t *pid = new Double_t[AliPID::kSPECIESN];
bdcfac30	98	AliDebug(3,"Fill particle lists");
	99
	100	//Get vertex for momentum calculation
	101	Double_t v[3] ; //vertex ;
	102	esd->GetVertex()->GetXYZ(v) ;
	103
	104	//########### PHOS ##############
	105
	106	Int_t begphos = 0;
	107	Int_t endphos = esd->GetNumberOfCaloClusters() ;
	108	Int_t indexPH = plPHOS->GetEntries() ;
	109	Int_t begem = 0;
	110	Int_t endem = endphos;
	111
	112	AliDebug(3,Form("First Calorimeter particle %d, last particle %d", begphos,endphos));
	113
	114	for (npar = begphos; npar < endphos; npar++) {//////////////PHOS track loop
	115	AliESDCaloCluster * clus = esd->GetCaloCluster(npar) ; // retrieve track from esd
	116
	117	//Check if it is PHOS cluster
	118	if(!clus->IsEMCAL())
	119	begem=npar+1;
	120	else
	121	continue;
	122
	123	AliDebug(4,Form("PHOS clusters: E %f, match %d", clus->E(),clus->GetTrackMatched()));
	124	if(clus->GetTrackMatched()==-1){
	125	//Create a TParticle to fill the particle list
	126	TLorentzVector momentum ;
	127	clus->GetMomentum(momentum, v);
	128	Double_t phi = momentum.Phi();
	129	if(phi<0) phi+=TMath::TwoPi() ;
	130	if(momentum.Pt() > fNeutralPtCut && TMath::Abs(momentum.Eta()) < fPHOSEtaCut &&
	131	phi > fPhiPHOSCut[0] && phi < fPhiPHOSCut[1] ) {
	132
	133	pid=clus->GetPid();
	134	Int_t pdg = 22;
4b707925	135
	136	if(IsPHOSPIDOn()){
	137	AliDebug(5,Form("E %1.2f; PID: ph %0.2f, pi0 %0.2f, el %0.2f, conv el %0.2f,pi %0.2f, k %0.2f, p %0.2f, k0 %0.2f, n %0.2f, mu %0.2f ",
	138	momentum.E(),pid[AliPID::kPhoton],pid[AliPID::kPi0],pid[AliPID::kElectron],pid[AliPID::kEleCon],pid[AliPID::kPion],
	139	pid[AliPID::kKaon],pid[AliPID::kProton], pid[AliPID::kKaon0],pid[AliPID::kNeutron], pid[AliPID::kMuon]));
	140
	141	Float_t wPhoton = fPHOSPhotonWeight;
	142	Float_t wPi0 = fPHOSPi0Weight;
	143
	144	if(fPHOSWeightFormula){
	145	wPhoton = fPHOSPhotonWeightFormula->Eval(momentum.E()) ;
	146	wPi0 = fPHOSPi0WeightFormula->Eval(momentum.E());
	147	}
	148
	149	if(pid[AliPID::kPhoton] > wPhoton)
	150	pdg = kPhoton ;
	151	else if(pid[AliPID::kPi0] > wPi0)
	152	pdg = kPi0 ;
	153	else if(pid[AliPID::kElectron] > fPHOSElectronWeight)
	154	pdg = kElectron ;
	155	else if(pid[AliPID::kEleCon] > fPHOSElectronWeight)
	156	pdg = kEleCon ;
	157	else if(pid[AliPID::kPion]+pid[AliPID::kKaon]+pid[AliPID::kProton] > fPHOSChargeWeight)
	158	pdg = kChargedHadron ;
	159	else if(pid[AliPID::kKaon0]+pid[AliPID::kNeutron] > fPHOSNeutralWeight)
	160	pdg = kNeutralHadron ;
	161
	162	else if(pid[AliPID::kElectron]+pid[AliPID::kEleCon]+pid[AliPID::kPion]+pid[AliPID::kKaon]+pid[AliPID::kProton] >
	163	pid[AliPID::kPhoton] + pid[AliPID::kPi0]+pid[AliPID::kKaon0]+pid[AliPID::kNeutron])
	164	pdg = kChargedUnknown ;
	165	else
	166	pdg = kNeutralUnknown ;
	167	//neutral cluster, unidentifed.
bdcfac30	168	}
bdcfac30	169
4b707925	170	if(pdg != kElectron && pdg != kEleCon && pdg !=kChargedHadron && pdg !=kChargedUnknown ){//keep only neutral particles in the array
	171	TParticle * particle = new TParticle(pdg, 1, -1, -1, -1, -1,
	172	momentum.Px(), momentum.Py(), momentum.Pz(), momentum.E(), v[0], v[1], v[2], 0);
	173
	174	AliDebug(4,Form("PHOS added: pdg %d, pt %f, phi %f, eta %f", pdg, particle->Pt(),particle->Phi(),particle->Eta()));
	175
	176	new((plPHOS)[indexPH++]) TParticle(particle) ;
	177	}
	178	else AliDebug(4,Form("PHOS charged cluster NOT added: pdg %d, pt %f, phi %f, eta %f\n",
	179	pdg, momentum.Pt(),momentum.Phi(),momentum.Eta()));
	180
bdcfac30	181	}//pt, eta, phi cut
	182	else AliDebug(4,"Particle not added");
	183	}//not charged
	184	}//cluster loop
	185
	186	//########### CTS (TPC+ITS) #####################
	187	Int_t begtpc = 0 ;
	188	Int_t endtpc = esd->GetNumberOfTracks() ;
	189	Int_t indexCh = plCTS->GetEntries() ;
	190	AliDebug(3,Form("First CTS particle %d, last particle %d", begtpc,endtpc));
	191
	192	for (npar = begtpc; npar < endtpc; npar++) {////////////// track loop
	193	AliESDtrack * track = esd->GetTrack(npar) ; // retrieve track from esd
	194
	195	//We want tracks fitted in the detectors:
	196	ULong_t status=AliESDtrack::kTPCrefit;
	197	status\|=AliESDtrack::kITSrefit;
	198
	199	//We want tracks whose PID bit is set:
	200	// ULong_t status =AliESDtrack::kITSpid;
	201	// status\|=AliESDtrack::kTPCpid;
	202
	203	if ( (track->GetStatus() & status) == status) {//Check if the bits we want are set
	204	// Do something with the tracks which were successfully
	205	// re-fitted
	206	Double_t en = 0; //track ->GetTPCsignal() ;
	207	Double_t mom[3];
	208	track->GetPxPyPz(mom) ;
	209	Double_t px = mom[0];
	210	Double_t py = mom[1];
	211	Double_t pz = mom[2]; //Check with TPC people if this is correct.
	212	Int_t pdg = 11; //Give any charged PDG code, in this case electron.
	213	//I just want to tag the particle as charged
	214	TParticle * particle = new TParticle(pdg, 1, -1, -1, -1, -1,
	215	px, py, pz, en, v[0], v[1], v[2], 0);
	216
	217	//TParticle * particle = new TParticle() ;
	218	//particle->SetMomentum(px,py,pz,en) ;
	219	if(particle->Pt() > fChargedPtCut && TMath::Abs(particle->Eta())<fCTSEtaCut)
	220	new((plCTS)[indexCh++]) TParticle(particle) ;
	221	}
	222	}
	223
	224	//################ EMCAL ##############
	225
	226	Int_t indexEM = plEMCAL->GetEntries() ;
bdcfac30	227
	228	for (npar = begem; npar < endem; npar++) {//////////////EMCAL track loop
	229	AliESDCaloCluster * clus = esd->GetCaloCluster(npar) ; // retrieve track from esd
	230	Int_t clustertype= clus->GetClusterType();
	231	AliDebug(4,Form("EMCAL clusters: E %f, match %d", clus->E(),clus->GetTrackMatched()));
	232
4b707925	233	if(clustertype == AliESDCaloCluster::kEMCALClusterv1 && clus->GetTrackMatched()==-1 ){
bdcfac30	234
	235	TLorentzVector momentum ;
	236	clus->GetMomentum(momentum, v);
	237	Double_t phi = momentum.Phi();
	238	if(phi<0) phi+=TMath::TwoPi() ;
	239	if(momentum.Pt() > fNeutralPtCut && TMath::Abs(momentum.Eta()) < fEMCALEtaCut &&
	240	phi > fPhiEMCALCut[0] && phi < fPhiEMCALCut[1] ) {
	241
	242	pid=clus->GetPid();
	243	Int_t pdg = 22;
bdcfac30	244
4b707925	245	if(IsEMCALPIDOn()){
	246	AliDebug(5,Form("E %1.2f; PID: ph %0.2f, pi0 %0.2f, el %0.2f, conv el %0.2f,pi %0.2f, k %0.2f, p %0.2f, k0 %0.2f, n %0.2f, mu %0.2f ",
	247	momentum.E(),pid[AliPID::kPhoton],pid[AliPID::kPi0],pid[AliPID::kElectron],pid[AliPID::kEleCon],pid[AliPID::kPion],
	248	pid[AliPID::kKaon],pid[AliPID::kProton], pid[AliPID::kKaon0],pid[AliPID::kNeutron], pid[AliPID::kMuon]));
	249
	250	if(pid[AliPID::kPhoton] > fEMCALPhotonWeight)
	251	pdg = kPhoton ;
	252	else if(pid[AliPID::kPi0] > fEMCALPi0Weight)
	253	pdg = kPi0 ;
	254	else if(pid[AliPID::kElectron] > fEMCALElectronWeight)
	255	pdg = kElectron ;
	256	else if(pid[AliPID::kEleCon] > fEMCALElectronWeight)
	257	pdg = kEleCon ;
	258	else if(pid[AliPID::kPion]+pid[AliPID::kKaon]+pid[AliPID::kProton] > fEMCALChargeWeight)
	259	pdg = kChargedHadron ;
	260	else if(pid[AliPID::kKaon0]+pid[AliPID::kNeutron] > fEMCALNeutralWeight)
	261	pdg = kNeutralHadron ;
	262	else if(pid[AliPID::kElectron]+pid[AliPID::kEleCon]+pid[AliPID::kPion]+pid[AliPID::kKaon]+pid[AliPID::kProton] >
	263	pid[AliPID::kPhoton] + pid[AliPID::kPi0]+pid[AliPID::kKaon0]+pid[AliPID::kNeutron])
	264	pdg = kChargedUnknown ;
	265	else
	266	pdg = kNeutralUnknown ;
	267	}
bdcfac30	268
4b707925	269	if(pdg != kElectron && pdg != kEleCon && pdg !=kChargedHadron && pdg !=kChargedUnknown){//keep only neutral particles in the array
	270
	271	TParticle * particle = new TParticle(pdg, 1, -1, -1, -1, -1,
	272	momentum.Px(), momentum.Py(), momentum.Pz(), momentum.E(), v[0], v[1], v[2], 0);
	273	AliDebug(4,Form("EMCAL cluster added: pdg %f, pt %f, phi %f, eta %f", pdg, particle->Pt(),particle->Phi(),particle->Eta()));
	274
	275	new((plEMCAL)[indexEM++]) TParticle(particle) ;
	276	}
	277	else AliDebug(4,Form("EMCAL charged cluster NOT added: pdg %d, pt %f, phi %f, eta %f",
	278	pdg, momentum.Pt(),momentum.Phi(),momentum.Eta()));
bdcfac30	279
bdcfac30	280	}//pt, phi, eta cut
	281	else AliDebug(4,"Particle not added");
	282	}//not charged, not pseudocluster
	283	}//Cluster loop
	284
4b707925	285	AliDebug(3,Form("Particle lists filled, PHOS clusters %d, EMCAL clusters %d, CTS tracks %d",
4b707925	286	indexPH,indexEM,indexCh));
bdcfac30	287
	288	}
	289