1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 // Class for EMCAL PID
17 // Implements the abstract base class AliHFEpidBase
18 // IsInitialized() does the PID decision with energy and
19 // momentum matching (e/p)
27 #include "AliESDInputHandler.h"
28 #include "AliAODpidUtil.h"
29 #include "AliAODTrack.h"
30 #include "AliESDpid.h"
32 #include "AliHFEdetPIDqa.h"
33 //#include "AliHFEpidEMCAL.h"
34 #include "AliHFEpidQAmanager.h"
36 //#include "AliEMCALRecoUtils.h"
37 //#include "AliEMCALGeometry.h"
38 #include "AliVCluster.h"
39 #include "AliVCaloCells.h"
40 #include "AliVEvent.h"
42 #include "AliEMCALPIDUtils.h"
44 #include "AliESDEvent.h"
45 #include "AliESDtrack.h"
46 //#include "AliEMCALTrack.h"
47 //#include "AliEMCALTracker.h"
48 #include "AliHFEpidEMCAL.h"
50 ClassImp(AliHFEpidEMCAL)
52 //___________________________________________________________________
53 AliHFEpidEMCAL::AliHFEpidEMCAL():
64 //___________________________________________________________________
65 AliHFEpidEMCAL::AliHFEpidEMCAL(const Char_t *name):
76 //___________________________________________________________________
77 AliHFEpidEMCAL::AliHFEpidEMCAL(const AliHFEpidEMCAL &c):
89 //___________________________________________________________________
90 AliHFEpidEMCAL &AliHFEpidEMCAL::operator=(const AliHFEpidEMCAL &ref){
92 // Assignment operator
101 //___________________________________________________________________
102 AliHFEpidEMCAL::~AliHFEpidEMCAL(){
106 if(fPID) delete fPID;
108 //___________________________________________________________________
109 void AliHFEpidEMCAL::Copy(TObject &ref) const {
111 // Performs the copying of the object
113 AliHFEpidEMCAL &target = dynamic_cast<AliHFEpidEMCAL &>(ref);
116 target.feopMax = feopMax;
117 target.feopMim = feopMim;
119 AliHFEpidBase::Copy(ref);
121 //___________________________________________________________________
122 Bool_t AliHFEpidEMCAL::InitializePID(){
124 // InitializePID: EMCAL experts have to implement code here
129 //___________________________________________________________________
130 Int_t AliHFEpidEMCAL::IsSelected(const AliHFEpidObject *track, AliHFEpidQAmanager */*pidqa*/) const
131 { //Function to a return a code indicating whether or not an electron candidate is selected
134 if((!fESDpid && track->IsESDanalysis()) || (!fAODpid && track->IsAODanalysis())) return 0;
135 AliDebug(2, "PID object available");
136 // EMCal not fESDpid (s.s Feb. 11)
139 AliHFEpidObject::AnalysisType_t anaType = track->IsESDanalysis() ? AliHFEpidObject::kESDanalysis : AliHFEpidObject::kAODanalysis;
140 if(!((dynamic_cast<const AliVTrack *>(track->GetRecTrack()))->GetStatus() & AliESDtrack::kEMCALpid)) return 0;
141 AliDebug(2, "Track Has EMCAL PID");
143 //if(pidqa) pidqa->ProcessTrack(track, AliHFEpid::kEMCALpid, AliHFEdetPIDqa::kBeforePID);
144 // not QA for EMCal, will be added (s.s Feb. 11)
146 Double_t eop = 0.;//MomentumEnergyMatch(track->GetRecTrack()); // get eop (What is GetRecTrack ?)
147 AliDebug(2, Form("Energy - Momentum Matching e/p : %f", eop));
149 if(eop>feopMim && eop<feopMax){
151 //if(pidqa) pidqa->ProcessTrack(track, AliHFEpid::kEMCALpid, AliHFEdetPIDqa::kAfterPID);
155 pdg = 211; // EMCal doesn't separate pi,k.p by e/p. return pion code as same as TRD
158 printf("eID %g ; %d \n",eop,pdg);
163 Double_t AliHFEpidEMCAL::MomentumEnergyMatch(const AliVParticle *track) const
164 { // Returns e/p if an electron is matched
168 Double_t matchclsE = 9999.9;
170 const AliESDtrack *esdtrack = dynamic_cast<const AliESDtrack *>(track);
171 AliESDEvent *evt = esdtrack->GetESDEvent();
172 Double_t magF = evt->GetMagneticField();
173 Double_t magSign = 1.0;
174 if(magF<0)magSign = -1.0;
175 //printf("magF ; %g ; %g \n", magF,magSign);
177 if (!TGeoGlobalMagField::Instance()->GetField()) {
178 printf("Loading field map...\n");
179 //AliMagF* field = new AliMagF("Maps","Maps", 1., 1., AliMagF::k5kG);
180 AliMagF* field = new AliMagF("Maps","Maps", magSign, magSign, AliMagF::k5kG); // for 10d
181 TGeoGlobalMagField::Instance()->SetField(field);
184 AliEMCALTrack *emctrack = new AliEMCALTrack(*esdtrack);
186 emctrack->GetBxByBz(fieldB);
187 //printf("%g %g %g \n", fieldB[0], fieldB[1], fieldB[2]);
189 for(Int_t icl=0; icl<evt->GetNumberOfCaloClusters(); icl++)
192 AliVCluster *cluster = (AliVCluster*) evt->GetCaloCluster(icl);
193 if(!cluster->IsEMCAL()) continue;
195 cluster->GetPosition(clsPos);
196 if(!emctrack->PropagateToGlobal(clsPos[0],clsPos[1],clsPos[2],0.,0.) ) continue;
197 emctrack->GetXYZ(trkPos);
199 TVector3 clsPosVec(clsPos[0],clsPos[1],clsPos[2]);
200 TVector3 trkPosVec(trkPos[0],trkPos[1],trkPos[2]);
202 Double_t delEmcphi = clsPosVec.Phi()-trkPosVec.Phi(); // track cluster matching
203 Double_t delEmceta = clsPosVec.Eta()-trkPosVec.Eta(); // track cluster matching
205 double rmatch = sqrt(pow(delEmcphi,2)+pow(delEmceta,2));
209 matchclsE = cluster->E();
214 double feop = -9999.9;
215 if(matchclsE<9999) feop = matchclsE/esdtrack->P();
217 // if(feop!=-9999.9)printf("%f\n",feop) ;