]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/AliAnaGammaHadron.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
replace ANALYSIS_NEW by ANALYSIS
[u/mrichter/AliRoot.git] / PWG4 / AliAnaGammaHadron.cxx
CommitLineData
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$
2a1d8a29 20 * Revision 1.1 2007/01/23 17:17:29 schutz
21 * New Gamma package
22 *
f9cea31c 23 *
24 */
25
26//_________________________________________________________________________
27// Class for the analysis of gamma-jet correlations
28// Basically it seaches for a prompt photon in the Calorimeters acceptance,
29// if so we construct a jet around the highest pt particle in the opposite
30// side in azimuth, inside the Central Tracking System (ITS+TPC) and
31// EMCAL acceptances (only when PHOS detects the gamma). First the leading
32// particle and then the jet have to fullfill several conditions
33// (energy, direction ..) to be accepted. Then the fragmentation function
34// of this jet is constructed
35// Class created from old AliPHOSGammaPion
36// (see AliRoot versions previous Release 4-09)
37//
38//*-- Author: Gustavo Conesa (LNF-INFN)
39//////////////////////////////////////////////////////////////////////////////
40
41
42// --- ROOT system ---
43
44#include <TFile.h>
45#include <TParticle.h>
46#include <TH2.h>
47
48#include "AliAnaGammaHadron.h"
49#include "AliESD.h"
50#include "AliESDtrack.h"
51#include "AliESDCaloCluster.h"
52#include "Riostream.h"
53#include "AliLog.h"
54
55ClassImp(AliAnaGammaHadron)
56
57//____________________________________________________________________________
58AliAnaGammaHadron::AliAnaGammaHadron(const char *name) :
59 AliAnaGammaDirect(name),
60 fPhiMaxCut(0.), fPhiMinCut(0.),
61 fInvMassMaxCut(0.), fInvMassMinCut(0.),
62 fMinPtPion(0),
63 fAngleMaxParam()
64{
65
66 // ctor
67 fAngleMaxParam.Set(4) ;
68 fAngleMaxParam.Reset(0.);
69
70 TList * list = gDirectory->GetListOfKeys() ;
71 TIter next(list) ;
72 TH2F * h = 0 ;
73 Int_t index ;
74 for (index = 0 ; index < list->GetSize()-1 ; index++) {
75 //-1 to avoid GammaPion Task
76 h = dynamic_cast<TH2F*>(gDirectory->Get(list->At(index)->GetName())) ;
77 fOutputContainer->Add(h) ;
78 }
79
80 // Input slot #0 works with an Ntuple
81 DefineInput(0, TChain::Class());
82 // Output slot #0 writes into a TH1 container
83 DefineOutput(0, TObjArray::Class()) ;
84
85}
86
87
88//____________________________________________________________________________
89AliAnaGammaHadron::AliAnaGammaHadron(const AliAnaGammaHadron & gj) :
90 AliAnaGammaDirect(gj),
91 fPhiMaxCut(gj.fPhiMaxCut), fPhiMinCut(gj.fPhiMinCut),
92 fInvMassMaxCut(gj.fInvMassMaxCut), fInvMassMinCut(gj.fInvMassMinCut),
93 fMinPtPion(gj.fMinPtPion),
94 fOutputContainer(0), fAngleMaxParam(gj.fAngleMaxParam)
95
96{
97 // cpy ctor
98 SetName (gj.GetName()) ;
99 SetTitle(gj.GetTitle()) ;
100
101}
102
103//____________________________________________________________________________
104AliAnaGammaHadron::~AliAnaGammaHadron()
105{
106 // Remove all pointers
107 fOutputContainer->Clear() ;
108 delete fOutputContainer ;
109
110 delete fhPhiCharged ;
111 delete fhPhiNeutral ;
112 delete fhEtaCharged ;
113 delete fhEtaNeutral ;
114 delete fhDeltaPhiGammaCharged ;
115 delete fhDeltaPhiGammaNeutral ;
116 delete fhDeltaEtaGammaCharged ;
117 delete fhDeltaEtaGammaNeutral ;
118
119 delete fhCorrelationGammaNeutral ;
120 delete fhCorrelationGammaCharged ;
121
122 delete fhAnglePairNoCut ;
123 delete fhAnglePairAzimuthCut ;
124 delete fhAnglePairOpeningAngleCut ;
125 delete fhAnglePairAllCut ;
126 delete fhInvMassPairNoCut ;
127 delete fhInvMassPairAzimuthCut ;
128 delete fhInvMassPairOpeningAngleCut ;
129 delete fhInvMassPairAllCut ;
130
131}
132
2a1d8a29 133//______________________________________________________________________________
134void AliAnaGammaHadron::ConnectInputData(const Option_t*)
135{
136 // Initialisation of branch container and histograms
137 AliAnaGammaDirect::ConnectInputData("");
138}
139
140//________________________________________________________________________
141void AliAnaGammaHadron::CreateOutputObjects()
142{
143
144 // Init parameteres and create histograms to be saved in output file and
145 // stores them in fOutputContainer
146 InitParameters();
147 AliAnaGammaDirect::CreateOutputObjects();
148
149 fOutputContainer = new TObjArray(100) ;
150
151 //Use histograms in AliAnaGammaDirect
152 TObjArray * outputContainer =GetOutputContainer();
153 for(Int_t i = 0; i < outputContainer->GetEntries(); i++ )
154 fOutputContainer->Add(outputContainer->At(i)) ;
155
156 fhPhiCharged = new TH2F
157 ("PhiCharged","#phi_{#pi^{#pm}} vs p_{T #gamma}",
158 120,0,120,120,0,7);
159 fhPhiCharged->SetYTitle("#phi_{#pi^{#pm}} (rad)");
160 fhPhiCharged->SetXTitle("p_{T #gamma} (GeV/c)");
161 fOutputContainer->Add(fhPhiCharged) ;
162
163 fhPhiNeutral = new TH2F
164 ("PhiNeutral","#phi_{#pi^{0}} vs p_{T #gamma}",
165 120,0,120,120,0,7);
166 fhPhiNeutral->SetYTitle("#phi_{#pi^{0}} (rad)");
167 fhPhiNeutral->SetXTitle("p_{T #gamma} (GeV/c)");
168 fOutputContainer->Add(fhPhiNeutral) ;
169
170 fhEtaCharged = new TH2F
171 ("EtaCharged","#eta_{#pi^{#pm}} vs p_{T #gamma}",
172 120,0,120,120,-1,1);
173 fhEtaCharged->SetYTitle("#eta_{#pi^{#pm}} (rad)");
174 fhEtaCharged->SetXTitle("p_{T #gamma} (GeV/c)");
175 fOutputContainer->Add(fhEtaCharged) ;
176
177 fhEtaNeutral = new TH2F
178 ("EtaNeutral","#eta_{#pi^{0}} vs p_{T #gamma}",
179 120,0,120,120,-1,1);
180 fhEtaNeutral->SetYTitle("#eta_{#pi^{0}} (rad)");
181 fhEtaNeutral->SetXTitle("p_{T #gamma} (GeV/c)");
182 fOutputContainer->Add(fhEtaNeutral) ;
f9cea31c 183
2a1d8a29 184 fhDeltaPhiGammaCharged = new TH2F
185 ("DeltaPhiGammaCharged","#phi_{#gamma} - #phi_{charged #pi} vs p_{T #gamma}",
186 200,0,120,200,0,6.4);
187 fhDeltaPhiGammaCharged->SetYTitle("#Delta #phi");
188 fhDeltaPhiGammaCharged->SetXTitle("p_{T #gamma} (GeV/c)");
189 fOutputContainer->Add(fhDeltaPhiGammaCharged) ;
190
191 fhDeltaEtaGammaCharged = new TH2F
192 ("DeltaEtaGammaCharged","#eta_{#gamma} - #eta_{#pi^{#pm}} vs p_{T #gamma}",
193 200,0,120,200,-2,2);
194 fhDeltaEtaGammaCharged->SetYTitle("#Delta #eta");
195 fhDeltaEtaGammaCharged->SetXTitle("p_{T #gamma} (GeV/c)");
196 fOutputContainer->Add(fhDeltaEtaGammaCharged) ;
197
198 fhDeltaPhiGammaNeutral = new TH2F
199 ("DeltaPhiGammaNeutral","#phi_{#gamma} - #phi_{#pi^{0}} vs p_{T #gamma}",
200 200,0,120,200,0,6.4);
201 fhDeltaPhiGammaNeutral->SetYTitle("#Delta #phi");
202 fhDeltaPhiGammaNeutral->SetXTitle("p_{T #gamma} (GeV/c)");
203 fOutputContainer->Add(fhDeltaPhiGammaNeutral) ;
204
205 fhDeltaEtaGammaNeutral = new TH2F
206 ("DeltaEtaGammaNeutral","#eta_{#gamma} - #eta_{#pi^{#pm}} vs p_{T #gamma}",
207 200,0,120,200,-2,2);
208 fhDeltaEtaGammaNeutral->SetYTitle("#Delta #eta");
209 fhDeltaEtaGammaNeutral->SetXTitle("p_{T #gamma} (GeV/c)");
210 fOutputContainer->Add(fhDeltaEtaGammaNeutral) ;
211
212 //
213 fhAnglePairAccepted = new TH2F
214 ("AnglePairAccepted",
215 "Angle between #pi^{0} #gamma pair vs p_{T #pi^{0}}, both #gamma in eta<0.7, inside window",
216 200,0,50,200,0,0.2);
217 fhAnglePairAccepted->SetYTitle("Angle (rad)");
218 fhAnglePairAccepted->SetXTitle("E_{ #pi^{0}} (GeV/c)");
219 fOutputContainer->Add(fhAnglePairAccepted) ;
220
221 fhAnglePairNoCut = new TH2F
222 ("AnglePairNoCut",
223 "Angle between all #gamma pair vs p_{T #pi^{0}}",200,0,50,200,0,0.2);
224 fhAnglePairNoCut->SetYTitle("Angle (rad)");
225 fhAnglePairNoCut->SetXTitle("E_{ #pi^{0}} (GeV/c)");
226 fOutputContainer->Add(fhAnglePairNoCut) ;
227
228 fhAnglePairAzimuthCut = new TH2F
229 ("AnglePairAzimuthCut",
230 "Angle between all #gamma pair that have a good phi and pt vs p_{T #pi^{0}}",
231 200,0,50,200,0,0.2);
232 fhAnglePairAzimuthCut->SetYTitle("Angle (rad)");
233 fhAnglePairAzimuthCut->SetXTitle("E_{ #pi^{0}} (GeV/c)");
234 fOutputContainer->Add(fhAnglePairAzimuthCut) ;
235
236 fhAnglePairOpeningAngleCut = new TH2F
237 ("AnglePairOpeningAngleCut",
238 "Angle between all #gamma pair (opening angle + azimuth cut) vs p_{T #pi^{0}}"
239 ,200,0,50,200,0,0.2);
240 fhAnglePairOpeningAngleCut->SetYTitle("Angle (rad)");
241 fhAnglePairOpeningAngleCut->SetXTitle("E_{ #pi^{0}} (GeV/c)");
242 fOutputContainer->Add(fhAnglePairOpeningAngleCut) ;
243
244 fhAnglePairAllCut = new TH2F
245 ("AnglePairAllCut",
246 "Angle between all #gamma pair (opening angle + inv mass cut+azimuth) vs p_{T #pi^{0}}"
247 ,200,0,50,200,0,0.2);
248 fhAnglePairAllCut->SetYTitle("Angle (rad)");
249 fhAnglePairAllCut->SetXTitle("E_{ #pi^{0}} (GeV/c)");
250 fOutputContainer->Add(fhAnglePairAllCut) ;
251
252
253 //
254 fhInvMassPairNoCut = new TH2F
255 ("InvMassPairNoCut","Invariant Mass of all #gamma pair vs p_{T #gamma}",
256 120,0,120,360,0,0.5);
257 fhInvMassPairNoCut->SetYTitle("Invariant Mass (GeV/c^{2})");
258 fhInvMassPairNoCut->SetXTitle("p_{T #gamma} (GeV/c)");
259 fOutputContainer->Add(fhInvMassPairNoCut) ;
260
261 fhInvMassPairAzimuthCut = new TH2F
262 ("InvMassPairAzimuthCut",
263 "Invariant Mass of #gamma pair (azimuth cuts) vs p_{T #gamma}",
264 120,0,120,360,0,0.5);
265 fhInvMassPairAzimuthCut->SetYTitle("Invariant Mass (GeV/c^{2})");
266 fhInvMassPairAzimuthCut->SetXTitle("p_{T #gamma} (GeV/c)");
267 fOutputContainer->Add(fhInvMassPairAzimuthCut) ;
268
269 fhInvMassPairOpeningAngleCut = new TH2F
270 ("InvMassPairOpeningAngleCut",
271 "Invariant Mass of #gamma pair (angle cut) vs p_{T #gamma}",
272 120,0,120,360,0,0.5);
273 fhInvMassPairOpeningAngleCut->SetYTitle("Invariant Mass (GeV/c^{2})");
274 fhInvMassPairOpeningAngleCut->SetXTitle("p_{T #gamma} (GeV/c)");
275 fOutputContainer->Add(fhInvMassPairOpeningAngleCut) ;
276
277 fhInvMassPairAllCut = new TH2F
278 ("InvMassPairAllCut",
279 "Invariant Mass of #gamma pair (opening angle+invmass cut+azimuth) vs p_{T #gamma}",
280 120,0,120,360,0,0.5);
281 fhInvMassPairAllCut->SetYTitle("Invariant Mass (GeV/c^{2})");
282 fhInvMassPairAllCut->SetXTitle("p_{T #gamma} (GeV/c)");
283 fOutputContainer->Add(fhInvMassPairAllCut) ;
284
285 //
286 fhCorrelationGammaCharged =
287 new TH2F("CorrelationGammaCharged","z_{#gamma #pi} = p_{T #pi^{#pm}} / p_{T #gamma}",
288 240,0.,120.,1000,0.,1.2);
289 fhCorrelationGammaCharged->SetYTitle("z_{#gamma #pi}");
290 fhCorrelationGammaCharged->SetXTitle("p_{T #gamma}");
291 fOutputContainer->Add(fhCorrelationGammaCharged) ;
292
293 fhCorrelationGammaNeutral =
294 new TH2F("CorrelationGammaNeutral","z_{#gamma #pi} = p_{T #pi^{0}} / p_{T #gamma}",
295 240,0.,120.,1000,0.,1.2);
296 fhCorrelationGammaNeutral->SetYTitle("z_{#gamma #pi}");
297 fhCorrelationGammaNeutral->SetXTitle("p_{T #gamma}");
298 fOutputContainer->Add(fhCorrelationGammaNeutral) ;
299
300}
f9cea31c 301
302//____________________________________________________________________________
303void AliAnaGammaHadron::Exec(Option_t *)
304{
305
306 // Processing of one event
2a1d8a29 307
f9cea31c 308 //Get ESDs
309 Long64_t entry = GetChain()->GetReadEntry() ;
310
311 if (!GetESD()) {
312 AliError("fESD is not connected to the input!") ;
313 return ;
314 }
315
316 if (GetPrintInfo())
317 AliInfo(Form("%s ----> Processing event # %lld", (dynamic_cast<TChain *>(GetChain()))->GetFile()->GetName(), entry)) ;
318
319 //CreateTLists with arrays of TParticles. Filled with particles only relevant for the analysis.
320
321 TClonesArray * particleList = new TClonesArray("TParticle",1000); // All particles refitted in CTS and detected in EMCAL (jet)
322 TClonesArray * plCTS = new TClonesArray("TParticle",1000); // All particles refitted in Central Tracking System (ITS+TPC)
323 TClonesArray * plNe = new TClonesArray("TParticle",1000); // All particles measured in Jet Calorimeter
324 TClonesArray * plCalo = new TClonesArray("TParticle",1000); // All particles measured in Prompt Gamma calorimeter
325
326
327 TParticle *pGamma = new TParticle(); //It will contain the kinematics of the found prompt gamma
328
329
330 Bool_t iIsInPHOSorEMCAL = kFALSE ; //To check if Gamma was in any calorimeter
331
332 //Fill lists with photons, neutral particles and charged particles
333 //look for the highest energy photon in the event inside fCalorimeter
334 //Fill particle lists
335 AliDebug(2, "Fill particle lists, get prompt gamma");
336
337 //Fill particle lists
338 if(GetCalorimeter() == "PHOS")
339 CreateParticleList(particleList, plCTS,plNe,plCalo);
340 else if(GetCalorimeter() == "EMCAL")
341 CreateParticleList(particleList, plCTS,plCalo,plNe);
342 else
343 AliError("No calorimeter selected");
344
345 //Search highest energy prompt gamma in calorimeter
346 GetPromptGamma(plCalo, plCTS, pGamma, iIsInPHOSorEMCAL) ;
347
348
349 AliDebug(1, Form("Is Gamma in %s? %d",GetCalorimeter().Data(),iIsInPHOSorEMCAL));
350 AliDebug(3,Form("Charged list entries %d, Neutral list entries %d, %s list entries %d",
351 plCTS->GetEntries(),plNe->GetEntries(), GetCalorimeter().Data(),plCalo->GetEntries()));
352
353 //If there is any prompt photon in fCalorimeter,
354 //search jet leading particle
355
356 if(iIsInPHOSorEMCAL){
357
358 if (GetPrintInfo())
359 AliInfo(Form("Prompt Gamma: pt %f, phi %f, eta %f", pGamma->Pt(),pGamma->Phi(),pGamma->Eta())) ;
360
361 AliDebug(2, "Make correlation");
362
363 //Search correlation
364 MakeGammaChargedCorrelation(plCTS, pGamma);
365 MakeGammaNeutralCorrelation(plNe, pGamma);
366
367 }//Gamma in Calo
368
369 AliDebug(2, "End of analysis, delete pointers");
370
371 particleList->Delete() ;
372 plCTS->Delete() ;
373 plNe->Delete() ;
374 plCalo->Delete() ;
375 pGamma->Delete();
376
377 delete plNe ;
378 delete plCalo ;
379 delete plCTS ;
380 delete particleList ;
381 // delete pGamma;
382
383 PostData(0, fOutputContainer);
384}
385
386
387//____________________________________________________________________________
388void AliAnaGammaHadron::MakeGammaChargedCorrelation(TClonesArray * pl, TParticle * pGamma) const
389{
390 //Search for the charged particle with highest with
391 //Phi=Phi_gamma-Pi and pT=0.1E_gamma
392 Double_t ptg = pGamma->Pt();
393 Double_t phig = pGamma->Phi();
394 Double_t pt = -100.;
395 Double_t rat = -100.;
396 Double_t phi = -100. ;
397
398 for(Int_t ipr = 0;ipr < pl->GetEntries() ; ipr ++ ){
399
400 TParticle * particle = dynamic_cast<TParticle *>(pl->At(ipr)) ;
401
402 pt = particle->Pt();
403 rat = pt/ptg ;
404 phi = particle->Phi() ;
405
406 AliDebug(3,Form("pt %f, phi %f, phi gamma %f. Cuts: delta phi min %f, max%f, pT min %f",pt,phi,phig,fPhiMinCut,fPhiMaxCut,fMinPtPion));
407
408 fhEtaCharged->Fill(ptg,particle->Eta());
409 fhPhiCharged->Fill(ptg,phi);
410 fhDeltaEtaGammaCharged->Fill(ptg,pGamma->Eta()-particle->Eta());
411 fhDeltaPhiGammaCharged->Fill(ptg,phig-phi);
412 //Selection within angular and energy limits
413 if(((phig-phi)> fPhiMinCut) && ((phig-phi)<fPhiMaxCut) && pt > fMinPtPion){
414 AliDebug(2,Form("Selected: pt %f, phi %f",pt,phi));
415 fhCorrelationGammaCharged->Fill(ptg,rat);
416 }
417 }//particle loop
418}
419
420//____________________________________________________________________________
421void AliAnaGammaHadron::MakeGammaNeutralCorrelation(TClonesArray * pl, TParticle * pGamma)
422{
423
424 //Search for the neutral pion with highest with
425 //Phi=Phi_gamma-Pi and pT=0.1E_gamma
426 Double_t pt = -100.;
427 Double_t rat = -100.;
428 Double_t phi = -100. ;
429 Double_t ptg = pGamma->Pt();
430 Double_t phig = pGamma->Phi();
431
432 TIter next(pl);
433 TParticle * particle = 0;
434
435 Int_t iPrimary = -1;
436 TLorentzVector gammai,gammaj;
437 Double_t angle = 0., e = 0., invmass = 0.;
438 Int_t ksPdg = 0;
439 Int_t jPrimary=-1;
440
441 while ( (particle = (TParticle*)next()) ) {
442 iPrimary++;
443 ksPdg = particle->GetPdgCode();
444
445 //2 gamma overlapped, found with PID
446 if(ksPdg == 111){
447 pt = particle->Pt();
448 rat = pt/ptg ;
449 phi = particle->Phi() ;
450 fhEtaCharged->Fill(ptg,particle->Eta());
451 fhPhiCharged->Fill(ptg,phi);
452 fhDeltaEtaGammaCharged->Fill(ptg,pGamma->Eta()-particle->Eta());
453 fhDeltaPhiGammaCharged->Fill(ptg,phig-phi);
454 //AliDebug(3,Form("pt %f, phi %f",pt,phi));
455 if (GetPrintInfo())
456 AliInfo(Form("pt %f, phi %f",pt,phi));
457 //Selection within angular and energy limits
458 if((pt> ptg)&& ((phig-phi)>fPhiMinCut)&&((phig-phi)<fPhiMaxCut)){
459 fhCorrelationGammaNeutral ->Fill(ptg,rat);
460 //AliDebug(2,Form("Selected: pt %f, phi %f",pt,phi));
461 if (GetPrintInfo())
462 AliInfo(Form("Selected: pt %f, phi %f",pt,phi));
463 }// cuts
464 }// pdg = 111
465
466 //Make invariant mass analysis
467 else if(ksPdg == 22){
468 //Search the photon companion in case it comes from a Pi0 decay
469 //Apply several cuts to select the good pair
470 particle->Momentum(gammai);
471 jPrimary=-1;
472 TIter next2(pl);
473 while ( (particle = (TParticle*)next2()) ) {
474 jPrimary++;
475 if(jPrimary>iPrimary){
476 ksPdg = particle->GetPdgCode();
477
478 if(ksPdg == 22){
479 particle->Momentum(gammaj);
480 //Info("GetLeadingPi0","Egammai %f, Egammaj %f",
481 //gammai.Pt(),gammaj.Pt());
482 pt = (gammai+gammaj).Pt();
483 phi = (gammai+gammaj).Phi();
484 if(phi < 0)
485 phi+=TMath::TwoPi();
486 rat = pt/ptg ;
487 invmass = (gammai+gammaj).M();
488 angle = gammaj.Angle(gammai.Vect());
489 e = (gammai+gammaj).E();
490 fhEtaNeutral->Fill(ptg,(gammai+gammaj).Eta());
491 fhPhiNeutral->Fill(ptg,phi);
492 fhDeltaEtaGammaNeutral->Fill(ptg,pGamma->Eta()-(gammai+gammaj).Eta());
493 fhDeltaPhiGammaNeutral->Fill(ptg,phig-phi);
494 // AliDebug(3,Form("pt %f, phi %f",pt,phi));
495 if (GetPrintInfo())
496 AliInfo(Form("pt %f, phi %f",pt,phi));
497
498 //Fill histograms with no cuts applied.
499 fhAnglePairNoCut->Fill(e,angle);
500 fhInvMassPairNoCut->Fill(ptg,invmass);
501
502 //First cut on the energy and azimuth of the pair
503 if((phig-phi) > fPhiMinCut && (phig-phi) < fPhiMaxCut
504 && pt > fMinPtPion){
505 fhAnglePairAzimuthCut ->Fill(e,angle);
506 fhInvMassPairAzimuthCut->Fill(ptg,invmass);
507 AliDebug(3,Form("1st cut: pt %f, phi %f",pt,phi));
508
509 //Second cut on the aperture of the pair
510 if(IsAngleInWindow(angle,e)){
511 fhAnglePairOpeningAngleCut ->Fill(e,angle);
512 fhInvMassPairOpeningAngleCut->Fill(ptg,invmass);
513 AliDebug(3,Form("2nd cut: pt %f, phi %f",pt,phi));
514
515 //Third cut on the invariant mass of the pair
516 if((invmass>fInvMassMinCut) && (invmass<fInvMassMaxCut)){
517 fhInvMassPairAllCut ->Fill(ptg,invmass);
518 fhAnglePairAllCut ->Fill(e,angle);
519 //Fill correlation histogram
520 fhCorrelationGammaNeutral ->Fill(ptg,rat);
521 //AliDebug(2,Form("Selected: pt %f, phi %f",pt,phi));
522 if (GetPrintInfo())
523 AliInfo(Form("Selected: pt %f, phi %f",pt,phi));
524 }//(invmass>0.125) && (invmass<0.145)
525 }//Opening angle cut
526 }//Azimuth and pt cut.
527 }//if pdg = 22
528 }//iprimary<jprimary
529 }//while
530 }// if pdg = 22
531 }//while
532}
533
534 //____________________________________________________________________________
2a1d8a29 535void AliAnaGammaHadron::InitParameters()
f9cea31c 536{
537 // Initialisation of branch container
2a1d8a29 538 //AliAnaGammaDirect::InitParameters();
f9cea31c 539
540 //Initialize the parameters of the analysis.
541 //fCalorimeter="PHOS";
542 fAngleMaxParam.Set(4) ;
543 fAngleMaxParam.AddAt(0.4,0);//={0.4,-0.25,0.025,-2e-4};
544 fAngleMaxParam.AddAt(-0.25,1) ;
545 fAngleMaxParam.AddAt(0.025,2) ;
546 fAngleMaxParam.AddAt(-2e-4,3) ;
547
548 //fPrintInfo = kTRUE;
549 fInvMassMaxCut = 0.16 ;
550 fInvMassMinCut = 0.11 ;
551 fPhiMaxCut = 4.5;
552 fPhiMinCut = 1.5 ;
553
554 fMinPtPion = 0. ;
555
556 //Fill particle lists when PID is ok
557 // fEMCALPID = kFALSE;
558 // fPHOSPID = kFALSE;
559
f9cea31c 560}
561
562//__________________________________________________________________________-
563Bool_t AliAnaGammaHadron::IsAngleInWindow(const Float_t angle,const Float_t e) {
564 //Check if the opening angle of the candidate pairs is inside
565 //our selection windowd
566
567 Bool_t result = kFALSE;
568 Double_t mpi0 = 0.1349766;
569 Double_t max = fAngleMaxParam.At(0)*TMath::Exp(fAngleMaxParam.At(1)*e)
570 +fAngleMaxParam.At(2)+fAngleMaxParam.At(3)*e;
571 Double_t arg = (e*e-2*mpi0*mpi0)/(e*e);
572 Double_t min = 100. ;
573 if(arg>0.)
574 min = TMath::ACos(arg);
575
576 if((angle<max)&&(angle>=min))
577 result = kTRUE;
578
579 return result;
580}
581
f9cea31c 582//____________________________________________________________________________
583void AliAnaGammaHadron::Print(const Option_t * opt) const
584{
585
586 //Print some relevant parameters set for the analysis
587 if(! opt)
588 return;
589
590 Info("Print", "%s %s", GetName(), GetTitle() ) ;
591 printf("pT Pion > %f\n", fMinPtPion) ;
592 printf("Phi Pion < %f\n", fPhiMaxCut) ;
593 printf("Phi Pion > %f\n", fPhiMinCut) ;
594 printf("M_pair < %f\n", fInvMassMaxCut) ;
595 printf("M_pair > %f\n", fInvMassMinCut) ;
596
597}
598
599//__________________________________________
600void AliAnaGammaHadron::Terminate(Option_t *)
601{
602 // The Terminate() function is the last function to be called during
603 // a query. It always runs on the client, it can be used to present
604 // the results graphically or save the results to file.
605
606
607}