]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/CaloTrackCorrelations/AliAnaGeneratorKine.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Dalitz: New histogram for photon effiVsRadius
[u/mrichter/AliRoot.git] / PWGGA / CaloTrackCorrelations / AliAnaGeneratorKine.cxx
CommitLineData
7b2086c3 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
16//_________________________________________________________________________
17// Do photon/pi0 analysis for isolation and correlation
18// at the generator level. Only for kine stack (ESDs)
19//
20//
21// -- Author: Gustavo Conesa (LPSC-CNRS-Grenoble)
22//////////////////////////////////////////////////////////////////////////////
23
24// --- ROOT system ---
25#include "TH2F.h"
26#include "TParticle.h"
27#include "TDatabasePDG.h"
28
29//---- ANALYSIS system ----
30#include "AliAnaGeneratorKine.h"
31#include "AliStack.h"
32#include "AliGenPythiaEventHeader.h"
33
34ClassImp(AliAnaGeneratorKine)
35
36
37//__________________________________________
38AliAnaGeneratorKine::AliAnaGeneratorKine() :
39AliAnaCaloTrackCorrBaseClass(),
7d409bf9 40fTriggerDetector(), fTriggerDetectorString(),
2b341b44 41fFidCutTrigger(0),
783b974c 42fMinChargedPt(0), fMinNeutralPt(0),
7b2086c3 43fStack(0),
44fParton2(0), fParton3(0),
7d409bf9 45fParton6(0), fParton7(0),
7b2086c3 46fJet6(), fJet7(),
7d409bf9 47fTrigger(), fLVTmp(),
7b2086c3 48fPtHard(0),
49fhPtHard(0), fhPtParton(0), fhPtJet(0),
50fhPtPartonPtHard(0), fhPtJetPtHard(0), fhPtJetPtParton(0),
51fhPtPhoton(0), fhPtPi0(0)
52{
53 //Default Ctor
54
55 //Initialize parameters
56 InitParameters();
57
58 for(Int_t i = 0; i < 4; i++)
59 {
dbb79a05 60 fhPtPhotonLeading[i] = fhPtPi0Leading[i] = 0;
c76fb00a 61 fhPtPhotonLeadingSumPt[i] = fhPtPi0LeadingSumPt[i] = 0;
dbb79a05 62 fhPtPhotonLeadingIsolated[i] = fhPtPi0LeadingIsolated[i] = 0;
63 for(Int_t j = 0; j < 2; j++)
64 {
65 fhZHardPhoton[j][i] = fhZHardPi0[j][i] = 0;
66 fhZHardPhotonIsolated[j][i] = fhZHardPi0Isolated[j][i] = 0;
67 fhZPartonPhoton[j][i] = fhZPartonPi0[j][i] = 0;
68 fhZPartonPhotonIsolated[j][i] = fhZPartonPi0Isolated[j][i] = 0;
69 fhZJetPhoton[j][i] = fhZJetPi0[j][i] = 0;
70 fhZJetPhotonIsolated[j][i] = fhZJetPi0Isolated[j][i] = 0;
71 fhXEPhoton[j][i] = fhXEPi0[j][i] = 0;
72 fhXEPhotonIsolated[j][i] = fhXEPi0Isolated[j][i] = 0;
73 fhXEUEPhoton[j][i] = fhXEUEPi0[j][i] = 0;
74 fhXEUEPhotonIsolated[j][i] = fhXEUEPi0Isolated[j][i] = 0;
75
76 fhPtPartonTypeNearPhoton[j][i] = fhPtPartonTypeNearPi0[j][i] = 0;
77 fhPtPartonTypeNearPhotonIsolated[j][i] = fhPtPartonTypeNearPi0Isolated[j][i] = 0;
78
79 fhPtPartonTypeAwayPhoton[j][i] = fhPtPartonTypeAwayPi0[j][i] = 0;
80 fhPtPartonTypeAwayPhotonIsolated[j][i] = fhPtPartonTypeAwayPi0Isolated[j][i] = 0;
81 }
7b2086c3 82 }
83
84}
85
b94e038e 86//___________________________________________________________________________
7d409bf9 87Bool_t AliAnaGeneratorKine::CorrelateWithPartonOrJet(Int_t indexTrig,
b94e038e 88 Int_t pdgTrig,
89 Bool_t leading[4],
90 Bool_t isolated[4],
2292cf03 91 Int_t & iparton )
7b2086c3 92{
93 //Correlate trigger with partons or jets, get z
94
8cc95539 95 if(GetDebug() > 1) printf("AliAnaGeneratorKine::CorrelateWithPartonOrJet() - Start \n");
96
7b2086c3 97 //Get the index of the mother
98 iparton = (fStack->Particle(indexTrig))->GetFirstMother();
99 TParticle * mother = fStack->Particle(iparton);
100 while (iparton > 7)
101 {
102 iparton = mother->GetFirstMother();
2292cf03 103 if(iparton < 0) { printf("AliAnaGeneratorKine::CorrelateWithPartonOrJet() - Negative index, skip event\n"); return kFALSE; }
7b2086c3 104 mother = fStack->Particle(iparton);
105 }
106
107 //printf("Mother is parton %d with pdg %d\n",imom,mother->GetPdgCode());
108
109 if(iparton < 6)
110 {
111 //printf("This particle is not from hard process - pdg %d, parton index %d\n",pdgTrig, iparton);
2292cf03 112 return kFALSE;
7b2086c3 113 }
114
7d409bf9 115 Float_t ptTrig = fTrigger.Pt();
7b2086c3 116 Float_t partonPt = fParton6->Pt();
117 Float_t jetPt = fJet6.Pt();
118 if(iparton==7)
119 {
120 partonPt = fParton6->Pt();
121 jetPt = fJet6.Pt();
122 }
123
23fa04c5 124 //Get id of parton in near and away side
125
126 Int_t away = -1;
127 Int_t near = -1;
128 Int_t nearPDG = -1;
129 Int_t awayPDG = -1;
130
131 //printf("parton 6 pdg = %d, parton 7 pdg = %d\n",fParton6->GetPdgCode(),fParton7->GetPdgCode());
132
133 if(iparton==6)
134 {
135 nearPDG = fParton6->GetPdgCode();
136 awayPDG = fParton7->GetPdgCode();
137 }
138 else
139 {
140 nearPDG = fParton7->GetPdgCode();
141 awayPDG = fParton6->GetPdgCode();
142 }
143
144 if (nearPDG == 22) near = 0;
145 else if(nearPDG == 21) near = 1;
146 else near = 2;
147
148 if (awayPDG == 22) away = 0;
149 else if(awayPDG == 21) away = 1;
150 else away = 2;
151
152 for( Int_t i = 0; i < 4; i++ )
153 {
154 if(pdgTrig==111)
155 {
dbb79a05 156
157 fhPtPartonTypeNearPi0[leading[i]][i]->Fill(ptTrig,near);
158 fhPtPartonTypeAwayPi0[leading[i]][i]->Fill(ptTrig,away);
159 if(isolated[i])
160 {
161 fhPtPartonTypeNearPi0Isolated[leading[i]][i]->Fill(ptTrig,near);
162 fhPtPartonTypeAwayPi0Isolated[leading[i]][i]->Fill(ptTrig,away);
23fa04c5 163 }
dbb79a05 164
23fa04c5 165 }// pi0
166 else if(pdgTrig==22)
167 {
dbb79a05 168 fhPtPartonTypeNearPhoton[leading[i]][i]->Fill(ptTrig,near);
169 fhPtPartonTypeAwayPhoton[leading[i]][i]->Fill(ptTrig,away);
170 if(isolated[i])
171 {
172 fhPtPartonTypeNearPhotonIsolated[leading[i]][i]->Fill(ptTrig,near);
173 fhPtPartonTypeAwayPhotonIsolated[leading[i]][i]->Fill(ptTrig,away);
23fa04c5 174 }
dbb79a05 175
23fa04c5 176 } // photon
dbb79a05 177 } // conditions loop
23fa04c5 178
179
180 // RATIOS
7b2086c3 181
182 Float_t zHard = ptTrig / fPtHard;
183 Float_t zPart = ptTrig / partonPt;
184 Float_t zJet = ptTrig / jetPt;
185
186 //if(zHard > 1 ) printf("*** Particle energy larger than pT hard z=%f\n",zHard);
187
12d80e80 188 //printf("Z: hard %2.2f, parton %2.2f, jet %2.2f\n",zHard,zPart,zJet);
7b2086c3 189
190 for( Int_t i = 0; i < 4; i++ )
191 {
192 if(pdgTrig==111)
193 {
dbb79a05 194 fhZHardPi0[leading[i]][i] ->Fill(ptTrig,zHard);
195 fhZPartonPi0[leading[i]][i]->Fill(ptTrig,zPart);
196 fhZJetPi0[leading[i]][i] ->Fill(ptTrig,zJet );
197
198 if(isolated[i])
199 {
200 fhZHardPi0Isolated[leading[i]][i] ->Fill(ptTrig,zHard);
201 fhZPartonPi0Isolated[leading[i]][i]->Fill(ptTrig,zPart);
202 fhZJetPi0Isolated[leading[i]][i] ->Fill(ptTrig,zJet);
7b2086c3 203 }
dbb79a05 204
7b2086c3 205 }// pi0
206 else if(pdgTrig==22)
207 {
dbb79a05 208
209 fhZHardPhoton[leading[i]][i] ->Fill(ptTrig,zHard);
210 fhZPartonPhoton[leading[i]][i]->Fill(ptTrig,zPart);
211 fhZJetPhoton[leading[i]][i] ->Fill(ptTrig,zJet );
212
213 if(isolated[i])
214 {
215 fhZHardPhotonIsolated[leading[i]][i] ->Fill(ptTrig,zHard);
216 fhZPartonPhotonIsolated[leading[i]][i]->Fill(ptTrig,zPart);
217 fhZJetPhotonIsolated[leading[i]][i] ->Fill(ptTrig,zJet);
218 }
219
7b2086c3 220 } // photon
dbb79a05 221 } // conditions loop
2292cf03 222
8cc95539 223 if(GetDebug() > 1) printf("AliAnaGeneratorKine::CorrelateWithPartonOrJet() - End TRUE \n");
224
2292cf03 225 return kTRUE;
dbb79a05 226}
7b2086c3 227
228
229//____________________________________________________
230TList * AliAnaGeneratorKine::GetCreateOutputObjects()
231{
232 // Create histograms to be saved in output file
233
234 TList * outputContainer = new TList() ;
235 outputContainer->SetName("GenKineHistos") ;
236
c76fb00a 237 Int_t nptbins = GetHistogramRanges()->GetHistoPtBins();
238 Float_t ptmax = GetHistogramRanges()->GetHistoPtMax();
239 Float_t ptmin = GetHistogramRanges()->GetHistoPtMin();
240
241 Int_t nptsumbins = GetHistogramRanges()->GetHistoNPtSumBins();
242 Float_t ptsummax = GetHistogramRanges()->GetHistoPtSumMax();
243 Float_t ptsummin = GetHistogramRanges()->GetHistoPtSumMin();
244
7b2086c3 245
246 fhPtHard = new TH1F("hPtHard"," pt hard for selected triggers",nptbins,ptmin,ptmax);
8cc95539 247 fhPtHard->SetXTitle("#it{p}_{T}^{hard} (GeV/#it{c})");
7b2086c3 248 outputContainer->Add(fhPtHard);
249
250 fhPtParton = new TH1F("hPtParton"," pt parton for selected triggers",nptbins,ptmin,ptmax);
8cc95539 251 fhPtParton->SetXTitle("#it{p}_{T}^{parton} (GeV/#it{c})");
7b2086c3 252 outputContainer->Add(fhPtParton);
253
254 fhPtJet = new TH1F("hPtJet"," pt jet for selected triggers",nptbins,ptmin,ptmax);
8cc95539 255 fhPtJet->SetXTitle("#it{p}_{T}^{jet} (GeV/#it{c})");
7b2086c3 256 outputContainer->Add(fhPtJet);
257
258 fhPtPartonPtHard = new TH2F("hPtPartonPtHard","parton pt / pt hard for selected triggers",nptbins,ptmin,ptmax,200,0,2);
8cc95539 259 fhPtPartonPtHard->SetXTitle("#it{p}_{T}^{hard} (GeV/#it{c})");
260 fhPtPartonPtHard->SetYTitle("#it{p}_{T}^{parton}/#it{p}_{T}^{hard}");
7b2086c3 261 outputContainer->Add(fhPtPartonPtHard);
262
263 fhPtJetPtHard = new TH2F("hPtJetPtHard","jet pt / pt hard for selected triggers",nptbins,ptmin,ptmax,200,0,2);
8cc95539 264 fhPtJetPtHard->SetXTitle("#it{p}_{T}^{hard} (GeV/#it{c})");
265 fhPtJetPtHard->SetYTitle("#it{p}_{T}^{jet}/#it{p}_{T}^{hard}");
7b2086c3 266 outputContainer->Add(fhPtJetPtHard);
267
268 fhPtJetPtParton = new TH2F("hPtJetPtParton","parton pt / pt hard for selected triggers",nptbins,ptmin,ptmax,200,0,2);
8cc95539 269 fhPtJetPtParton->SetXTitle("#it{p}_{T}^{hard} (GeV/#it{c})");
270 fhPtJetPtParton->SetYTitle("#it{p}_{T}^{jet}/#it{p}_{T}^{parton}");
7b2086c3 271 outputContainer->Add(fhPtJetPtParton);
272
12d80e80 273 fhPtPhoton = new TH1F("hPtPhoton","Input #gamma",nptbins,ptmin,ptmax);
8cc95539 274 fhPtPhoton->SetXTitle("#it{p}_{T} (GeV/#it{c})");
7b2086c3 275 outputContainer->Add(fhPtPhoton);
276
12d80e80 277 fhPtPi0 = new TH1F("hPtPi0","Input #pi^{0}",nptbins,ptmin,ptmax);
8cc95539 278 fhPtPi0->SetXTitle("#it{p}_{T} (GeV/#it{c})");
7b2086c3 279 outputContainer->Add(fhPtPi0);
280
dbb79a05 281 TString name [] = {"","_EMC","_Photon","_EMC_Photon"};
12d80e80 282 TString title [] = {"",", neutral in EMCal",", neutral only #gamma-like",", neutral in EMCal and only #gamma-like"};
dbb79a05 283 TString leading[] = {"NotLeading","Leading"};
284
7b2086c3 285 for(Int_t i = 0; i < 4; i++)
286 {
287
288 // Pt
289
290 fhPtPhotonLeading[i] = new TH1F(Form("hPtPhotonLeading%s",name[i].Data()),
12d80e80 291 Form("#gamma: Leading of all particles%s",title[i].Data()),
dbb79a05 292 nptbins,ptmin,ptmax);
8cc95539 293 fhPtPhotonLeading[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
7b2086c3 294 outputContainer->Add(fhPtPhotonLeading[i]);
295
296 fhPtPi0Leading[i] = new TH1F(Form("hPtPi0Leading%s",name[i].Data()),
12d80e80 297 Form("#pi^{0}: Leading of all particles%s",title[i].Data()),
dbb79a05 298 nptbins,ptmin,ptmax);
8cc95539 299 fhPtPi0Leading[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
dbb79a05 300 outputContainer->Add(fhPtPi0Leading[i]);
7b2086c3 301
302 fhPtPhotonLeadingIsolated[i] = new TH1F(Form("hPtPhotonLeadingIsolated%s",name[i].Data()),
12d80e80 303 Form("#gamma: Leading of all particles%s, isolated",title[i].Data()),
dbb79a05 304 nptbins,ptmin,ptmax);
8cc95539 305 fhPtPhotonLeadingIsolated[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
7b2086c3 306 outputContainer->Add(fhPtPhotonLeadingIsolated[i]);
307
308 fhPtPi0LeadingIsolated[i] = new TH1F(Form("hPtPi0LeadingIsolated%s",name[i].Data()),
12d80e80 309 Form("#pi^{0}: Leading of all particles%s, isolated",title[i].Data()),
dbb79a05 310 nptbins,ptmin,ptmax);
8cc95539 311 fhPtPi0LeadingIsolated[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
dbb79a05 312 outputContainer->Add(fhPtPi0LeadingIsolated[i]);
7b2086c3 313
c76fb00a 314 fhPtPhotonLeadingSumPt[i] = new TH2F(Form("hPtPhotonLeadingSumPt%s",name[i].Data()),
12d80e80 315 Form("#gamma: Leading of all particles%s",title[i].Data()),
c76fb00a 316 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
8cc95539 317 fhPtPhotonLeadingSumPt[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
318 fhPtPhotonLeadingSumPt[i]->SetYTitle("#Sigma #it{p}_{T} (GeV/#it{c})");
c76fb00a 319 outputContainer->Add(fhPtPhotonLeadingSumPt[i]);
320
321 fhPtPi0LeadingSumPt[i] = new TH2F(Form("hPtPi0LeadingSumPt%s",name[i].Data()),
12d80e80 322 Form("#pi^{0}: Leading of all particles%s",title[i].Data()),
c76fb00a 323 nptbins,ptmin,ptmax,nptsumbins,ptsummin,ptsummax);
8cc95539 324 fhPtPi0LeadingSumPt[i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
325 fhPtPi0LeadingSumPt[i]->SetYTitle("#Sigma #it{p}_{T} (GeV/#it{c})");
c76fb00a 326 outputContainer->Add(fhPtPi0LeadingSumPt[i]);
327
328
dbb79a05 329 // Leading or not loop
330 for(Int_t j = 0; j < 2; j++)
331 {
332 // Near side parton
333
334 fhPtPartonTypeNearPhoton[j][i] = new TH2F(Form("hPtPartonTypeNearPhoton%s%s",leading[j].Data(),name[i].Data()),
12d80e80 335 Form("#gamma: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 336 nptbins,ptmin,ptmax,3,0,3);
8cc95539 337 fhPtPartonTypeNearPhoton[j][i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
dbb79a05 338 fhPtPartonTypeNearPhoton[j][i]->SetYTitle("Parton type");
339 fhPtPartonTypeNearPhoton[j][i]->GetYaxis()->SetBinLabel(1,"#gamma");
340 fhPtPartonTypeNearPhoton[j][i]->GetYaxis()->SetBinLabel(2,"g");
341 fhPtPartonTypeNearPhoton[j][i]->GetYaxis()->SetBinLabel(3,"q");
342 outputContainer->Add(fhPtPartonTypeNearPhoton[j][i]);
343
344 fhPtPartonTypeNearPi0[j][i] = new TH2F(Form("hPtPartonTypeNearPi0%s%s",leading[j].Data(),name[i].Data()),
12d80e80 345 Form("#pi^{0}: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 346 nptbins,ptmin,ptmax,3,0,3);
8cc95539 347 fhPtPartonTypeNearPi0[j][i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
dbb79a05 348 fhPtPartonTypeNearPi0[j][i]->SetYTitle("Parton type");
349 fhPtPartonTypeNearPi0[j][i]->GetYaxis()->SetBinLabel(1,"#gamma");
350 fhPtPartonTypeNearPi0[j][i]->GetYaxis()->SetBinLabel(2,"g");
351 fhPtPartonTypeNearPi0[j][i]->GetYaxis()->SetBinLabel(3,"q");
352 outputContainer->Add(fhPtPartonTypeNearPi0[j][i]);
353
354 fhPtPartonTypeNearPhotonIsolated[j][i] = new TH2F(Form("hPtPartonTypeNearPhoton%sIsolated%s",leading[j].Data(),name[i].Data()),
12d80e80 355 Form("#gamma: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 356 nptbins,ptmin,ptmax,3,0,3);
8cc95539 357 fhPtPartonTypeNearPhotonIsolated[j][i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
dbb79a05 358 fhPtPartonTypeNearPhotonIsolated[j][i]->SetYTitle("Parton type");
359 fhPtPartonTypeNearPhotonIsolated[j][i]->GetYaxis()->SetBinLabel(1,"#gamma");
360 fhPtPartonTypeNearPhotonIsolated[j][i]->GetYaxis()->SetBinLabel(2,"g");
361 fhPtPartonTypeNearPhotonIsolated[j][i]->GetYaxis()->SetBinLabel(3,"q");
362 outputContainer->Add(fhPtPartonTypeNearPhotonIsolated[j][i]);
363
364 fhPtPartonTypeNearPi0Isolated[j][i] = new TH2F(Form("hPtPartonTypeNearPi0%sIsolated%s",leading[j].Data(),name[i].Data()),
12d80e80 365 Form("#pi^{0}: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 366 nptbins,ptmin,ptmax,3,0,3);
8cc95539 367 fhPtPartonTypeNearPi0Isolated[j][i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
dbb79a05 368 fhPtPartonTypeNearPi0Isolated[j][i]->SetYTitle("Parton type");
369 fhPtPartonTypeNearPi0Isolated[j][i]->GetYaxis()->SetBinLabel(1,"#gamma");
370 fhPtPartonTypeNearPi0Isolated[j][i]->GetYaxis()->SetBinLabel(2,"g");
371 fhPtPartonTypeNearPi0Isolated[j][i]->GetYaxis()->SetBinLabel(3,"q");
372 outputContainer->Add(fhPtPartonTypeNearPi0Isolated[j][i]);
373
374
375 // Away side parton
376
377 fhPtPartonTypeAwayPhoton[j][i] = new TH2F(Form("hPtPartonTypeAwayPhoton%s%s",leading[j].Data(),name[i].Data()),
12d80e80 378 Form("#gamma: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 379 nptbins,ptmin,ptmax,3,0,3);
8cc95539 380 fhPtPartonTypeAwayPhoton[j][i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
dbb79a05 381 fhPtPartonTypeAwayPhoton[j][i]->SetYTitle("Parton type");
382 fhPtPartonTypeAwayPhoton[j][i]->GetYaxis()->SetBinLabel(1,"#gamma");
383 fhPtPartonTypeAwayPhoton[j][i]->GetYaxis()->SetBinLabel(2,"g");
384 fhPtPartonTypeAwayPhoton[j][i]->GetYaxis()->SetBinLabel(3,"q");
385 outputContainer->Add(fhPtPartonTypeAwayPhoton[j][i]);
386
387 fhPtPartonTypeAwayPi0[j][i] = new TH2F(Form("hPtPartonTypeAwayPi0%s%s",leading[j].Data(),name[i].Data()),
12d80e80 388 Form("#pi^{0}: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 389 nptbins,ptmin,ptmax,3,0,3);
8cc95539 390 fhPtPartonTypeAwayPi0[j][i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
dbb79a05 391 fhPtPartonTypeAwayPi0[j][i]->SetYTitle("Parton type");
392 fhPtPartonTypeAwayPi0[j][i]->GetYaxis()->SetBinLabel(1,"#gamma");
393 fhPtPartonTypeAwayPi0[j][i]->GetYaxis()->SetBinLabel(2,"g");
394 fhPtPartonTypeAwayPi0[j][i]->GetYaxis()->SetBinLabel(3,"q");
395 outputContainer->Add(fhPtPartonTypeAwayPi0[j][i]);
396
397 fhPtPartonTypeAwayPhotonIsolated[j][i] = new TH2F(Form("hPtPartonTypeAwayPhoton%sIsolated%s",leading[j].Data(),name[i].Data()),
12d80e80 398 Form("#gamma: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 399 nptbins,ptmin,ptmax,3,0,3);
8cc95539 400 fhPtPartonTypeAwayPhotonIsolated[j][i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
dbb79a05 401 fhPtPartonTypeAwayPhotonIsolated[j][i]->SetYTitle("Parton type");
402 fhPtPartonTypeAwayPhotonIsolated[j][i]->GetYaxis()->SetBinLabel(1,"#gamma");
403 fhPtPartonTypeAwayPhotonIsolated[j][i]->GetYaxis()->SetBinLabel(2,"g");
404 fhPtPartonTypeAwayPhotonIsolated[j][i]->GetYaxis()->SetBinLabel(3,"q");
405 outputContainer->Add(fhPtPartonTypeAwayPhotonIsolated[j][i]);
406
407 fhPtPartonTypeAwayPi0Isolated[j][i] = new TH2F(Form("hPtPartonTypeAwayPi0%sIsolated%s",leading[j].Data(),name[i].Data()),
12d80e80 408 Form("#pi^{0}: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 409 nptbins,ptmin,ptmax,3,0,3);
8cc95539 410 fhPtPartonTypeAwayPi0Isolated[j][i]->SetXTitle("#it{p}_{T} (GeV/#it{c})");
dbb79a05 411 fhPtPartonTypeAwayPi0Isolated[j][i]->SetYTitle("Parton type");
412 fhPtPartonTypeAwayPi0Isolated[j][i]->GetYaxis()->SetBinLabel(1,"#gamma");
413 fhPtPartonTypeAwayPi0Isolated[j][i]->GetYaxis()->SetBinLabel(2,"g");
414 fhPtPartonTypeAwayPi0Isolated[j][i]->GetYaxis()->SetBinLabel(3,"q");
415 outputContainer->Add(fhPtPartonTypeAwayPi0Isolated[j][i]);
416
417 // zHard
418
419 fhZHardPhoton[j][i] = new TH2F(Form("hZHardPhoton%s%s",leading[j].Data(),name[i].Data()),
8cc95539 420 Form("#it{z}_{Hard} of #gamma: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 421 nptbins,ptmin,ptmax,200,0,2);
8cc95539 422 fhZHardPhoton[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
423 fhZHardPhoton[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 424 outputContainer->Add(fhZHardPhoton[j][i]);
425
426 fhZHardPi0[j][i] = new TH2F(Form("hZHardPi0%s%s",leading[j].Data(),name[i].Data()),
8cc95539 427 Form("#it{z}_{Hard} of #pi^{0}: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 428 nptbins,ptmin,ptmax,200,0,2);
8cc95539 429 fhZHardPi0[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
430 fhZHardPi0[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 431 outputContainer->Add(fhZHardPi0[j][i]);
432
433 fhZHardPhotonIsolated[j][i] = new TH2F(Form("hZHardPhoton%sIsolated%s",leading[j].Data(),name[i].Data()),
8cc95539 434 Form("#it{z}_{Hard} of #gamma: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 435 nptbins,ptmin,ptmax,200,0,2);
8cc95539 436 fhZHardPhotonIsolated[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
437 fhZHardPhotonIsolated[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 438 outputContainer->Add(fhZHardPhotonIsolated[j][i]);
439
440 fhZHardPi0Isolated[j][i] = new TH2F(Form("hZHardPi0%sIsolated%s",leading[j].Data(),name[i].Data()),
8cc95539 441 Form("#it{z}_{Hard} of #pi^{0}: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 442 nptbins,ptmin,ptmax,200,0,2);
8cc95539 443 fhZHardPi0Isolated[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
444 fhZHardPi0Isolated[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 445 outputContainer->Add(fhZHardPi0Isolated[j][i]);
446
447 // zHard
448
449 fhZPartonPhoton[j][i] = new TH2F(Form("hZPartonPhoton%s%s",leading[j].Data(),name[i].Data()),
8cc95539 450 Form("#it{z}_{Parton} of #gamma: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 451 nptbins,ptmin,ptmax,200,0,2);
8cc95539 452 fhZPartonPhoton[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
453 fhZPartonPhoton[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 454 outputContainer->Add(fhZPartonPhoton[j][i]);
455
456 fhZPartonPi0[j][i] = new TH2F(Form("hZPartonPi0%s%s",leading[j].Data(),name[i].Data()),
8cc95539 457 Form("#it{z}_{Parton} of #pi^{0}: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 458 nptbins,ptmin,ptmax,200,0,2);
8cc95539 459 fhZPartonPi0[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
460 fhZPartonPi0[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 461 outputContainer->Add(fhZPartonPi0[j][i]);
462
463 fhZPartonPhotonIsolated[j][i] = new TH2F(Form("hZPartonPhoton%sIsolated%s",leading[j].Data(),name[i].Data()),
8cc95539 464 Form("#it{z}_{Parton} of #gamma: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 465 nptbins,ptmin,ptmax,200,0,2);
8cc95539 466 fhZPartonPhotonIsolated[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
467 fhZPartonPhotonIsolated[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 468 outputContainer->Add(fhZPartonPhotonIsolated[j][i]);
469
470 fhZPartonPi0Isolated[j][i] = new TH2F(Form("hZPartonPi0%sIsolated%s",leading[j].Data(),name[i].Data()),
8cc95539 471 Form("#it{z}_{Parton} of #pi^{0}: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 472 nptbins,ptmin,ptmax,200,0,2);
8cc95539 473 fhZPartonPi0Isolated[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
474 fhZPartonPi0Isolated[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 475 outputContainer->Add(fhZPartonPi0Isolated[j][i]);
476
477
478 // zJet
479
480 fhZJetPhoton[j][i] = new TH2F(Form("hZJetPhoton%s%s",leading[j].Data(),name[i].Data()),
8cc95539 481 Form("#it{z}_{Jet} of #gamma: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 482 nptbins,ptmin,ptmax,200,0,2);
8cc95539 483 fhZJetPhoton[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
484 fhZJetPhoton[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 485 outputContainer->Add(fhZJetPhoton[j][i]);
486
487 fhZJetPi0[j][i] = new TH2F(Form("hZJetPi0%s%s",leading[j].Data(),name[i].Data()),
8cc95539 488 Form("#it{z}_{Jet} of #pi^{0}: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 489 nptbins,ptmin,ptmax,200,0,2);
8cc95539 490 fhZJetPi0[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
491 fhZJetPi0[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 492 outputContainer->Add(fhZJetPi0[j][i]);
493
494 fhZJetPhotonIsolated[j][i] = new TH2F(Form("hZJetPhoton%sIsolated%s",leading[j].Data(),name[i].Data()),
8cc95539 495 Form("#it{z}_{Jet} of #gamma: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 496 nptbins,ptmin,ptmax,200,0,2);
8cc95539 497 fhZJetPhotonIsolated[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
498 fhZJetPhotonIsolated[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 499 outputContainer->Add(fhZJetPhotonIsolated[j][i]);
500
501 fhZJetPi0Isolated[j][i] = new TH2F(Form("hZJetPi0%sIsolated%s",leading[j].Data(),name[i].Data()),
8cc95539 502 Form("#it{z}_{Jet} of #pi^{0}: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 503 nptbins,ptmin,ptmax,200,0,2);
8cc95539 504 fhZJetPi0Isolated[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
505 fhZJetPi0Isolated[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 506 outputContainer->Add(fhZJetPi0Isolated[j][i]);
507
508
509 // XE
510
511 fhXEPhoton[j][i] = new TH2F(Form("hXEPhoton%s%s",leading[j].Data(),name[i].Data()),
8cc95539 512 Form("#it{z}_{Jet} of #gamma: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 513 nptbins,ptmin,ptmax,200,0,2);
8cc95539 514 fhXEPhoton[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
515 fhXEPhoton[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 516 outputContainer->Add(fhXEPhoton[j][i]);
517
518 fhXEPi0[j][i] = new TH2F(Form("hXEPi0%s%s",leading[j].Data(),name[i].Data()),
8cc95539 519 Form("#it{z}_{Jet} of #pi^{0}: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 520 nptbins,ptmin,ptmax,200,0,2);
8cc95539 521 fhXEPi0[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
522 fhXEPi0[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 523 outputContainer->Add(fhXEPi0[j][i]);
524
525 fhXEPhotonIsolated[j][i] = new TH2F(Form("hXEPhoton%sIsolated%s",leading[j].Data(),name[i].Data()),
8cc95539 526 Form("#it{z}_{Jet} of #gamma: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 527 nptbins,ptmin,ptmax,200,0,2);
8cc95539 528 fhXEPhotonIsolated[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
529 fhXEPhotonIsolated[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 530 outputContainer->Add(fhXEPhotonIsolated[j][i]);
531
532 fhXEPi0Isolated[j][i] = new TH2F(Form("hXEPi0%sIsolated%s",leading[j].Data(),name[i].Data()),
8cc95539 533 Form("#it{z}_{Jet} of #pi^{0}: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 534 nptbins,ptmin,ptmax,200,0,2);
8cc95539 535 fhXEPi0Isolated[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
536 fhXEPi0Isolated[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 537 outputContainer->Add(fhXEPi0Isolated[j][i]);
538
539
540 // XE from UE
541
542 fhXEUEPhoton[j][i] = new TH2F(Form("hXEUEPhoton%s%s",leading[j].Data(),name[i].Data()),
8cc95539 543 Form("#it{z}_{Jet} of #gamma: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 544 nptbins,ptmin,ptmax,200,0,2);
8cc95539 545 fhXEUEPhoton[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
546 fhXEUEPhoton[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 547 outputContainer->Add(fhXEUEPhoton[j][i]);
548
549 fhXEUEPi0[j][i] = new TH2F(Form("hXEUEPi0%s%s",leading[j].Data(),name[i].Data()),
8cc95539 550 Form("#it{z}_{Jet} of #pi^{0}: %s of all particles%s",leading[j].Data(),title[i].Data()),
dbb79a05 551 nptbins,ptmin,ptmax,200,0,2);
8cc95539 552 fhXEUEPi0[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
553 fhXEUEPi0[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 554 outputContainer->Add(fhXEUEPi0[j][i]);
555
556 fhXEUEPhotonIsolated[j][i] = new TH2F(Form("hXEUEPhoton%sIsolated%s",leading[j].Data(),name[i].Data()),
8cc95539 557 Form("#it{z}_{Jet} of #gamma: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
dbb79a05 558 nptbins,ptmin,ptmax,200,0,2);
8cc95539 559 fhXEUEPhotonIsolated[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
560 fhXEUEPhotonIsolated[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 561 outputContainer->Add(fhXEUEPhotonIsolated[j][i]);
562
563 fhXEUEPi0Isolated[j][i] = new TH2F(Form("hXEUEPi0%sIsolated%s",leading[j].Data(),name[i].Data()),
8cc95539 564 Form("#it{z}_{Jet} of #pi^{0}: %s of all particles%s, isolated",leading[j].Data(),title[i].Data()),
7b2086c3 565 nptbins,ptmin,ptmax,200,0,2);
8cc95539 566 fhXEUEPi0Isolated[j][i]->SetYTitle("#it{p}_{T}^{particle}/#it{p}_{T}^{hard}");
567 fhXEUEPi0Isolated[j][i]->SetXTitle("#it{p}_{T}^{particle} (GeV/#it{c})");
dbb79a05 568 outputContainer->Add(fhXEUEPi0Isolated[j][i]);
569 }
7b2086c3 570 }
571
572 return outputContainer;
573
574}
575
576//____________________________________________
577void AliAnaGeneratorKine::GetPartonsAndJets()
578{
579 // Fill data members with partons,jets and generated pt hard
580
8cc95539 581 if(GetDebug() > 1) printf("AliAnaGeneratorKine::GetPartonsAndJets() - Start \n");
582
7b2086c3 583 fStack = GetMCStack() ;
584
585 if(!fStack)
586 {
587 printf("AliAnaGeneratorKine::MakeAnalysisFillHistograms() - No Stack available, STOP\n");
588 abort();
589 }
590
591 fParton2 = fStack->Particle(2) ;
592 fParton3 = fStack->Particle(3) ;
593 fParton6 = fStack->Particle(6) ;
594 fParton7 = fStack->Particle(7) ;
595
596 Float_t p6phi = fParton6->Phi();
597 if(p6phi < 0) p6phi +=TMath::TwoPi();
598 Float_t p7phi = fParton7->Phi();
599 if(p7phi < 0) p7phi +=TMath::TwoPi();
600
601 //printf("parton6: pt %2.2f, eta %2.2f, phi %2.2f with pdg %d\n",fParton6->Pt(),fParton6->Eta(),p6phi, fParton6->GetPdgCode());
602 //printf("parton7: pt %2.2f, eta %2.2f, phi %2.2f with pdg %d\n",fParton7->Pt(),fParton7->Eta(),p7phi, fParton7->GetPdgCode());
603
604 // Get the jets, only for pythia
605 if(!strcmp(GetReader()->GetGenEventHeader()->ClassName(), "AliGenPythiaEventHeader"))
606 {
607 AliGenPythiaEventHeader* pygeh= (AliGenPythiaEventHeader*) GetReader()->GetGenEventHeader();
608
609 fPtHard = pygeh->GetPtHard();
610
611 //printf("pt Hard %2.2f\n",fPtHard);
612
613 const Int_t nTriggerJets = pygeh->NTriggerJets();
614
615 Float_t tmpjet[]={0,0,0,0};
616
617 // select the closest jet to parton
618 Float_t jet7R = 100;
619 Float_t jet6R = 100;
620
621 for(Int_t ijet = 0; ijet< nTriggerJets; ijet++)
622 {
623 pygeh->TriggerJet(ijet, tmpjet);
624
7d409bf9 625 fLVTmp.SetPxPyPzE(tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3]);
626 Float_t jphi = fLVTmp.Phi();
7b2086c3 627 if(jphi < 0) jphi +=TMath::TwoPi();
628
7d409bf9 629 Double_t radius6 = GetIsolationCut()->Radius(fParton6->Eta(), p6phi, fLVTmp.Eta() , jphi) ;
630 Double_t radius7 = GetIsolationCut()->Radius(fParton7->Eta(), p7phi, fLVTmp.Eta() , jphi) ;
7b2086c3 631
632 //printf("jet %d: pt %2.2f, eta %2.2f, phi %2.2f, r6 %2.2f, r7 %2.2f\n",ijet,jet.Pt(),jet.Eta(),jphi,radius6, radius7);
633
634 if (radius6 < jet6R)
635 {
636 jet6R = radius6;
7d409bf9 637 fJet6 = fLVTmp;
7b2086c3 638
639 }
640 if (radius7 < jet7R)
641 {
642 jet7R = radius7;
7d409bf9 643 fJet7 = fLVTmp;
7b2086c3 644 }
645
646 } // jet loop
647
648 //printf("jet6: pt %2.2f, eta %2.2f, phi %2.2f\n",fJet6.Pt(),fJet6.Eta(),fJet6.Phi());
649 //printf("jet7: pt %2.2f, eta %2.2f, phi %2.2f\n",fJet7.Pt(),fJet7.Eta(),fJet6.Phi());
650
651 } // pythia header
652
653 fhPtHard ->Fill(fPtHard);
654 fhPtJet ->Fill(fJet6.Pt());
655 fhPtJet ->Fill(fJet7.Pt());
656 fhPtParton ->Fill(fParton6->Pt());
657 fhPtParton ->Fill(fParton7->Pt());
658
8cc95539 659 fhPtPartonPtHard->Fill(fPtHard, fParton6->Pt()/fPtHard);
660 fhPtPartonPtHard->Fill(fPtHard, fParton7->Pt()/fPtHard);
661 fhPtJetPtHard ->Fill(fPtHard, fJet6.Pt()/fPtHard);
662 fhPtJetPtHard ->Fill(fPtHard, fJet7.Pt()/fPtHard);
663 fhPtJetPtParton ->Fill(fPtHard, fJet6.Pt()/fParton6->Pt());
664 fhPtJetPtParton ->Fill(fPtHard, fJet7.Pt()/fParton7->Pt());
665
666 if(GetDebug() > 1) printf("AliAnaGeneratorKine::GetPartonsAndJets() - End \n");
667
7b2086c3 668}
669
b94e038e 670//_____________________________________________________
7d409bf9 671void AliAnaGeneratorKine::GetXE(Int_t indexTrig,
b94e038e 672 Int_t pdgTrig,
673 Bool_t leading[4],
674 Bool_t isolated[4],
675 Int_t iparton)
7b2086c3 676{
677
678 // Calculate the real XE and the UE XE
679
8cc95539 680 if(GetDebug() > 1) printf("AliAnaGeneratorKine::GetXE() - Start \n");
681
7d409bf9 682 Float_t ptTrig = fTrigger.Pt();
683 Float_t phiTrig = fTrigger.Phi();
7b2086c3 684 if(phiTrig < 0 ) phiTrig += TMath::TwoPi();
685
686 //Loop on primaries, start from position 8, no partons
687 for(Int_t ipr = 8; ipr < fStack->GetNprimary(); ipr ++ )
688 {
7b2086c3 689 TParticle * particle = fStack->Particle(ipr) ;
690
691 if(ipr==indexTrig) continue;
692
693 Int_t pdg = particle->GetPdgCode();
694 Int_t status = particle->GetStatusCode();
695
696 // Compare trigger with final state particles
94386973 697 if( status != 1 ) continue ;
7b2086c3 698
699 Double_t charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
700
94386973 701 if(charge == 0 ) continue; // construct xe only with charged
7b2086c3 702
703 Float_t pt = particle->Pt();
7b2086c3 704 Float_t phi = particle->Phi();
705 if(phi < 0 ) phi += TMath::TwoPi();
706
783b974c 707 if( pt < fMinChargedPt) continue ;
7b2086c3 708
7d409bf9 709 particle->Momentum(fLVTmp);
710 Bool_t inTPC = GetFiducialCut()->IsInFiducialCut(fLVTmp.Eta(),fLVTmp.Phi(),kCTS) ;
94386973 711
712 if(!inTPC) continue;
7b2086c3 713
714 Float_t xe = -pt/ptTrig*TMath::Cos(phi-phiTrig);
715
94386973 716 // ---------------------------------------------------
717 // Get the index of the mother, get from what parton
7b2086c3 718 Int_t ipartonAway = particle->GetFirstMother();
23fa04c5 719 if(ipartonAway < 0) return;
94386973 720
7b2086c3 721 TParticle * mother = fStack->Particle(ipartonAway);
722 while (ipartonAway > 7)
723 {
724 ipartonAway = mother->GetFirstMother();
23fa04c5 725 if(ipartonAway < 0) break;
7b2086c3 726 mother = fStack->Particle(ipartonAway);
727 }
728
94386973 729 //-----------------------------------------
730 // Get XE of particles belonging to the jet
731 // on the opposite side of the trigger
7b2086c3 732 if((ipartonAway==6 || ipartonAway==7) && iparton!=ipartonAway)
733 {
7b2086c3 734 for( Int_t i = 0; i < 4; i++ )
735 {
736 if(pdgTrig==111)
737 {
dbb79a05 738 fhXEPi0[leading[i]][i] ->Fill(ptTrig,xe);
739
740 if(isolated[i])
741 {
742 fhXEPi0Isolated[leading[i]][i] ->Fill(ptTrig,xe);
7b2086c3 743 }
dbb79a05 744
7b2086c3 745 }// pi0
746 else if(pdgTrig==22)
747 {
dbb79a05 748 fhXEPhoton[leading[i]][i] ->Fill(ptTrig,xe);
749
750 if(isolated[i])
751 {
752 fhXEPhotonIsolated[leading[i]][i] ->Fill(ptTrig,xe);
753 }
754
7b2086c3 755 } // photon
dbb79a05 756 } // conditions loop
7b2086c3 757 } // Away side
758
94386973 759 //----------------------------------------------------------
760 // Get the XE from particles not attached to any of the jets
dbb79a05 761 if(ipartonAway!=6 && ipartonAway!=7)
7b2086c3 762 {
7b2086c3 763 for( Int_t i = 0; i < 4; i++ )
764 {
765 if(pdgTrig==111)
766 {
dbb79a05 767 fhXEUEPi0[leading[i]][i] ->Fill(ptTrig,xe);
768
769 if(isolated[i])
770 {
771 fhXEUEPi0Isolated[leading[i]][i] ->Fill(ptTrig,xe);
7b2086c3 772 }
dbb79a05 773
7b2086c3 774 }// pi0
775 else if(pdgTrig==22)
776 {
dbb79a05 777 fhXEUEPhoton[leading[i]][i] ->Fill(ptTrig,xe);
778
779 if(isolated[i])
780 {
781 fhXEUEPhotonIsolated[leading[i]][i] ->Fill(ptTrig,xe);
782 }
783
7b2086c3 784 } // photon
785 } // conditions loop
786 } // Away side
787
788 } // primary loop
789
8cc95539 790 if(GetDebug() > 1) printf("AliAnaGeneratorKine::GetPartonsAndJets() - End \n");
7b2086c3 791
792}
793
794//________________________________________
795void AliAnaGeneratorKine::InitParameters()
796{
797
798 //Initialize the parameters of the analysis.
799 AddToHistogramsName("AnaGenKine_");
800
1290eee4 801 fTriggerDetector = kEMCAL;
783b974c 802
803 fMinChargedPt = 0.2;
804 fMinNeutralPt = 0.3;
805
7b2086c3 806}
807
b94e038e 808//_____________________________________________________________________
7d409bf9 809void AliAnaGeneratorKine::IsLeadingAndIsolated(Int_t indexTrig,
b94e038e 810 Int_t pdgTrig,
7b2086c3 811 Bool_t leading[4],
812 Bool_t isolated[4])
813{
dbb79a05 814 // Check if the trigger is the leading particle and if it is isolated
7b2086c3 815 // In case of neutral particles check all neutral or neutral in EMCAL acceptance
816
8cc95539 817 if(GetDebug() > 1) printf("AliAnaGeneratorKine::GetIsLeadingAndIsolated() - Start \n");
818
7b2086c3 819 Float_t ptMaxCharged = 0; // all charged
820 Float_t ptMaxNeutral = 0; // all neutral
821 Float_t ptMaxNeutEMCAL = 0; // for neutral, select them in EMCAL acceptance
822 Float_t ptMaxNeutPhot = 0; // for neutral, take only photons
823 Float_t ptMaxNeutEMCALPhot = 0; // for neutral, take only photons in EMCAL acceptance
824
825 leading[0] = 0;
826 leading[1] = 0;
827 leading[2] = 0;
828 leading[3] = 0;
829
830 isolated[0] = 0;
831 isolated[1] = 0;
832 isolated[2] = 0;
833 isolated[3] = 0;
834
7d409bf9 835 Float_t ptTrig = fTrigger.Pt();
836 Float_t etaTrig = fTrigger.Eta();
837 Float_t phiTrig = fTrigger.Phi();
7b2086c3 838 if(phiTrig < 0 ) phiTrig += TMath::TwoPi();
839
840 // Minimum track or cluster energy
7b2086c3 841
842 //Isolation cuts
c76fb00a 843 Float_t ptThresIC = GetIsolationCut()->GetPtThreshold();
844 Float_t sumThresIC = GetIsolationCut()->GetPtThreshold();
783b974c 845 Float_t rThresIC = GetIsolationCut()->GetConeSize();
c76fb00a 846 Float_t isoMethod = GetIsolationCut()->GetICMethod();
847
848 //Counters
7b2086c3 849 Int_t nICTrack = 0;
850 Int_t nICNeutral = 0;
851 Int_t nICNeutEMCAL = 0;
852 Int_t nICNeutPhot = 0;
853 Int_t nICNeutEMCALPhot = 0;
854
c76fb00a 855 // Sum of pT
856 Float_t sumNePt = 0;
857 Float_t sumNePtPhot = 0;
858 Float_t sumNePtEMC = 0;
859 Float_t sumNePtEMCPhot = 0;
860 Float_t sumChPt = 0;
861
7b2086c3 862 //Loop on primaries, start from position 8, no partons
863 for(Int_t ipr = 8; ipr < fStack->GetNprimary(); ipr ++ )
864 {
865 if(ipr == indexTrig) continue;
866 TParticle * particle = fStack->Particle(ipr) ;
867
868 Int_t imother= particle->GetFirstMother();
869 //printf("Leading ipr %d - mother %d\n",ipr, imother);
870
871 if(imother==indexTrig) continue ;
872
873 Int_t pdg = particle->GetPdgCode();
874 Int_t status = particle->GetStatusCode();
875
876 // Compare trigger with final state particles
877 if( status != 1) continue ;
cecdfdb3 878
879 // Select all particles in at least the TPC acceptance
7d409bf9 880 Bool_t inTPC = GetFiducialCut()->IsInFiducialCut(fTrigger.Eta(),fTrigger.Phi(),kCTS) ;
cecdfdb3 881 if(!inTPC) continue;
7b2086c3 882
883 Float_t pt = particle->Pt();
884 Float_t eta = particle->Eta();
885 Float_t phi = particle->Phi();
886 if(phi < 0 ) phi += TMath::TwoPi();
887
7b2086c3 888 Double_t charge = TDatabasePDG::Instance()->GetParticle(pdg)->Charge();
889
890 //Isolation
891 Double_t radius = GetIsolationCut()->Radius(etaTrig, phiTrig, eta , phi) ;
892
893 if(charge==0)
894 {
783b974c 895 if(pt < fMinNeutralPt) continue ;
7b2086c3 896
c76fb00a 897 if( ptMaxNeutral < pt ) ptMaxNeutral = pt;
898
899 if( radius < rThresIC )
900 {
901 if( pt > ptThresIC ) nICNeutral++ ;
902 sumNePt+= pt;
903 }
904
7b2086c3 905 Bool_t phPDG = kFALSE;
906 if(pdg==22 || pdg==111) phPDG = kTRUE;
907
908 //if(pt > ptTrig) printf(" --- pdg %d, phPDG %d pT %2.2f, pTtrig %2.2f, eta %2.2f, phi %2.2f\n",pdg,phPDG,pt,ptTrig,particle->Eta(), particle->Phi()*TMath::RadToDeg());
909 if(phPDG)
910 {
c76fb00a 911 if( ptMaxNeutPhot < pt) ptMaxNeutPhot = pt;
912
913 if( radius < rThresIC )
914 {
915 if(pt > ptThresIC) nICNeutPhot++ ;
916 sumNePtPhot += pt;
917 }
7b2086c3 918 }
919
229bed7a 920 //Calorimeter acceptance
7d409bf9 921 Bool_t inCalo = GetFiducialCut()->IsInFiducialCut(fTrigger.Eta(),fTrigger.Phi(),GetCalorimeter()) ;
229bed7a 922 if(!inCalo) continue;
923
c76fb00a 924 if( ptMaxNeutEMCAL < pt ) ptMaxNeutEMCAL = pt;
925 if( radius < rThresIC )
926 {
927 if( pt > ptThresIC ) nICNeutEMCAL++ ;
928 sumNePtEMC += pt;
929 }
229bed7a 930
931 if(phPDG)
7b2086c3 932 {
229bed7a 933 if( ptMaxNeutEMCALPhot < pt ) ptMaxNeutEMCALPhot = pt;
c76fb00a 934 if( radius < rThresIC )
935 {
936 if (pt > ptThresIC) nICNeutEMCALPhot++ ;
937 sumNePtEMCPhot += pt;
938 }
7b2086c3 939 }
940 }
229bed7a 941 else
7b2086c3 942 {
783b974c 943 if( pt < fMinChargedPt) continue ;
c76fb00a 944
7b2086c3 945 if( ptMaxCharged < pt ) ptMaxCharged = pt;
946
c76fb00a 947 if( radius < rThresIC )
7b2086c3 948 {
949 //printf("UE track? pTtrig %2.2f, pt %2.2f, etaTrig %2.2f, eta %2.2f, phiTrig %2.2f, phi %2.2f, radius %2.2f\n",
950 // ptTrig, pt,etaTrig, eta, phiTrig*TMath::RadToDeg(), phi*TMath::RadToDeg(),radius);
c76fb00a 951 if( pt > ptThresIC ) nICTrack++ ;
952 sumChPt += pt;
7b2086c3 953 }
954 }
955
956 } // particle loop
957
958 // Leding decision
959 if(ptTrig > ptMaxCharged)
960 {
961 //printf("pt charged %2.2f, pt neutral %2.2f, pt neutral emcal %2.2f, pt photon %2.2f, pt photon emcal %2.2f\n",
962 // ptMaxCharged, ptMaxNeutral, ptMaxNeutEMCAL,ptMaxNeutPhot, ptMaxNeutEMCALPhot);
963 if(ptTrig > ptMaxNeutral ) leading[0] = kTRUE ;
964 if(ptTrig > ptMaxNeutEMCAL ) leading[1] = kTRUE ;
965 if(ptTrig > ptMaxNeutPhot ) leading[2] = kTRUE ;
966 if(ptTrig > ptMaxNeutEMCALPhot) leading[3] = kTRUE ;
967 }
968
12d80e80 969 //printf("N in cone over threshold: tracks %d, neutral %d, neutral emcal %d, photon %d, photon emcal %d\n",
7b2086c3 970 // nICTrack, nICNeutral ,nICNeutEMCAL,nICNeutPhot, nICNeutEMCALPhot);
971
c76fb00a 972 //------------------
7b2086c3 973 // Isolation decision
c76fb00a 974 if( isoMethod == AliIsolationCut::kPtThresIC )
975 {
976 if( nICTrack == 0 )
977 {
978 if(nICNeutral == 0 ) isolated[0] = kTRUE ;
979 if(nICNeutEMCAL == 0 ) isolated[1] = kTRUE ;
980 if(nICNeutPhot == 0 ) isolated[2] = kTRUE ;
981 if(nICNeutEMCALPhot == 0 ) isolated[3] = kTRUE ;
982 }
983 }
984 else if( isoMethod == AliIsolationCut::kSumPtIC )
7b2086c3 985 {
c76fb00a 986 if(sumChPt + sumNePt < sumThresIC ) isolated[0] = kTRUE ;
987 if(sumChPt + sumNePtEMC < sumThresIC ) isolated[1] = kTRUE ;
988 if(sumChPt + sumNePtPhot < sumThresIC ) isolated[2] = kTRUE ;
989 if(sumChPt + sumNePtEMCPhot < sumThresIC ) isolated[3] = kTRUE ;
7b2086c3 990 }
991
c76fb00a 992
993 //----------------------------------------------------
7b2086c3 994 // Fill histograms if conditions apply for all 4 cases
995 for( Int_t i = 0; i < 4; i++ )
996 {
997 if(pdgTrig==111)
998 {
999 if(leading[i])
1000 {
1001 fhPtPi0Leading[i]->Fill(ptTrig);
c76fb00a 1002
1003 if (i == 0) fhPtPi0LeadingSumPt[i]->Fill(ptTrig, sumChPt + sumNePt);
1004 else if(i == 1) fhPtPi0LeadingSumPt[i]->Fill(ptTrig, sumChPt + sumNePtEMC);
1005 else if(i == 2) fhPtPi0LeadingSumPt[i]->Fill(ptTrig, sumChPt + sumNePtPhot);
1006 else if(i == 3) fhPtPi0LeadingSumPt[i]->Fill(ptTrig, sumChPt + sumNePtEMCPhot);
1007
7b2086c3 1008 if(isolated[i]) fhPtPi0LeadingIsolated[i]->Fill(ptTrig);
1009 }
1010 }// pi0
1011 else if(pdgTrig==22)
1012 {
1013 if(leading[i])
1014 {
1015 fhPtPhotonLeading[i]->Fill(ptTrig);
c76fb00a 1016
1017 if (i == 0) fhPtPhotonLeadingSumPt[i]->Fill(ptTrig, sumChPt + sumNePt);
1018 else if(i == 1) fhPtPhotonLeadingSumPt[i]->Fill(ptTrig, sumChPt + sumNePtEMC);
1019 else if(i == 2) fhPtPhotonLeadingSumPt[i]->Fill(ptTrig, sumChPt + sumNePtPhot);
1020 else if(i == 3) fhPtPhotonLeadingSumPt[i]->Fill(ptTrig, sumChPt + sumNePtEMCPhot);
1021
7b2086c3 1022 if(isolated[i]) fhPtPhotonLeadingIsolated[i]->Fill(ptTrig);
1023 }
1024 } // photon
1025 } // conditions loop
8cc95539 1026
1027 if(GetDebug() > 1) printf("AliAnaGeneratorKine::IsLeadingAndIsolated() - End \n");
1028
7b2086c3 1029}
1030
1031//_____________________________________________________
1032void AliAnaGeneratorKine::MakeAnalysisFillHistograms()
1033{
1034 //Particle-Parton Correlation Analysis, fill histograms
1035
8cc95539 1036 if(GetDebug() > 1) printf("AliAnaGeneratorKine::MakeAnalysisFillHistograms() - Start \n");
7b2086c3 1037
1038 GetPartonsAndJets();
1039
1040 for(Int_t ipr = 0; ipr < fStack->GetNprimary(); ipr ++ )
1041 {
1042 TParticle * particle = fStack->Particle(ipr) ;
1043
1044 Int_t pdgTrig = particle->GetPdgCode();
1045 Int_t statusTrig = particle->GetStatusCode();
1046 Int_t imother = particle->GetFirstMother();
1047 Float_t ptTrig = particle->Pt();
1048
1049 // Select final state photons (prompt, fragmentation) or pi0s
1050
1051 //Check the origin of the photon, accept if prompt or initial/final state radiation
1052 if(pdgTrig == 22 && statusTrig == 1)
1053 {
1054 if(imother > 8) continue;
1055 }
1056 // If not photon, trigger on pi0
1057 else if(pdgTrig != 111) continue;
1058
1059 // Acceptance and kinematical cuts
7d409bf9 1060 if( ptTrig < GetMinPt() ) continue ;
7b2086c3 1061
8cc95539 1062 // Recover the kinematics:
7d409bf9 1063 particle->Momentum(fTrigger);
8cc95539 1064
7d409bf9 1065 Bool_t in = GetFiducialCutForTrigger()->IsInFiducialCut(fTrigger.Eta(),fTrigger.Phi(),fTriggerDetector) ;
229bed7a 1066 if(! in ) continue ;
1067
8cc95539 1068 if( GetDebug() > 2) printf("Select trigger particle %d: pdg %d status %d, mother index %d, pT %2.2f, eta %2.2f, phi %2.2f \n",
1069 ipr, pdgTrig, statusTrig, imother, ptTrig, particle->Eta(), particle->Phi()*TMath::RadToDeg());
7b2086c3 1070
1071// if(pdgTrig==111)
1072// {
1073// printf("\t pi0 daughters %d, %d\n", particle->GetDaughter(0), particle->GetDaughter(1));
1074// }
7b2086c3 1075
1076 if (pdgTrig==22 ) fhPtPhoton->Fill(ptTrig);
1077 else if(pdgTrig==111) fhPtPi0 ->Fill(ptTrig);
1078
1079 // Check if it is leading
1080 Bool_t leading[4] ;
1081 Bool_t isolated[4] ;
1082
7d409bf9 1083 IsLeadingAndIsolated(ipr, pdgTrig, leading, isolated);
7b2086c3 1084
1085 Int_t iparton = -1;
7d409bf9 1086 Int_t ok = CorrelateWithPartonOrJet(ipr, pdgTrig, leading, isolated, iparton);
2292cf03 1087 if(!ok) continue;
1088
7d409bf9 1089 GetXE(ipr,pdgTrig,leading,isolated,iparton) ;
7b2086c3 1090
1091 }
1092
1093 if(GetDebug() > 1) printf("AliAnaGeneratorKine::MakeAnalysisFillHistograms() - End fill histograms \n");
1094
7d409bf9 1095}
1096
1097//_________________________________________________________
1098void AliAnaGeneratorKine::SetTriggerDetector(TString & det)
1099{
1100 // Set the detrimeter for the analysis
1101
1102 fTriggerDetectorString = det;
1103
1104 if (det=="EMCAL") fTriggerDetector = kEMCAL;
1105 else if(det=="PHOS" ) fTriggerDetector = kPHOS;
1106 else if(det=="CTS") fTriggerDetector = kCTS;
1107 else if(det=="DCAL") fTriggerDetector = kDCAL;
1108 else if(det.Contains("DCAL") && det.Contains("PHOS")) fTriggerDetector = kDCALPHOS;
1109 else AliFatal(Form("Detector < %s > not known!", det.Data()));
1110
1111}
1112
1113//_____________________________________________________
1114void AliAnaGeneratorKine::SetTriggerDetector(Int_t det)
1115{
1116 // Set the detrimeter for the analysis
1117
1118 fTriggerDetector = det;
1119
1120 if (det==kEMCAL) fTriggerDetectorString = "EMCAL";
1121 else if(det==kPHOS ) fTriggerDetectorString = "PHOS";
1122 else if(det==kCTS) fTriggerDetectorString = "CTS";
1123 else if(det==kDCAL) fTriggerDetectorString = "DCAL";
1124 else if(det==kDCALPHOS) fTriggerDetectorString = "DCAL_PHOS";
1125 else AliFatal(Form("Detector < %d > not known!", det));
1126
1127}
1128