]>
Commit | Line | Data |
---|---|---|
f9cea31c | 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 | /* $Id$ */ | |
16 | ||
17 | /* History of cvs commits: | |
18 | * | |
19 | * $Log$ | |
84e10d55 | 20 | * Revision 1.3 2007/03/08 10:24:32 schutz |
21 | * Coding convention | |
22 | * | |
463ee300 | 23 | * Revision 1.2 2007/02/09 18:40:40 schutz |
24 | * B\bNew version from Gustavo | |
25 | * | |
2a1d8a29 | 26 | * Revision 1.1 2007/01/23 17:17:29 schutz |
27 | * New Gamma package | |
28 | * | |
f9cea31c | 29 | * |
30 | */ | |
31 | ||
32 | //_________________________________________________________________________ | |
33 | // Class for the analysis of gamma correlations (gamma-jet, | |
34 | // gamma-hadron and isolation cut. | |
35 | // This class contains 3 main methods: one to fill lists of particles (ESDs) comming | |
36 | // from the CTS (ITS+TPC) and the calorimeters; the other one tags a candidate | |
37 | // cluster as isolated; the last one search in the | |
38 | // corresponing calorimeter for the highest energy cluster, identified it as | |
39 | // prompt photon; | |
40 | // | |
41 | // Class created from old AliPHOSGammaJet | |
42 | // (see AliRoot versions previous Release 4-09) | |
43 | // | |
44 | //*-- Author: Gustavo Conesa (LNF-INFN) | |
45 | ////////////////////////////////////////////////////////////////////////////// | |
46 | ||
47 | ||
48 | // --- ROOT system --- | |
49 | ||
50 | #include <TFile.h> | |
51 | #include <TParticle.h> | |
52 | #include <TH2.h> | |
2a1d8a29 | 53 | #include <TChain.h> |
f9cea31c | 54 | #include "AliAnaGammaDirect.h" |
55 | #include "AliESD.h" | |
56 | #include "AliESDtrack.h" | |
57 | #include "AliESDCaloCluster.h" | |
58 | #include "Riostream.h" | |
59 | #include "AliLog.h" | |
60 | ||
61 | ClassImp(AliAnaGammaDirect) | |
62 | ||
63 | //____________________________________________________________________________ | |
64 | AliAnaGammaDirect::AliAnaGammaDirect(const char *name) : | |
65 | AliAnalysisTask(name,""), fChain(0), fESD(0), | |
66 | fOutputContainer(new TObjArray(100)), | |
67 | fPrintInfo(0), fMinGammaPt(0.), | |
68 | fCalorimeter(""), fEMCALPID(0),fPHOSPID(0), | |
84e10d55 | 69 | fEMCALPhotonWeight(0.), fEMCALPi0Weight(0.), fPHOSPhotonWeight(0.), |
f9cea31c | 70 | fConeSize(0.),fPtThreshold(0.),fPtSumThreshold(0), |
463ee300 | 71 | fMakeICMethod(0),fhNGamma(0),fhPhiGamma(0),fhEtaGamma(0) |
f9cea31c | 72 | { |
463ee300 | 73 | //Ctor |
74 | ||
75 | //Initialize parameters | |
76 | InitParameters(); | |
77 | ||
f9cea31c | 78 | // Input slot #0 works with an Ntuple |
79 | DefineInput(0, TChain::Class()); | |
80 | // Output slot #0 writes into a TH1 container | |
81 | DefineOutput(0, TObjArray::Class()) ; | |
82 | ||
83 | } | |
84 | ||
85 | ||
86 | //____________________________________________________________________________ | |
87 | AliAnaGammaDirect::AliAnaGammaDirect(const AliAnaGammaDirect & g) : | |
88 | AliAnalysisTask(g), fChain(g.fChain), fESD(g.fESD), | |
89 | fOutputContainer(g. fOutputContainer), fPrintInfo(g.fPrintInfo), | |
90 | fMinGammaPt(g.fMinGammaPt), fCalorimeter(g.fCalorimeter), | |
91 | fEMCALPID(g.fEMCALPID),fPHOSPID(g.fPHOSPID), | |
84e10d55 | 92 | fEMCALPhotonWeight(g.fEMCALPhotonWeight), |
93 | fEMCALPi0Weight(g.fEMCALPi0Weight), | |
94 | fPHOSPhotonWeight(g.fPHOSPhotonWeight), | |
f9cea31c | 95 | fConeSize(g.fConeSize), |
96 | fPtThreshold(g.fPtThreshold), | |
97 | fPtSumThreshold(g.fPtSumThreshold), | |
463ee300 | 98 | fMakeICMethod(g.fMakeICMethod), |
99 | fhNGamma(g.fhNGamma),fhPhiGamma(g.fhPhiGamma),fhEtaGamma(g.fhEtaGamma) | |
f9cea31c | 100 | { |
101 | // cpy ctor | |
102 | SetName (g.GetName()) ; | |
103 | SetTitle(g.GetTitle()) ; | |
f9cea31c | 104 | } |
105 | ||
463ee300 | 106 | //_________________________________________________________________________ |
107 | AliAnaGammaDirect & AliAnaGammaDirect::operator = (const AliAnaGammaDirect & source) | |
108 | { | |
109 | // assignment operator | |
110 | ||
111 | if(&source == this) return *this; | |
112 | ||
113 | fChain = source.fChain ; | |
114 | fESD = source.fESD ; | |
115 | fOutputContainer = source.fOutputContainer ; | |
116 | fPrintInfo = source.fPrintInfo ; | |
117 | fMinGammaPt = source.fMinGammaPt ; | |
118 | fCalorimeter = source. fCalorimeter ; | |
119 | fEMCALPID = source.fEMCALPID ; | |
120 | fPHOSPID = source.fPHOSPID ; | |
84e10d55 | 121 | fEMCALPhotonWeight = source. fEMCALPhotonWeight ; |
122 | fEMCALPi0Weight = source.fEMCALPi0Weight ; | |
123 | fPHOSPhotonWeight = source.fPHOSPhotonWeight ; | |
463ee300 | 124 | fConeSize = source.fConeSize ; |
125 | fPtThreshold = source.fPtThreshold ; | |
126 | fPtSumThreshold = source.fPtSumThreshold ; | |
127 | fMakeICMethod = source.fMakeICMethod ; | |
128 | fhNGamma = source.fhNGamma ; | |
129 | fhPhiGamma = source.fhPhiGamma ; | |
130 | fhEtaGamma = source.fhEtaGamma ; | |
131 | ||
84e10d55 | 132 | |
463ee300 | 133 | return *this; |
134 | ||
135 | } | |
136 | ||
f9cea31c | 137 | //____________________________________________________________________________ |
138 | AliAnaGammaDirect::~AliAnaGammaDirect() | |
139 | { | |
140 | // Remove all pointers | |
141 | fOutputContainer->Clear() ; | |
142 | delete fOutputContainer ; | |
143 | delete fhNGamma ; | |
144 | delete fhPhiGamma ; | |
145 | delete fhEtaGamma ; | |
2a1d8a29 | 146 | |
147 | } | |
148 | ||
149 | //______________________________________________________________________________ | |
150 | void AliAnaGammaDirect::ConnectInputData(const Option_t*) | |
151 | { | |
152 | // Initialisation of branch container and histograms | |
153 | ||
154 | AliInfo(Form("*** Initialization of %s", GetName())) ; | |
155 | ||
156 | // Get input data | |
157 | fChain = dynamic_cast<TChain *>(GetInputData(0)) ; | |
158 | if (!fChain) { | |
159 | AliError(Form("Input 0 for %s not found\n", GetName())); | |
160 | return ; | |
161 | } | |
162 | ||
163 | // One should first check if the branch address was taken by some other task | |
164 | char ** address = (char **)GetBranchAddress(0, "ESD"); | |
165 | if (address) { | |
166 | fESD = (AliESD*)(*address); | |
167 | } else { | |
168 | fESD = new AliESD(); | |
169 | SetBranchAddress(0, "ESD", &fESD); | |
170 | } | |
171 | ||
172 | } | |
173 | ||
174 | //________________________________________________________________________ | |
175 | void AliAnaGammaDirect::CreateOutputObjects() | |
176 | { | |
177 | ||
463ee300 | 178 | // Create histograms to be saved in output file and |
179 | // store them in fOutputContainer | |
2a1d8a29 | 180 | |
181 | fOutputContainer = new TObjArray(3) ; | |
182 | ||
183 | //Histograms of highest gamma identified in Event | |
184 | fhNGamma = new TH1F("NGamma","Number of #gamma over PHOS",240,0,120); | |
185 | fhNGamma->SetYTitle("N"); | |
186 | fhNGamma->SetXTitle("p_{T #gamma}(GeV/c)"); | |
187 | fOutputContainer->Add(fhNGamma) ; | |
188 | ||
189 | fhPhiGamma = new TH2F | |
190 | ("PhiGamma","#phi_{#gamma}",200,0,120,200,0,7); | |
191 | fhPhiGamma->SetYTitle("#phi"); | |
192 | fhPhiGamma->SetXTitle("p_{T #gamma} (GeV/c)"); | |
193 | fOutputContainer->Add(fhPhiGamma) ; | |
194 | ||
195 | fhEtaGamma = new TH2F | |
196 | ("EtaGamma","#phi_{#gamma}",200,0,120,200,-0.8,0.8); | |
197 | fhEtaGamma->SetYTitle("#eta"); | |
198 | fhEtaGamma->SetXTitle("p_{T #gamma} (GeV/c)"); | |
199 | fOutputContainer->Add(fhEtaGamma) ; | |
f9cea31c | 200 | |
201 | } | |
202 | ||
203 | //____________________________________________________________________________ | |
204 | void AliAnaGammaDirect::CreateParticleList(TClonesArray * pl, | |
205 | TClonesArray * plCTS, | |
206 | TClonesArray * plEMCAL, | |
207 | TClonesArray * plPHOS){ | |
208 | ||
209 | //Create a list of particles from the ESD. These particles have been measured | |
210 | //by the Central Tracking system (TPC+ITS), PHOS and EMCAL | |
211 | ||
212 | Int_t index = pl->GetEntries() ; | |
213 | Int_t npar = 0 ; | |
214 | Float_t *pid = new Float_t[AliPID::kSPECIESN]; | |
215 | AliDebug(3,"Fill particle lists"); | |
216 | if(fPrintInfo) | |
217 | AliInfo(Form("fCalorimeter %s",fCalorimeter.Data())); | |
218 | ||
219 | Double_t v[3] ; //vertex ; | |
220 | fESD->GetVertex()->GetXYZ(v) ; | |
221 | ||
222 | //########### PHOS ############## | |
223 | if(fCalorimeter == "PHOS"){ | |
224 | ||
225 | Int_t begphos = fESD->GetFirstPHOSCluster(); | |
226 | Int_t endphos = fESD->GetFirstPHOSCluster() + | |
227 | fESD->GetNumberOfPHOSClusters() ; | |
228 | Int_t indexNePHOS = plPHOS->GetEntries() ; | |
229 | AliDebug(3,Form("First PHOS particle %d, last particle %d", begphos,endphos)); | |
230 | ||
231 | if(fCalorimeter == "PHOS"){ | |
232 | for (npar = begphos; npar < endphos; npar++) {//////////////PHOS track loop | |
233 | AliESDCaloCluster * clus = fESD->GetCaloCluster(npar) ; // retrieve track from esd | |
234 | ||
235 | //Create a TParticle to fill the particle list | |
463ee300 | 236 | TLorentzVector momentum ; |
237 | clus->GetMomentum(momentum); | |
f9cea31c | 238 | TParticle * particle = new TParticle() ; |
463ee300 | 239 | //particle->SetMomentum(px,py,pz,en) ; |
240 | particle->SetMomentum(momentum) ; | |
241 | ||
f9cea31c | 242 | AliDebug(4,Form("PHOS clusters: pt %f, phi %f, eta %f", particle->Pt(),particle->Phi(),particle->Eta())); |
243 | ||
244 | //Select only photons | |
245 | ||
246 | pid=clus->GetPid(); | |
247 | //cout<<"pid "<<pid[AliPID::kPhoton]<<endl ; | |
248 | if( !fPHOSPID) | |
249 | new((*plPHOS)[indexNePHOS++]) TParticle(*particle) ; | |
84e10d55 | 250 | else if( pid[AliPID::kPhoton] > fPHOSPhotonWeight) |
f9cea31c | 251 | new((*plPHOS)[indexNePHOS++]) TParticle(*particle) ; |
252 | } | |
253 | } | |
254 | } | |
255 | ||
463ee300 | 256 | //########## ##################### |
257 | //Prepare bool array for EMCAL track matching | |
258 | ||
259 | // step 1, set the flag in a way that it rejects all not-V1 clusters | |
260 | // but at this point all V1 clusters are accepted | |
261 | ||
262 | Int_t begem = fESD->GetFirstEMCALCluster(); | |
263 | Int_t endem = fESD->GetFirstEMCALCluster() + | |
264 | fESD->GetNumberOfEMCALClusters() ; | |
265 | ||
84e10d55 | 266 | // if(endem < begem+12) |
267 | // AliError("Number of pseudoclusters smaller than 12"); | |
463ee300 | 268 | Bool_t *useCluster = new Bool_t[endem+1]; |
269 | ||
84e10d55 | 270 | // for (npar = 0; npar < endem; npar++){ |
271 | // if(npar < begem+12) | |
272 | // useCluster[npar] =kFALSE; //EMCAL Pseudoclusters and PHOS clusters | |
273 | // else | |
274 | // useCluster[npar] =kTRUE; //EMCAL clusters | |
275 | // } | |
276 | for (npar = 0; npar < endem; npar++) | |
277 | useCluster[npar] =kFALSE; //EMCAL Pseudoclusters and clusters | |
463ee300 | 278 | |
f9cea31c | 279 | //########### CTS (TPC+ITS) ##################### |
280 | Int_t begtpc = 0 ; | |
281 | Int_t endtpc = fESD->GetNumberOfTracks() ; | |
282 | Int_t indexCh = plCTS->GetEntries() ; | |
283 | AliDebug(3,Form("First CTS particle %d, last particle %d", begtpc,endtpc)); | |
463ee300 | 284 | |
285 | Int_t iemcalMatch = -1 ; | |
f9cea31c | 286 | for (npar = begtpc; npar < endtpc; npar++) {////////////// track loop |
287 | AliESDtrack * track = fESD->GetTrack(npar) ; // retrieve track from esd | |
463ee300 | 288 | |
289 | // step 2 for EMCAL matching, change the flag for all matched clusters found in tracks | |
290 | iemcalMatch = track->GetEMCALcluster(); | |
84e10d55 | 291 | if(iemcalMatch > 0) useCluster[iemcalMatch] = kTRUE; // reject matched cluster |
463ee300 | 292 | |
f9cea31c | 293 | //We want tracks fitted in the detectors: |
294 | ULong_t status=AliESDtrack::kTPCrefit; | |
295 | status|=AliESDtrack::kITSrefit; | |
296 | ||
297 | //We want tracks whose PID bit is set: | |
298 | // ULong_t status =AliESDtrack::kITSpid; | |
299 | // status|=AliESDtrack::kTPCpid; | |
300 | ||
301 | if ( (track->GetStatus() & status) == status) {//Check if the bits we want are set | |
302 | // Do something with the tracks which were successfully | |
303 | // re-fitted | |
304 | Double_t en = 0; //track ->GetTPCsignal() ; | |
305 | Double_t mom[3]; | |
306 | track->GetPxPyPz(mom) ; | |
307 | Double_t px = mom[0]; | |
308 | Double_t py = mom[1]; | |
309 | Double_t pz = mom[2]; //Check with TPC people if this is correct. | |
310 | Int_t pdg = 11; //Give any charged PDG code, in this case electron. | |
311 | //I just want to tag the particle as charged | |
312 | TParticle * particle = new TParticle(pdg, 1, -1, -1, -1, -1, | |
463ee300 | 313 | px, py, pz, en, v[0], v[1], v[2], 0); |
f9cea31c | 314 | |
315 | //TParticle * particle = new TParticle() ; | |
316 | //particle->SetMomentum(px,py,pz,en) ; | |
317 | ||
318 | new((*plCTS)[indexCh++]) TParticle(*particle) ; | |
319 | new((*pl)[index++]) TParticle(*particle) ; | |
320 | } | |
321 | } | |
322 | ||
323 | //################ EMCAL ############## | |
324 | ||
f9cea31c | 325 | Int_t indexNe = plEMCAL->GetEntries() ; |
326 | ||
327 | AliDebug(3,Form("First EMCAL particle %d, last particle %d",begem,endem)); | |
328 | ||
329 | for (npar = begem; npar < endem; npar++) {//////////////EMCAL track loop | |
330 | AliESDCaloCluster * clus = fESD->GetCaloCluster(npar) ; // retrieve track from esd | |
331 | Int_t clustertype= clus->GetClusterType(); | |
84e10d55 | 332 | if(clustertype == AliESDCaloCluster::kClusterv1 && !useCluster[npar] ){ |
463ee300 | 333 | TLorentzVector momentum ; |
334 | clus->GetMomentum(momentum); | |
335 | TParticle * particle = new TParticle() ; | |
336 | //particle->SetMomentum(px,py,pz,en) ; | |
337 | particle->SetMomentum(momentum) ; | |
84e10d55 | 338 | cout<<"GOOD EMCAL "<<particle->Pt()<<endl; |
f9cea31c | 339 | pid=clus->GetPid(); |
340 | if(fCalorimeter == "EMCAL") | |
341 | { | |
342 | TParticle * particle = new TParticle() ; | |
463ee300 | 343 | //particle->SetMomentum(px,py,pz,en) ; |
f9cea31c | 344 | AliDebug(4,Form("EMCAL clusters: pt %f, phi %f, eta %f", particle->Pt(),particle->Phi(),particle->Eta())); |
345 | if(!fEMCALPID) //Only identified particles | |
346 | new((*plEMCAL)[indexNe++]) TParticle(*particle) ; | |
84e10d55 | 347 | else if(pid[AliPID::kPhoton] > fEMCALPhotonWeight) |
f9cea31c | 348 | new((*plEMCAL)[indexNe++]) TParticle(*particle) ; |
349 | } | |
350 | else | |
351 | { | |
463ee300 | 352 | Int_t pdg = 0; |
353 | if(fEMCALPID) | |
354 | { | |
84e10d55 | 355 | if( pid[AliPID::kPhoton] > fEMCALPhotonWeight) |
463ee300 | 356 | pdg = 22; |
84e10d55 | 357 | else if( pid[AliPID::kPi0] > fEMCALPi0Weight) |
463ee300 | 358 | pdg = 111; |
359 | } | |
360 | else | |
361 | pdg = 22; //No PID, assume all photons | |
362 | ||
363 | TParticle * particle = new TParticle(pdg, 1, -1, -1, -1, -1, | |
364 | momentum.Px(), momentum.Py(), momentum.Pz(), momentum.E(), v[0], v[1], v[2], 0); | |
365 | AliDebug(4,Form("EMCAL clusters: pt %f, phi %f, eta %f", particle->Pt(),particle->Phi(),particle->Eta())); | |
366 | ||
367 | new((*plEMCAL)[indexNe++]) TParticle(*particle) ; | |
368 | new((*pl)[index++]) TParticle(*particle) ; | |
f9cea31c | 369 | } |
370 | } | |
371 | } | |
463ee300 | 372 | |
f9cea31c | 373 | AliDebug(3,"Particle lists filled"); |
374 | ||
375 | } | |
376 | ||
377 | ||
378 | ||
379 | //____________________________________________________________________________ | |
380 | void AliAnaGammaDirect::Exec(Option_t *) | |
381 | { | |
382 | ||
383 | // Processing of one event | |
2a1d8a29 | 384 | |
f9cea31c | 385 | //Get ESDs |
386 | Long64_t entry = fChain->GetReadEntry() ; | |
387 | ||
388 | if (!fESD) { | |
389 | AliError("fESD is not connected to the input!") ; | |
390 | return ; | |
391 | } | |
392 | ||
393 | if (fPrintInfo) | |
394 | AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(fChain))->GetFile()->GetName(), entry)) ; | |
395 | ||
396 | //CreateTLists with arrays of TParticles. Filled with particles only relevant for the analysis. | |
397 | ||
398 | TClonesArray * particleList = new TClonesArray("TParticle",1000); // All particles refitted in CTS and detected in EMCAL (jet) | |
399 | TClonesArray * plCTS = new TClonesArray("TParticle",1000); // All particles refitted in Central Tracking System (ITS+TPC) | |
400 | TClonesArray * plNe = new TClonesArray("TParticle",1000); // All particles measured in Jet Calorimeter (EMCAL) | |
401 | TClonesArray * plPHOS = new TClonesArray("TParticle",1000); // All particles measured in PHOS as Gamma calorimeter | |
402 | TClonesArray * plEMCAL = new TClonesArray("TParticle",1000); // All particles measured in EMCAL as Gamma calorimeter | |
403 | ||
404 | TParticle *pGamma = new TParticle(); //It will contain the kinematics of the found prompt gamma | |
2a1d8a29 | 405 | |
f9cea31c | 406 | //Fill lists with photons, neutral particles and charged particles |
407 | //look for the highest energy photon in the event inside fCalorimeter | |
2a1d8a29 | 408 | |
f9cea31c | 409 | AliDebug(2, "Fill particle lists, get prompt gamma"); |
2a1d8a29 | 410 | |
f9cea31c | 411 | //Fill particle lists |
412 | CreateParticleList(particleList, plCTS,plEMCAL,plPHOS); | |
2a1d8a29 | 413 | |
f9cea31c | 414 | if(fCalorimeter == "PHOS") |
415 | plNe = plPHOS; | |
416 | if(fCalorimeter == "EMCAL") | |
417 | plNe = plEMCAL; | |
2a1d8a29 | 418 | |
419 | ||
420 | Bool_t iIsInPHOSorEMCAL = kFALSE ; //To check if Gamma was in any calorimeter | |
421 | //Search highest energy prompt gamma in calorimeter | |
422 | GetPromptGamma(plNe, plCTS, pGamma, iIsInPHOSorEMCAL) ; | |
423 | ||
424 | AliDebug(1, Form("Is Gamma in %s? %d",fCalorimeter.Data(),iIsInPHOSorEMCAL)); | |
425 | ||
426 | //If there is any photon candidate in fCalorimeter | |
427 | if(iIsInPHOSorEMCAL){ | |
428 | if (fPrintInfo) | |
429 | AliInfo(Form("Prompt Gamma: pt %f, phi %f, eta %f", pGamma->Pt(),pGamma->Phi(),pGamma->Eta())) ; | |
f9cea31c | 430 | |
2a1d8a29 | 431 | }//Gamma in Calo |
f9cea31c | 432 | |
433 | AliDebug(2, "End of analysis, delete pointers"); | |
434 | ||
435 | particleList->Delete() ; | |
436 | plCTS->Delete() ; | |
437 | plNe->Delete() ; | |
2a1d8a29 | 438 | plEMCAL->Delete() ; |
f9cea31c | 439 | plPHOS->Delete() ; |
f9cea31c | 440 | pGamma->Delete(); |
2a1d8a29 | 441 | |
442 | PostData(0, fOutputContainer); | |
f9cea31c | 443 | |
2a1d8a29 | 444 | delete plNe ; |
f9cea31c | 445 | delete plCTS ; |
2a1d8a29 | 446 | //delete plPHOS ; |
447 | //delete plEMCAL ; | |
f9cea31c | 448 | delete particleList ; |
2a1d8a29 | 449 | |
f9cea31c | 450 | // delete pGamma; |
451 | ||
f9cea31c | 452 | } |
453 | ||
454 | ||
455 | //____________________________________________________________________________ | |
456 | void AliAnaGammaDirect::GetPromptGamma(TClonesArray * pl, TClonesArray * plCTS, TParticle *pGamma, Bool_t &Is) const | |
457 | { | |
458 | //Search for the prompt photon in Calorimeter with pt > fMinGammaPt | |
459 | ||
460 | Double_t pt = 0; | |
461 | Int_t index = -1; | |
462 | for(Int_t ipr = 0;ipr < pl->GetEntries() ; ipr ++ ){ | |
463 | TParticle * particle = dynamic_cast<TParticle *>(pl->At(ipr)) ; | |
464 | ||
465 | if((particle->Pt() > fMinGammaPt) && (particle->Pt() > pt)){ | |
466 | index = ipr ; | |
467 | pt = particle->Pt(); | |
468 | pGamma->SetMomentum(particle->Px(),particle->Py(),particle->Pz(),particle->Energy()); | |
469 | AliDebug(4,Form("Cluster in calo: pt %f, phi %f, eta %f", pGamma->Pt(),pGamma->Phi(),pGamma->Eta())) ; | |
470 | Is = kTRUE; | |
471 | } | |
472 | } | |
473 | ||
474 | //Do Isolation? | |
475 | if(fMakeICMethod && Is) | |
476 | { | |
477 | Float_t coneptsum = 0 ; | |
478 | Bool_t icPtThres = kFALSE; | |
479 | Bool_t icPtSum = kFALSE; | |
480 | MakeIsolationCut(plCTS,pl, pGamma, index, | |
481 | icPtThres, icPtSum,coneptsum); | |
482 | if(fMakeICMethod == 1) //Pt thres method | |
483 | Is = icPtThres ; | |
484 | if(fMakeICMethod == 2) //Pt cone sum method | |
485 | Is = icPtSum ; | |
486 | } | |
487 | ||
488 | if(Is){ | |
489 | AliDebug(3,Form("Cluster with p_{T} larger than %f found in calorimeter ", fMinGammaPt)) ; | |
490 | AliDebug(3,Form("Gamma: pt %f, phi %f, eta %f", pGamma->Pt(),pGamma->Phi(),pGamma->Eta())) ; | |
491 | //Fill prompt gamma histograms | |
492 | fhNGamma ->Fill(pGamma->Pt()); | |
493 | fhPhiGamma->Fill( pGamma->Pt(),pGamma->Phi()); | |
494 | fhEtaGamma->Fill(pGamma->Pt(),pGamma->Eta()); | |
495 | } | |
496 | else | |
497 | AliDebug(1,Form("NO Cluster with pT larger than %f found in calorimeter ", fMinGammaPt)) ; | |
498 | } | |
499 | ||
500 | //____________________________________________________________________________ | |
2a1d8a29 | 501 | void AliAnaGammaDirect::InitParameters() |
f9cea31c | 502 | { |
f9cea31c | 503 | |
2a1d8a29 | 504 | //Initialize the parameters of the analysis. |
f9cea31c | 505 | fCalorimeter="PHOS"; |
506 | fPrintInfo = kTRUE; | |
2a1d8a29 | 507 | fMinGammaPt = 5. ; |
f9cea31c | 508 | |
509 | //Fill particle lists when PID is ok | |
510 | fEMCALPID = kFALSE; | |
511 | fPHOSPID = kFALSE; | |
84e10d55 | 512 | fEMCALPhotonWeight = 0.5 ; |
513 | fEMCALPi0Weight = 0.5 ; | |
514 | fPHOSPhotonWeight = 0.8 ; | |
f9cea31c | 515 | fConeSize = 0.2 ; |
516 | fPtThreshold = 2.0; | |
517 | fPtSumThreshold = 1.; | |
518 | ||
2a1d8a29 | 519 | fMakeICMethod = 1; // 0 don't isolate, 1 pt thresh method, 2 cone pt sum method |
f9cea31c | 520 | } |
521 | ||
522 | //__________________________________________________________________ | |
523 | void AliAnaGammaDirect::MakeIsolationCut(TClonesArray * plCTS, | |
524 | TClonesArray * plNe, | |
525 | TParticle * pCandidate, | |
526 | Int_t index, | |
527 | Bool_t &icmpt, Bool_t &icms, | |
528 | Float_t &coneptsum) const | |
529 | { | |
530 | //Search in cone around a candidate particle if it is isolated | |
531 | Float_t phiC = pCandidate->Phi() ; | |
532 | Float_t etaC = pCandidate->Eta() ; | |
533 | Float_t pt = -100. ; | |
534 | Float_t eta = -100. ; | |
535 | Float_t phi = -100. ; | |
536 | Float_t rad = -100 ; | |
537 | Int_t n = 0 ; | |
538 | TParticle * particle = new TParticle(); | |
539 | ||
540 | coneptsum = 0; | |
541 | icmpt = kFALSE; | |
542 | icms = kFALSE; | |
543 | ||
544 | //Check charged particles in cone. | |
545 | for(Int_t ipr = 0;ipr < plCTS->GetEntries() ; ipr ++ ){ | |
546 | particle = dynamic_cast<TParticle *>(plCTS->At(ipr)) ; | |
547 | pt = particle->Pt(); | |
548 | eta = particle->Eta(); | |
549 | phi = particle->Phi() ; | |
550 | ||
551 | //Check if there is any particle inside cone with pt larger than fPtThreshold | |
552 | rad = TMath::Sqrt(TMath::Power((eta-etaC),2) + | |
553 | TMath::Power((phi-phiC),2)); | |
554 | if(rad<fConeSize){ | |
555 | AliDebug(3,Form("charged in cone pt %f, phi %f, eta %f, R %f ",pt,phi,eta,rad)); | |
556 | coneptsum+=pt; | |
557 | if(pt > fPtThreshold ) n++; | |
558 | } | |
559 | }// charged particle loop | |
560 | ||
561 | //Check neutral particles in cone. | |
562 | for(Int_t ipr = 0;ipr < plNe->GetEntries() ; ipr ++ ){ | |
563 | if(ipr != index){//Do not count the candidate | |
564 | particle = dynamic_cast<TParticle *>(plNe->At(ipr)) ; | |
565 | pt = particle->Pt(); | |
566 | eta = particle->Eta(); | |
567 | phi = particle->Phi() ; | |
568 | ||
569 | //Check if there is any particle inside cone with pt larger than fPtThreshold | |
570 | rad = TMath::Sqrt(TMath::Power((eta-etaC),2) + | |
571 | TMath::Power((phi-phiC),2)); | |
572 | if(rad<fConeSize){ | |
573 | AliDebug(3,Form("charged in cone pt %f, phi %f, eta %f, R %f ",pt,phi,eta,rad)); | |
574 | coneptsum+=pt; | |
575 | if(pt > fPtThreshold ) n++; | |
576 | } | |
577 | } | |
578 | }// neutral particle loop | |
579 | ||
580 | if(n == 0) | |
581 | icmpt = kTRUE ; | |
582 | if(coneptsum < fPtSumThreshold) | |
583 | icms = kTRUE ; | |
584 | ||
585 | } | |
586 | ||
f9cea31c | 587 | void AliAnaGammaDirect::Print(const Option_t * opt) const |
588 | { | |
589 | ||
590 | //Print some relevant parameters set for the analysis | |
591 | if(! opt) | |
592 | return; | |
593 | ||
594 | Info("Print", "%s %s", GetName(), GetTitle() ) ; | |
2a1d8a29 | 595 | printf("IC method = %d\n", fMakeICMethod) ; |
f9cea31c | 596 | printf("Cone Size = %f\n", fConeSize) ; |
597 | printf("pT threshold = %f\n", fPtThreshold) ; | |
598 | printf("pT sum threshold = %f\n", fPtSumThreshold) ; | |
2a1d8a29 | 599 | printf("Min Gamma pT = %f\n", fMinGammaPt) ; |
f9cea31c | 600 | printf("Calorimeter = %s\n", fCalorimeter.Data()) ; |
601 | } | |
602 | ||
603 | void AliAnaGammaDirect::Terminate(Option_t *) | |
604 | { | |
605 | // The Terminate() function is the last function to be called during | |
606 | // a query. It always runs on the client, it can be used to present | |
607 | // the results graphically or save the results to file. | |
608 | ||
609 | ||
610 | } |