]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG2/FLOW/oldcode/AliFlowAnalyser.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Moved old reacton plane code to /oldcode
[u/mrichter/AliRoot.git] / PWG2 / FLOW / oldcode / AliFlowAnalyser.cxx
CommitLineData
9777bfcb 1//////////////////////////////////////////////////////////////////////
2//
448e8856 3// $Id: AliFlowAnalyser.cxx 18618 2007-05-16 15:38:22Z snelling $
9777bfcb 4//
5// Author: Emanuele Simili
6//
7//////////////////////////////////////////////////////////////////////
8//_____________________________________________________________
9//
10// Description:
11// the AliFlowAnalyser provides the methods to perform an Event
12// plane flow analysis over AliFlowEvents.
13// - The method Init() generates the histograms.
4e566f2f 14// - The method Analyze(AliFlowEvent*) calculates/extracts flow observables,
9777bfcb 15// fills some histograms and performs the E.P. analysis.
16// - The method Resolution() calculates the resolution correction and
17// applies it to observed v_n values.
18// - The method Finish() saves the histograms and closes the file.
19// Weights calculation has been moved to the class AliFlowWeighter.
20//
21// The most important features of the Flow Analysis are:
22// - Reaction Plane calculation: for different harmonics and selections
23// with user defined cuts. (The calculation of the R.P. is a method in
24// the AliFlowEvent class).
25// - Correlation Analysis: particles are related to the R.P. and Vn
26// coefficient are extracted. The AliFlowAnalysisMaker also removes
27// auto-correlations of each particle with respect to the event plane,
28// by subtracting track by track from the Q vector before the correlation
29// study. Observable Flow is plotted.
30// - Differential Flow: for each harmonic and selection, flow values
31// (Vn) are plotted versus rapidity/pseudorapidity and pT.
32// - Resulution Estimate: event plane resolution is estimated via the
33// sub-event method (res ~ sqrt(cos(psi1-psi2))), and flow values are
34// then corrected for it. "True" Flow is plotted.
35// - Integrated Flow: for the doubly-integrated (pT & eta) Vn values
36// the error on the event plane resolution is folded in to the error.
37// - Phi and Baiesian weights calculation (from event observables).
38//
39// The AliFlowAnalysisMaker Class is adapted from the original
40// StFlowAnalysisMaker, succesfully used to study flow in the
41// STAR experiment at RICH.
42// Original Authors: Raimond Snellings & Art Poskanzer
43//
44//////////////////////////////////////////////////////////////////////
45
46#ifndef ALIFLOWANALYSER_CXX
47#define ALIFLOWANALYSER_CXX
48
49
50// ROOT things
51#include <TROOT.h>
52#include <TMath.h>
53#include <TFile.h>
54#include <TTree.h>
55#include <TChain.h>
56#include <TString.h>
57#include <TObject.h>
92016a03 58#include <TClonesArray.h>
9777bfcb 59#include <TOrdCollection.h>
60#include <TParticle.h>
61#include <TParticlePDG.h>
62#include <TDatabasePDG.h>
63#include <TH1.h>
64#include <TH2.h>
65#include <TH3.h>
66#include <TProfile.h>
67#include <TProfile2D.h>
68#include <TVector2.h>
69
70// Flow things
71#include "AliFlowEvent.h"
72#include "AliFlowTrack.h"
73#include "AliFlowV0.h"
74#include "AliFlowConstants.h"
75#include "AliFlowSelection.h"
76#include "AliFlowAnalyser.h"
77
78// ANSI things
79#include <math.h>
80#include <stdio.h>
81#include <stdlib.h>
82#include <iostream>
83using namespace std; //required for resolving the 'cout' symbol
84
da4d7f17 85ClassImp(AliFlowAnalyser)
9777bfcb 86//-----------------------------------------------------------------------
4e566f2f 87AliFlowAnalyser::AliFlowAnalyser(const AliFlowSelection* flowSelect):
7eb9830d 88 fTrackLoop(kTRUE),
89 fV0loop(kTRUE),
90 fShuffle(kFALSE),
91 fV1Ep1Ep2(kFALSE),
92 fEtaSub(kFALSE),
93 fReadPhiWgt(kFALSE),
94 fBayWgt(kFALSE),
95 fRePid(kFALSE),
96 fCustomRespFunc(kFALSE),
da5aa0a0 97 fSub(0),
7eb9830d 98 fPtWgt(kFALSE),
99 fEtaWgt(kFALSE),
100 fOnePhiWgt(kTRUE),
4e566f2f 101 fHistFile(0x0), fPhiWgtFile(0x0),
7eb9830d 102 fHistFileName("flowAnalysPlot.root"),
103 fEventNumber(0),
104 fTrackNumber(0),
105 fV0Number(0),
106 fNumberOfTracks(0),
107 fNumberOfV0s(0),
108 fPidId(0),
109 fSelParts(0),
110 fSelV0s(0),
da5aa0a0 111 fTotalNumberOfEvents(0),
112 fTotalNumberOfTracks(0),
113 fTotalNumberOfV0s(0),
7eb9830d 114 fFlowEvent(0x0),
115 fFlowTrack(0x0),
116 fFlowV0(0x0),
117 fFlowSelect(0x0),
118 fFlowTracks(0x0),
119 fFlowV0s(0x0),
120 fPhiBins(AliFlowConstants::kPhiBins),
121 fPhiMin(0.),
122 fPhiMax(2*TMath::Pi()),
123 fLabel(0),
124 fEtaMin(0.),
125 fEtaMax(0.),
126 fPtMin(0.),
127 fPtMax(0.),
128 fPtMaxPart(0),
129 fEtaBins(0),
130 fPtBins(0),
131 fPtBinsPart(0),
132 fMaxLabel(0),
133 fPhiWgtHistList(0x0),
134 fVnResHistList(0x0),
135 fHistTrigger(0),
136 fHistMult(0),
137 fHistV0Mult(0),
138 fHistOrigMult(0),
139 fHistMultOverOrig(0),
140 fHistMultEta(0),
141 fHistCent(0),
142 fHistVertexZ(0),
143 fHistVertexXY2D(0),
144 fHistEnergyZDC(0),
145 fHistPartZDC(0),
146 fHistPidMult(0),
147 fHistBayPidMult(0),
148 fHistEtaSym(0),
149 fHistEtaSymPart(0),
150 fHistEtaSymVerZ2D(0),
151 fHistEtaSymVerZ2DPart(0),
152 fHistMultPart(0),
153 fHistV0MultPart(0),
154 fHistBayPidMultPart(0),
155 fHistMultPartUnit(0),
156 fHistPtot(0),
157 fHistPt(0),
158 fHistCharge(0),
159 fHistDcaGlobal(0),
160 fHistDca(0),
161 fHistTransDca(0),
162 fHistChi2(0),
163 fHistLenght(0),
164 fHistInvMass(0),
165 fHistFitOverMax(0),
166 fHistPhiPtCon(0),
167 fHistPhiPtUnc(0),
168 fHistPtPhiPos(0),
169 fHistPtPhiNeg(0),
170 fHistAllEtaPtPhi3D(0),
171 fHistCosPhi(0),
172 fHistPidPt(0),
173 fHistPhi(0),
174 fHistPhiCons(0),
175 fHistYieldAll2D(0),
176 fHistYieldCon2D(0),
177 fHistYieldUnc2D(0),
178 fHistConsEtaPtPhi3D(0),
179 fHistGlobEtaPtPhi3D(0),
180 fHistUncEtaPtPhi3D(0),
181 fHistChi2ITS(0),
182 fHistChi2normITS(0),
183 fHistFitPtsITS(0),
184 fHistMaxPtsITS(0),
185 fHistMeanDedxPos2DITS(0),
186 fHistMeanDedxNeg2DITS(0),
187 fHistChi2TPC(0),
188 fHistChi2normTPC(0),
189 fHistFitPtsTPC(0),
190 fHistMaxPtsTPC(0),
191 fHistFitOverMaxTPC(0),
192 fHistMeanDedxPos2D(0),
193 fHistMeanDedxNeg2D(0),
194 fHistChi2TRD(0),
195 fHistChi2normTRD(0),
196 fHistFitPtsTRD(0),
197 fHistMaxPtsTRD(0),
198 fHistMeanDedxPos2DTRD(0),
199 fHistMeanDedxNeg2DTRD(0),
200 fHistChi2TOF(0),
201 fHistChi2normTOF(0),
202 fHistFitPtsTOF(0),
203 fHistMaxPtsTOF(0),
204 fHistMeanDedxPos2DTOF(0),
205 fHistMeanDedxNeg2DTOF(0),
206 fHistMeanTPCPiPlus(0),
207 fHistMeanTPCPiMinus(0),
208 fHistMeanTPCProton(0),
209 fHistMeanTPCPbar(0),
210 fHistMeanTPCKplus(0),
211 fHistMeanTPCKminus(0),
212 fHistMeanTPCDeuteron(0),
213 fHistMeanTPCAntiDeuteron(0),
214 fHistMeanTPCPositron(0),
215 fHistMeanTPCElectron(0),
216 fHistMeanTPCMuonPlus(0),
217 fHistMeanTPCMuonMinus(0),
218 fHistMeanITSPiPlus(0),
219 fHistMeanITSPiMinus(0),
220 fHistMeanITSProton(0),
221 fHistMeanITSPbar(0),
222 fHistMeanITSKplus(0),
223 fHistMeanITSKminus(0),
224 fHistMeanITSDeuteron(0),
225 fHistMeanITSAntiDeuteron(0),
226 fHistMeanITSPositron(0),
227 fHistMeanITSElectron(0),
228 fHistMeanITSMuonPlus(0),
229 fHistMeanITSMuonMinus(0),
230 fHistMeanTOFPiPlus(0),
231 fHistMeanTOFPiMinus(0),
232 fHistMeanTOFProton(0),
233 fHistMeanTOFPbar(0),
234 fHistMeanTOFKplus(0),
235 fHistMeanTOFKminus(0),
236 fHistMeanTOFDeuteron(0),
237 fHistMeanTOFAntiDeuteron(0),
238 fHistMeanTOFPositron(0),
239 fHistMeanTOFElectron(0),
240 fHistMeanTOFMuonPlus(0),
241 fHistMeanTOFMuonMinus(0),
242 fHistMeanTRDPiPlus(0),
243 fHistMeanTRDPiMinus(0),
244 fHistMeanTRDProton(0),
245 fHistMeanTRDPbar(0),
246 fHistMeanTRDKplus(0),
247 fHistMeanTRDKminus(0),
248 fHistMeanTRDDeuteron(0),
249 fHistMeanTRDAntiDeuteron(0),
250 fHistMeanTRDPositron(0),
251 fHistMeanTRDElectron(0),
252 fHistMeanTRDMuonPlus(0),
253 fHistMeanTRDMuonMinus(0),
254 fHistPidPiPlus(0),
255 fHistPidPiMinus(0),
256 fHistPidProton(0),
257 fHistPidAntiProton(0),
258 fHistPidKplus(0),
259 fHistPidKminus(0),
260 fHistPidDeuteron(0),
261 fHistPidAntiDeuteron(0),
262 fHistPidElectron(0),
263 fHistPidPositron(0),
264 fHistPidMuonMinus(0),
265 fHistPidMuonPlus(0),
266 fHistPidPiPlusPart(0),
267 fHistPidPiMinusPart(0),
268 fHistPidProtonPart(0),
269 fHistPidAntiProtonPart(0),
270 fHistPidKplusPart(0),
271 fHistPidKminusPart(0),
272 fHistPidDeuteronPart(0),
273 fHistPidAntiDeuteronPart(0),
274 fHistPidElectronPart(0),
275 fHistPidPositronPart(0),
276 fHistPidMuonMinusPart(0),
277 fHistPidMuonPlusPart(0),
278 fHistBinEta(0),
279 fHistBinPt(0),
280 fHistEtaPtPhi3DPart(0),
281 fHistYieldPart2D(0),
282 fHistDcaGlobalPart(0),
283 fHistInvMassPart(0),
284 fHistEtaPtPhi3DOut(0),
285 fHistYieldOut2D(0),
286 fHistDcaGlobalOut(0),
287 fHistInvMassOut(0),
288 fHistMeanDedxPos3DPart(0),
289 fHistMeanDedxNeg3DPart(0),
290 fHistMeanDedxPos3DPartITS(0),
291 fHistMeanDedxNeg3DPartITS(0),
292 fHistV0Mass(0),
293 fHistV0EtaPtPhi3D(0),
294 fHistV0YieldAll2D(0),
295 fHistV0PYall2D(0),
296 fHistV0Dca(0),
297 fHistV0Chi2(0),
298 fHistV0Lenght(0),
299 fHistV0Sigma(0),
300 fHistV0CosPhi(0),
301 fHistV0MassPtSlices(0),
302 fHistV0BinEta(0),
303 fHistV0BinPt(0),
304 fHistV0sbBinEta(0),
305 fHistV0sbBinPt(0),
306 fHistV0MassWin(0),
307 fHistV0EtaPtPhi3DPart(0),
308 fHistV0YieldPart2D(0),
309 fHistV0DcaPart(0),
310 fHistV0LenghtPart(0),
311 fHistV0sbMassSide(0),
312 fHistV0sbEtaPtPhi3DPart(0),
313 fHistV0sbYieldPart2D(0),
314 fHistV0sbDcaPart(0),
315 fHistV0sbLenghtPart(0)
9777bfcb 316{
317 // default constructor (selection given or default selection)
318
319 if(flowSelect) { fFlowSelect = new AliFlowSelection(*flowSelect) ; }
320 else { fFlowSelect = new AliFlowSelection() ; }
9777bfcb 321}
322//-----------------------------------------------------------------------
323AliFlowAnalyser::~AliFlowAnalyser()
324{
325 // default distructor (no actions)
326}
327//-----------------------------------------------------------------------
328Bool_t AliFlowAnalyser::Init()
329{
330 // sets some defaults for the analysis
331
332 cout << "* FlowAnalysis * - Init()" << endl ; cout << endl ;
333
334 // Open output files (->plots)
335 fHistFile = new TFile(fHistFileName.Data(), "RECREATE");
336 fHistFile->cd() ;
337
338 // counters and pointers to 0
339 fEventNumber = 0 ;
340 fTrackNumber = 0 ;
341 fV0Number = 0 ;
342 fPidId = -1 ;
343 //fNumberOfEvents = 0 ;
344 fNumberOfTracks = 0 ;
345 fNumberOfV0s = 0 ;
346 fFlowEvent = 0 ;
347 fFlowTrack = 0 ;
348 fFlowV0 = 0 ;
349 fFlowTracks = 0 ;
350 fFlowV0s = 0 ;
92016a03 351 fSelParts = 0 ;
352 fSelV0s = 0 ;
da5aa0a0 353 for(Int_t ii=0;ii<3;ii++) { fVertex[ii] = 0 ; }
354
355 fTotalNumberOfEvents = 0 ;
356 fTotalNumberOfTracks = 0 ;
357 fTotalNumberOfV0s = 0 ;
9777bfcb 358
359 // Check for Reading weights: if weight file is not there, wgt arrays are not used (phiwgt & bayesian)
360 if(!fPhiWgtFile)
361 {
362 fReadPhiWgt = kFALSE ;
363 fBayWgt = kFALSE ;
364 cout << " No wgt file. Phi & Bayesian weights will not be used ! " << endl ;
365 }
366
367 // Histogram settings
327288af 368 fLabel = "Pseudorapidity" ; if(strlen(fFlowSelect->PidPart()) != 0) { fLabel = "Rapidity"; }
369 fPtMaxPart = AliFlowConstants::fgPtMaxPart ; if(fFlowSelect->PtMaxPart()) { fPtMaxPart = fFlowSelect->PtMaxPart() ; }
370 fPtBinsPart = AliFlowConstants::kPtBinsPart ; if(fFlowSelect->PtBinsPart()) { fPtBinsPart = fFlowSelect->PtBinsPart() ; }
371 fPtMin = AliFlowConstants::fgPtMin ;
372 fPtMax = AliFlowConstants::fgPtMax ;
373 fEtaMin = AliFlowConstants::fgEtaMin ;
374 fEtaMax = AliFlowConstants::fgEtaMax ;
375 fPtBins = AliFlowConstants::kPtBins ;
376 fEtaBins = AliFlowConstants::kEtaBins ;
9777bfcb 377 // -
327288af 378 SetPtRangevEta(1.,0.); // (AliFlowConstants::fgPtMin , AliFlowConstants::fgPtMax) ;
379 SetEtaRangevPt(1.,0.); // (AliFlowConstants::fgEtaMin , AliFlowConstants::fgEtaMax) ;
9777bfcb 380 // -
327288af 381 float triggerMin = -1.5;
382 float triggerMax = 10.5;
383 float chargeMin = -2.5;
384 float chargeMax = 2.5;
385 float dcaMin = 0.0;
386 float dcaMax = 10.0; // 0.5;
387 float glDcaMax = 50.0;
388 float chi2Min = 0.;
389 float chi2Max = 500.;
390 float chi2normMax = 50.;
391 float chi2MaxC = 30.;
392 float chi2MaxI = 60.;
393 float fitPtsMin = -0.5;
394 float fitPtsMaxTpc = 160.5;
395 float fitPtsMaxIts = 6.5;
396 float fitPtsMaxTrd = 130.5;
397 float fitPtsMaxTof = 1.5;
398 float fitOverMaxMin = 0.;
399 float fitOverMaxMax = 1.1;
400 float multOverOrigMin = 0.;
401 float multOverOrigMax = 1.;
402 float vertexZMin = -10.;
403 float vertexZMax = 10.;
404 float vertexXYMin = -0.1;
405 float vertexXYMax = 0.1;
406 float etaSymMinPart = -1.;
407 float etaSymMaxPart = 1.;
408 float etaSymMin = -2.;
409 float etaSymMax = 2.;
410 float psiMin = 0.;
411 float psiMax = 2*TMath::Pi() ;
412 float multMin = 0. ;
413 float multMax =25000. ;
414 float multV0 = 5000. ;
415 float qMin = 0. ;
416 float qMax = 20. ;
417 float pidMin = 0. ;
418 float pidMax = 1. ;
419 float centMin = -0.5 ;
420 float centMax = 9.5 ;
421 float logpMin = -2.5 ;
422 float logpMax = 2.5 ;
423 float pMin = 0. ;
da5aa0a0 424 float pMax = 10. ;
327288af 425 float dEdxMax = 1000. ;
426 float dEdxMaxTPC = 1500. ;
427 float tofmin = 10000. ;
428 float tofmax = 50000. ;
429 float trdmax = 2000. ;
430 float lgMin = -1000. ;
431 float lgMax = 1000. ;
432 float lgMinV0 = 0. ;
433 float lgMaxV0 = 50. ;
434 float massMin = -0.1 ;
435 float massMax = 2.1 ;
436 float zdcpartMin = -0.5 ;
437 float zdcpartMax = 1000.5 ;
438 float zdceMin = 0. ;
439 float zdceMax = 10000. ;
9777bfcb 440 // -
327288af 441 enum { kTriggerBins = 12,
442 kChargeBins = 5,
443 kDcaBins = 1000,
444 kChi2Bins = 150,
445 kFitPtsBinsTpc = 161,
446 kFitPtsBinsIts = 7,
447 kFitPtsBinsTrd = 131,
448 kFitPtsBinsTof = 2,
449 kFitOverMaxBins = 55,
450 kMultOverOrigBins =100,
451 kVertexZBins = 200,
452 kVertexXYBins = 1000,
453 kEtaSymBins = 200,
454 kPhi3DBins = 60,
455 kPsiBins = 36,
456 kMultBins = 250,
457 kPidBins = 100,
458 kCentBins = 10,
459 kDedxBins = 1000,
460 kTofBins = 1000,
461 kTrdBins = 100,
da5aa0a0 462 kMomenBins = 100,
327288af 463 kQbins = 50,
464 kLgBins = 500,
465 kMassBins = 2200,
466 kZDCpartBins = 1001,
467 kZDCeBins = 200
9777bfcb 468 } ;
469
470 // Histograms Booking ...
471 // Trigger
327288af 472 fHistTrigger = new TH1F("Flow_Trigger", "Flow_Trigger", kTriggerBins, triggerMin, triggerMax);
9777bfcb 473 fHistTrigger->SetXTitle("Trigger: 1 minbias, 2 central, 3 laser, 10 other");
474 fHistTrigger->SetYTitle("Counts");
475
476 // Charge
327288af 477 fHistCharge = new TH1F("Flow_Charge", "Flow_Charge", kChargeBins, chargeMin, chargeMax);
9777bfcb 478 fHistCharge->SetXTitle("Charge");
479 fHistCharge->SetYTitle("Counts");
480
481 // Reconstructed Momentum (constrained, if so)
482 fHistPtot = new TH1F("Flow_totalP","Flow_totalP", fPtBins, fPtMin, fPtMax);
483 fHistPtot->Sumw2();
484 fHistPtot->SetXTitle("P (GeV/c)");
485 fHistPtot->SetYTitle("Counts");
486
487 // Reconstructed pT (constrained, if so)
327288af 488 fHistPt = new TH1F("Flow_Pt","Flow_Pt", kMomenBins, pMin, pMax);
9777bfcb 489 fHistPt->Sumw2();
490 fHistPt->SetXTitle("Pt (GeV/c)");
491 fHistPt->SetYTitle("Counts");
492
493 // Reconstructed P.id vs pT
494 fHistPidPt = new TH2F("Flow_PidPt","Flow_PidPt", 12, 0., 12., fPtBins, fPtMin, fPtMax);
495 fHistPidPt->Sumw2();
496 fHistPidPt->SetXTitle("e-, e+, mu-, mu+, pi-, pi+, K-, K+, pr-, pr+, d-, d+");
497 fHistPidPt->SetYTitle("Pt (GeV/c)");
498
499 // Distance of closest approach - Transverse, Unsigned
327288af 500 fHistDca = new TH1F("Flow_TransDCA", "Flow_TransverseDCA", kDcaBins, dcaMin, dcaMax);
9777bfcb 501 fHistDca->Sumw2();
502 fHistDca->SetXTitle("|Track's Signed dca (cm)|");
503 fHistDca->SetYTitle("Counts");
504
505 // Distance of closest approach - Transverse
327288af 506 fHistTransDca = new TH1F("Flow_TransDCA_Signed", "Flow_TransverseDCA_Signed", kDcaBins, -dcaMax, dcaMax);
9777bfcb 507 fHistTransDca->Sumw2();
508 fHistTransDca->SetXTitle("Track's Transverse dca (cm)");
509 fHistTransDca->SetYTitle("Counts");
510
511 // Distance of closest approach - 3d
327288af 512 fHistDcaGlobal = new TH1F("Flow_3dDcaGlobal", "Flow_3dDcaGlobal", kDcaBins, dcaMin, glDcaMax);
9777bfcb 513 fHistDcaGlobal->Sumw2();
514 fHistDcaGlobal->SetXTitle("|3d Global Track's dca to Vertex (cm)|");
515 fHistDcaGlobal->SetYTitle("Counts");
516
517 // Distance of closest approach for particles correlated with the event plane - 3d
327288af 518 fHistDcaGlobalPart = new TH1F("Flow_3dDcaGlobalPart", "Flow_3dDcaGlobalPart", kDcaBins, dcaMin, glDcaMax);
9777bfcb 519 fHistDcaGlobalPart->Sumw2();
520 fHistDcaGlobalPart->SetXTitle("|3d Global Track's dca to Vertex (cm)|");
521 fHistDcaGlobalPart->SetYTitle("Counts");
522
523 // Distance of closest approach for particles NOT SELECTED for correlation with the event plane - 3d
327288af 524 fHistDcaGlobalOut = new TH1F("Flow_3dDcaGlobalOut", "Flow_3dDcaGlobalOut", kDcaBins, dcaMin, glDcaMax);
9777bfcb 525 fHistDcaGlobalOut->Sumw2();
526 fHistDcaGlobalOut->SetXTitle("|3d Global Track's dca to Vertex (cm)|");
527 fHistDcaGlobalOut->SetYTitle("Counts");
528
529 // Yield Pt-Phi for constrained tracks
530 fHistPhiPtCon = new TH2D("Flow_PtPhi_Con", "Flow_PtPhi_Con", fPhiBins, fPhiMin, fPhiMax, fPtBins, fPtMin, fPtMax);
531 fHistPhiPtCon->Sumw2();
532 fHistPhiPtCon->SetXTitle("Phi");
533 fHistPhiPtCon->SetYTitle("Pt (GeV/c)");
534
535 // Yield Pt-Phi for UNconstrained tracks
536 fHistPhiPtUnc = new TH2D("Flow_PtPhi_Unc", "Flow_PtPhi_Unc", fPhiBins, fPhiMin, fPhiMax, fPtBins, fPtMin, fPtMax);
537 fHistPhiPtUnc->Sumw2();
538 fHistPhiPtUnc->SetXTitle("Phi");
539 fHistPhiPtUnc->SetYTitle("Pt (GeV/c)");
540
541 // Chi2
542 // at the main vertex
327288af 543 fHistChi2 = new TH1F("Flow_Chi2", "Flow_Chi2", kChi2Bins, chi2Min, chi2MaxC);
9777bfcb 544 fHistChi2->SetXTitle("Chi square at Main Vertex");
545 fHistChi2->SetYTitle("Counts");
546 // TPC
327288af 547 fHistChi2TPC = new TH1F("Flow_Chi2_TPC", "Flow_Chi2_TPC", kChi2Bins, chi2Min, chi2Max);
9777bfcb 548 fHistChi2TPC->SetXTitle("Chi square for TPC");
549 fHistChi2TPC->SetYTitle("Counts");
550 // ITS
327288af 551 fHistChi2ITS = new TH1F("Flow_Chi2_ITS", "Flow_Chi2_ITS", kChi2Bins, chi2Min, chi2MaxI);
9777bfcb 552 fHistChi2ITS->SetXTitle("Chi square for ITS");
553 fHistChi2ITS->SetYTitle("Counts");
554 // TRD
327288af 555 fHistChi2TRD = new TH1F("Flow_Chi2_TRD", "Flow_Chi2_TRD", kChi2Bins, chi2Min, chi2Max);
9777bfcb 556 fHistChi2TRD->SetXTitle("Chi square for TRD");
557 fHistChi2TRD->SetYTitle("Counts");
558 // TOF
327288af 559 fHistChi2TOF = new TH1F("Flow_Chi2_TOF", "Flow_Chi2_TOF", kChi2Bins, chi2Min, chi2Max);
9777bfcb 560 fHistChi2TOF->SetXTitle("Chi square for TOF");
561 fHistChi2TOF->SetYTitle("Counts");
562
563 // Chi2 normalized
564 // TPC
327288af 565 fHistChi2normTPC = new TH1F("Flow_Chi2norm_TPC", "Flow_Chi2norm_TPC", kChi2Bins, chi2Min, chi2normMax);
9777bfcb 566 fHistChi2normTPC->SetXTitle("Normalized #Chi^{2} for TPC");
567 fHistChi2normTPC->SetYTitle("Counts");
568 // ITS
327288af 569 fHistChi2normITS = new TH1F("Flow_Chi2norm_ITS", "Flow_Chi2norm_ITS", kChi2Bins, chi2Min, chi2normMax);
9777bfcb 570 fHistChi2normITS->SetXTitle("Normalized #Chi^{2} for ITS");
571 fHistChi2normITS->SetYTitle("Counts");
572 // TRD
327288af 573 fHistChi2normTRD = new TH1F("Flow_Chi2norm_TRD", "Flow_Chi2norm_TRD", kChi2Bins, chi2Min, chi2normMax);
9777bfcb 574 fHistChi2normTRD->SetXTitle("Normalized #Chi^{2} for TRD");
575 fHistChi2normTRD->SetYTitle("Counts");
576 // TOF
327288af 577 fHistChi2normTOF = new TH1F("Flow_Chi2norm_TOF", "Flow_Chi2norm_TOF", kChi2Bins, chi2Min, chi2normMax);
9777bfcb 578 fHistChi2normTOF->SetXTitle("Normalized #Chi^{2} for TOF");
579 fHistChi2normTOF->SetYTitle("Counts");
580
581 // FitPts
582 // TPC
327288af 583 fHistFitPtsTPC = new TH1F("Flow_FitPts_TPC", "Flow_FitPts_TPC", kFitPtsBinsTpc, fitPtsMin, fitPtsMaxTpc);
9777bfcb 584 fHistFitPtsTPC->SetXTitle("Fit Points");
585 fHistFitPtsTPC->SetYTitle("Counts");
586 // ITS
327288af 587 fHistFitPtsITS = new TH1F("Flow_HitPts_ITS", "Flow_HitPts_ITS", kFitPtsBinsIts, fitPtsMin, fitPtsMaxIts);
9777bfcb 588 fHistFitPtsITS->SetXTitle("Fit Points");
589 fHistFitPtsITS->SetYTitle("Counts");
590 // TRD
327288af 591 fHistFitPtsTRD = new TH1F("Flow_FitPts_TRD", "Flow_FitPts_TRD", kFitPtsBinsTrd, fitPtsMin, fitPtsMaxTrd);
9777bfcb 592 fHistFitPtsTRD->SetXTitle("Fit Points");
593 fHistFitPtsTRD->SetYTitle("Counts");
594 // TOF
327288af 595 fHistFitPtsTOF = new TH1F("Flow_HitPts_TOF", "Flow_HitPts_TOF", kFitPtsBinsTof, fitPtsMin, fitPtsMaxTof);
9777bfcb 596 fHistFitPtsTOF->SetXTitle("Fit Points");
597 fHistFitPtsTOF->SetYTitle("Counts");
598
599 // MaxPts
600 // TPC
327288af 601 fHistMaxPtsTPC = new TH1F("Flow_MaxPts_TPC", "Flow_MaxPts_TPC", kFitPtsBinsTpc, fitPtsMin, fitPtsMaxTpc);
9777bfcb 602 fHistMaxPtsTPC->SetXTitle("Max Points");
603 fHistMaxPtsTPC->SetYTitle("Counts");
604 // ITS
327288af 605 fHistMaxPtsITS = new TH1F("Flow_MaxPts_ITS", "Flow_MaxPts_ITS", kFitPtsBinsIts, fitPtsMin, fitPtsMaxIts);
9777bfcb 606 fHistMaxPtsITS->SetXTitle("Max Points");
607 fHistMaxPtsITS->SetYTitle("Counts");
608 // TRD
327288af 609 fHistMaxPtsTRD = new TH1F("Flow_MaxPts_TRD", "Flow_MaxPts_TRD", kFitPtsBinsTrd, fitPtsMin, fitPtsMaxTrd);
9777bfcb 610 fHistMaxPtsTRD->SetXTitle("Max Points");
611 fHistMaxPtsTRD->SetYTitle("Counts");
612 // TOF
327288af 613 fHistMaxPtsTOF = new TH1F("Flow_MaxPts_TOF", "Flow_MaxPts_TOF", kFitPtsBinsTof, fitPtsMin, fitPtsMaxTof);
9777bfcb 614 fHistMaxPtsTOF->SetXTitle("Max Points");
615 fHistMaxPtsTOF->SetYTitle("Counts");
616
617 // FitOverMax
618 // Tpc
327288af 619 fHistFitOverMaxTPC = new TH1F("Flow_FitOverMax_TPC", "Flow_FitOverMax_TPC", kFitOverMaxBins, fitOverMaxMin, fitOverMaxMax);
9777bfcb 620 fHistFitOverMaxTPC->SetXTitle("(Fit Points - 1) / Max Points");
621 fHistFitOverMaxTPC->SetYTitle("Counts");
622 // All
327288af 623 fHistFitOverMax = new TH1F("Flow_FitOverMax", "Flow_FitOverMax", kFitOverMaxBins, fitOverMaxMin, fitOverMaxMax);
9777bfcb 624 fHistFitOverMax->SetXTitle("(Fit Points - 1) / Max Points");
625 fHistFitOverMax->SetYTitle("Counts");
626
627 // lenght
327288af 628 fHistLenght = new TH1F("Flow_TrackLenght", "Flow_TrackLenght", kLgBins, lgMin, lgMax);
9777bfcb 629 fHistLenght->SetXTitle("Lenght of the Track (cm)");
630 fHistLenght->SetYTitle("Counts");
631
632 // OrigMult
327288af 633 fHistOrigMult = new TH1F("Flow_OrigMult", "Flow_OrigMult", kMultBins, multMin, multMax);
9777bfcb 634 fHistOrigMult->SetXTitle("Original Mult");
635 fHistOrigMult->SetYTitle("Counts");
636
637 // MultEta
327288af 638 fHistMultEta = new TH1F("Flow_MultEta", "Flow_MultEta", kMultBins, multMin, multMax);
9777bfcb 639 fHistMultEta->SetXTitle("Mult for Centrality");
640 fHistMultEta->SetYTitle("Counts");
641
642 // Mult
327288af 643 fHistMult = new TH1F("Flow_Mult", "Flow_Mult", kMultBins, multMin, multMax);
9777bfcb 644 fHistMult->SetXTitle("Mult");
645 fHistMult->SetYTitle("Counts");
646
647 // V0s multiplicity
327288af 648 fHistV0Mult = new TH1F("FlowV0_Mult","FlowV0_Mult", kMultBins, multMin, multV0);
9777bfcb 649 fHistV0Mult->SetXTitle("V0s Multiplicity");
650 fHistV0Mult->SetYTitle("Counts");
651
652 // MultOverOrig
327288af 653 fHistMultOverOrig = new TH1F("Flow_MultOverOrig", "Flow_MultOverOrig", kMultOverOrigBins, multOverOrigMin, multOverOrigMax);
9777bfcb 654 fHistMultOverOrig->SetXTitle("Mult / Orig. Mult");
655 fHistMultOverOrig->SetYTitle("Counts");
656
657 // Mult correlated with the event planes
327288af 658 fHistMultPart = new TH1F("Flow_MultPart", "Flow_MultPart", 2*kMultBins, multMin, multMax);
9777bfcb 659 fHistMultPart->SetXTitle("Mult of Correlated Particles");
660 fHistMultPart->SetYTitle("Counts");
661
662 // Mult correlated with the event planes in 1 unit rapidity
327288af 663 fHistMultPartUnit = new TH1F("Flow_MultPartUnit", "Flow_MultPartUnit", 2*kMultBins, multMin, multMax/2);
9777bfcb 664 fHistMultPartUnit->SetXTitle("Mult of Correlated Particles (-0.5 < eta < 0.5)");
665 fHistMultPartUnit->SetYTitle("Counts");
666
667 // Mult of V0s correlated with the event planes
327288af 668 fHistV0MultPart = new TH1F("FlowV0_MultPart", "FlowV0_MultPart", kMultBins, multMin, multV0);
9777bfcb 669 fHistV0MultPart->SetXTitle("Mult of Correlated V0s");
670 fHistV0MultPart->SetYTitle("Counts");
671
672 // VertexZ
327288af 673 fHistVertexZ = new TH1F("Flow_VertexZ", "Flow_VertexZ", kVertexZBins, vertexZMin, vertexZMax);
9777bfcb 674 fHistVertexZ->SetXTitle("Vertex Z (cm)");
675 fHistVertexZ->SetYTitle("Counts");
676
677 // VertexXY
327288af 678 fHistVertexXY2D = new TH2F("Flow_VertexXY2D", "Flow_VertexXY2D", kVertexXYBins, vertexXYMin, vertexXYMax, kVertexXYBins, vertexXYMin, vertexXYMax);
9777bfcb 679 fHistVertexXY2D->SetXTitle("Vertex X (cm)");
680 fHistVertexXY2D->SetYTitle("Vertex Y (cm)");
681
682 // EtaSym vs. Vertex Z Tpc
327288af 683 fHistEtaSymVerZ2D = new TH2F("Flow_EtaSymVerZ2D", "Flow_EtaSymVerZ2D", kVertexZBins, vertexZMin, vertexZMax, kEtaSymBins, etaSymMin, etaSymMax);
9777bfcb 684 fHistEtaSymVerZ2D->SetXTitle("Vertex Z (cm)");
685 fHistEtaSymVerZ2D->SetYTitle("Eta Symmetry TPC");
686
687 // EtaSym Tpc
327288af 688 fHistEtaSym = new TH1F("Flow_EtaSym_TPC", "Flow_EtaSym_TPC", kEtaSymBins, etaSymMin, etaSymMax);
9777bfcb 689 fHistEtaSym->SetXTitle("Eta Symmetry Ratio TPC");
690 fHistEtaSym->SetYTitle("Counts");
691
692 // EtaSym vs. Vertex Z Tpc - correlation analysis
327288af 693 fHistEtaSymVerZ2DPart = new TH2F("Flow_EtaSymVerZ2D_part", "Flow_EtaSymVerZ2D_part", kVertexZBins, vertexZMin, vertexZMax, kEtaSymBins, etaSymMinPart, etaSymMaxPart);
9777bfcb 694 fHistEtaSymVerZ2DPart->SetXTitle("Vertex Z (cm)");
695 fHistEtaSymVerZ2DPart->SetYTitle("Eta Symmetry TPC");
696
697 // EtaSym Tpc - correlation analysis
327288af 698 fHistEtaSymPart = new TH1F("Flow_EtaSym_TPC_part", "Flow_EtaSym_TPC_part", kEtaSymBins, etaSymMinPart, etaSymMaxPart);
9777bfcb 699 fHistEtaSymPart->SetXTitle("Eta Symmetry Ratio TPC");
700 fHistEtaSymPart->SetYTitle("Counts");
701
702 // phi , whatever
703 fHistPhi = new TH1F("Flow_xPhi", "Flow_xPhi (all)", fPhiBins, fPhiMin, fPhiMax);
704 fHistPhi->SetXTitle("Phi (rad)");
705 fHistPhi->SetYTitle("Counts");
706
707 // phi constrained
708 fHistPhiCons = new TH1F("Flow_cPhi", "Flow_cPhi", fPhiBins, fPhiMin, fPhiMax);
709 fHistPhiCons->SetXTitle("cPhi (rad)");
710 fHistPhiCons->SetYTitle("Counts");
711
712 // EtaPtPhi , whatever
327288af 713 fHistAllEtaPtPhi3D = new TH3F("Flow_EtaPtPhi3Dall", "Flow_EtaPtPhi3Dall (whatever)", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax, kPhi3DBins, fPhiMin, fPhiMax);
9777bfcb 714 fHistAllEtaPtPhi3D->SetXTitle("Eta");
715 fHistAllEtaPtPhi3D->SetYTitle("Pt (GeV/c)");
716 fHistAllEtaPtPhi3D->SetZTitle("Phi (rad)");
717
718 // Constrained EtaPtPhi
327288af 719 fHistConsEtaPtPhi3D = new TH3F("Flow_consEtaPtPhi3D", "Flow_consEtaPtPhi3D (constrainable)", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax, kPhi3DBins, fPhiMin, fPhiMax);
9777bfcb 720 fHistConsEtaPtPhi3D->SetXTitle("cEta");
721 fHistConsEtaPtPhi3D->SetYTitle("cPt (GeV/c)");
722 fHistConsEtaPtPhi3D->SetZTitle("cPhi (rad)");
723
724 // Global EtaPtPhi
327288af 725 fHistGlobEtaPtPhi3D = new TH3F("Flow_globEtaPtPhi3D", "Flow_globEtaPtPhi3D (constrainable)", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax, kPhi3DBins, fPhiMin, fPhiMax);
9777bfcb 726 fHistGlobEtaPtPhi3D->SetXTitle("gEta");
727 fHistGlobEtaPtPhi3D->SetYTitle("gPt (GeV/c)");
728 fHistGlobEtaPtPhi3D->SetZTitle("gPhi (rad)");
729
730 // UnConstrained EtaPtPhi
327288af 731 fHistUncEtaPtPhi3D = new TH3F("Flow_uncEtaPtPhi3D", "Flow_uncEtaPtPhi3D (un-constrainable)", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax, kPhi3DBins, fPhiMin, fPhiMax);
9777bfcb 732 fHistUncEtaPtPhi3D->SetXTitle("gEta");
733 fHistUncEtaPtPhi3D->SetYTitle("gPt (GeV/c)");
734 fHistUncEtaPtPhi3D->SetZTitle("gPhi (rad)");
735
736 // EtaPtPhi for particles correlated with the event plane
327288af 737 fHistEtaPtPhi3DPart = new TH3F("Flow_EtaPtPhi3Dpart", "Flow_EtaPtPhi3Dpart (selected part)", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax, kPhi3DBins, fPhiMin, fPhiMax);
9777bfcb 738 fHistEtaPtPhi3DPart->SetXTitle("Eta");
739 fHistEtaPtPhi3DPart->SetYTitle("Pt (GeV/c)");
740 fHistEtaPtPhi3DPart->SetZTitle("Phi (rad)");
741
742 // EtaPtPhi for particles NOT SELECTED for correlation with the event plane
327288af 743 fHistEtaPtPhi3DOut = new TH3F("Flow_EtaPtPhi3Dout", "Flow_EtaPtPhi3Dout (NOT selected part)", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax, kPhi3DBins, fPhiMin, fPhiMax);
9777bfcb 744 fHistEtaPtPhi3DOut->SetXTitle("Eta");
745 fHistEtaPtPhi3DOut->SetYTitle("Pt (GeV/c)");
746 fHistEtaPtPhi3DOut->SetZTitle("Phi (rad)");
747
748 // Yield Pt-Phi for all positive
749 fHistPtPhiPos = new TH2D("Flow_PtPhi_Plus", "Flow_PtPhi_Plus", fPhiBins, fPhiMin, fPhiMax, 50, fPtMin, fPtMax);
750 fHistPtPhiPos->Sumw2();
751 fHistPtPhiPos->SetXTitle("Phi");
752 fHistPtPhiPos->SetYTitle("Pt (GeV/c)");
753
754 // Yield Pt-Phi for all negative
755 fHistPtPhiNeg = new TH2D("Flow_PtPhi_Minus", "Flow_PtPhi_Minus", fPhiBins, fPhiMin, fPhiMax, 50, fPtMin, fPtMax);
756 fHistPtPhiNeg->Sumw2();
757 fHistPtPhiNeg->SetXTitle("Phi");
758 fHistPtPhiNeg->SetYTitle("Pt (GeV/c)");
759
760 // Yield for all particles
761 fHistYieldAll2D = new TH2D("Flow_YieldAll2D", "Flow_YieldAll2D", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax);
762 fHistYieldAll2D->Sumw2();
763 fHistYieldAll2D->SetXTitle("Pseudorapidty");
764 fHistYieldAll2D->SetYTitle("Pt (GeV/c)");
765
766 // Yield for constrainable tracks
767 fHistYieldCon2D = new TH2D("Flow_YieldCons2D", "Flow_YieldCons2D", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax);
768 fHistYieldCon2D->Sumw2();
769 fHistYieldCon2D->SetXTitle("Pseudorapidty");
770 fHistYieldCon2D->SetYTitle("Pt (GeV/c)");
771
772 // Yield for un-constrainable tracks
773 fHistYieldUnc2D = new TH2D("Flow_YieldUnc2D", "Flow_YieldUnc2D", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax);
774 fHistYieldUnc2D->Sumw2();
775 fHistYieldUnc2D->SetXTitle("Pseudorapidty");
776 fHistYieldUnc2D->SetYTitle("Pt (GeV/c)");
777
778 // Yield for particles correlated with the event plane
779 fHistYieldPart2D = new TH2D("Flow_YieldPart2D", "Flow_YieldPart2D (selected part)", fEtaBins, fEtaMin, fEtaMax, fPtBinsPart, fPtMin, fPtMaxPart);
780 fHistYieldPart2D->Sumw2();
781 fHistYieldPart2D->SetXTitle((char*)fLabel.Data());
782 fHistYieldPart2D->SetYTitle("Pt (GeV/c)");
783
784 // Yield for particles NOT SELECTED for correlation with the event plane
785 fHistYieldOut2D = new TH2D("Flow_YieldOut2D", "Flow_YieldOut2D (NOT selected part)", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax);
786 fHistYieldOut2D->Sumw2();
787 fHistYieldOut2D->SetXTitle("Pseudorapidty");
788 fHistYieldOut2D->SetYTitle("Pt (GeV/c)");
789
790 // invariant Mass for all particles (from TOF)
327288af 791 fHistInvMass = new TH1F("Flow_InvMass", "Flow_InvMass (tof)", kMassBins, massMin, massMax);
9777bfcb 792 fHistInvMass->SetXTitle("Invariant Mass (GeV)");
793 fHistInvMass->SetYTitle("Counts");
794
795 // invariant Mass for particles correlated with the event plane (from TOF)
327288af 796 fHistInvMassPart = new TH1F("Flow_InvMassPart", "Flow_InvMassPart (tof)", kMassBins, massMin, massMax);
9777bfcb 797 fHistInvMassPart->SetXTitle("Invariant Mass (GeV)");
798 fHistInvMassPart->SetYTitle("Counts");
799
800 // invariant Mass for particles NOT SELECTED for correlation with the event plane (from TOF)
327288af 801 fHistInvMassOut = new TH1F("Flow_InvMassOut", "Flow_InvMassOut (tof)", kMassBins, massMin, massMax);
9777bfcb 802 fHistInvMassOut->SetXTitle("Invariant Mass (GeV)");
803 fHistInvMassOut->SetYTitle("Counts");
804
805 // Mean Eta in each bin for particles correlated with the event plane
806 fHistBinEta = new TProfile("Flow_Bin_Eta", "Flow_Bin_Eta_part (selected part)", fEtaBins, fEtaMin, fEtaMax, fEtaMin, fEtaMax, "");
807 fHistBinEta->SetXTitle((char*)fLabel.Data());
808 fHistBinEta->SetYTitle((char*)fLabel.Data());
809
810 // Mean Pt in each bin for particles correlated with the event plane
811 fHistBinPt = new TProfile("Flow_Bin_Pt", "Flow_Bin_Pt_part (selected part)", fPtBinsPart, fPtMin, fPtMaxPart, fPtMin, fPtMaxPart, "");
812 fHistBinPt->SetXTitle("Pt (GeV/c)");
813 fHistBinPt->SetYTitle("<Pt> (GeV/c)");
814
815 // cos(n*phiLab)
327288af 816 fHistCosPhi = new TProfile("Flow_CosPhiLab", "Flow_CosPhiLab", AliFlowConstants::kHars, 0.5, (float)(AliFlowConstants::kHars) + 0.5, -100., 100., "");
9777bfcb 817 fHistCosPhi->SetXTitle("Harmonic");
818 fHistCosPhi->SetYTitle("<cos(n*PhiLab)> (%)");
819
820 // PID pi+
327288af 821 fHistPidPiPlus = new TH1F("Flow_PidPiPlus", "Flow_PidPiPlus", kPidBins, pidMin, pidMax);
9777bfcb 822 fHistPidPiPlus->SetXTitle("ALICE P.Id.");
823 fHistPidPiPlus->SetYTitle("Counts");
824
825 // PID pi-
327288af 826 fHistPidPiMinus = new TH1F("Flow_PidPiMinus", "Flow_PidPiMinus", kPidBins, pidMin, pidMax);
9777bfcb 827 fHistPidPiMinus->SetXTitle("ALICE P.Id.");
828 fHistPidPiMinus->SetYTitle("Counts");
829
830 // PID proton
327288af 831 fHistPidProton = new TH1F("Flow_PidProton", "Flow_PidProton", kPidBins, pidMin, pidMax);
9777bfcb 832 fHistPidProton->SetXTitle("ALICE P.Id.");
833 fHistPidProton->SetYTitle("Counts");
834
835 // PID anti proton
327288af 836 fHistPidAntiProton = new TH1F("Flow_PidAntiProton", "Flow_PidAntiProton", kPidBins, pidMin, pidMax);
9777bfcb 837 fHistPidAntiProton->SetXTitle("ALICE P.Id.");
838 fHistPidAntiProton->SetYTitle("Counts");
839
840 // PID Kplus
327288af 841 fHistPidKplus = new TH1F("Flow_PidKplus", "Flow_PidKplus", kPidBins, pidMin, pidMax);
9777bfcb 842 fHistPidKplus->SetXTitle("ALICE P.Id.");
843 fHistPidKplus->SetYTitle("Counts");
844
845 // PID Kminus
327288af 846 fHistPidKminus = new TH1F("Flow_PidKminus", "Flow_PidKminus", kPidBins, pidMin, pidMax);
9777bfcb 847 fHistPidKminus->SetXTitle("ALICE P.Id.");
848 fHistPidKminus->SetYTitle("Counts");
849
850 // PID deuteron
327288af 851 fHistPidDeuteron = new TH1F("Flow_PidDeuteron", "Flow_PidDeuteron", kPidBins, pidMin, pidMax);
9777bfcb 852 fHistPidDeuteron->SetXTitle("ALICE P.Id.");
853 fHistPidDeuteron->SetYTitle("Counts");
854
855 // PID anti deuteron
327288af 856 fHistPidAntiDeuteron = new TH1F("Flow_PidAntiDeuteron", "Flow_PidAntiDeuteron", kPidBins, pidMin, pidMax);
9777bfcb 857 fHistPidAntiDeuteron->SetXTitle("ALICE P.Id.");
858 fHistPidAntiDeuteron->SetYTitle("Counts");
859
860 // PID electron
327288af 861 fHistPidElectron = new TH1F("Flow_PidElectron", "Flow_PidElectron", kPidBins, pidMin, pidMax);
9777bfcb 862 fHistPidElectron->SetXTitle("ALICE P.Id.");
863 fHistPidElectron->SetYTitle("Counts");
864
865 // PID positron
327288af 866 fHistPidPositron = new TH1F("Flow_PidPositron", "Flow_PidPositron", kPidBins, pidMin, pidMax);
9777bfcb 867 fHistPidPositron->SetXTitle("ALICE P.Id.");
868 fHistPidPositron->SetYTitle("Counts");
869
870 // PID Muon+
327288af 871 fHistPidMuonPlus = new TH1F("Flow_PidMuonPlus", "Flow_PidMuonPlus", kPidBins, pidMin, pidMax);
9777bfcb 872 fHistPidMuonPlus->SetXTitle("ALICE P.Id.");
873 fHistPidMuonPlus->SetYTitle("Counts");
874
875 // PID Muon-
327288af 876 fHistPidMuonMinus = new TH1F("Flow_PidMuonMinus", "Flow_PidMuonMinus", kPidBins, pidMin, pidMax);
9777bfcb 877 fHistPidMuonMinus->SetXTitle("ALICE P.Id.");
878 fHistPidMuonMinus->SetYTitle("Counts");
879
880 // PID pi+ selected
327288af 881 fHistPidPiPlusPart = new TH1F("Flow_PidPiPlusPart", "Flow_PidPiPlusPart", kPidBins, pidMin, pidMax);
9777bfcb 882 fHistPidPiPlusPart->SetXTitle("ALICE P.Id. (pid = pi+)");
883 fHistPidPiPlusPart->SetYTitle("Counts");
884
885 // PID pi- selected
327288af 886 fHistPidPiMinusPart = new TH1F("Flow_PidPiMinusPart", "Flow_PidPiMinusPart", kPidBins, pidMin, pidMax);
9777bfcb 887 fHistPidPiMinusPart->SetXTitle("ALICE P.Id. (pid = pi-)");
888 fHistPidPiMinusPart->SetYTitle("Counts");
889
890 // PID proton selected
327288af 891 fHistPidProtonPart = new TH1F("Flow_PidProtonPart", "Flow_PidProtonPart", kPidBins, pidMin, pidMax);
9777bfcb 892 fHistPidProtonPart->SetXTitle("ALICE P.Id. (pid = p)");
893 fHistPidProtonPart->SetYTitle("Counts");
894
895 // PID anti proton selected
327288af 896 fHistPidAntiProtonPart = new TH1F("Flow_PidAntiProtonPart", "Flow_PidAntiProtonPart", kPidBins, pidMin, pidMax);
9777bfcb 897 fHistPidAntiProtonPart->SetXTitle("ALICE P.Id. (pid = p-)");
898 fHistPidAntiProtonPart->SetYTitle("Counts");
899
900 // PID Kplus selected
327288af 901 fHistPidKplusPart = new TH1F("Flow_PidKplusPart", "Flow_PidKplusPart", kPidBins, pidMin, pidMax);
9777bfcb 902 fHistPidKplusPart->SetXTitle("ALICE P.Id. (pid = K+)");
903 fHistPidKplusPart->SetYTitle("Counts");
904
905 // PID Kminus selected
327288af 906 fHistPidKminusPart = new TH1F("Flow_PidKminusPart", "Flow_PidKminusPart", kPidBins, pidMin, pidMax);
9777bfcb 907 fHistPidKminusPart->SetXTitle("ALICE P.Id. (pid = K-)");
908 fHistPidKminusPart->SetYTitle("Counts");
909
910 // PID deuteron selected
327288af 911 fHistPidDeuteronPart = new TH1F("Flow_PidDeuteronPart", "Flow_PidDeuteronPart", kPidBins, pidMin, pidMax);
9777bfcb 912 fHistPidDeuteronPart->SetXTitle("ALICE P.Id. (pid = d)");
913 fHistPidDeuteronPart->SetYTitle("Counts");
914
915 // PID anti deuteron selected
327288af 916 fHistPidAntiDeuteronPart = new TH1F("Flow_PidAntiDeuteronPart", "Flow_PidAntiDeuteronPart", kPidBins, pidMin, pidMax);
9777bfcb 917 fHistPidAntiDeuteronPart->SetXTitle("ALICE P.Id. (pid = d--)");
918 fHistPidAntiDeuteronPart->SetYTitle("Counts");
919
920 // PID electron selected
327288af 921 fHistPidElectronPart = new TH1F("Flow_PidElectronPart", "Flow_PidElectronPart", kPidBins, pidMin, pidMax);
9777bfcb 922 fHistPidElectronPart->SetXTitle("ALICE P.Id. (pid = e-)");
923 fHistPidElectronPart->SetYTitle("Counts");
924
925 // PID positron selected
327288af 926 fHistPidPositronPart = new TH1F("Flow_PidPositronPart", "Flow_PidPositronPart", kPidBins, pidMin, pidMax);
9777bfcb 927 fHistPidPositronPart->SetXTitle("ALICE P.Id. (pid = e+)");
928 fHistPidPositronPart->SetYTitle("Counts");
929
930 // PID Muon+ selected
327288af 931 fHistPidMuonPlusPart = new TH1F("Flow_PidMuonPlusPart", "Flow_PidMuonPlusPart", kPidBins, pidMin, pidMax);
9777bfcb 932 fHistPidMuonPlusPart->SetXTitle("ALICE P.Id. (pid = mu+)");
933 fHistPidMuonPlusPart->SetYTitle("Counts");
934
935 // PID Muon- selected
327288af 936 fHistPidMuonMinusPart = new TH1F("Flow_PidMuonMinusPart", "Flow_PidMuonMinusPart", kPidBins, pidMin, pidMax);
9777bfcb 937 fHistPidMuonMinusPart->SetXTitle("ALICE P.Id. (pid = mu-)");
938 fHistPidMuonMinusPart->SetYTitle("Counts");
939
940 // PID multiplicities (all)
941 fHistPidMult = new TProfile("Flow_PidMult", "Flow_PidMult", 15, 0.5, 15.5, "");
942 fHistPidMult->SetXTitle("all, h+, h-, pi+, pi-, pr+, pr-, K+, K-, d+, d-, e-, e+, mu-, mu+");
943 fHistPidMult->SetYTitle("Multiplicity");
944
945 // PID for all tracks
327288af 946 fHistBayPidMult = new TH1F("Flow_BayPidMult","Flow_BayPidMult",AliFlowConstants::kPid,-0.5,((float)AliFlowConstants::kPid-0.5));
9777bfcb 947 fHistBayPidMult->Sumw2() ;
948 fHistBayPidMult->SetXTitle("e+/- , mu+/- , pi+/- , K+/- , p+/- , d+/- ");
949 fHistBayPidMult->SetYTitle("Counts");
950
951 // PID for particles correlated with the event plane
327288af 952 fHistBayPidMultPart = new TH1F("Flow_BayPidMult_Part","Flow_BayPidMult_Part",AliFlowConstants::kPid,-0.5,((float)AliFlowConstants::kPid-0.5));
9777bfcb 953 fHistBayPidMultPart->Sumw2() ;
954 fHistBayPidMultPart->SetXTitle("e+/- , mu+/- , pi+/- , K+/- , p+/- , d+/- ");
955 fHistBayPidMultPart->SetYTitle("Counts");
956
957 // Centrality
327288af 958 fHistCent = new TH1F("Flow_Cent", "Flow_Cent", kCentBins, centMin, centMax);
9777bfcb 959 fHistCent->SetXTitle("Centrality Bin");
960 fHistCent->SetYTitle("Counts");
961
962 // Deposited Energy in ZDC
327288af 963 fHistEnergyZDC = new TH2F("Flow_ZDC_E", "Flow_ZDC_E", 3, 0., 3., kZDCeBins, zdceMin, zdceMax);
9777bfcb 964 fHistEnergyZDC->SetXTitle("neutron , proton , e.m.");
965 fHistEnergyZDC->SetYTitle("ZDC energy");
966
967 // Part. versus Energy in ZDC
327288af 968 fHistPartZDC = new TH1F("Flow_ZDC_Participants", "Flow_ZDC_Participants", kZDCpartBins, zdcpartMin, zdcpartMax);
9777bfcb 969 fHistPartZDC->SetXTitle("ZDC part");
970 fHistPartZDC->SetYTitle("Counts");
971
972 // MeanDedxPos TPC
327288af 973 fHistMeanDedxPos2D = new TH2F("Flow_MeanDedxPos2D_TPC","Flow_MeanDedxPos2D_TPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 974 fHistMeanDedxPos2D->SetXTitle("log(momentum) (GeV/c)");
975 fHistMeanDedxPos2D->SetYTitle("mean dEdx (+)");
976
977 // MeanDedxNeg TPC
327288af 978 fHistMeanDedxNeg2D = new TH2F("Flow_MeanDedxNeg2D_TPC","Flow_MeanDedxNeg2D_TPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 979 fHistMeanDedxNeg2D->SetXTitle("log(momentum) (GeV/c)");
980 fHistMeanDedxNeg2D->SetYTitle("mean dEdx (-)");
981
982 // MeanDedxPos ITS
327288af 983 fHistMeanDedxPos2DITS = new TH2F("Flow_MeanDedxPos2D_ITS","Flow_MeanDedxPos2D_ITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 984 fHistMeanDedxPos2DITS->SetXTitle("log(momentum) (GeV/c)");
985 fHistMeanDedxPos2DITS->SetYTitle("mean dEdx (+)");
986
987 // MeanDedxNeg ITS
327288af 988 fHistMeanDedxNeg2DITS = new TH2F("Flow_MeanDedxNeg2D_ITS","Flow_MeanDedxNeg2D_ITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 989 fHistMeanDedxNeg2DITS->SetXTitle("log(momentum) (GeV/c)");
990 fHistMeanDedxNeg2DITS->SetYTitle("mean dEdx (-)");
991
992 // MeanDedxPos TPC 3D selected Part
327288af 993 fHistMeanDedxPos3DPart = new TH3F("Flow_MeanDedxPos3DPart_TPC","Flow_MeanDedxPos3DPart_TPC", kMomenBins, logpMin, logpMax, kDedxBins, 0., dEdxMaxTPC, AliFlowConstants::kPid, 0., AliFlowConstants::kPid);
9777bfcb 994 fHistMeanDedxPos3DPart->SetXTitle("log(momentum) (GeV/c)");
995 fHistMeanDedxPos3DPart->SetYTitle("mean dEdx (+)");
996 fHistMeanDedxPos3DPart->SetZTitle("e , mu , pi , k , p , d");
997
998 // MeanDedxNeg TPC 3D selected Part
327288af 999 fHistMeanDedxNeg3DPart = new TH3F("Flow_MeanDedxNeg3DPart_TPC","Flow_MeanDedxNeg3DPart_TPC", kMomenBins, logpMin, logpMax, kDedxBins, 0., dEdxMaxTPC, AliFlowConstants::kPid, 0., AliFlowConstants::kPid);
9777bfcb 1000 fHistMeanDedxNeg3DPart->SetXTitle("log(momentum) (GeV/c)");
1001 fHistMeanDedxNeg3DPart->SetYTitle("mean dEdx (-)");
1002 fHistMeanDedxNeg3DPart->SetZTitle("e , mu , pi , k , p , d");
1003
1004 // MeanDedxPos ITS 3D selected Part
327288af 1005 fHistMeanDedxPos3DPartITS = new TH3F("Flow_MeanDedxPos3DPart_ITS","Flow_MeanDedxPos3DPart_ITS", kMomenBins, logpMin, logpMax, kDedxBins, 0., dEdxMax, AliFlowConstants::kPid, 0., AliFlowConstants::kPid);
9777bfcb 1006 fHistMeanDedxPos3DPartITS->SetXTitle("log(momentum) (GeV/c)");
1007 fHistMeanDedxPos3DPartITS->SetYTitle("mean dEdx (+)");
1008 fHistMeanDedxPos3DPartITS->SetZTitle("e , mu , pi , k , p , d");
1009
1010 // MeanDedxNeg ITS 3D selected Part
327288af 1011 fHistMeanDedxNeg3DPartITS = new TH3F("Flow_MeanDedxNeg3DPart_ITS","Flow_MeanDedxNeg3DPart_ITS", kMomenBins, logpMin, logpMax, kDedxBins, 0., dEdxMax, AliFlowConstants::kPid, 0., AliFlowConstants::kPid);
9777bfcb 1012 fHistMeanDedxNeg3DPartITS->SetXTitle("log(momentum) (GeV/c)");
1013 fHistMeanDedxNeg3DPartITS->SetYTitle("mean dEdx (-)");
1014 fHistMeanDedxNeg3DPartITS->SetZTitle("e , mu , pi , k , p , d");
1015
1016 // MeanSignalPos TRD
327288af 1017 fHistMeanDedxPos2DTRD = new TH2F("Flow_MeanSignalPos2D_TRD","Flow_MeanSignalPos2D_TRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1018 fHistMeanDedxPos2DTRD->SetXTitle("log(momentum) (GeV/c)");
1019 fHistMeanDedxPos2DTRD->SetYTitle("mean TRD (+)");
1020
1021 // MeanSignalNeg TRD
327288af 1022 fHistMeanDedxNeg2DTRD = new TH2F("Flow_MeanSignalNeg2D_TRD","Flow_MeanSignalNeg2D_TRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1023 fHistMeanDedxNeg2DTRD->SetXTitle("log(momentum) (GeV/c)");
1024 fHistMeanDedxNeg2DTRD->SetYTitle("mean TRD (-)");
1025
1026 // MeanTimePos TOF
327288af 1027 fHistMeanDedxPos2DTOF = new TH2F("Flow_MeanTimePos2D_TOF","Flow_MeanTimePos2D_TOF", kMomenBins, logpMin, logpMax, kTofBins, tofmin, tofmax);
9777bfcb 1028 fHistMeanDedxPos2DTOF->SetXTitle("log(momentum) (GeV/c)");
1029 fHistMeanDedxPos2DTOF->SetYTitle("mean Time of Flight (+)");
1030
1031 // MeanTimeNeg TOF
327288af 1032 fHistMeanDedxNeg2DTOF = new TH2F("Flow_MeanTimeNeg2D_TOF","Flow_MeanTimeNeg2D_TOF", kMomenBins, logpMin, logpMax, kTofBins, tofmin, tofmax);
9777bfcb 1033 fHistMeanDedxNeg2DTOF->SetXTitle("log(momentum) (GeV/c)");
1034 fHistMeanDedxNeg2DTOF->SetYTitle("mean Time of Flight (-)");
1035
1036 // Detector response for each particle { ...
1037 // MeanDedx PiPlus - TPC
327288af 1038 fHistMeanTPCPiPlus = new TH2F("Flow_MeanDedxPiPlusTPC", "Flow_MeanDedxPiPlusTPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 1039 fHistMeanTPCPiPlus->SetXTitle("Log10(p) (GeV/c)");
1040 fHistMeanTPCPiPlus->SetYTitle("mean dEdx (pi+)");
1041 // MeanDedx PiMinus
327288af 1042 fHistMeanTPCPiMinus = new TH2F("Flow_MeanDedxPiMinusTPC", "Flow_MeanDedxPiMinusTPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 1043 fHistMeanTPCPiMinus->SetXTitle("Log10(p) (GeV/c)");
1044 fHistMeanTPCPiMinus->SetYTitle("mean dEdx (pi-)");
1045 // MeanDedx Proton
327288af 1046 fHistMeanTPCProton = new TH2F("Flow_MeanDedxProtonTPC", "Flow_MeanDedxProtonTPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 1047 fHistMeanTPCProton->SetXTitle("Log10(p) (GeV/c)");
1048 fHistMeanTPCProton->SetYTitle("mean dEdx (pr+)");
1049 // MeanDedx Pbar
327288af 1050 fHistMeanTPCPbar = new TH2F("Flow_MeanDedxPbarTPC", "Flow_MeanDedxPbarTPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 1051 fHistMeanTPCPbar->SetXTitle("Log10(p) (GeV/c)");
1052 fHistMeanTPCPbar->SetYTitle("mean dEdx (pr-)");
1053 // MeanDedx Kplus
327288af 1054 fHistMeanTPCKplus = new TH2F("Flow_MeanDedxKplusTPC", "Flow_MeanDedxKplusTPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 1055 fHistMeanTPCKplus->SetXTitle("Log10(p) (GeV/c)");
1056 fHistMeanTPCKplus->SetYTitle("mean dEdx (K+)");
1057 // MeanDedx Kminus
327288af 1058 fHistMeanTPCKminus = new TH2F("Flow_MeanDedxKminusTPC", "Flow_MeanDedxKminusTPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 1059 fHistMeanTPCKminus->SetXTitle("Log10(p) (GeV/c)");
1060 fHistMeanTPCKminus->SetYTitle("mean dEdx (K-)");
1061 // MeanDedx Deuteron
327288af 1062 fHistMeanTPCDeuteron = new TH2F("Flow_MeanDedxDeuteronTPC", "Flow_MeanDedxDeuteronTPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 1063 fHistMeanTPCDeuteron->SetXTitle("Log10(p) (GeV/c)");
1064 fHistMeanTPCDeuteron->SetYTitle("mean dEdx (d+)");
1065 // MeanDedx AntiDeuteron
327288af 1066 fHistMeanTPCAntiDeuteron = new TH2F("Flow_MeanDedxAntiDeuteronTPC", "Flow_MeanDedxAntiDeuteronTPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 1067 fHistMeanTPCAntiDeuteron->SetXTitle("Log10(p) (GeV/c)");
1068 fHistMeanTPCAntiDeuteron->SetYTitle("mean dEdx (d-)");
1069 // MeanDedxElectron
327288af 1070 fHistMeanTPCElectron = new TH2F("Flow_MeanDedxElectronTPC", "Flow_MeanDedxElectronTPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 1071 fHistMeanTPCElectron->SetXTitle("Log10(p) (GeV/c)");
1072 fHistMeanTPCElectron->SetYTitle("mean dEdx (e-)");
1073 // MeanDedx Positron
327288af 1074 fHistMeanTPCPositron = new TH2F("Flow_MeanDedxPositronTPC", "Flow_MeanDedxPositronTPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 1075 fHistMeanTPCPositron->SetXTitle("Log10(p) (GeV/c)");
1076 fHistMeanTPCPositron->SetYTitle("mean dEdx (e+)");
1077 // MeanDedx MuonPlus
327288af 1078 fHistMeanTPCMuonPlus = new TH2F("Flow_MeanDedxMuonPlusTPC", "Flow_MeanDedxMuonPlusTPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 1079 fHistMeanTPCMuonPlus->SetXTitle("Log10(p) (GeV/c)");
1080 fHistMeanTPCMuonPlus->SetYTitle("mean dEdx (mu+)");
1081 // MeanDedx MuonMinus
327288af 1082 fHistMeanTPCMuonMinus = new TH2F("Flow_MeanDedxMuonMinusTPC", "Flow_MeanDedxMuonMinusTPC", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMaxTPC);
9777bfcb 1083 fHistMeanTPCMuonMinus->SetXTitle("Log10(p) (GeV/c)");
1084 fHistMeanTPCMuonMinus->SetYTitle("mean dEdx (mu-)");
1085
1086 // MeanDedx PiPlus - ITS
327288af 1087 fHistMeanITSPiPlus = new TH2F("Flow_MeanDedxPiPlusITS", "Flow_MeanDedxPiPlusITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 1088 fHistMeanITSPiPlus->SetXTitle("Log10(p) (GeV/c)");
1089 fHistMeanITSPiPlus->SetYTitle("mean dEdx (pi+)");
1090 // MeanDedx PiMinus
327288af 1091 fHistMeanITSPiMinus = new TH2F("Flow_MeanDedxPiMinusITS", "Flow_MeanDedxPiMinusITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 1092 fHistMeanITSPiMinus->SetXTitle("Log10(p) (GeV/c)");
1093 fHistMeanITSPiMinus->SetYTitle("mean dEdx (pi-)");
1094 // MeanDedx Proton
327288af 1095 fHistMeanITSProton = new TH2F("Flow_MeanDedxProtonITS", "Flow_MeanDedxProtonITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 1096 fHistMeanITSProton->SetXTitle("Log10(p) (GeV/c)");
1097 fHistMeanITSProton->SetYTitle("mean dEdx (pr+)");
1098 // MeanDedx Pbar
327288af 1099 fHistMeanITSPbar = new TH2F("Flow_MeanDedxPbarITS", "Flow_MeanDedxPbarITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 1100 fHistMeanITSPbar->SetXTitle("Log10(p) (GeV/c)");
1101 fHistMeanITSPbar->SetYTitle("mean dEdx (pr-)");
1102 // MeanDedx Kplus
327288af 1103 fHistMeanITSKplus = new TH2F("Flow_MeanDedxKplusITS", "Flow_MeanDedxKplusITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 1104 fHistMeanITSKplus->SetXTitle("Log10(p) (GeV/c)");
1105 fHistMeanITSKplus->SetYTitle("mean dEdx (K+)");
1106 // MeanDedx Kminus
327288af 1107 fHistMeanITSKminus = new TH2F("Flow_MeanDedxKminusITS", "Flow_MeanDedxKminusITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 1108 fHistMeanITSKminus->SetXTitle("Log10(p) (GeV/c)");
1109 fHistMeanITSKminus->SetYTitle("mean dEdx (K-)");
1110 // MeanDedx Deuteron
327288af 1111 fHistMeanITSDeuteron = new TH2F("Flow_MeanDedxDeuteronITS", "Flow_MeanDedxDeuteronITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 1112 fHistMeanITSDeuteron->SetXTitle("Log10(p) (GeV/c)");
1113 fHistMeanITSDeuteron->SetYTitle("mean dEdx (d+)");
1114 // MeanDedx AntiDeuteron
327288af 1115 fHistMeanITSAntiDeuteron = new TH2F("Flow_MeanDedxAntiDeuteronITS", "Flow_MeanDedxAntiDeuteronITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 1116 fHistMeanITSAntiDeuteron->SetXTitle("Log10(p) (GeV/c)");
1117 fHistMeanITSAntiDeuteron->SetYTitle("mean dEdx (d-)");
1118 // MeanDedx Electron
327288af 1119 fHistMeanITSElectron = new TH2F("Flow_MeanDedxElectronITS", "Flow_MeanDedxElectronITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 1120 fHistMeanITSElectron->SetXTitle("Log10(p) (GeV/c)");
1121 fHistMeanITSElectron->SetYTitle("mean dEdx (e-)");
1122 // MeanDedx Positron
327288af 1123 fHistMeanITSPositron = new TH2F("Flow_MeanDedxPositronITS", "Flow_MeanDedxPositronITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 1124 fHistMeanITSPositron->SetXTitle("Log10(p) (GeV/c)");
1125 fHistMeanITSPositron->SetYTitle("mean dEdx (e+)");
1126 // MeanDedx MuonPlus
327288af 1127 fHistMeanITSMuonPlus = new TH2F("Flow_MeanDedxMuonPlusITS", "Flow_MeanDedxMuonPlusITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 1128 fHistMeanITSMuonPlus->SetXTitle("Log10(p) (GeV/c)");
1129 fHistMeanITSMuonPlus->SetYTitle("mean dEdx (mu+)");
1130 // MeanDedx MuonMinus
327288af 1131 fHistMeanITSMuonMinus = new TH2F("Flow_MeanDedxMuonMinusITS", "Flow_MeanDedxMuonMinusITS", kMomenBins, logpMin, logpMax, kDedxBins, 0, dEdxMax);
9777bfcb 1132 fHistMeanITSMuonMinus->SetXTitle("Log10(p) (GeV/c)");
1133 fHistMeanITSMuonMinus->SetYTitle("mean dEdx (mu-)");
1134
1135 // MeanTrd PiPlus - TRD
327288af 1136 fHistMeanTRDPiPlus = new TH2F("Flow_MeanTrdPiPlusTRD", "Flow_MeanTrdPiPlusTRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1137 fHistMeanTRDPiPlus->SetXTitle("Log10(p) (GeV/c)");
1138 fHistMeanTRDPiPlus->SetYTitle("mean t.r.[] (pi+)");
1139 // MeanTrd PiMinus
327288af 1140 fHistMeanTRDPiMinus = new TH2F("Flow_MeanTrdPiMinusTRD", "Flow_MeanTrdPiMinusTRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1141 fHistMeanTRDPiMinus->SetXTitle("Log10(p) (GeV/c)");
1142 fHistMeanTRDPiMinus->SetYTitle("mean t.r.[] (pi-)");
1143 // MeanTrd Proton
327288af 1144 fHistMeanTRDProton = new TH2F("Flow_MeanTrdProtonTRD", "Flow_MeanTrdProtonTRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1145 fHistMeanTRDProton->SetXTitle("Log10(p) (GeV/c)");
1146 fHistMeanTRDProton->SetYTitle("mean t.r.[] (pr+)");
1147 // MeanTrd Pbar
327288af 1148 fHistMeanTRDPbar = new TH2F("Flow_MeanTrdPbarTRD", "Flow_MeanTrdPbarTRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1149 fHistMeanTRDPbar->SetXTitle("Log10(p) (GeV/c)");
1150 fHistMeanTRDPbar->SetYTitle("mean t.r.[] (pr-)");
1151 // MeanTrd Kplus
327288af 1152 fHistMeanTRDKplus = new TH2F("Flow_MeanTrdKplusTRD", "Flow_MeanTrdKplusTRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1153 fHistMeanTRDKplus->SetXTitle("Log10(p) (GeV/c)");
1154 fHistMeanTRDKplus->SetYTitle("mean t.r.[] (K+)");
1155 // MeanTrd Kminus
327288af 1156 fHistMeanTRDKminus = new TH2F("Flow_MeanTrdKminusTRD", "Flow_MeanTrdKminusTRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1157 fHistMeanTRDKminus->SetXTitle("Log10(p) (GeV/c)");
1158 fHistMeanTRDKminus->SetYTitle("mean t.r.[] (K-)");
1159 // MeanTrd Deuteron
327288af 1160 fHistMeanTRDDeuteron = new TH2F("Flow_MeanTrdDeuteronTRD", "Flow_MeanTrdDeuteronTRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1161 fHistMeanTRDDeuteron->SetXTitle("Log10(p) (GeV/c)");
1162 fHistMeanTRDDeuteron->SetYTitle("mean t.r.[] (d+)");
1163 // MeanTrd AntiDeuteron
327288af 1164 fHistMeanTRDAntiDeuteron = new TH2F("Flow_MeanTrdAntiDeuteronTRD", "Flow_MeanTrdAntiDeuteronTRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1165 fHistMeanTRDAntiDeuteron->SetXTitle("Log10(p) (GeV/c)");
1166 fHistMeanTRDAntiDeuteron->SetYTitle("mean t.r.[] (d-)");
1167 // MeanTrd Electron
327288af 1168 fHistMeanTRDElectron = new TH2F("Flow_MeanTrdElectronTRD", "Flow_MeanTrdElectronTRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1169 fHistMeanTRDElectron->SetXTitle("Log10(p) (GeV/c)");
1170 fHistMeanTRDElectron->SetYTitle("mean t.r.[] (e-)");
1171 // MeanTrd Positron
327288af 1172 fHistMeanTRDPositron = new TH2F("Flow_MeanTrdPositronTRD", "Flow_MeanTrdPositronTRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1173 fHistMeanTRDPositron->SetXTitle("Log10(p) (GeV/c)");
1174 fHistMeanTRDPositron->SetYTitle("mean t.r.[] (e+)");
1175 // MeanTrd MuonPlus
327288af 1176 fHistMeanTRDMuonPlus = new TH2F("Flow_MeanTrdMuonPlusTRD", "Flow_MeanTrdMuonPlusTRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1177 fHistMeanTRDMuonPlus->SetXTitle("Log10(p) (GeV/c)");
1178 fHistMeanTRDMuonPlus->SetYTitle("mean t.r.[] (mu+)");
1179 // MeanTrd MuonMinus
327288af 1180 fHistMeanTRDMuonMinus = new TH2F("Flow_MeanTrdMuonMinusTRD", "Flow_MeanTrdMuonMinusTRD", kMomenBins, logpMin, logpMax, kTrdBins, 0, trdmax);
9777bfcb 1181 fHistMeanTRDMuonMinus->SetXTitle("Log10(p) (GeV/c)");
1182 fHistMeanTRDMuonMinus->SetYTitle("mean t.r.[] (mu-)");
1183
1184 // T.O.F. PiPlus - TOF
327288af 1185 fHistMeanTOFPiPlus = new TH2F("Flow_MeanTofPiPlusTOF", "Flow_MeanTofPiPlusTOF", kMassBins, massMin, massMax, kTofBins, tofmin, tofmax);
9777bfcb 1186 fHistMeanTOFPiPlus->SetXTitle("invariant mass (GeV)");
1187 fHistMeanTOFPiPlus->SetYTitle("mean t.o.f.[psec] (pi+)");
1188 // MeanTof PiMinus
327288af 1189 fHistMeanTOFPiMinus = new TH2F("Flow_MeanTofPiMinusTOF", "Flow_MeanTofPiMinusTOF", kMassBins, massMin, massMax, kTofBins, tofmin, tofmax);
9777bfcb 1190 fHistMeanTOFPiMinus->SetXTitle("invariant mass (GeV)");
1191 fHistMeanTOFPiMinus->SetYTitle("mean t.o.f.[psec] (pi-)");
1192 // MeanTof Proton
327288af 1193 fHistMeanTOFProton = new TH2F("Flow_MeanTofProtonTOF", "Flow_MeanTofProtonTOF", kMassBins, massMin, massMax, kTofBins, tofmin, tofmax);
9777bfcb 1194 fHistMeanTOFProton->SetXTitle("invariant mass (GeV)");
1195 fHistMeanTOFProton->SetYTitle("mean t.o.f.[psec] (pr+)");
1196 // Mean TofPbar
327288af 1197 fHistMeanTOFPbar = new TH2F("Flow_MeanTofPbarTOF", "Flow_MeanTofPbarTOF", kMassBins, massMin, massMax, kTofBins, tofmin, tofmax);
9777bfcb 1198 fHistMeanTOFPbar->SetXTitle("invariant mass (GeV)");
1199 fHistMeanTOFPbar->SetYTitle("mean t.o.f.[psec] (pr-)");
1200 // mean t.o.f.[psec]Kplus
327288af 1201 fHistMeanTOFKplus = new TH2F("Flow_MeanTofKplusTOF", "Flow_MeanTofKplusTOF", kMassBins, massMin, massMax, kTofBins, tofmin, tofmax);
9777bfcb 1202 fHistMeanTOFKplus->SetXTitle("invariant mass (GeV)");
1203 fHistMeanTOFKplus->SetYTitle("mean t.o.f.[psec] (K+)");
1204 // mean t.o.f.[psec]Kminus
327288af 1205 fHistMeanTOFKminus = new TH2F("Flow_MeanTofKminusTOF", "Flow_MeanTofKminusTOF", kMassBins, massMin, massMax, kTofBins, tofmin, tofmax);
9777bfcb 1206 fHistMeanTOFKminus->SetXTitle("invariant mass (GeV)");
1207 fHistMeanTOFKminus->SetYTitle("mean t.o.f.[psec] (K-)");
1208 // MeanTof Deuteron
327288af 1209 fHistMeanTOFDeuteron = new TH2F("Flow_MeanTofDeuteronTOF", "Flow_MeanTofDeuteronTOF", kMassBins, massMin, massMax, kTofBins, tofmin, tofmax);
9777bfcb 1210 fHistMeanTOFDeuteron->SetXTitle("invariant mass (GeV)");
1211 fHistMeanTOFDeuteron->SetYTitle("mean t.o.f.[psec] (d+)");
1212 // MeanTof AntiDeuteron
327288af 1213 fHistMeanTOFAntiDeuteron = new TH2F("Flow_MeanTofAntiDeuteronTOF", "Flow_MeanTofAntiDeuteronTOF", kMassBins, massMin, massMax, kTofBins, tofmin, tofmax);
9777bfcb 1214 fHistMeanTOFAntiDeuteron->SetXTitle("invariant mass (GeV)");
1215 fHistMeanTOFAntiDeuteron->SetYTitle("mean t.o.f.[psec] (d-)");
1216 // MeanTof Electron
327288af 1217 fHistMeanTOFElectron = new TH2F("Flow_MeanTofElectronTOF", "Flow_MeanTofElectronTOF", kMassBins, massMin, massMax, kTofBins, tofmin, tofmax);
9777bfcb 1218 fHistMeanTOFElectron->SetXTitle("invariant mass (GeV)");
1219 fHistMeanTOFElectron->SetYTitle("mean t.o.f.[psec] (e-)");
1220 // MeanTof Positron
327288af 1221 fHistMeanTOFPositron = new TH2F("Flow_MeanTofPositronTOF", "Flow_MeanTofPositronTOF", kMassBins, massMin, massMax, kTofBins, tofmin, tofmax);
9777bfcb 1222 fHistMeanTOFPositron->SetXTitle("invariant mass (GeV)");
1223 fHistMeanTOFPositron->SetYTitle("mean t.o.f.[psec] (e+)");
1224 // MeanTof MuonPlus
327288af 1225 fHistMeanTOFMuonPlus = new TH2F("Flow_MeanTofMuonPlusTOF", "Flow_MeanTofMuonPlusTOF", kMassBins, massMin, massMax, kTofBins, tofmin, tofmax);
9777bfcb 1226 fHistMeanTOFMuonPlus->SetXTitle("invariant mass (GeV)");
1227 fHistMeanTOFMuonPlus->SetYTitle("mean t.o.f.[psec] (mu+)");
1228 // MeanTof MuonMinus
327288af 1229 fHistMeanTOFMuonMinus = new TH2F("Flow_MeanTofMuonMinusTOF", "Flow_MeanTofMuonMinusTOF", kMassBins, massMin, massMax, kTofBins, tofmin, tofmax);
9777bfcb 1230 fHistMeanTOFMuonMinus->SetXTitle("invariant mass (GeV)");
1231 fHistMeanTOFMuonMinus->SetYTitle("mean t.o.f.[psec] (mu-)");
1232 // ... }
1233
1234 TString* histTitle;
327288af 1235 for (int n = 0; n < AliFlowConstants::kSubs; n++) // for sub-events
9777bfcb 1236 {
327288af 1237 for (int k = 0; k < AliFlowConstants::kSels; k++)
9777bfcb 1238 {
327288af 1239 for (int j = 0; j < AliFlowConstants::kHars; j++)
9777bfcb 1240 {
1241 float order = (float)(j + 1);
327288af 1242 int i = AliFlowConstants::kSubs * k + n ;
9777bfcb 1243
1244 // event planes
1245 histTitle = new TString("Flow_Psi_Sub");
1246 *histTitle += n+1;
1247 histTitle->Append("_Sel");
1248 *histTitle += k+1;
1249 histTitle->Append("_Har");
1250 *histTitle += j+1;
327288af 1251 fHistSub[i].fHistSubHar[j].fHistPsiSubs = new TH1F(histTitle->Data(),histTitle->Data(), kPsiBins, psiMin, (psiMax / order));
9777bfcb 1252 fHistSub[i].fHistSubHar[j].fHistPsiSubs->SetXTitle("Event Plane Angle (rad)");
1253 fHistSub[i].fHistSubHar[j].fHistPsiSubs->SetYTitle("Counts");
1254 delete histTitle;
1255 }
1256 }
1257 }
1258
1259 if(fV0loop) // All V0s (if there, if flag on)
1260 {
1261 // Mass
327288af 1262 fHistV0Mass = new TH1F("FlowV0_InvMass", "FlowV0_InvMass", kMassBins, massMin, massMax);
9777bfcb 1263 fHistV0Mass->SetXTitle("Invariant Mass (GeV)");
1264 fHistV0Mass->SetYTitle("Counts");
1265 // Distance of closest approach
327288af 1266 fHistV0Dca = new TH1F("FlowV0_Dca", "FlowV0_Dca", kDcaBins, dcaMin, glDcaMax);
9777bfcb 1267 fHistV0Dca->SetXTitle("dca between tracks (cm)");
1268 fHistV0Dca->SetYTitle("Counts");
1269 // lenght
327288af 1270 fHistV0Lenght = new TH1F("FlowV0_Lenght", "FlowV0_Lenght", kLgBins, lgMinV0, lgMaxV0);
9777bfcb 1271 fHistV0Lenght->SetXTitle("Distance of V0s (cm)");
1272 fHistV0Lenght->SetYTitle("Counts");
1273 // Sigma for all particles
327288af 1274 fHistV0Sigma = new TH1F("FlowV0_Sigma", "FlowV0_Sigma", kLgBins, lgMinV0, lgMaxV0 );
9777bfcb 1275 fHistV0Sigma->SetXTitle("Sigma");
1276 fHistV0Sigma->SetYTitle("Counts");
1277 // Chi2
327288af 1278 fHistV0Chi2 = new TH1F("FlowV0_Chi2", "FlowV0_Chi2", kChi2Bins, chi2Min, chi2MaxC);
9777bfcb 1279 fHistV0Chi2->SetXTitle("Chi square at Main Vertex");
1280 fHistV0Chi2->SetYTitle("Counts");
1281 // EtaPtPhi
327288af 1282 fHistV0EtaPtPhi3D = new TH3F("FlowV0_EtaPtPhi3D", "FlowV0_EtaPtPhi3D", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax, kPhi3DBins, fPhiMin, fPhiMax);
9777bfcb 1283 fHistV0EtaPtPhi3D->SetXTitle("Eta");
1284 fHistV0EtaPtPhi3D->SetYTitle("Pt (GeV/c)");
1285 fHistV0EtaPtPhi3D->SetZTitle("Phi (rad)");
1286 // Yield for all v0s
da5aa0a0 1287 fHistV0YieldAll2D = new TH2D("FlowV0_PtEta2Dall", "FlowV0_PtEta2Dall", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax);
9777bfcb 1288 fHistV0YieldAll2D->Sumw2();
1289 fHistV0YieldAll2D->SetXTitle("Pseudorapidty");
1290 fHistV0YieldAll2D->SetYTitle("Pt (GeV/c)");
1291 // Mass slices on pT
327288af 1292 fHistV0MassPtSlices = new TH2D("FlowV0_MassPtSlices", "FlowV0_MassPtSlices", kMassBins, massMin, massMax, fPtBins, fPtMin, fPtMax);
9777bfcb 1293 fHistV0MassPtSlices->Sumw2();
1294 fHistV0MassPtSlices->SetXTitle("Invariant Mass (GeV)");
1295 fHistV0MassPtSlices->SetYTitle("Pt (GeV/c)");
1296 // Yield v0s (total P and Rapidity)
1297 fHistV0PYall2D = new TH2D("FlowV0_PYall2D", "FlowV0_PYall2D", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax);
1298 fHistV0PYall2D->Sumw2();
1299 fHistV0PYall2D->SetXTitle("Rapidty");
1300 fHistV0PYall2D->SetYTitle("P (GeV/c)");
1301
1302 // Selected V0s ...
1303 // Yield
da5aa0a0 1304 fHistV0YieldPart2D = new TH2D("FlowV0_PtEta2Dpart", "FlowV0_PtEta2Dpart", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax);
9777bfcb 1305 fHistV0YieldPart2D->Sumw2();
1306 fHistV0YieldPart2D->SetXTitle("Pseudorapidty");
1307 fHistV0YieldPart2D->SetYTitle("Pt (GeV/c)");
1308 // Mass Window
327288af 1309 fHistV0MassWin = new TH1F("FlowV0_MassWinPart", "FlowV0_MassWinPart", kMassBins, massMin, massMax);
9777bfcb 1310 fHistV0MassWin->SetXTitle("Invariant Mass (GeV)");
1311 fHistV0MassWin->SetYTitle("Counts");
1312 // EtaPtPhi
327288af 1313 fHistV0EtaPtPhi3DPart = new TH3F("FlowV0_EtaPtPhi3Dpart", "FlowV0_EtaPtPhi3Dpart", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax, kPhi3DBins, fPhiMin, fPhiMax);
9777bfcb 1314 fHistV0EtaPtPhi3DPart->SetXTitle("Eta");
1315 fHistV0EtaPtPhi3DPart->SetYTitle("Pt (GeV/c)");
1316 fHistV0EtaPtPhi3DPart->SetZTitle("Phi (rad)");
1317 // Distance of closest approach
327288af 1318 fHistV0DcaPart = new TH1F("FlowV0_DcaPart", "FlowV0_DcaPart", kDcaBins, dcaMin, dcaMax);
9777bfcb 1319 fHistV0DcaPart->Sumw2();
1320 fHistV0DcaPart->SetXTitle("dca between tracks (cm)");
1321 fHistV0DcaPart->SetYTitle("Counts");
1322 // lenght
327288af 1323 fHistV0LenghtPart = new TH1F("FlowV0_LenghtPart", "FlowV0_LenghtPart", kLgBins, lgMinV0, lgMaxV0);
9777bfcb 1324 fHistV0LenghtPart->SetXTitle("Distance of V0s (cm)");
1325 fHistV0LenghtPart->SetYTitle("Counts");
1326 // SideBand Mass (sidebands)
327288af 1327 fHistV0sbMassSide = new TH1F("FlowV0sb_MassWinSideBands", "FlowV0sb_MassWinSideBands", kMassBins, massMin, massMax);
9777bfcb 1328 fHistV0sbMassSide->SetXTitle("Invariant Mass (GeV)");
1329 fHistV0sbMassSide->SetYTitle("Counts");
1330 // EtaPtPhi (sidebands)
327288af 1331 fHistV0sbEtaPtPhi3DPart = new TH3F("FlowV0sb_EtaPtPhi3D", "FlowV0sb_EtaPtPhi3D", fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax, kPhi3DBins, fPhiMin, fPhiMax);
9777bfcb 1332 fHistV0sbEtaPtPhi3DPart->SetXTitle("Eta");
1333 fHistV0sbEtaPtPhi3DPart->SetYTitle("Pt (GeV/c)");
1334 fHistV0sbEtaPtPhi3DPart->SetZTitle("Phi (rad)");
1335 // Distance of closest approach (sidebands)
327288af 1336 fHistV0sbDcaPart = new TH1F("FlowV0sb_Dca", "FlowV0sb_Dca", kDcaBins, dcaMin, dcaMax);
9777bfcb 1337 fHistV0sbDcaPart->Sumw2();
1338 fHistV0sbDcaPart->SetXTitle("dca between tracks (cm)");
1339 fHistV0sbDcaPart->SetYTitle("Counts");
1340 // lenght (sidebands)
327288af 1341 fHistV0sbLenghtPart = new TH1F("FlowV0sb_Lenght", "FlowV0sb_Lenght", kLgBins, lgMinV0, lgMaxV0);
9777bfcb 1342 fHistV0sbLenghtPart->SetXTitle("Distance of V0s (cm)");
1343 fHistV0sbLenghtPart->SetYTitle("Counts");
1344
1345 // Mean Eta in each bin
1346 fHistV0BinEta = new TProfile("FlowV0_Bin_Eta", "FlowV0_Bin_Eta", fEtaBins, fEtaMin, fEtaMax, fEtaMin, fEtaMax, "");
1347 fHistV0BinEta->SetXTitle((char*)fLabel.Data());
1348 fHistV0BinEta->SetYTitle("<Eta>");
1349 // Mean Pt in each bin
1350 fHistV0BinPt = new TProfile("FlowV0_Bin_Pt", "FlowV0_Bin_Pt", fPtBinsPart, fPtMin, fPtMaxPart, fPtMin, fPtMaxPart, "");
1351 fHistV0BinPt->SetXTitle("Pt (GeV/c)");
1352 fHistV0BinPt->SetYTitle("<Pt> (GeV/c)");
1353 // Mean Eta in each bin (sidebands)
1354 fHistV0sbBinEta = new TProfile("FlowV0sb_Bin_Eta", "FlowV0sb_Bin_Eta", fEtaBins, fEtaMin, fEtaMax, fEtaMin, fEtaMax, "");
1355 fHistV0sbBinEta->SetXTitle((char*)fLabel.Data());
1356 fHistV0sbBinEta->SetYTitle("<Eta>");
1357 // Mean Pt in each bin (sidebands)
1358 fHistV0sbBinPt = new TProfile("FlowV0sb_Bin_Pt", "FlowV0sb_Bin_Pt", fPtBinsPart, fPtMin, fPtMaxPart, fPtMin, fPtMaxPart, "");
1359 fHistV0sbBinPt->SetXTitle("Pt (GeV/c)");
1360 fHistV0sbBinPt->SetYTitle("<Pt> (GeV/c)");
1361 }
1362
327288af 1363 for (int k = 0; k < AliFlowConstants::kSels; k++) // for each selection
9777bfcb 1364 {
1365 // cos(n*delta_Psi)
1366 histTitle = new TString("Flow_Cos_Sel");
1367 *histTitle += k+1;
327288af 1368 fHistFull[k].fHistCos = new TProfile(histTitle->Data(), histTitle->Data(), AliFlowConstants::kHars, 0.5, (float)(AliFlowConstants::kHars) + 0.5, -1., 1., "");
9777bfcb 1369 fHistFull[k].fHistCos->SetXTitle("Harmonic");
1370 fHistFull[k].fHistCos->SetYTitle("<cos(n*delta_Psi)>");
1371 delete histTitle;
1372
1373 // resolution
1374 histTitle = new TString("Flow_Res_Sel");
1375 *histTitle += k+1;
327288af 1376 fHistFull[k].fHistRes = new TH1F(histTitle->Data(), histTitle->Data(), AliFlowConstants::kHars, 0.5, (float)(AliFlowConstants::kHars) + 0.5);
9777bfcb 1377 fHistFull[k].fHistRes->SetXTitle("Harmonic");
1378 fHistFull[k].fHistRes->SetYTitle("Resolution");
1379 delete histTitle;
1380
1381 // vObs
1382 histTitle = new TString("Flow_vObs_Sel");
1383 *histTitle += k+1;
327288af 1384 fHistFull[k].fHistvObs = new TProfile(histTitle->Data(), histTitle->Data(), AliFlowConstants::kHars, 0.5, (float)(AliFlowConstants::kHars) + 0.5, -100., 100., "");
9777bfcb 1385 fHistFull[k].fHistvObs->SetXTitle("Harmonic");
1386 fHistFull[k].fHistvObs->SetYTitle("vObs (%)");
1387 delete histTitle;
1388
1389 // vObs V0
1390 histTitle = new TString("FlowV0_vObs_Sel");
1391 *histTitle += k+1;
327288af 1392 fHistFull[k].fHistV0vObs = new TProfile(histTitle->Data(), histTitle->Data(), AliFlowConstants::kHars, 0.5, (float)(AliFlowConstants::kHars) + 0.5, -100., 100., "");
9777bfcb 1393 fHistFull[k].fHistV0vObs->SetXTitle("Harmonic");
1394 fHistFull[k].fHistV0vObs->SetYTitle("vObs (%)");
1395 delete histTitle;
1396
1397 // vObs V0 sideband SX
1398 histTitle = new TString("FlowV0sb_vObs_sx_Sel");
1399 *histTitle += k+1;
327288af 1400 fHistFull[k].fHistV0sbvObsSx = new TProfile(histTitle->Data(), histTitle->Data(), AliFlowConstants::kHars, 0.5, (float)(AliFlowConstants::kHars) + 0.5, -100., 100., "");
9777bfcb 1401 fHistFull[k].fHistV0sbvObsSx->SetXTitle("Harmonic");
1402 fHistFull[k].fHistV0sbvObsSx->SetYTitle("vObs (%)");
1403 delete histTitle;
1404
1405 // vObs V0 sideband DX
1406 histTitle = new TString("FlowV0sb_vObs_dx_Sel");
1407 *histTitle += k+1;
327288af 1408 fHistFull[k].fHistV0sbvObsDx = new TProfile(histTitle->Data(), histTitle->Data(), AliFlowConstants::kHars, 0.5, (float)(AliFlowConstants::kHars) + 0.5, -100., 100., "");
9777bfcb 1409 fHistFull[k].fHistV0sbvObsDx->SetXTitle("Harmonic");
1410 fHistFull[k].fHistV0sbvObsDx->SetYTitle("vObs (%)");
1411 delete histTitle;
1412
1413 // PID for tracks used in R.P.
1414 histTitle = new TString("Flow_BayPidMult_Sel");
1415 *histTitle += k+1;
327288af 1416 fHistFull[k].fHistBayPidMult = new TH1F(histTitle->Data(), histTitle->Data(),AliFlowConstants::kPid,-0.5,((float)AliFlowConstants::kPid-0.5));
9777bfcb 1417 fHistFull[k].fHistBayPidMult->Sumw2() ;
1418 fHistFull[k].fHistBayPidMult->SetXTitle("e+/- , mu+/- , pi+/- , K+/- , p+/- , d+/- ");
1419 fHistFull[k].fHistBayPidMult->SetYTitle("Counts");
1420 delete histTitle;
1421
327288af 1422 for (int j = 0; j < AliFlowConstants::kHars; j++) // for each harmonic
9777bfcb 1423 {
1424 float order = (float)(j+1);
1425
1426 // multiplicity
1427 histTitle = new TString("Flow_Mul_Sel");
1428 *histTitle += k+1;
1429 histTitle->Append("_Har");
1430 *histTitle += j+1;
327288af 1431 fHistFull[k].fHistFullHar[j].fHistMult = new TH1F(histTitle->Data(),histTitle->Data(), kMultBins, multMin, multMax);
9777bfcb 1432 fHistFull[k].fHistFullHar[j].fHistMult->SetXTitle("Multiplicity");
1433 fHistFull[k].fHistFullHar[j].fHistMult->SetYTitle("Counts");
1434 delete histTitle;
1435
1436 // event plane
1437 histTitle = new TString("Flow_Psi_Sel");
1438 *histTitle += k+1;
1439 histTitle->Append("_Har");
1440 *histTitle += j+1;
327288af 1441 fHistFull[k].fHistFullHar[j].fHistPsi = new TH1F(histTitle->Data(), histTitle->Data(), kPsiBins, psiMin, psiMax / order);
9777bfcb 1442 fHistFull[k].fHistFullHar[j].fHistPsi->SetXTitle("Event Plane Angle (rad)");
1443 fHistFull[k].fHistFullHar[j].fHistPsi->SetYTitle("Counts");
1444 delete histTitle;
1445
1446 // event plane difference of two selections
1447 histTitle = new TString("Flow_Psi_Diff_Sel");
1448 *histTitle += k+1;
1449 histTitle->Append("_Har");
1450 *histTitle += j+1;
1451 if (k == 0 )
1452 {
327288af 1453 Int_t myOrder = 1;
1454 if (j == 1) { myOrder = 2 ; }
1455 fHistFull[k].fHistFullHar[j].fHistPsiDiff = new TH1F(histTitle->Data(), histTitle->Data(), kPsiBins, -psiMax/myOrder/2., psiMax/myOrder/2.);
9777bfcb 1456 }
1457 else
1458 {
327288af 1459 fHistFull[k].fHistFullHar[j].fHistPsiDiff = new TH1F(histTitle->Data(), histTitle->Data(), kPsiBins, -psiMax/2., psiMax/2.);
9777bfcb 1460 }
1461 if (k == 0)
1462 {
1463 if (j == 0)
1464 {
1465 fHistFull[k].fHistFullHar[j].fHistPsiDiff->SetXTitle("#Psi_{1,Sel1} - #Psi_{1,Sel2}(rad)");
1466 }
1467 else if (j == 1)
1468 {
1469 fHistFull[k].fHistFullHar[j].fHistPsiDiff->SetXTitle("#Psi_{2,Sel1} - #Psi_{2,Sel2}(rad)");
1470 }
1471 }
1472 else if (k == 1)
1473 {
1474 if (j == 0)
1475 {
1476 fHistFull[k].fHistFullHar[j].fHistPsiDiff->SetXTitle("#Psi_{1,Sel1} - #Psi_{2,Sel2}(rad)");
1477 }
1478 else if (j == 1)
1479 {
1480 fHistFull[k].fHistFullHar[j].fHistPsiDiff->SetXTitle("#Psi_{1,Sel1} - #Psi_{2,Sel1}(rad)");
1481 }
1482 }
1483 fHistFull[k].fHistFullHar[j].fHistPsiDiff->SetYTitle("Counts");
1484 delete histTitle;
1485
1486 // correlation of sub-event planes
1487 histTitle = new TString("Flow_Psi_Sub_Corr_Sel");
1488 *histTitle += k+1;
1489 histTitle->Append("_Har");
1490 *histTitle += j+1;
327288af 1491 fHistFull[k].fHistFullHar[j].fHistPsiSubCorr = new TH1F(histTitle->Data(), histTitle->Data(), kPsiBins, psiMin, psiMax / order);
9777bfcb 1492 fHistFull[k].fHistFullHar[j].fHistPsiSubCorr->Sumw2();
1493 fHistFull[k].fHistFullHar[j].fHistPsiSubCorr->SetXTitle("Sub-Event Correlation (rad)");
1494 fHistFull[k].fHistFullHar[j].fHistPsiSubCorr->SetYTitle("Counts");
1495 delete histTitle;
1496
1497 // correlation of sub-event planes of different order
1498 histTitle = new TString("Flow_Psi_Sub_Corr_Diff_Sel");
1499 *histTitle += k+1;
1500 histTitle->Append("_Har");
1501 *histTitle += j+1;
327288af 1502 fHistFull[k].fHistFullHar[j].fHistPsiSubCorrDiff = new TH1F(histTitle->Data(), histTitle->Data(), kPsiBins, psiMin, psiMax / (order+1.));
9777bfcb 1503 fHistFull[k].fHistFullHar[j].fHistPsiSubCorrDiff->Sumw2();
1504 fHistFull[k].fHistFullHar[j].fHistPsiSubCorrDiff->SetXTitle("Sub-Event Correlation (rad)");
1505 fHistFull[k].fHistFullHar[j].fHistPsiSubCorrDiff->SetYTitle("Counts");
1506 delete histTitle;
1507
1508 // q
1509 histTitle = new TString("Flow_NormQ_Sel");
1510 *histTitle += k+1;
1511 histTitle->Append("_Har");
1512 *histTitle += j+1;
327288af 1513 fHistFull[k].fHistFullHar[j].fHistQnorm = new TH1F(histTitle->Data(), histTitle->Data(), kQbins, qMin, qMax);
9777bfcb 1514 fHistFull[k].fHistFullHar[j].fHistQnorm->Sumw2();
1515 fHistFull[k].fHistFullHar[j].fHistQnorm->SetXTitle("q = |Q|/sqrt(Mult)");
1516 fHistFull[k].fHistFullHar[j].fHistQnorm->SetYTitle("Counts");
1517 delete histTitle;
1518
1519 // particle-plane azimuthal correlation
1520 histTitle = new TString("Flow_Phi_Corr_Sel");
1521 *histTitle += k+1;
1522 histTitle->Append("_Har");
1523 *histTitle += j+1;
1524 fHistFull[k].fHistFullHar[j].fHistPhiCorr = new TH1F(histTitle->Data(), histTitle->Data(), fPhiBins, fPhiMin, fPhiMax / order);
1525 fHistFull[k].fHistFullHar[j].fHistPhiCorr->Sumw2();
1526 fHistFull[k].fHistFullHar[j].fHistPhiCorr->SetXTitle("Particle-Plane Correlation (rad)");
1527 fHistFull[k].fHistFullHar[j].fHistPhiCorr->SetYTitle("Counts");
1528 delete histTitle;
1529
1530 // neutral particle-plane azimuthal correlation
1531 histTitle = new TString("FlowV0_Phi_Corr_Sel");
1532 *histTitle += k+1;
1533 histTitle->Append("_Har");
1534 *histTitle += j+1;
1535 fHistFull[k].fHistFullHar[j].fHistV0PhiCorr = new TH1F(histTitle->Data(), histTitle->Data(), fPhiBins, fPhiMin, fPhiMax / order);
1536 fHistFull[k].fHistFullHar[j].fHistV0PhiCorr->Sumw2();
1537 fHistFull[k].fHistFullHar[j].fHistV0PhiCorr->SetXTitle("V0-Plane Correlation (rad)");
1538 fHistFull[k].fHistFullHar[j].fHistV0PhiCorr->SetYTitle("Counts");
1539 delete histTitle;
1540
1541 // neutral sidebands-plane azimuthal correlation
1542 histTitle = new TString("FlowV0sb_Phi_Corr_Sel");
1543 *histTitle += k+1;
1544 histTitle->Append("_Har");
1545 *histTitle += j+1;
1546 fHistFull[k].fHistFullHar[j].fHistV0sbPhiCorr = new TH1F(histTitle->Data(), histTitle->Data(), fPhiBins, fPhiMin, fPhiMax / order);
1547 fHistFull[k].fHistFullHar[j].fHistV0sbPhiCorr->Sumw2();
1548 fHistFull[k].fHistFullHar[j].fHistV0sbPhiCorr->SetXTitle("V0sideBands-Plane Correlation (rad)");
1549 fHistFull[k].fHistFullHar[j].fHistV0sbPhiCorr->SetYTitle("Counts");
1550 delete histTitle;
1551
1552 // Yield(pt)
1553 histTitle = new TString("Flow_YieldPt_Sel");
1554 *histTitle += k+1;
1555 histTitle->Append("_Har");
1556 *histTitle += j+1;
1557 fHistFull[k].fHistFullHar[j].fHistYieldPt = new TH1F(histTitle->Data(), histTitle->Data(), fPtBins, fPtMin, fPtMax);
1558 fHistFull[k].fHistFullHar[j].fHistYieldPt->Sumw2();
1559 fHistFull[k].fHistFullHar[j].fHistYieldPt->SetXTitle("Pt (GeV/c)");
1560 fHistFull[k].fHistFullHar[j].fHistYieldPt->SetYTitle("Yield");
1561 delete histTitle;
1562
1563 // EtaPtPhi
1564 histTitle = new TString("Flow_EtaPtPhi3D_Sel");
1565 *histTitle += k+1;
1566 histTitle->Append("_Har");
1567 *histTitle += j+1;
327288af 1568 fHistFull[k].fHistFullHar[j].fHistEtaPtPhi3D = new TH3F(histTitle->Data(), histTitle->Data(), fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax, kPhi3DBins, fPhiMin, fPhiMax);
9777bfcb 1569 fHistFull[k].fHistFullHar[j].fHistEtaPtPhi3D->SetXTitle("Eta");
1570 fHistFull[k].fHistFullHar[j].fHistEtaPtPhi3D->SetYTitle("Pt (GeV/c)");
1571 fHistFull[k].fHistFullHar[j].fHistEtaPtPhi3D->SetZTitle("Phi (rad)");
1572 delete histTitle;
1573
1574 // Yield(eta,pt)
1575 histTitle = new TString("Flow_Yield2D_Sel");
1576 *histTitle += k+1;
1577 histTitle->Append("_Har");
1578 *histTitle += j+1;
1579 fHistFull[k].fHistFullHar[j].fHistYield2D = new TH2D(histTitle->Data(), histTitle->Data(), fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax);
1580 fHistFull[k].fHistFullHar[j].fHistYield2D->Sumw2();
1581 fHistFull[k].fHistFullHar[j].fHistYield2D->SetXTitle("Pseudorapidty");
1582 fHistFull[k].fHistFullHar[j].fHistYield2D->SetYTitle("Pt (GeV/c)");
1583 delete histTitle;
1584
1585 // Dca - 3D
1586 histTitle = new TString("Flow_3dDca_Sel") ;
1587 *histTitle += k+1;
1588 histTitle->Append("_Har");
1589 *histTitle += j+1;
327288af 1590 fHistFull[k].fHistFullHar[j].fHistDcaGlob = new TH1F(histTitle->Data(), histTitle->Data(), kDcaBins, dcaMin, glDcaMax);
9777bfcb 1591 fHistFull[k].fHistFullHar[j].fHistDcaGlob->Sumw2();
1592 fHistFull[k].fHistFullHar[j].fHistDcaGlob->SetXTitle("|3d Global Track's dca to Vertex (cm)|");
1593 delete histTitle;
1594
1595 // Yield(pt) - excluded from R.P.
1596 histTitle = new TString("Flow_YieldPt_outSel");
1597 *histTitle += k+1;
1598 histTitle->Append("_Har");
1599 *histTitle += j+1;
1600 fHistFull[k].fHistFullHar[j].fHistYieldPtout = new TH1F(histTitle->Data(), histTitle->Data(), fPtBins, fPtMin, fPtMax);
1601 fHistFull[k].fHistFullHar[j].fHistYieldPtout->Sumw2();
1602 fHistFull[k].fHistFullHar[j].fHistYieldPtout->SetXTitle("Pt (GeV/c)");
1603 fHistFull[k].fHistFullHar[j].fHistYieldPtout->SetYTitle("Yield");
1604 delete histTitle;
1605
1606 // EtaPtPhi - excluded from R.P.
1607 histTitle = new TString("Flow_EtaPtPhi3D_outSel");
1608 *histTitle += k+1;
1609 histTitle->Append("_Har");
1610 *histTitle += j+1;
327288af 1611 fHistFull[k].fHistFullHar[j].fHistEtaPtPhi3Dout = new TH3F(histTitle->Data(), histTitle->Data(), fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax, kPhi3DBins, fPhiMin, fPhiMax);
9777bfcb 1612 fHistFull[k].fHistFullHar[j].fHistEtaPtPhi3Dout->SetXTitle("Eta");
1613 fHistFull[k].fHistFullHar[j].fHistEtaPtPhi3Dout->SetYTitle("Pt (GeV/c)");
1614 fHistFull[k].fHistFullHar[j].fHistEtaPtPhi3Dout->SetZTitle("Phi (rad)");
1615 delete histTitle;
1616
1617 // Yield(eta,pt) - excluded from R.P.
1618 histTitle = new TString("Flow_Yield2D_outSel");
1619 *histTitle += k+1;
1620 histTitle->Append("_Har");
1621 *histTitle += j+1;
1622 fHistFull[k].fHistFullHar[j].fHistYield2Dout = new TH2D(histTitle->Data(), histTitle->Data(), fEtaBins, fEtaMin, fEtaMax, fPtBins, fPtMin, fPtMax);
1623 fHistFull[k].fHistFullHar[j].fHistYield2Dout->Sumw2();
1624 fHistFull[k].fHistFullHar[j].fHistYield2Dout->SetXTitle("Pseudorapidty");
1625 fHistFull[k].fHistFullHar[j].fHistYield2Dout->SetYTitle("Pt (GeV/c)");
1626 delete histTitle;
1627
1628 // Dca - 3D - excluded from R.P.
1629 histTitle = new TString("Flow_3dDca_outSel") ;
1630 *histTitle += k+1;
1631 histTitle->Append("_Har");
1632 *histTitle += j+1;
327288af 1633 fHistFull[k].fHistFullHar[j].fHistDcaGlobout = new TH1F(histTitle->Data(), histTitle->Data(), kDcaBins, dcaMin, glDcaMax);
9777bfcb 1634 fHistFull[k].fHistFullHar[j].fHistDcaGlobout->Sumw2();
1635 fHistFull[k].fHistFullHar[j].fHistDcaGlobout->SetXTitle("|3d Global Track's dca to Vertex (cm)|");
1636 delete histTitle;
1637
1638 // Flow observed - v_obs,Pt,Eta
1639 histTitle = new TString("Flow_vObs2D_Sel");
1640 *histTitle += k+1;
1641 histTitle->Append("_Har");
1642 *histTitle += j+1;
1643 fHistFull[k].fHistFullHar[j].fHistvObs2D = new TProfile2D(histTitle->Data(), histTitle->Data(), fEtaBins, fEtaMin, fEtaMax, fPtBinsPart, fPtMin, fPtMaxPart, -100., 100., "");
1644 fHistFull[k].fHistFullHar[j].fHistvObs2D->SetXTitle((char*)fLabel.Data());
1645 fHistFull[k].fHistFullHar[j].fHistvObs2D->SetYTitle("Pt (GeV/c)");
1646 delete histTitle;
1647
1648 // v_obs,Eta
1649 histTitle = new TString("Flow_vObsEta_Sel");
1650 *histTitle += k+1;
1651 histTitle->Append("_Har");
1652 *histTitle += j+1;
1653 fHistFull[k].fHistFullHar[j].fHistvObsEta = new TProfile(histTitle->Data(), histTitle->Data(), fEtaBins, fEtaMin, fEtaMax, -100., 100., "");
1654 fHistFull[k].fHistFullHar[j].fHistvObsEta->SetXTitle((char*)fLabel.Data());
1655 fHistFull[k].fHistFullHar[j].fHistvObsEta->SetYTitle("v (%)");
1656 delete histTitle;
1657
1658 // v_obs,Pt
1659 histTitle = new TString("Flow_vObsPt_Sel");
1660 *histTitle += k+1;
1661 histTitle->Append("_Har");
1662 *histTitle += j+1;
1663 fHistFull[k].fHistFullHar[j].fHistvObsPt = new TProfile(histTitle->Data(), histTitle->Data(), fPtBinsPart, fPtMin, fPtMaxPart, -100., 100., "");
1664 fHistFull[k].fHistFullHar[j].fHistvObsPt->SetXTitle("Pt (GeV/c)");
1665 fHistFull[k].fHistFullHar[j].fHistvObsPt->SetYTitle("v (%)");
1666 delete histTitle;
1667
1668 // neutral Flow observed - Pt,Eta
1669 histTitle = new TString("FlowV0_vObs2D_Sel");
1670 *histTitle += k+1;
1671 histTitle->Append("_Har");
1672 *histTitle += j+1;
1673 fHistFull[k].fHistFullHar[j].fHistV0vObs2D = new TProfile2D(histTitle->Data(), histTitle->Data(), fEtaBins, fEtaMin, fEtaMax, fPtBinsPart, fPtMin, fPtMaxPart, -100., 100., "");
1674 fHistFull[k].fHistFullHar[j].fHistV0vObs2D->SetXTitle((char*)fLabel.Data());
1675 fHistFull[k].fHistFullHar[j].fHistV0vObs2D->SetYTitle("Pt (GeV/c)");
1676 delete histTitle;
1677
1678 // neutral Flow observed - Eta
1679 histTitle = new TString("FlowV0_vObsEta_Sel");
1680 *histTitle += k+1;
1681 histTitle->Append("_Har");
1682 *histTitle += j+1;
1683 fHistFull[k].fHistFullHar[j].fHistV0vObsEta = new TProfile(histTitle->Data(), histTitle->Data(), fEtaBins, fEtaMin, fEtaMax, -100., 100., "");
1684 fHistFull[k].fHistFullHar[j].fHistV0vObsEta->SetXTitle((char*)fLabel.Data());
1685 fHistFull[k].fHistFullHar[j].fHistV0vObsEta->SetYTitle("v (%)");
1686 delete histTitle;
1687
1688 // neutral Flow observed - Pt
1689 histTitle = new TString("FlowV0_vObsPt_Sel");
1690 *histTitle += k+1;
1691 histTitle->Append("_Har");
1692 *histTitle += j+1;
1693 fHistFull[k].fHistFullHar[j].fHistV0vObsPt = new TProfile(histTitle->Data(), histTitle->Data(), fPtBinsPart, fPtMin, fPtMaxPart, -100., 100., "");
1694 fHistFull[k].fHistFullHar[j].fHistV0vObsPt->SetXTitle("Pt (GeV/c)");
1695 fHistFull[k].fHistFullHar[j].fHistV0vObsPt->SetYTitle("v (%)");
1696 delete histTitle;
1697
1698 // neutral sidebands Flow observed - Pt,Eta
1699 histTitle = new TString("FlowV0sb_vObs2D_Sel");
1700 *histTitle += k+1;
1701 histTitle->Append("_Har");
1702 *histTitle += j+1;
1703 fHistFull[k].fHistFullHar[j].fHistV0sbvObs2D = new TProfile2D(histTitle->Data(), histTitle->Data(), fEtaBins, fEtaMin, fEtaMax, fPtBinsPart, fPtMin, fPtMaxPart, -100., 100., "");
1704 fHistFull[k].fHistFullHar[j].fHistV0sbvObs2D->SetXTitle((char*)fLabel.Data());
1705 fHistFull[k].fHistFullHar[j].fHistV0sbvObs2D->SetYTitle("Pt (GeV/c)");
1706 delete histTitle;
1707
1708 // neutral sidebands Flow observed - Eta
1709 histTitle = new TString("FlowV0sb_vObsEta_Sel");
1710 *histTitle += k+1;
1711 histTitle->Append("_Har");
1712 *histTitle += j+1;
1713 fHistFull[k].fHistFullHar[j].fHistV0sbvObsEta = new TProfile(histTitle->Data(), histTitle->Data(), fEtaBins, fEtaMin, fEtaMax, -100., 100., "");
1714 fHistFull[k].fHistFullHar[j].fHistV0sbvObsEta->SetXTitle((char*)fLabel.Data());
1715 fHistFull[k].fHistFullHar[j].fHistV0sbvObsEta->SetYTitle("v (%)");
1716 delete histTitle;
1717
1718 // neutral sidebands Flow observed - Pt
1719 histTitle = new TString("FlowV0sb_vObsPt_Sel");
1720 *histTitle += k+1;
1721 histTitle->Append("_Har");
1722 *histTitle += j+1;
1723 fHistFull[k].fHistFullHar[j].fHistV0sbvObsPt = new TProfile(histTitle->Data(), histTitle->Data(), fPtBinsPart, fPtMin, fPtMaxPart, -100., 100., "");
1724 fHistFull[k].fHistFullHar[j].fHistV0sbvObsPt->SetXTitle("Pt (GeV/c)");
1725 fHistFull[k].fHistFullHar[j].fHistV0sbvObsPt->SetYTitle("v (%)");
1726 delete histTitle;
1727
1728 // SX neutral sidebands Flow observed - Eta
1729 histTitle = new TString("FlowV0sb_vObsEta_sx_Sel");
1730 *histTitle += k+1;
1731 histTitle->Append("_Har");
1732 *histTitle += j+1;
1733 fHistFull[k].fHistFullHar[j].fHistV0sbvObsEtaSx = new TProfile(histTitle->Data(), histTitle->Data(), fEtaBins, fEtaMin, fEtaMax, -100., 100., "");
1734 fHistFull[k].fHistFullHar[j].fHistV0sbvObsEtaSx->SetXTitle((char*)fLabel.Data());
1735 fHistFull[k].fHistFullHar[j].fHistV0sbvObsEtaSx->SetYTitle("v (%)");
1736 delete histTitle;
1737
1738 // SX neutral sidebands Flow observed - Pt
1739 histTitle = new TString("FlowV0sb_vObsPt_sx_Sel");
1740 *histTitle += k+1;
1741 histTitle->Append("_Har");
1742 *histTitle += j+1;
1743 fHistFull[k].fHistFullHar[j].fHistV0sbvObsPtSx = new TProfile(histTitle->Data(), histTitle->Data(), fPtBinsPart, fPtMin, fPtMaxPart, -100., 100., "");
1744 fHistFull[k].fHistFullHar[j].fHistV0sbvObsPtSx->SetXTitle("Pt (GeV/c)");
1745 fHistFull[k].fHistFullHar[j].fHistV0sbvObsPtSx->SetYTitle("v (%)");
1746 delete histTitle;
1747
1748 // DX neutral sidebands Flow observed - Eta
1749 histTitle = new TString("FlowV0sb_vObsEta_dx_Sel");
1750 *histTitle += k+1;
1751 histTitle->Append("_Har");
1752 *histTitle += j+1;
1753 fHistFull[k].fHistFullHar[j].fHistV0sbvObsEtaDx = new TProfile(histTitle->Data(), histTitle->Data(), fEtaBins, fEtaMin, fEtaMax, -100., 100., "");
1754 fHistFull[k].fHistFullHar[j].fHistV0sbvObsEtaDx->SetXTitle((char*)fLabel.Data());
1755 fHistFull[k].fHistFullHar[j].fHistV0sbvObsEtaDx->SetYTitle("v (%)");
1756 delete histTitle;
1757
1758 // DX neutral sidebands Flow observed - Pt
1759 histTitle = new TString("FlowV0sb_vObsPt_dx_Sel");
1760 *histTitle += k+1;
1761 histTitle->Append("_Har");
1762 *histTitle += j+1;
1763 fHistFull[k].fHistFullHar[j].fHistV0sbvObsPtDx = new TProfile(histTitle->Data(), histTitle->Data(), fPtBinsPart, fPtMin, fPtMaxPart, -100., 100., "");
1764 fHistFull[k].fHistFullHar[j].fHistV0sbvObsPtDx->SetXTitle("Pt (GeV/c)");
1765 fHistFull[k].fHistFullHar[j].fHistV0sbvObsPtDx->SetYTitle("v (%)");
1766 delete histTitle;
1767
1768 // Phi lab
1769 // Tpc (plus)
1770 histTitle = new TString("Flow_Phi_TPCplus_Sel");
1771 *histTitle += k+1;
1772 histTitle->Append("_Har");
1773 *histTitle += j+1;
1774 fHistFull[k].fHistFullHar[j].fHistPhiPlus = new TH1D(histTitle->Data(), histTitle->Data(), fPhiBins, fPhiMin, fPhiMax);
1775 fHistFull[k].fHistFullHar[j].fHistPhiPlus->SetXTitle("Azimuthal Angles (rad)");
1776 fHistFull[k].fHistFullHar[j].fHistPhiPlus->SetYTitle("Counts");
1777 delete histTitle;
1778
1779 // Tpc (minus)
1780 histTitle = new TString("Flow_Phi_TPCminus_Sel");
1781 *histTitle += k+1;
1782 histTitle->Append("_Har");
1783 *histTitle += j+1;
1784 fHistFull[k].fHistFullHar[j].fHistPhiMinus = new TH1D(histTitle->Data(), histTitle->Data(), fPhiBins, fPhiMin, fPhiMax);
1785 fHistFull[k].fHistFullHar[j].fHistPhiMinus->SetXTitle("Azimuthal Angles (rad)");
1786 fHistFull[k].fHistFullHar[j].fHistPhiMinus->SetYTitle("Counts");
1787 delete histTitle;
1788
1789 // Tpc (cross)
1790 histTitle = new TString("Flow_Phi_TPCcross_Sel");
1791 *histTitle += k+1;
1792 histTitle->Append("_Har");
1793 *histTitle += j+1;
1794 fHistFull[k].fHistFullHar[j].fHistPhiAll = new TH1D(histTitle->Data(), histTitle->Data(), fPhiBins, fPhiMin, fPhiMax);
1795 fHistFull[k].fHistFullHar[j].fHistPhiAll->SetXTitle("Azimuthal Angles (rad)");
1796 fHistFull[k].fHistFullHar[j].fHistPhiAll->SetYTitle("Counts");
1797 delete histTitle;
1798
1799 // Tpc
1800 histTitle = new TString("Flow_Phi_TPC_Sel");
1801 *histTitle += k+1;
1802 histTitle->Append("_Har");
1803 *histTitle += j+1;
1804 fHistFull[k].fHistFullHar[j].fHistPhi = new TH1D(histTitle->Data(), histTitle->Data(), fPhiBins, fPhiMin, fPhiMax);
1805 fHistFull[k].fHistFullHar[j].fHistPhi->SetXTitle("Azimuthal Angles (rad)");
1806 fHistFull[k].fHistFullHar[j].fHistPhi->SetYTitle("Counts");
1807 delete histTitle;
1808
1809 // Phi lab flattened
1810 // Tpc (Plus)
1811 histTitle = new TString("Flow_Phi_Flat_TPCplus_Sel");
1812 *histTitle += k+1;
1813 histTitle->Append("_Har");
1814 *histTitle += j+1;
1815 fHistFull[k].fHistFullHar[j].fHistPhiFlatPlus = new TH1D(histTitle->Data(), histTitle->Data(), fPhiBins, fPhiMin, fPhiMax);
1816 fHistFull[k].fHistFullHar[j].fHistPhiFlatPlus->SetXTitle("Azimuthal Angles (rad)");
1817 fHistFull[k].fHistFullHar[j].fHistPhiFlatPlus->SetYTitle("Counts");
1818 delete histTitle;
1819
1820 // Tpc (Minus)
1821 histTitle = new TString("Flow_Phi_Flat_TPCminus_Sel");
1822 *histTitle += k+1;
1823 histTitle->Append("_Har");
1824 *histTitle += j+1;
1825 fHistFull[k].fHistFullHar[j].fHistPhiFlatMinus = new TH1D(histTitle->Data(), histTitle->Data(), fPhiBins, fPhiMin, fPhiMax);
1826 fHistFull[k].fHistFullHar[j].fHistPhiFlatMinus->SetXTitle("Azimuthal Angles (rad)");
1827 fHistFull[k].fHistFullHar[j].fHistPhiFlatMinus->SetYTitle("Counts");
1828 delete histTitle;
1829
1830 // Tpc (cross)
1831 histTitle = new TString("Flow_Phi_Flat_TPCcross_Sel");
1832 *histTitle += k+1;
1833 histTitle->Append("_Har");
1834 *histTitle += j+1;
1835 fHistFull[k].fHistFullHar[j].fHistPhiFlatAll = new TH1D(histTitle->Data(), histTitle->Data(), fPhiBins, fPhiMin, fPhiMax);
1836 fHistFull[k].fHistFullHar[j].fHistPhiFlatAll->SetXTitle("Azimuthal Angles (rad)");
1837 fHistFull[k].fHistFullHar[j].fHistPhiFlatAll->SetYTitle("Counts");
1838 delete histTitle;
1839
1840 // Tpc
1841 histTitle = new TString("Flow_Phi_Flat_TPC_Sel");
1842 *histTitle += k+1;
1843 histTitle->Append("_Har");
1844 *histTitle += j+1;
1845 fHistFull[k].fHistFullHar[j].fHistPhiFlat = new TH1D(histTitle->Data(), histTitle->Data(), fPhiBins, fPhiMin, fPhiMax);
1846 fHistFull[k].fHistFullHar[j].fHistPhiFlat->SetXTitle("Azimuthal Angles (rad)");
1847 fHistFull[k].fHistFullHar[j].fHistPhiFlat->SetYTitle("Counts");
1848 delete histTitle;
1849 }
1850 }
1851 cout << " Init() - Histograms booked" << endl ;
1852
1853 return kTRUE ;
1854}
1855//-----------------------------------------------------------------------
1856
1857//-----------------------------------------------------------------------
1858Bool_t AliFlowAnalyser::Finish()
1859{
1860 // Close the analysis and saves the histograms on the histFile .
1861
1862 cout << "* FlowAnalysis * - Finish()" << endl ; cout << endl ;
1863
1864 // Write all histograms
1865 fHistFile->cd() ;
1866 fHistFile->Write() ;
1867
1868 // Write Resolution corrected histograms
1869 if(fVnResHistList)
1870 {
1871 fVnResHistList->Write();
1872 delete fVnResHistList ;
1873 }
1874 else { cout << " E.P. resolution has not been calculated. No v_n histograms!" << endl ; }
1875
1876 // Write PhiWgt histograms
1877 if(fPhiWgtHistList)
1878 {
1879 fPhiWgtHistList->Write();
1880 delete fPhiWgtHistList ;
1881 }
1882
1883 fFlowSelect->Write();
1884 // delete fFlowSelect ;
1885
1886 fHistFile->Close() ;
1887
da5aa0a0 1888 cout << endl ;
1889 cout << " Finish() - tot. " << fTotalNumberOfEvents << " have been analyzed (with " << fTotalNumberOfTracks << " tracks and " << fTotalNumberOfV0s << " v0s) ." << endl ;
1890 cout << endl ;
9777bfcb 1891 cout << " Finish() - Histograms saved : " << fHistFileName.Data() << endl ; cout << endl ;
1892
1893 return kTRUE ;
1894}
1895//-----------------------------------------------------------------------
1896// ###
1897//----------------------------------------------------------------------
327288af 1898Float_t AliFlowAnalyser::GetRunBayesian(Int_t nPid, Int_t selN) const
9777bfcb 1899{
1900 // Returns the normalized particle abundance of "e","mu","pi","k","p","d"
1901 // in all the analysed events (in selection selN).
1902 // Call at the end of the analysis.
1903
327288af 1904 if(selN>AliFlowConstants::kSels) { selN = 0 ; }
9777bfcb 1905 Double_t totCount = (fHistFull[selN].fHistBayPidMult)->GetSumOfWeights() ;
1906 if(totCount) { return (fHistFull[selN].fHistBayPidMult->GetBinContent(nPid+1) / totCount) ; }
1907 else { return 1. ; }
1908}
1909//-----------------------------------------------------------------------
327288af 1910void AliFlowAnalyser::PrintRunBayesian(Int_t selN) const
9777bfcb 1911{
1912 // Prints the normalized particle abundance of all the analysed events
1913 // (in selection selN).
1914
327288af 1915 if(selN>AliFlowConstants::kSels) { selN = 0 ; }
1916 Char_t* names[AliFlowConstants::kPid] = {"e","mu","pi","k","p","d"} ;
9777bfcb 1917 Double_t bayes = 0. ;
1918 cout << " selN = " << selN << " particles normalized abundance : " ;
327288af 1919 for(int i=0;i<AliFlowConstants::kPid;i++)
9777bfcb 1920 {
1921 bayes = GetRunBayesian(i, selN) ;
1922 cout << bayes << "_" << names[i] << " ; " ;
1923 }
1924 cout << endl ;
1925
1926 return ;
1927}
1928//-----------------------------------------------------------------------
1929void AliFlowAnalyser::FillWgtArrays(TFile* wgtFile)
1930{
1931 // Loads PhiWeights & Bayesian particles' abundance from file (default: flowPhiWgt.hist.root).
1932 // Weights are stored in a static TH1D* data member, ready to be plugged into the AliFlowEvent.
1933 // The plugging is done by the method ::FillEvtPhiWgt() (if wgt file is there).
1934
1935 fPhiWgtFile = wgtFile ;
1936
1937 TString* histTitle ;
327288af 1938 TH1D* hTPCall ; TH1D* hTPCplus ; TH1D* hTPCminus ; TH1D* hTPCcross ;
1939 TH1D* hPIDbay ;
1940 for(int k=0;k<AliFlowConstants::kSels;k++)
9777bfcb 1941 {
327288af 1942 for(int j=0;j<AliFlowConstants::kHars;j++)
9777bfcb 1943 {
1944 // Tpc (plus)
1945 histTitle = new TString("Flow_Phi_Weight_TPCplus_Sel");
1946 *histTitle += k+1;
1947 histTitle->Append("_Har");
1948 *histTitle += j+1;
327288af 1949 hTPCplus = (TH1D*)fPhiWgtFile->Get(histTitle->Data());
9777bfcb 1950 delete histTitle;
1951 // Tpc (minus)
1952 histTitle = new TString("Flow_Phi_Weight_TPCminus_Sel");
1953 *histTitle += k+1;
1954 histTitle->Append("_Har");
1955 *histTitle += j+1;
327288af 1956 hTPCminus = (TH1D*)fPhiWgtFile->Get(histTitle->Data());
9777bfcb 1957 delete histTitle;
1958 // Tpc (cross)
1959 histTitle = new TString("Flow_Phi_Weight_TPCcross_Sel");
1960 *histTitle += k+1;
1961 histTitle->Append("_Har");
1962 *histTitle += j+1;
327288af 1963 hTPCcross = (TH1D*)fPhiWgtFile->Get(histTitle->Data());
9777bfcb 1964 delete histTitle;
1965
1966 // Tpc
1967 histTitle = new TString("Flow_Phi_Weight_TPC_Sel");
1968 *histTitle += k+1;
1969 histTitle->Append("_Har");
1970 *histTitle += j+1;
327288af 1971 hTPCall = (TH1D*)fPhiWgtFile->Get(histTitle->Data());
9777bfcb 1972 delete histTitle;
1973
1974 for(int n=0;n<fPhiBins;n++)
1975 {
327288af 1976 fPhiWgtPlus[k][j][n] = hTPCplus->GetBinContent(n+1) ;
1977 fPhiWgtMinus[k][j][n] = hTPCminus->GetBinContent(n+1) ;
1978 fPhiWgtCross[k][j][n] = hTPCcross->GetBinContent(n+1) ;
1979 fPhiWgt[k][j][n] = hTPCall->GetBinContent(n+1) ;
9777bfcb 1980 // cout << " Weights: " << fPhiWgt[k][j][n] << " ; " << fPhiWgtPlus[k][j][n] << " | " << fPhiWgtMinus[k][j][n] << " | " << fPhiWgtCross[k][j][n] << endl ;
1981 }
1982 }
1983
1984 // Bayesian weights
1985 histTitle = new TString("Flow_BayPidMult_Sel");
1986 *histTitle += k+1;
327288af 1987 hPIDbay = (TH1D*)fPhiWgtFile->Get(histTitle->Data());
9777bfcb 1988 delete histTitle;
327288af 1989 Double_t totCount = hPIDbay->GetSumOfWeights() ;
1990 for (int n=0;n<AliFlowConstants::kPid;n++)
9777bfcb 1991 {
327288af 1992 if(totCount) { fBayesianWgt[k][n] = hPIDbay->GetBinContent(n+1) / totCount ; }
9777bfcb 1993 else { fBayesianWgt[k][n] = 1. ; }
1994 // cout << " Bayesian Weights (" << n << ") : " << fBayesianWgt[k][n] << endl ;
1995 }
1996 }
1997
327288af 1998 delete hTPCall ; delete hTPCplus ; delete hTPCminus ; delete hTPCcross ;
1999 delete hPIDbay ;
9777bfcb 2000
2001 return ;
2002}
2003//-----------------------------------------------------------------------
2004void AliFlowAnalyser::FillEvtPhiWgt(AliFlowEvent* fFlowEvent)
2005{
2006 // Plugs phi weights into the static dwgt data member of the AliFlowEvent class.
327288af 2007 // Weights are given in special AliFlowConstants::PhiWgt_t arrays (see AliFlowConstants),
9777bfcb 2008 // which are read from the wgt histograms by the method FillWgtArrays(...).
2009
2010 fFlowEvent->SetPhiWeight(fPhiWgt);
2011 fFlowEvent->SetPhiWeightPlus(fPhiWgtPlus);
2012 fFlowEvent->SetPhiWeightMinus(fPhiWgtMinus);
2013 fFlowEvent->SetPhiWeightCross(fPhiWgtCross);
2014
2015 return ;
2016}
2017//-----------------------------------------------------------------------
2018void AliFlowAnalyser::FillBayesianWgt(AliFlowEvent* fFlowEvent)
2019{
2020 // Plugs Bayesian particle abundance into the current AliFlowEvent.
2021 // A bayesian vector should be used for the PId of any different selection
2022 // (different sets of cuts -> different particle abundance), but for now
2023 // just the Selection n.0 (with no cuts) is used .
2024 // (AliFlowEvent::mBayesianCs[6] is a 1-dimensional array, change that first!).
2025
327288af 2026 Double_t bayes[AliFlowConstants::kPid] ;
9777bfcb 2027 Double_t bayCheck = 0. ;
327288af 2028 for (int n=0;n<AliFlowConstants::kPid;n++)
9777bfcb 2029 {
2030 bayes[n] = fBayesianWgt[0][n] ;
2031 bayCheck += bayes[n] ;
2032 // cout << "Bayesian V[" << n << "] = " << fBayesianWgt[0][n] << endl ;
2033 }
2034 if(bayCheck) { fFlowEvent->SetBayesian(bayes) ; fRePid = kTRUE ; }
2035 else { cout << "An empty bayesian vector is stored !!! - Bayesian weights = {1,1,1,1,1,1} " << endl ; }
2036
2037 return ;
2038}
2039//-----------------------------------------------------------------------
2040void AliFlowAnalyser::Weightening()
2041{
2042 // Calculates weights, and fills PhiWgt histograms .
2043 // This is called at the end of the event analysis.
2044
2045 cout << " AliFlowAnalyser::Weightening() " << endl ; cout << endl ;
2046
2047 // PhiWgt histogram collection
327288af 2048 fPhiWgtHistList = new TOrdCollection(4*AliFlowConstants::kSels*AliFlowConstants::kHars) ;
9777bfcb 2049
2050 // Creates PhiWgt Histograms
2051 TString* histTitle ;
327288af 2052 for(int k = 0; k < AliFlowConstants::kSels; k++)
9777bfcb 2053 {
327288af 2054 for(int j = 0; j < AliFlowConstants::kHars; j++)
9777bfcb 2055 {
2056 // Tpc (plus)
2057 histTitle = new TString("Flow_Phi_Weight_TPCplus_Sel");
2058 *histTitle += k+1;
2059 histTitle->Append("_Har");
2060 *histTitle += j+1;
2061 fHistFull[k].fHistFullHar[j].fHistPhiWgtPlus = new TH1D(histTitle->Data(),histTitle->Data(), fPhiBins, fPhiMin, fPhiMax);
2062 fHistFull[k].fHistFullHar[j].fHistPhiWgtPlus->Sumw2();
2063 fHistFull[k].fHistFullHar[j].fHistPhiWgtPlus->SetXTitle("Azimuthal Angles (rad)");
2064 fHistFull[k].fHistFullHar[j].fHistPhiWgtPlus->SetYTitle("PhiWgt");
2065 delete histTitle;
2066 // Tpc (minus)
2067 histTitle = new TString("Flow_Phi_Weight_TPCminus_Sel");
2068 *histTitle += k+1;
2069 histTitle->Append("_Har");
2070 *histTitle += j+1;
2071 fHistFull[k].fHistFullHar[j].fHistPhiWgtMinus = new TH1D(histTitle->Data(),histTitle->Data(), fPhiBins, fPhiMin, fPhiMax);
2072 fHistFull[k].fHistFullHar[j].fHistPhiWgtMinus->Sumw2();
2073 fHistFull[k].fHistFullHar[j].fHistPhiWgtMinus->SetXTitle("Azimuthal Angles (rad)");
2074 fHistFull[k].fHistFullHar[j].fHistPhiWgtMinus->SetYTitle("PhiWgt");
2075 delete histTitle;
2076 // Tpc (cross)
2077 histTitle = new TString("Flow_Phi_Weight_TPCcross_Sel");
2078 *histTitle += k+1;
2079 histTitle->Append("_Har");
2080 *histTitle += j+1;
2081 fHistFull[k].fHistFullHar[j].fHistPhiWgtAll = new TH1D(histTitle->Data(),histTitle->Data(), fPhiBins, fPhiMin, fPhiMax);
2082 fHistFull[k].fHistFullHar[j].fHistPhiWgtAll->Sumw2();
2083 fHistFull[k].fHistFullHar[j].fHistPhiWgtAll->SetXTitle("Azimuthal Angles (rad)");
2084 fHistFull[k].fHistFullHar[j].fHistPhiWgtAll->SetYTitle("PhiWgt");
2085 delete histTitle;
2086 // Tpc
2087 histTitle = new TString("Flow_Phi_Weight_TPC_Sel");
2088 *histTitle += k+1;
2089 histTitle->Append("_Har");
2090 *histTitle += j+1;
2091 fHistFull[k].fHistFullHar[j].fHistPhiWgt = new TH1D(histTitle->Data(),histTitle->Data(), fPhiBins, fPhiMin, fPhiMax);
2092 fHistFull[k].fHistFullHar[j].fHistPhiWgt->Sumw2();
2093 fHistFull[k].fHistFullHar[j].fHistPhiWgt->SetXTitle("Azimuthal Angles (rad)");
2094 fHistFull[k].fHistFullHar[j].fHistPhiWgt->SetYTitle("PhiWgt");
2095 delete histTitle;
2096
2097 // Calculate PhiWgt
2098 double meanPlus = fHistFull[k].fHistFullHar[j].fHistPhiPlus->Integral() / (double)fPhiBins ;
2099 double meanMinus = fHistFull[k].fHistFullHar[j].fHistPhiMinus->Integral() / (double)fPhiBins ;
2100 double meanCross = fHistFull[k].fHistFullHar[j].fHistPhiAll->Integral() / (double)fPhiBins ;
2101 double meanTPC = fHistFull[k].fHistFullHar[j].fHistPhi->Integral() / (double)fPhiBins ;
2102
2103 // Tpc
2104 for (int i=0;i<fPhiBins;i++)
2105 {
2106 fHistFull[k].fHistFullHar[j].fHistPhiWgtPlus->SetBinContent(i+1,meanPlus);
2107 fHistFull[k].fHistFullHar[j].fHistPhiWgtPlus->SetBinError(i+1, 0.);
2108 fHistFull[k].fHistFullHar[j].fHistPhiWgtMinus->SetBinContent(i+1,meanMinus);
2109 fHistFull[k].fHistFullHar[j].fHistPhiWgtMinus->SetBinError(i+1, 0.);
2110 fHistFull[k].fHistFullHar[j].fHistPhiWgtAll->SetBinContent(i+1,meanCross);
2111 fHistFull[k].fHistFullHar[j].fHistPhiWgtAll->SetBinError(i+1, 0.);
2112 fHistFull[k].fHistFullHar[j].fHistPhiWgt->SetBinContent(i+1,meanTPC);
2113 fHistFull[k].fHistFullHar[j].fHistPhiWgt->SetBinError(i+1, 0.);
2114 }
2115
2116 if(meanTPC==0) { cout << " Sel." << k << " , Har." << j << " : empty phi histogram ! " << endl ; }
2117 else
2118 {
2119 fHistFull[k].fHistFullHar[j].fHistPhiWgtPlus->Divide(fHistFull[k].fHistFullHar[j].fHistPhiPlus);
2120 fHistFull[k].fHistFullHar[j].fHistPhiWgtMinus->Divide(fHistFull[k].fHistFullHar[j].fHistPhiMinus);
2121 fHistFull[k].fHistFullHar[j].fHistPhiWgtAll->Divide(fHistFull[k].fHistFullHar[j].fHistPhiAll);
2122 fHistFull[k].fHistFullHar[j].fHistPhiWgt->Divide(fHistFull[k].fHistFullHar[j].fHistPhi);
2123 }
2124
2125 fPhiWgtHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistPhiWgtPlus);
2126 fPhiWgtHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistPhiWgtMinus);
2127 fPhiWgtHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistPhiWgtAll);
2128 fPhiWgtHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistPhiWgt);
2129 }
2130 }
2131
2132 return ;
2133}
2134//-----------------------------------------------------------------------
2135// ###
2136//-----------------------------------------------------------------------
4e566f2f 2137Bool_t AliFlowAnalyser::Analyze(AliFlowEvent* flowEvent)
9777bfcb 2138{
2139 // Runs the analysis on the AliFlowEvent (* fFlowEvent).
2140 // This method can be inserted in a loop over a collection of
2141 // AliFlowEvents or for on-fly analysis on AliESDs if used toghether
2142 // with the AliFlowMaker.
2143
2144 cout << " AliFlowAnalyser::Analyze(" << fFlowEvent << " ) - " << fEventNumber << endl ;
2145 if(!flowEvent) { return kFALSE ; }
2146 else { fFlowEvent = flowEvent ; }
2147
2148 if(fFlowSelect->Select(fFlowEvent)) // event selected - here below the ANALYSIS FLAGS are setted -
2149 {
2150 cout << " * 1 . Load event (track & v0s) and set flags . " << endl ;
2151 fFlowTracks = fFlowEvent->TrackCollection() ;
2152 fNumberOfTracks = fFlowTracks->GetEntries() ;
2153 fFlowV0s = fFlowEvent->V0Collection() ;
2154 fNumberOfV0s = fFlowV0s->GetEntries() ;
2155 cout << " event ID = " << fFlowEvent->EventID() << " : found " << fNumberOfTracks << " AliFlowTracks, and " << fNumberOfV0s << " AliFlowV0s . " << endl ;
2156
2157 if(fReadPhiWgt) { FillEvtPhiWgt(fFlowEvent) ; } // phi and bayesian weights are filled previous to the loop (FillWgtArrays(TFile*))
2158 else { fFlowEvent->SetNoWgt() ; } // phi weights can be used or not , this plugs them into the event
2159
4e566f2f 2160 fFlowEvent->SetCustomRespFunc(fCustomRespFunc) ; // a custom "detector response function" is used for assigning P.Id
2161
9777bfcb 2162 if(fBayWgt) { FillBayesianWgt(fFlowEvent) ; } // bayesian weights can be used or not , this plugs them into the event
2163 if(fRePid) { fFlowEvent->SetPids() ; } // re-calculate all p.id. hypotesis with the (new) bayesian array
2164
2165 if(fOnePhiWgt) { fFlowEvent->SetOnePhiWgt() ; } // one phi-wgt histogram
2166 else { fFlowEvent->SetFirstLastPhiWgt() ; } // three phi-wgt histogram
2167 if(fPtWgt) { fFlowEvent->SetPtWgt(); ; } // pT as a weight
2168 if(fEtaWgt) { fFlowEvent->SetEtaWgt() ; } // eta as a weight
2169
2170 if(fShuffle) { fFlowEvent->RandomShuffle() ; } // tracks re-shuffling
2171
2172 fFlowEvent->SetSelections(fFlowSelect) ; // does the selection of tracks for r.p. calculation (sets flags in AliFlowTrack)
da5aa0a0 2173
2174 if(fSub == 1) { fFlowEvent->SetEtaSubs() ; } // setting for the subevents (eta, random, charged)
2175 else if(fSub == -1) { fFlowEvent->SetChrSubs() ; } //
2176 else { fFlowEvent->SetRndSubs() ; } // ...
2177
2178 fEtaSub = fFlowEvent->EtaSubs() ;
2179
9777bfcb 2180 fFlowEvent->MakeSubEvents() ; // makes the subevent, eta or random basing on the previous flag
2181
2182 cout << " * 2 . Calculating event quantities all in one shoot . " << endl ;
2183 fFlowEvent->MakeAll() ;
2184
2185 if(FillFromFlowEvent(fFlowEvent)) // calculates event quantities
2186 {
2187 cout << " * 3 . Event Histograms and Particles loop . " << endl ;
2188 FillEventHistograms(fFlowEvent); // fill histograms from AliFlowEvents
2189 if(fTrackLoop) { FillParticleHistograms(fFlowTracks) ; } // fill histograms from AliFlowTracks
2190 if(fV0loop) { FillV0Histograms(fFlowV0s) ; } // fill histograms from AliFlowV0s
2191 //FillLabels() ; // fill the histogram of MC labels (from the simulation)
2192 }
2193 else
2194 {
2195 cout << " * 3 . Event psi = 0 - Skipping! " << endl ;
2196 return kFALSE ;
2197 }
2198 }
2199 else
2200 {
2201 cout << " * 0 . Event " << fEventNumber << " (event ID = " << fFlowEvent->EventID() << ") discarded . " << endl ;
92016a03 2202 // delete fFlowEvent ; fFlowEvent = 0 ;
9777bfcb 2203 return kFALSE ;
2204 }
2205 fEventNumber++ ;
2206
2207 return kTRUE ;
2208}
2209//-----------------------------------------------------------------------
2210Bool_t AliFlowAnalyser::FillFromFlowEvent(AliFlowEvent* fFlowEvent)
2211{
2212 // gets event quantities, returns kFALSE if Q vector is always 0
2213
2214 Int_t selCheck = 0 ;
327288af 2215 for(Int_t k = 0; k < AliFlowConstants::kSels; k++)
9777bfcb 2216 {
2217 fFlowSelect->SetSelection(k) ;
327288af 2218 for(Int_t j = 0; j < AliFlowConstants::kHars; j++)
9777bfcb 2219 {
2220 fFlowSelect->SetHarmonic(j) ;
327288af 2221 for(Int_t n = 0; n < AliFlowConstants::kSubs; n++)
9777bfcb 2222 {
2223 fFlowSelect->SetSubevent(n) ;
2224 fPsiSub[n][k][j] = fFlowEvent->Psi(fFlowSelect) ; // sub-event quantities
2225 fMultSub[n][k][j] = fFlowEvent->Mult(fFlowSelect) ;
2226 }
2227 fFlowSelect->SetSubevent(-1);
2228 fQ[k][j] = fFlowEvent->Q(fFlowSelect) ; // full event quantities
2229 fPsi[k][j] = fFlowEvent->Psi(fFlowSelect) ;
2230 fQnorm[k][j] = fFlowEvent->NormQ(fFlowSelect).Mod() ; // was: fFlowEvent->q(fFlowSelect) ; // but the normalization was bad (no pT,eta weight)
2231 fMult[k][j] = fFlowEvent->Mult(fFlowSelect) ;
2232 selCheck += fMult[k][j] ;
2233 }
2234 }
2235 if(!selCheck) { return kFALSE ; } // if there are no particles in the selection -> skip the event
2236
2237 return kTRUE ;
2238}
2239//-----------------------------------------------------------------------
2240void AliFlowAnalyser::FillEventHistograms(AliFlowEvent* fFlowEvent)
2241{
2242 // fill event histograms
2243
2244 cout << " Fill Event Histograms ... " << endl ;
2245
2246 float trigger = (float)fFlowEvent->L0TriggerWord() ;
2247 fHistTrigger->Fill(trigger);
2248
2249 // no selections: OrigMult, Centrality, Mult, MultOverOrig, VertexZ, VertexXY
2250 int origMult = fFlowEvent->OrigMult();
2251 fHistOrigMult->Fill((float)origMult);
2252 fHistMultEta->Fill((float)fFlowEvent->MultEta());
2253
2254 int cent = fFlowEvent->Centrality();
2255 fHistCent->Fill((float)cent);
2256
2257 fHistMult->Fill((float)fNumberOfTracks) ;
2258 fHistV0Mult->Fill((float)fNumberOfV0s) ;
2259 if(origMult) { fHistMultOverOrig->Fill((float)fNumberOfTracks/(float)origMult) ; }
2260
2261 fFlowEvent->VertexPos(fVertex) ;
2262 fHistVertexZ->Fill(fVertex[2]) ;
2263 fHistVertexXY2D->Fill(fVertex[0],fVertex[1]) ;
2264
2265 // ZDC info
2266 fHistPartZDC->Fill(fFlowEvent->ZDCpart()) ;
2267 for(int ii=0;ii<3;ii++) { fHistEnergyZDC->Fill(ii,fFlowEvent->ZDCenergy(ii)) ; }
2268
2269 // sub-event Psi_Subs
327288af 2270 for(int k = 0; k < AliFlowConstants::kSels; k++)
9777bfcb 2271 {
327288af 2272 for(int j = 0; j < AliFlowConstants::kHars; j++)
9777bfcb 2273 {
327288af 2274 for(int n = 0; n < AliFlowConstants::kSubs; n++)
9777bfcb 2275 {
327288af 2276 int iii = AliFlowConstants::kSubs * k + n ; //cout << " " << k << j << n << " , " << iii << endl ;
9777bfcb 2277 fHistSub[iii].fHistSubHar[j].fHistPsiSubs->Fill(fPsiSub[n][k][j]) ;
2278 }
2279 }
2280 }
2281
2282 // full event Psi, PsiSubCorr, PsiSubCorrDiff, cos, mult, q
327288af 2283 for(int k = 0; k < AliFlowConstants::kSels; k++)
9777bfcb 2284 {
327288af 2285 for(int j = 0; j < AliFlowConstants::kHars; j++)
9777bfcb 2286 {
2287 float order = (float)(j+1);
2288 fHistFull[k].fHistFullHar[j].fHistPsi->Fill(fPsi[k][j]);
2289 if(k<2 && j<2)
2290 {
2291 if(k==0)
2292 {
2293 float psi1 = fPsi[0][j] ;
2294 float psi2 = fPsi[1][j] ;
2295 float diff = psi1 - psi2 ;
2296 if(diff < -TMath::Pi()/(j+1)) { diff += 2*TMath::Pi()/(j+1) ; }
2297 else if(diff > +TMath::Pi()/(j+1)) { diff -= 2*TMath::Pi()/(j+1) ; }
2298 fHistFull[k].fHistFullHar[j].fHistPsiDiff->Fill(diff) ; // k=0
2299 }
2300 else if(k==1)
2301 {
51a9b7d8 2302 float psi1 = 0. ; float psi2 = 0. ;
9777bfcb 2303 if (j==0) { psi1 = fPsi[0][0] ; psi2 = fPsi[1][1] ; }
2304 else if(j==1) { psi1 = fPsi[0][0] ; psi2 = fPsi[0][1] ; }
2305 float diff = psi1 - psi2 ;
2306 diff = (TMath::Abs(diff) > TMath::Pi()) ? ((diff > 0.) ? -(2*TMath::Pi()-diff) : -(diff+2*TMath::Pi())) : diff ;
2307 fHistFull[k].fHistFullHar[j].fHistPsiDiff->Fill(diff) ; // k=1
2308 }
2309 }
2310
2311 if(fPsiSub[0][k][j] != 0. && fPsiSub[1][k][j] != 0.)
2312 {
2313 float psiSubCorr; // this is: delta_Psi
2314 if(fV1Ep1Ep2 == kFALSE || order != 1)
2315 {
2316 psiSubCorr = fPsiSub[0][k][j] - fPsiSub[1][k][j];
2317 }
2318 else // i.e. (fV1Ep1Ep2 == kTRUE && order == 1)
2319 {
2320 psiSubCorr = fPsiSub[0][k][0] + fPsiSub[1][k][0] - 2*fPsi[k][1];
2321 }
2322 fHistFull[k].fHistCos->Fill(order, (float)cos(order * psiSubCorr)) ;
2323 if(psiSubCorr < 0.) { psiSubCorr += 2*TMath::Pi()/order ; }
2324 if(psiSubCorr > 2*TMath::Pi()/order) { psiSubCorr -= 2*TMath::Pi()/order ; } // for v1Ep1Ep2 which gives -2*TMath::Pi() < psiSubCorr < 2*2*TMath::Pi()
2325 fHistFull[k].fHistFullHar[j].fHistPsiSubCorr->Fill(psiSubCorr);
2326 }
2327
327288af 2328 if(j < AliFlowConstants::kHars - 1) // subevents of different harmonics
9777bfcb 2329 {
51a9b7d8 2330 int j1 = 0 ; int j2 = 0 ;
9777bfcb 2331 float psiSubCorrDiff;
2332 if(j==0) { j1 = 1, j2 = 2 ; }
2333 else if(j==1) { j1 = 1, j2 = 3 ; }
2334 else if(j==2) { j1 = 2, j2 = 4 ; }
2335 psiSubCorrDiff = fmod((double)fPsiSub[0][k][j1-1],2*TMath::Pi()/(double)j2)-fmod((double)fPsiSub[1][k][j2-1],2*TMath::Pi()/(double)j2) ;
2336 if(psiSubCorrDiff < 0.) { psiSubCorrDiff += 2*TMath::Pi()/(float)j2 ; }
2337 fHistFull[k].fHistFullHar[j].fHistPsiSubCorrDiff->Fill(psiSubCorrDiff) ;
2338 psiSubCorrDiff = fmod((double)fPsiSub[0][k][j2-1],2*TMath::Pi()/(double)j2)-fmod((double)fPsiSub[1][k][j1-1],2*TMath::Pi()/(double)j2) ;
2339 if(psiSubCorrDiff < 0.) { psiSubCorrDiff += 2*TMath::Pi()/(float)j2 ; }
2340 fHistFull[k].fHistFullHar[j].fHistPsiSubCorrDiff->Fill(psiSubCorrDiff) ;
2341 }
2342
2343 fHistFull[k].fHistFullHar[j].fHistMult->Fill((float)fMult[k][j]) ;
2344 fHistFull[k].fHistFullHar[j].fHistQnorm->Fill(fQnorm[k][j]) ;
2345 }
da5aa0a0 2346 }
2347 fTotalNumberOfEvents++ ;
9777bfcb 2348
2349 return ;
2350}
2351//-----------------------------------------------------------------------
92016a03 2352void AliFlowAnalyser::FillParticleHistograms(TClonesArray* fFlowTracks)
9777bfcb 2353{
2354 // fills tracks histograms
2355
2356 cout << " Tracks Loop . " << endl ;
2357
2358 float corrMultUnit = 0. ;
2359 float corrMultN = 0. ;
2360 float etaSymPosTpcN = 0. ;
2361 float etaSymNegTpcN = 0. ;
2362 float etaSymPosTpcNpart = 0. ;
2363 float etaSymNegTpcNpart = 0. ;
2364 float hPlusN = 0. ;
2365 float hMinusN = 0. ;
2366 float piPlusN = 0. ;
2367 float piMinusN = 0. ;
2368 float protonN = 0. ;
2369 float pbarN = 0. ;
2370 float kMinusN = 0. ;
2371 float kPlusN = 0. ;
2372 float deuteronN = 0. ;
2373 float dbarN = 0. ;
2374 float electronN = 0. ;
2375 float positronN = 0. ;
2376 float muonMinusN = 0. ;
2377 float muonPlusN = 0. ;
2378
92016a03 2379 fSelParts = 0 ;
2380 fSelV0s = 0 ;
9777bfcb 2381 for(fTrackNumber=0;fTrackNumber<fNumberOfTracks;fTrackNumber++)
2382 {
92016a03 2383 //fFlowTrack = (AliFlowTrack*)fFlowTracks->At(fTrackNumber) ;
2384
2385 fFlowTrack = (AliFlowTrack*)(fFlowTracks->AddrAt(fTrackNumber)) ;
2386 //cout << "Track n. " << fTrackNumber << endl ; // fFlowTrack->Dump() ;
9777bfcb 2387
2388 bool constrainable = fFlowTrack->IsConstrainable() ;
2389 // int label = fFlowTrack->Label() ;
2390 float dcaGlobal = TMath::Abs(fFlowTrack->Dca()) ;
2391 float dcaSigned = fFlowTrack->TransDcaSigned() ;
2392 float dcaTrans = TMath::Abs(dcaSigned) ;
2393 float eta = fFlowTrack->Eta() ;
2394 float phi = fFlowTrack->Phi() ;
2395 float pt = fFlowTrack->Pt() ;
2396 float etaGlob = fFlowTrack->EtaGlobal() ;
2397 float phiGlob = fFlowTrack->PhiGlobal() ;
2398 float ptGlob = fFlowTrack->PtGlobal() ;
2399 float totalp = fFlowTrack->P() ;
2400 // float logp = TMath::Log10(totalp) ;
2401 // float zFirstPoint = fFlowTrack->ZFirstPoint() ;
2402 // float zLastPoint = fFlowTrack->ZLastPoint() ;
2403 float lenght = fFlowTrack->TrackLength() ;
2404 int charge = fFlowTrack->Charge() ;
2405 float chi2 = fFlowTrack->Chi2() ;
2406 int fitPtsTPC = (int)((float)fFlowTrack->FitPtsTPC()) ;
2407 int maxPtsTPC = fFlowTrack->MaxPtsTPC() ;
2408 float chi2TPC = fFlowTrack->Chi2TPC() ;
2409 int fitPtsITS = fFlowTrack->FitPtsITS() ;
2410 int maxPtsITS = fFlowTrack->MaxPtsITS() ;
2411 float chi2ITS = fFlowTrack->Chi2ITS() ;
2412 int fitPtsTRD = fFlowTrack->NhitsTRD() ;
2413 int maxPtsTRD = fFlowTrack->MaxPtsTRD() ;
2414 float chi2TRD = fFlowTrack->Chi2TRD() ;
2415 int fitPtsTOF = fFlowTrack->NhitsTOF() ;
2416 int maxPtsTOF = fFlowTrack->MaxPtsTOF() ;
2417 float chi2TOF = fFlowTrack->Chi2TOF() ;
2418 float dEdx = fFlowTrack->DedxTPC() ;
2419 float its = fFlowTrack->DedxITS() ;
2420 float trd = fFlowTrack->SigTRD() ;
2421 float tof = fFlowTrack->TofTOF() ;
2422 float lpTPC = 0 ; if(fFlowTrack->PatTPC()>0) { lpTPC = TMath::Log10(fFlowTrack->PatTPC()) ; }
2423 float lpITS = 0 ; if(fFlowTrack->PatITS()>0) { lpITS = TMath::Log10(fFlowTrack->PatITS()) ; }
2424 float lpTRD = 0 ; if(fFlowTrack->PatTRD()>0) { lpTRD = TMath::Log10(fFlowTrack->PatTRD()) ; }
2425 float lpTOF = 0 ; if(fFlowTrack->PatTOF()>0) { lpTOF = TMath::Log10(fFlowTrack->PatTOF()) ; }
2426 float invMass = fFlowTrack->InvMass() ;
2427 Char_t pid[10]="0" ; strcpy(pid,fFlowTrack->Pid()) ;
2428 fPidId = -1 ; // assigned later
2429
2430 // no selections: Charge, Dca, Chi2, FitPts, MaxPts, FitOverMax, PID
2431 fHistCharge->Fill((float)charge);
2432 fHistDcaGlobal->Fill(dcaGlobal);
2433 fHistDca->Fill(dcaTrans) ;
2434 fHistTransDca->Fill(dcaSigned);
2435 fHistChi2->Fill(chi2);
2436
9777bfcb 2437 // - here TPC (chi2 & nHits)
2438 fHistChi2TPC->Fill(chi2TPC);
92016a03 2439 if(fitPtsTPC>0) { fHistChi2normTPC->Fill(chi2TPC/((float)fitPtsTPC)) ; }
9777bfcb 2440 fHistFitPtsTPC->Fill((float)fitPtsTPC);
2441 fHistMaxPtsTPC->Fill((float)maxPtsTPC);
92016a03 2442 if(maxPtsTPC>0) { fHistFitOverMaxTPC->Fill((float)(fitPtsTPC)/(float)maxPtsTPC) ; }
2443
2444 // - here ITS (chi2 & nHits)
2445 fHistChi2ITS->Fill(chi2ITS);
2446 if(fitPtsITS>0) { fHistChi2normITS->Fill(chi2ITS/((float)fitPtsITS)) ; }
2447 fHistFitPtsITS->Fill((float)fitPtsITS);
2448 fHistMaxPtsITS->Fill((float)maxPtsITS);
9777bfcb 2449
92016a03 2450 // - here TRD (chi2 & nHits)
9777bfcb 2451 fHistChi2TRD->Fill(chi2TRD);
92016a03 2452 if(fitPtsTRD>0) { fHistChi2normTRD->Fill(chi2TRD/((float)fitPtsTRD)) ; }
9777bfcb 2453 fHistFitPtsTRD->Fill((float)fitPtsTRD);
2454 fHistMaxPtsTRD->Fill((float)maxPtsTRD);
2455
2456 // - here TOF (chi2 & nHits)
2457 fHistChi2TOF->Fill(chi2TOF);
92016a03 2458 if(fitPtsTOF>0) { fHistChi2normTOF->Fill(chi2TOF/((float)fitPtsTOF)) ; }
9777bfcb 2459 fHistFitPtsTOF->Fill((float)fitPtsTOF);
2460 fHistMaxPtsTOF->Fill((float)maxPtsTOF);
2461
2462 // fit over max (all)
2463 int maxPts = maxPtsITS + maxPtsTPC + maxPtsTRD + maxPtsTOF ;
2464 int fitPts = fitPtsITS + fitPtsTPC + fitPtsTRD + fitPtsTOF ;
92016a03 2465 if(maxPts>0) { fHistFitOverMax->Fill((float)(fitPts)/(float)maxPts) ; }
9777bfcb 2466
2467 // lenght
2468 fHistLenght->Fill(lenght) ;
2469 fHistInvMass->Fill(invMass) ;
2470
2471 // PID histograms & multiplicity count (for bayesian histogram)
2472 if(charge == 1)
2473 {
2474 hPlusN++ ;
2475 fHistMeanDedxPos2D->Fill(lpTPC, dEdx) ;
2476 fHistMeanDedxPos2DITS->Fill(lpITS, its) ;
2477 fHistMeanDedxPos2DTRD->Fill(lpTRD, trd) ;
2478 fHistMeanDedxPos2DTOF->Fill(lpTOF, tof) ;
2479 //
2480 float positron = fFlowTrack->ElectronPositronProb() ;
2481 fHistPidPositron->Fill(positron) ;
2482 if(strcmp(pid, "e+") == 0)
2483 {
2484 fPidId = 0 ; positronN++ ; fHistPidPt->Fill(1.5,pt) ;
2485 fHistMeanTPCPositron->Fill(lpTPC, dEdx) ;
2486 fHistMeanITSPositron->Fill(lpITS, its);
2487 fHistMeanTRDPositron->Fill(lpTRD, trd);
2488 fHistMeanTOFPositron->Fill(invMass, tof);
2489 fHistPidPositronPart->Fill(positron) ;
2490 }
2491 float muonPlus = fFlowTrack->MuonPlusMinusProb() ;
2492 fHistPidMuonPlus->Fill(muonPlus) ;
2493 if(strcmp(pid, "mu+") == 0)
2494 {
2495 fPidId = 1 ; muonPlusN++ ; fHistPidPt->Fill(3.5,pt) ;
2496 fHistMeanTPCMuonPlus->Fill(lpTPC, dEdx) ;
2497 fHistMeanITSMuonPlus->Fill(lpITS, its);
2498 fHistMeanTRDMuonPlus->Fill(lpTRD, trd);
2499 fHistMeanTOFMuonPlus->Fill(invMass, tof);
2500 fHistPidMuonPlusPart->Fill(muonPlus) ;
2501 }
2502 float piPlus = fFlowTrack->PionPlusMinusProb() ;
2503 fHistPidPiPlus->Fill(piPlus) ;
2504 if(strcmp(pid, "pi+") == 0)
2505 {
2506 fPidId = 2 ; piPlusN++ ; fHistPidPt->Fill(5.5,pt) ;
2507 fHistMeanTPCPiPlus->Fill(lpTPC, dEdx) ;
2508 fHistMeanITSPiPlus->Fill(lpITS, its);
2509 fHistMeanTRDPiPlus->Fill(lpTRD, trd);
2510 fHistMeanTOFPiPlus->Fill(invMass, tof);
2511 fHistPidPiPlusPart->Fill(piPlus) ;
2512 }
2513 float kplus = fFlowTrack->KaonPlusMinusProb() ;
2514 fHistPidKplus->Fill(kplus) ;
2515 if(strcmp(pid, "k+") == 0)
2516 {
2517 fPidId = 3 ; kPlusN++ ; fHistPidPt->Fill(7.5,pt) ;
2518 fHistMeanTPCKplus->Fill(lpTPC, dEdx) ;
2519 fHistMeanITSKplus->Fill(lpITS, its);
2520 fHistMeanTRDKplus->Fill(lpTRD, trd);
2521 fHistMeanTOFKplus->Fill(invMass, tof);
2522 fHistPidKplusPart->Fill(kplus) ;
2523 }
2524 float proton = fFlowTrack->ProtonPbarProb() ;
2525 fHistPidProton->Fill(proton) ;
2526 if(strcmp(pid, "pr+") == 0)
2527 {
2528 fPidId = 4 ; protonN++ ; fHistPidPt->Fill(9.5,pt) ;
2529 fHistMeanTPCProton->Fill(lpTPC, dEdx) ;
2530 fHistMeanITSProton->Fill(lpITS, its);
2531 fHistMeanTRDProton->Fill(lpTRD, trd);
2532 fHistMeanTOFProton->Fill(invMass, tof);
2533 fHistPidProtonPart->Fill(proton) ;
2534 }
2535 float deuteron = fFlowTrack->DeuteriumAntiDeuteriumProb() ;
2536 fHistPidDeuteron->Fill(deuteron) ;
2537 if(strcmp(pid, "d+") == 0)
2538 {
2539 fPidId = 6 ; deuteronN++ ; fHistPidPt->Fill(11.5,pt) ;
2540 fHistMeanTPCDeuteron->Fill(lpTPC, dEdx) ;
2541 fHistMeanITSDeuteron->Fill(lpITS, its);
2542 fHistMeanTRDDeuteron->Fill(lpTRD, trd);
2543 fHistMeanTOFDeuteron->Fill(invMass, tof);
2544 fHistPidDeuteronPart->Fill(deuteron) ;
2545 }
92016a03 2546 }
9777bfcb 2547 else if(charge == -1)
2548 {
2549 hMinusN++ ;
2550 fHistMeanDedxNeg2D->Fill(lpTPC, dEdx) ;
2551 fHistMeanDedxNeg2DITS->Fill(lpITS, its) ;
2552 fHistMeanDedxNeg2DTRD->Fill(lpTRD, trd) ;
2553 fHistMeanDedxNeg2DTOF->Fill(lpTOF, tof) ;
2554 //
2555 float electron = fFlowTrack->ElectronPositronProb() ;
2556 fHistPidElectron->Fill(electron);
2557 if(strcmp(pid, "e-") == 0)
2558 {
2559 fPidId = 0 ; electronN++ ; fHistPidPt->Fill(0.5,pt) ;
2560 fHistMeanTPCElectron->Fill(lpTPC, dEdx);
2561 fHistMeanITSElectron->Fill(lpITS, its);
2562 fHistMeanTRDElectron->Fill(lpTRD, trd);
2563 fHistMeanTOFElectron->Fill(invMass, tof);
2564 fHistPidElectronPart->Fill(electron);
2565 }
2566 float muonMinus = fFlowTrack->MuonPlusMinusProb() ;
2567 fHistPidMuonMinus->Fill(muonMinus) ;
2568 if(strcmp(pid, "mu-") == 0)
2569 {
2570 fPidId = 1 ; muonMinusN++ ; fHistPidPt->Fill(2.5,pt) ;
2571 fHistMeanTPCMuonMinus->Fill(lpTPC, dEdx) ;
2572 fHistMeanITSMuonMinus->Fill(lpITS, its);
2573 fHistMeanTRDMuonMinus->Fill(lpTRD, trd);
2574 fHistMeanTOFMuonMinus->Fill(invMass, tof);
2575 fHistPidMuonMinusPart->Fill(muonMinus) ;
2576 }
2577 float piMinus = fFlowTrack->PionPlusMinusProb() ;
2578 fHistPidPiMinus->Fill(piMinus) ;
2579 if(strcmp(pid, "pi-") == 0)
2580 {
2581 fPidId = 2 ; piMinusN++ ; fHistPidPt->Fill(4.5,pt) ;
2582 fHistMeanTPCPiMinus->Fill(lpTPC, dEdx);
2583 fHistMeanITSPiMinus->Fill(lpITS, its);
2584 fHistMeanTRDPiMinus->Fill(lpTRD, trd);
2585 fHistMeanTOFPiMinus->Fill(invMass, tof);
2586 fHistPidPiMinusPart->Fill(piMinus);
2587 }
2588 float kminus = fFlowTrack->KaonPlusMinusProb() ;
2589 fHistPidKminus->Fill(kminus);
2590 if(strcmp(pid, "k-") == 0)
2591 {
2592 fPidId = 3 ; kMinusN++ ; fHistPidPt->Fill(6.5,pt) ;
2593 fHistMeanTPCKminus->Fill(lpTPC, dEdx);
2594 fHistMeanITSKminus->Fill(lpITS, its);
2595 fHistMeanTRDKminus->Fill(lpTRD, trd);
2596 fHistMeanTOFKminus->Fill(invMass, tof);
2597 fHistPidKminusPart->Fill(kminus);
2598 }
2599 float antiproton = fFlowTrack->ProtonPbarProb() ;
2600 fHistPidAntiProton->Fill(antiproton);
2601 if(strcmp(pid, "pr-") == 0)
2602 {
2603 fPidId = 4 ; pbarN++ ; fHistPidPt->Fill(8.5,pt) ;
2604 fHistMeanTPCPbar->Fill(lpTPC, dEdx);
2605 fHistMeanITSPbar->Fill(lpITS, its);
2606 fHistMeanTRDPbar->Fill(lpTRD, trd);
2607 fHistMeanTOFPbar->Fill(invMass, tof);
2608 fHistPidAntiProtonPart->Fill(antiproton);
2609 }
2610 float antideuteron = fFlowTrack->DeuteriumAntiDeuteriumProb() ;
2611 fHistPidAntiDeuteron->Fill(antideuteron);
2612 if(strcmp(pid, "d-") == 0)
2613 {
2614 fPidId = 6 ; dbarN++ ; fHistPidPt->Fill(10.5,pt) ;
2615 fHistMeanTPCAntiDeuteron->Fill(lpTPC, dEdx);
2616 fHistMeanITSAntiDeuteron->Fill(lpITS, its);
2617 fHistMeanTRDAntiDeuteron->Fill(lpTRD, trd);
2618 fHistMeanTOFAntiDeuteron->Fill(invMass, tof);
2619 fHistPidAntiDeuteronPart->Fill(antideuteron);
2620 }
2621 }
2622
2623 // Yield3D, Yield2D, Eta, Pt, Phi, bayP.Id.
2624 fHistPtot->Fill(totalp) ;
2625 fHistPt->Fill(pt) ;
2626 fHistPhi->Fill(phi);
2627 fHistAllEtaPtPhi3D->Fill(eta, pt, phi) ;
2628 fHistYieldAll2D->Fill(eta, pt) ;
2629 fHistBayPidMult->Fill(fPidId) ;
2630 if(constrainable)
2631 {
2632 fHistPhiCons->Fill(phi);
2633 fHistPhiPtCon->Fill(phi, pt);
2634 fHistYieldCon2D->Fill(eta, pt) ;
2635 fHistConsEtaPtPhi3D->Fill(eta, pt, phi) ;
2636 fHistGlobEtaPtPhi3D->Fill(etaGlob, ptGlob, phiGlob) ;
2637 }
2638 else
2639 {
2640 fHistYieldUnc2D->Fill(etaGlob, ptGlob) ;
2641 fHistUncEtaPtPhi3D->Fill(etaGlob, ptGlob, phiGlob) ;
2642 fHistPhiPtUnc->Fill(phiGlob, ptGlob) ;
2643 }
2644 if(fFlowTrack->Charge()>0) { fHistPtPhiPos->Fill(phi, pt); }
2645 else if(fFlowTrack->Charge()<0) { fHistPtPhiNeg->Fill(phi, pt); }
2646
2647 // fills selected part histograms
2648 if(fFlowSelect->SelectPart(fFlowTrack))
2649 {
92016a03 2650 fSelParts++ ;
2651
9777bfcb 2652 if(strlen(fFlowSelect->PidPart()) != 0)
2653 {
2654 float rapidity = fFlowTrack->Y();
2655 fHistBinEta->Fill(rapidity, rapidity);
2656 fHistYieldPart2D->Fill(rapidity, pt);
2657 }
2658 else
2659 {
2660 fHistBinEta->Fill(eta, eta) ;
2661 fHistYieldPart2D->Fill(eta, pt) ;
2662 }
2663 fHistBayPidMultPart->Fill(fPidId) ;
2664 fHistBinPt->Fill(pt, pt) ;
2665 fHistEtaPtPhi3DPart->Fill(eta,pt,phi) ;
2666 fHistDcaGlobalPart->Fill(dcaGlobal) ;
2667 fHistInvMassPart->Fill(invMass) ;
2668 if(charge == 1)
2669 {
2670 fHistMeanDedxPos3DPart->Fill(lpITS, dEdx, fPidId) ;
2671 fHistMeanDedxPos3DPartITS->Fill(lpITS, its, fPidId) ;
2672 }
2673 else if(charge == -1)
2674 {
2675 fHistMeanDedxNeg3DPart->Fill(lpITS, dEdx, fPidId) ;
2676 fHistMeanDedxNeg3DPartITS->Fill(lpITS, its, fPidId) ;
2677 }
2678 if(eta > 0.) { etaSymPosTpcNpart++ ; } // for fHistEtaSymPart
2679 else { etaSymNegTpcNpart++ ; }
2680 }
2681 else
2682 {
2683 fHistEtaPtPhi3DOut->Fill(eta,pt,phi) ;
2684 fHistYieldOut2D->Fill(eta, pt) ;
2685 fHistDcaGlobalOut->Fill(dcaGlobal) ;
2686 fHistInvMassOut->Fill(invMass) ;
2687 }
2688
2689 // cos(n*phiLab)
327288af 2690 for(int j = 0; j < AliFlowConstants::kHars; j++)
9777bfcb 2691 {
2692 bool oddHar = (j+1) % 2 ;
2693 float order = (float)(j+1) ;
2694 float vIn = 100 * cos((double)order * phi) ;
2695 if(eta < 0 && oddHar) { vIn *= -1 ; }
2696 fHistCosPhi->Fill(order, vIn);
2697 }
2698
2699 //For Eta symmetry TPC
2700 if(fFlowTrack->FitPtsTPC())
2701 {
2702 if(eta > 0.) { etaSymPosTpcN++ ; } // for fHistEtaSym
2703 else { etaSymNegTpcN++ ; }
2704 }
2705
2706 // not to call it twice ...
2707 int nEtaS = HarmonicsLoop(fFlowTrack) ;
2708
2709 //For Correlated Multiplicity
2710 corrMultN += (float)nEtaS ;
2711
2712 //For Correlated Multiplicity in 1 unit rapidity
2713 if(TMath::Abs(eta) <= 0.5) { corrMultUnit += (float)nEtaS ; }
2714
2715 //delete pointer
2716 fFlowTrack = 0 ;
2717 } // end of tracks loop
2718
2719 // EtaSym
2720 float etaSymTpc = 0 ;
2721 if(etaSymPosTpcN || etaSymNegTpcN) { etaSymTpc = (etaSymPosTpcN - etaSymNegTpcN) / (etaSymPosTpcN + etaSymNegTpcN); }
2722 float etaSymTpcPart = 0 ;
2723 if(etaSymPosTpcNpart || etaSymNegTpcNpart) { etaSymTpcPart = (etaSymPosTpcNpart - etaSymNegTpcNpart) / (etaSymPosTpcNpart + etaSymNegTpcNpart) ; }
2724 Float_t vertexZ = fVertex[2] ;
2725 // all
2726 fHistEtaSym->Fill(etaSymTpc);
2727 fHistEtaSymVerZ2D->Fill(vertexZ , etaSymTpc);
2728 // selected
2729 fHistEtaSymPart->Fill(etaSymTpc);
2730 fHistEtaSymVerZ2DPart->Fill(vertexZ , etaSymTpcPart);
2731
2732 // PID multiplicities
2733 float totalMult = (float)fFlowTracks->GetEntries() ;
2734 fHistPidMult->Fill(1., totalMult);
2735 fHistPidMult->Fill(2., hPlusN);
2736 fHistPidMult->Fill(3., hMinusN);
2737 fHistPidMult->Fill(4., piPlusN);
2738 fHistPidMult->Fill(5., piMinusN);
2739 fHistPidMult->Fill(6., protonN);
2740 fHistPidMult->Fill(7., pbarN);
2741 fHistPidMult->Fill(8., kPlusN);
2742 fHistPidMult->Fill(9., kMinusN);
2743 fHistPidMult->Fill(10., deuteronN);
2744 fHistPidMult->Fill(11., dbarN);
2745 fHistPidMult->Fill(12., electronN);
2746 fHistPidMult->Fill(13., positronN);
2747 fHistPidMult->Fill(14., muonMinusN);
2748 fHistPidMult->Fill(15., muonPlusN);
2749
2750 // Multiplicity of particles correlated with the event planes
327288af 2751 corrMultN /= (float)(AliFlowConstants::kHars * AliFlowConstants::kSels) ;
92016a03 2752 if(fSelParts == corrMultN) { fHistMultPart->Fill(corrMultN) ; } // crosscheck
2753 else { fHistMultPart->Fill(-1) ; }
9777bfcb 2754 // ...in one unit rapidity
327288af 2755 corrMultUnit /= (float)(AliFlowConstants::kHars * AliFlowConstants::kSels) ;
9777bfcb 2756 fHistMultPartUnit->Fill(corrMultUnit) ;
2757
da5aa0a0 2758 fTotalNumberOfTracks += fNumberOfTracks ;
2759
9777bfcb 2760 return ;
2761}
2762//-----------------------------------------------------------------------
2763Int_t AliFlowAnalyser::HarmonicsLoop(AliFlowTrack* fFlowTrack)
2764{
2765 // HarmonicsLoop, does the correlation analysis, fills histograms
2766 // (internally called by ::FillParticleHistograms(..) loop) .
2767
2768 float eta = fFlowTrack->Eta() ;
2769 float phi = fFlowTrack->Phi() ;
2770 float pt = fFlowTrack->Pt() ;
2771 float zFirstPoint = fFlowTrack->ZFirstPoint() ;
2772 // float zLastPoint = fFlowTrack->ZLastPoint() ;
2773
2774 // Looping over Selections and Harmonics
2775 int corrMultN = 0 ;
327288af 2776 for (int k = 0; k < AliFlowConstants::kSels; k++)
9777bfcb 2777 {
2778 fFlowSelect->SetSelection(k) ;
327288af 2779 for (int j = 0; j < AliFlowConstants::kHars; j++)
9777bfcb 2780 {
2781 bool oddHar = (j+1) % 2;
2782 fFlowSelect->SetHarmonic(j);
2783 double order = (double)(j+1);
327288af 2784 float psii = 0. ; float psi2 = 0. ;
da5aa0a0 2785 if(fEtaSub) // particles with the opposite subevent
9777bfcb 2786 {
327288af 2787 if(eta > 0) { psii = fPsiSub[1][k][j] ; } //check
2788 else { psii = fPsiSub[0][k][j] ; }
9777bfcb 2789 }
2790 else if(order > 3. && !oddHar)
2791 {
327288af 2792 psii = fPsi[k][1]; // 2nd harmomic event plane
2793 if(psii > 2*TMath::Pi()/order) { psii -= 2*TMath::Pi()/order ; }
2794 if(psii > 2*TMath::Pi()/order) { psii -= 2*TMath::Pi()/order ; }
9777bfcb 2795 }
2796 else // random subevents
2797 {
327288af 2798 psii = fPsi[k][j] ;
9777bfcb 2799 }
2800
2801 if(fFlowSelect->Select(fFlowTrack)) // Get detID
2802 {
2803 Bool_t kTpcPlus = kFALSE ;
2804 Bool_t kTpcMinus = kFALSE ;
2805 Bool_t kTpcAll = kFALSE ;
2806
2807 fHistFull[k].fHistFullHar[j].fHistYieldPt->Fill(pt);
2808 fHistFull[k].fHistFullHar[j].fHistEtaPtPhi3D->Fill(eta, pt, phi);
2809 fHistFull[k].fHistFullHar[j].fHistYield2D->Fill(eta, pt);
2810 fHistFull[k].fHistFullHar[j].fHistDcaGlob->Fill(TMath::Abs(fFlowTrack->Dca()));
2811
2812 // Set Tpc (+ and -)
2813 if(fFlowTrack->FitPtsTPC()) //OR*: AliESDtrack:: "TBits& GetTPCClusterMap()" or "Int_t GetTPCclusters(Int_t* idx)" ...
2814 {
2815 if(zFirstPoint >= 0. && eta > 0.) { kTpcPlus = kTRUE ; }
2816 else if(zFirstPoint <= 0. && eta < 0.) { kTpcMinus = kTRUE ; }
2817 else { kTpcAll = kTRUE ; }
2818 }
2819
2820 // PID Multiplicities (particle for R.P.) - done just one time for each selection
2821 if(j==0) { fHistFull[k].fHistBayPidMult->Fill(fPidId) ; }
2822
2823 // Calculate weights for filling histograms
2824 float wt = 1. ; // TMath::Abs(fFlowEvent->Weight(k, j, fFlowTrack)) ;
2825
2826 // Fill histograms with selections
2827 if(kTpcPlus) { fHistFull[k].fHistFullHar[j].fHistPhiPlus->Fill(phi,wt) ; }
2828 else if(kTpcMinus) { fHistFull[k].fHistFullHar[j].fHistPhiMinus->Fill(phi,wt) ; }
2829 else if(kTpcAll) { fHistFull[k].fHistFullHar[j].fHistPhiAll->Fill(phi,wt) ; }
2830 fHistFull[k].fHistFullHar[j].fHistPhi->Fill(phi,wt) ;
2831
2832 // Get phiWgt from file
2833 double phiWgt;
2834 if(order > 3. && !oddHar) { phiWgt = fFlowEvent->PhiWeight(k, 1, fFlowTrack) ; }
2835 else { phiWgt = fFlowEvent->PhiWeight(k, j, fFlowTrack) ; }
2836 if(oddHar && eta<0.) { phiWgt /= -1. ; } // only for flat hists
2837 // cout << " Weight [" << k << "][" << j << "] (" << fFlowTrack->GetName() << ") = " << phiWgt << " or " << wt << " for phi = " << phi << endl ; // chk
2838
2839 // Fill Flat histograms
2840 if(kTpcPlus) { fHistFull[k].fHistFullHar[j].fHistPhiFlatPlus->Fill(phi, phiWgt) ; }
2841 else if(kTpcMinus) { fHistFull[k].fHistFullHar[j].fHistPhiFlatMinus->Fill(phi, phiWgt) ; }
2842 else if(kTpcAll) { fHistFull[k].fHistFullHar[j].fHistPhiFlatAll->Fill(phi, phiWgt) ; }
2843 fHistFull[k].fHistFullHar[j].fHistPhiFlat->Fill(phi,phiWgt) ;
2844
2845 if(oddHar && eta<0.) { phiWgt *= -1. ; } // restore value
2846
2847 // Remove autocorrelations
327288af 2848 TVector2 qi;
da5aa0a0 2849 if(!fEtaSub) // random subevents (or other)
9777bfcb 2850 {
2851 if(order > 3. && !oddHar) // 2nd harmonic event plane
2852 {
327288af 2853 qi.Set(phiWgt * cos(phi * 2), phiWgt * sin(phi * 2));
2854 TVector2 mQi = fQ[k][1] - qi;
2855 psii = mQi.Phi() / 2;
2856 if(psii < 0.) { psii += TMath::Pi() ; }
9777bfcb 2857 }
2858 else
2859 {
327288af 2860 qi.Set(phiWgt * cos(phi * order), phiWgt * sin(phi * order));
2861 TVector2 mQi = fQ[k][j] - qi;
2862 psii = mQi.Phi() / order;
2863 if(psii < 0.) { psii += 2*TMath::Pi()/order ; }
9777bfcb 2864 }
2865 }
2866
2867 // Remove autocorrelations of the second order 'particles' which are used for v1{EP1,EP2}.
2868 if (fV1Ep1Ep2 == kTRUE && order == 1)
2869 {
2870 AliFlowSelection usedForPsi2 = *fFlowSelect ;
2871 usedForPsi2.SetHarmonic(1);
2872 if(usedForPsi2.Select(fFlowTrack)) // particle was used for Psi2
2873 {
327288af 2874 qi.Set(phiWgt * cos(phi * 2), phiWgt * sin(phi * 2));
2875 TVector2 mQi = fQ[k][1] - qi;
2876 psi2 = mQi.Phi() / 2;
2877 if(psi2 < 0.) { psi2 += TMath::Pi() ; }
9777bfcb 2878 }
2879 else // particle was not used for Psi2
2880 {
327288af 2881 psi2 = fPsi[k][1];
9777bfcb 2882 }
2883 }
2884 }
2885 else
2886 {
2887 fHistFull[k].fHistFullHar[j].fHistYieldPtout->Fill(pt);
2888 fHistFull[k].fHistFullHar[j].fHistEtaPtPhi3Dout->Fill(eta, pt, phi);
2889 fHistFull[k].fHistFullHar[j].fHistYield2Dout->Fill(eta, pt);
2890 fHistFull[k].fHistFullHar[j].fHistDcaGlobout->Fill(TMath::Abs(fFlowTrack->Dca()));
2891 }
2892
2893 // Caculate v for all particles selected for correlation analysis
2894 if(fFlowSelect->SelectPart(fFlowTrack))
2895 {
2896 corrMultN++;
2897
2898 float v;
2899 if (fV1Ep1Ep2 == kFALSE || order != 1)
2900 {
327288af 2901 v = 100 * cos(order * (phi - psii)) ;
9777bfcb 2902 }
2903 else // i.e. (fV1Ep1Ep2 == kTRUE && order == 1)
2904 {
327288af 2905 v = 100 * cos(phi + psii - 2*psi2) ;
9777bfcb 2906 }
2907
2908 float vFlip = v;
2909 if(eta < 0 && oddHar) { vFlip *= -1 ; }
2910 if(strlen(fFlowSelect->PidPart()) != 0) // pid, fill rapidity
2911 {
2912 float rapidity = fFlowTrack->Y();
2913 fHistFull[k].fHistFullHar[j].fHistvObs2D->Fill(rapidity, pt, v);
2914
2915 if(fPtRangevEta[1] > fPtRangevEta[0]) // cut is used
2916 {
2917 if(pt < fPtRangevEta[1] && pt >= fPtRangevEta[0]) // check cut range, fill if in range
2918 {
2919 fHistFull[k].fHistFullHar[j].fHistvObsEta->Fill(rapidity, v);
2920 }
2921 }
2922 else // cut is not used, fill in any case
2923 {
2924 fHistFull[k].fHistFullHar[j].fHistvObsEta->Fill(rapidity, v);
2925 }
2926 }
2927 else // no pid, fill eta
2928 {
2929 fHistFull[k].fHistFullHar[j].fHistvObs2D->Fill(eta, pt, v);
2930
2931 if(fPtRangevEta[1] > fPtRangevEta[0]) // cut is used
2932 {
2933 if(pt < fPtRangevEta[1] && pt >= fPtRangevEta[0]) // check cut range, fill if in range
2934 {
2935 fHistFull[k].fHistFullHar[j].fHistvObsEta->Fill(eta, v);
2936 }
2937 }
2938 else // cut is not used, fill in any case
2939 {
2940 fHistFull[k].fHistFullHar[j].fHistvObsEta->Fill(eta, v);
2941 }
2942 }
2943
2944 if(fEtaRangevPt[1] > fEtaRangevPt[0]) // cut is used
2945 {
2946 if(TMath::Abs(eta) < fEtaRangevPt[1] && TMath::Abs(eta) >= fEtaRangevPt[0]) // check cut range, fill if in range
2947 {
2948 fHistFull[k].fHistFullHar[j].fHistvObsPt->Fill(pt, vFlip); // for odd harmonis /-/
2949 }
2950 }
2951 else // cut is not used, fill in any case
2952 {
2953 fHistFull[k].fHistFullHar[j].fHistvObsPt->Fill(pt, vFlip);
2954 }
2955
2956 // v_
2957 Bool_t etaPtNoCut = kTRUE;
2958 if(fPtRangevEta[1] > fPtRangevEta[0] && (pt < fPtRangevEta[0] || pt >= fPtRangevEta[1]))
2959 {
2960 etaPtNoCut = kFALSE;
2961 }
2962 if(fEtaRangevPt[1] > fEtaRangevPt[0] && (TMath::Abs(eta) < fEtaRangevPt[0] || TMath::Abs(eta) >= fEtaRangevPt[1]))
2963 {
2964 etaPtNoCut = kFALSE;
2965 }
2966 if(etaPtNoCut) { fHistFull[k].fHistvObs->Fill(order, vFlip) ; }
2967
2968 // Correlation of Phi of selected particles with Psi
327288af 2969 float phii = phi;
9777bfcb 2970 if(eta < 0 && oddHar)
2971 {
327288af 2972 phii += TMath::Pi() ; // backward particle and odd harmonic
2973 if(phii > 2*TMath::Pi()) { phii -= 2*TMath::Pi() ; }
9777bfcb 2974 }
327288af 2975 float dPhi = phii - psii;
9777bfcb 2976 if(dPhi < 0.) { dPhi += 2*TMath::Pi() ; }
2977 fHistFull[k].fHistFullHar[j].fHistPhiCorr->Fill(fmod((double)dPhi, 2*TMath::Pi() / order));
2978 }
2979 }
2980 }
2981
2982 return corrMultN ;
2983}
2984//-----------------------------------------------------------------------
92016a03 2985void AliFlowAnalyser::FillV0Histograms(TClonesArray* fFlowV0s)
9777bfcb 2986{
2987 // v0s histograms
2988
2989 cout << " V0s Loop . " << endl ;
2990
92016a03 2991 fSelV0s = 0 ;
9777bfcb 2992 for(fV0Number=0;fV0Number<fNumberOfV0s;fV0Number++)
2993 {
92016a03 2994 //fFlowV0 = (AliFlowV0*)fFlowV0s->At(fV0Number) ;
2995
2996 fFlowV0 = (AliFlowV0*)(fFlowV0s->AddrAt(fV0Number)) ;
2997 // cout << " looping v0 n. " << fV0Number << " * " << fFlowV0 << endl ;
9777bfcb 2998
2999 float mass = fFlowV0->Mass() ;
3000 float eta = fFlowV0->Eta() ;
3001 float rapidity = fFlowV0->Y() ;
3002 float phi = fFlowV0->Phi() ;
3003 float pt = fFlowV0->Pt() ;
3004 float totalp = fFlowV0->P() ;
3005 //int charge = fFlowV0->Charge() ;
3006 float dca = fFlowV0->Dca() ;
3007 float lenght = fFlowV0->V0Lenght() ;
3008 float sigma = fFlowV0->Sigma() ;
3009 float chi2 = fFlowV0->Chi2() ;
3010 Char_t pid[10] ; strcpy(pid, fFlowV0->Pid()) ;
3011 int labelPlus = -1 ;
3012 int labelMinus = -1 ;
da5aa0a0 3013 AliFlowTrack* daughterPlus = (AliFlowTrack*)(fFlowTracks->AddrAt(fFlowV0->DaughterP())) ;
3014 AliFlowTrack* daughterMinus = (AliFlowTrack*)(fFlowTracks->AddrAt(fFlowV0->DaughterN())) ; ;
9777bfcb 3015 if(daughterPlus) { labelPlus = daughterPlus->Label() ; }
3016 if(daughterMinus) { labelMinus = daughterMinus->Label() ; }
3017
3018 fHistV0Mass->Fill(mass) ;
3019 fHistV0EtaPtPhi3D->Fill(eta, pt, phi) ;
3020 fHistV0YieldAll2D->Fill(eta, pt) ;
3021 fHistV0Dca->Fill(dca);
3022 fHistV0Chi2->Fill(chi2);
3023 fHistV0Lenght->Fill(lenght);
3024 fHistV0Sigma->Fill(sigma);
3025 fHistV0MassPtSlices->Fill(mass,pt);
3026 fHistV0PYall2D->Fill(rapidity, totalp) ;
3027
3028 if(fFlowSelect->SelectPart(fFlowV0))
3029 {
3030 bool inWin = fFlowSelect->SelectV0Part(fFlowV0) ;
3031 bool sx = fFlowSelect->SelectV0sxSide(fFlowV0) ;
3032 bool dx = fFlowSelect->SelectV0dxSide(fFlowV0) ;
92016a03 3033 fSelV0s++ ;
9777bfcb 3034
3035 fHistV0YieldPart2D->Fill(eta, pt) ;
3036 if(inWin)
3037 {
3038 fHistV0EtaPtPhi3DPart->Fill(eta, pt, phi) ;
3039 fHistV0LenghtPart->Fill(lenght);
3040 fHistV0DcaPart->Fill(dca);
3041 fHistV0MassWin->Fill(mass) ;
3042 fHistV0BinEta->Fill(eta, eta);
3043 fHistV0BinPt->Fill(pt, pt);
3044 }
3045 else
3046 {
3047 fHistV0sbEtaPtPhi3DPart->Fill(eta, pt, phi) ;
3048 fHistV0sbLenghtPart->Fill(lenght);
3049 fHistV0sbDcaPart->Fill(dca);
3050 fHistV0sbMassSide->Fill(mass) ;
3051 fHistV0sbBinEta->Fill(eta, eta);
3052 fHistV0sbBinPt->Fill(pt, pt);
3053 }
3054
327288af 3055 for(int k = 0; k < AliFlowConstants::kSels; k++) // sort of HarmonicsLoop - selection number used
9777bfcb 3056 {
3057 fFlowSelect->SetSelection(k) ;
327288af 3058 for(Int_t j=0;j<AliFlowConstants::kHars;j++)
9777bfcb 3059 {
3060 Bool_t oddHar = (j+1) % 2 ;
3061 Float_t order = (Float_t)(j+1) ;
3062 fFlowSelect->SetHarmonic(j);
3063
3064 // Remove autocorrelations
327288af 3065 Float_t psii, psi2 ;
3066 TVector2 q1, q2 ;
9777bfcb 3067 Float_t phiDaughter1 = 0. ; Float_t phiDaughter2 = 0. ;
3068 Double_t phiWgt1 = 0. ; Double_t phiWgt2 = 0. ;
3069 // -
3070 if(daughterPlus)
3071 { // 1
3072 fFlowTrack = daughterPlus ;
3073 phiDaughter1 = fFlowTrack->Phi() ;
3074 if(fFlowSelect->Select(fFlowTrack)) // Get phiWgt from file
3075 {
3076 if(order > 3. && !oddHar) { phiWgt1 = fFlowEvent->PhiWeight(k, 1, fFlowTrack) ; }
3077 else { phiWgt1 = fFlowEvent->PhiWeight(k, j, fFlowTrack) ; }
3078 }
3079 }
3080 if(daughterMinus)
3081 { // 2
3082 fFlowTrack = daughterMinus ;
3083 phiDaughter2 = fFlowTrack->Phi() ;
3084 if(fFlowSelect->Select(fFlowTrack)) // Get phiWgt from file
3085 {
3086 if(order > 3. && !oddHar) { phiWgt2 = fFlowEvent->PhiWeight(k, 1, fFlowTrack) ; }
3087 else { phiWgt2 = fFlowEvent->PhiWeight(k, j, fFlowTrack) ; }
3088 }
3089 }
3090
327288af 3091 // psi2
3092 q1.Set(phiWgt1 * cos(phiDaughter1 * 2), phiWgt1 * sin(phiDaughter1 * 2));
3093 q2.Set(phiWgt2 * cos(phiDaughter2 * 2), phiWgt2 * sin(phiDaughter2 * 2));
3094 TVector2 mQi = fQ[k][1] ; mQi -= q1 ; mQi -= q2 ;
3095 psi2 = mQi.Phi() / 2 ;
3096 if(psi2 < 0.) { psi2 += TMath::Pi() ; }
3097 if(psi2 > TMath::Pi()) { psi2 -= TMath::Pi() ; }
9777bfcb 3098
327288af 3099 // psii
3100 if(order > 3. && !oddHar) { psii = psi2 ; }
9777bfcb 3101 else
3102 {
327288af 3103 q1.Set(phiWgt1 * cos(phiDaughter1 * order), phiWgt1 * sin(phiDaughter1 * order));
3104 q2.Set(phiWgt2 * cos(phiDaughter2 * order), phiWgt2 * sin(phiDaughter2 * order));
3105 TVector2 mQi = fQ[k][j] ; mQi -= q1 ; mQi -= q2 ;
3106 psii = mQi.Phi()/order ;
3107 if(psii < 0.) { psii += 2*TMath::Pi()/order ; }
3108 if(psii > 2*TMath::Pi()/order) { psii -= 2*TMath::Pi()/order ; }
9777bfcb 3109 }
3110
3111 // Caculate v for all V0s selected for correlation analysis
3112 float v ;
327288af 3113 if(fV1Ep1Ep2 == kFALSE || order != 1) { v = 100 * cos(order * (phi - psii)) ; }
3114 else { v = 100 * cos(phi + psii - 2*psi2) ; }
9777bfcb 3115 float vFlip = v ; if(eta < 0 && oddHar) { vFlip *= -1 ; }
3116
3117 // invariant mass windows & sidebands
3118 if(inWin) { fHistFull[k].fHistFullHar[j].fHistV0vObs2D->Fill(eta, pt, v) ; }
3119 else { fHistFull[k].fHistFullHar[j].fHistV0sbvObs2D->Fill(eta, pt, v) ; }
3120
3121 if(fPtRangevEta[1] > fPtRangevEta[0]) // cut is used
3122 {
3123 if(pt < fPtRangevEta[1] && pt >= fPtRangevEta[0]) // check cut range, fill if in range
3124 {
3125 if(inWin) { fHistFull[k].fHistFullHar[j].fHistV0vObsEta->Fill(eta, v) ; }
3126 else
3127 {
3128 fHistFull[k].fHistFullHar[j].fHistV0sbvObsEta->Fill(eta, v) ;
3129 if(sx) { fHistFull[k].fHistFullHar[j].fHistV0sbvObsEtaSx->Fill(eta, v) ; }
3130 else if(dx) { fHistFull[k].fHistFullHar[j].fHistV0sbvObsEtaDx->Fill(eta, v) ; }
3131 }
3132 }
3133 }
3134 else // cut is not used, fill in any case
3135 {
3136 if(inWin) { fHistFull[k].fHistFullHar[j].fHistV0vObsEta->Fill(eta, v) ; }
3137 else
3138 {
3139 fHistFull[k].fHistFullHar[j].fHistV0sbvObsEta->Fill(eta, v) ;
3140 if(sx) { fHistFull[k].fHistFullHar[j].fHistV0sbvObsEtaSx->Fill(eta, v) ; }
3141 else if(dx) { fHistFull[k].fHistFullHar[j].fHistV0sbvObsEtaDx->Fill(eta, v) ; }
3142 }
3143 }
3144 if(fEtaRangevPt[1] > fEtaRangevPt[0]) // cut is used
3145 {
3146 if(TMath::Abs(eta) < fEtaRangevPt[1] && TMath::Abs(eta) >= fEtaRangevPt[0]) // check cut range, fill if in range
3147 {
3148 if(inWin) { fHistFull[k].fHistFullHar[j].fHistV0vObsPt->Fill(pt, vFlip) ; } // for odd harmonis /-/
3149 else
3150 {
3151 fHistFull[k].fHistFullHar[j].fHistV0sbvObsPt->Fill(pt, vFlip) ;
3152 if(sx) { fHistFull[k].fHistFullHar[j].fHistV0sbvObsPtSx->Fill(eta, v) ; }
3153 else if(dx) { fHistFull[k].fHistFullHar[j].fHistV0sbvObsPtDx->Fill(eta, v) ; }
3154 }
3155 }
3156 }
3157 else // cut is not used, fill in any case
3158 {
3159 if(inWin) { fHistFull[k].fHistFullHar[j].fHistV0vObsPt->Fill(pt, vFlip) ; }
3160 else
3161 {
3162 fHistFull[k].fHistFullHar[j].fHistV0sbvObsPt->Fill(pt, vFlip) ;
3163 if(sx) { fHistFull[k].fHistFullHar[j].fHistV0sbvObsPtSx->Fill(eta, v) ; }
3164 else if(dx) { fHistFull[k].fHistFullHar[j].fHistV0sbvObsPtDx->Fill(eta, v) ; }
3165 }
3166 }
3167 // v_
3168 Bool_t etaPtNoCut = kTRUE;
3169 if(fPtRangevEta[1] > fPtRangevEta[0] && (pt < fPtRangevEta[0] || pt >= fPtRangevEta[1]))
3170 {
3171 etaPtNoCut = kFALSE;
3172 }
3173 if(fEtaRangevPt[1] > fEtaRangevPt[0] && (TMath::Abs(eta) < fEtaRangevPt[0] || TMath::Abs(eta) >= fEtaRangevPt[1]))
3174 {
3175 etaPtNoCut = kFALSE;
3176 }
3177 if(etaPtNoCut)
3178 {
3179 if(inWin) { fHistFull[k].fHistV0vObs->Fill(order, vFlip) ; }
3180 else
3181 {
3182 if(sx) { fHistFull[k].fHistV0sbvObsSx->Fill(order, vFlip) ; }
3183 else if(dx) { fHistFull[k].fHistV0sbvObsDx->Fill(order, vFlip) ; }
3184 }
3185 }
3186
3187 // Correlation of Phi of selected v0s with Psi
327288af 3188 float phii = phi;
9777bfcb 3189 if(eta < 0 && oddHar)
3190 {
327288af 3191 phii += TMath::Pi() ; // backward particle and odd harmonic
3192 if(phii > 2*TMath::Pi()) { phii -= 2*TMath::Pi() ; }
9777bfcb 3193 }
327288af 3194 float dPhi = phii - psii;
9777bfcb 3195 if(dPhi < 0.) { dPhi += 2*TMath::Pi() ; }
3196
3197 if(inWin) { fHistFull[k].fHistFullHar[j].fHistV0PhiCorr->Fill(fmod((double)dPhi, 2*TMath::Pi() / order)) ; }
3198 else { fHistFull[k].fHistFullHar[j].fHistV0sbPhiCorr->Fill(fmod((double)dPhi, 2*TMath::Pi() / order)) ; }
3199 }
3200 }
3201 }
3202 //delete fFlowV0 ; fFlowV0 = 0 ;
3203 //delete daughterPlus ; daughterPlus = 0 ;
3204 //delete daughterMinus ; daughterMinus = 0 ;
3205 }
92016a03 3206 fHistV0MultPart->Fill(fSelV0s) ;
da5aa0a0 3207 fTotalNumberOfV0s += fNumberOfV0s ;
9777bfcb 3208
3209 return ;
3210}
3211//-----------------------------------------------------------------------
3212// void AliFlowAnalyser::FillLabels()
3213// {
3214// // Reads the tracks label (i.e. the link from ESD tracking to KineTree)
3215// // and ...
3216//
3217// cout << " Fill Labels . " << endl ;
3218//
327288af 3219// fLabHist ;
9777bfcb 3220// return ;
3221// }
3222//-----------------------------------------------------------------------
327288af 3223void AliFlowAnalyser::PrintEventQuantities() const
9777bfcb 3224{
3225 // prints event by event calculated quantities
3226
3227 cout << endl ; cout << " # " << fEventNumber << " - Event quantities : " << endl ; cout << endl ;
3228
3229 cout << "fQ[k][j] = " ;
327288af 3230 for(int k=0;k<AliFlowConstants::kSels;k++)
9777bfcb 3231 {
327288af 3232 for(int j=0;j<AliFlowConstants::kHars;j++)
9777bfcb 3233 {
3234 cout << (Float_t)fQ[k][j].X() << "," << (Float_t)fQ[k][j].Y() << " ; " ;
3235 }
3236 cout << endl ;
3237 cout << " " ;
3238 }
3239 cout << endl ; cout << "fPsi[k][j] = " ;
327288af 3240 for(int k=0;k<AliFlowConstants::kSels;k++)
9777bfcb 3241 {
327288af 3242 for(int j=0;j<AliFlowConstants::kHars;j++) { Float_t aaa = (Float_t)fPsi[k][j] ; cout << aaa << " , " ; }
9777bfcb 3243 cout << endl ;
3244 cout << " " ;
3245 }
3246 cout << endl ; cout << "fQnorm[k][j] = " ;
327288af 3247 for(int k=0;k<AliFlowConstants::kSels;k++)
9777bfcb 3248 {
327288af 3249 for(int j=0;j<AliFlowConstants::kHars;j++) { Float_t aaa = (Float_t)fQnorm[k][j] ; cout << aaa << " , " ; }
9777bfcb 3250 cout << endl ;
3251 cout << " " ;
3252 }
3253 cout << endl ; cout << "fMult[k][j] = " ;
327288af 3254 for(int k=0;k<AliFlowConstants::kSels;k++)
9777bfcb 3255 {
327288af 3256 for(int j=0;j<AliFlowConstants::kHars;j++) { Float_t aaa = (Float_t)fMult[k][j] ; cout << aaa << " , " ; }
9777bfcb 3257 cout << endl ;
3258 cout << " " ;
3259 }
3260 cout << endl ; cout << "fMultSub[n][k][j] = " ;
327288af 3261 for(int n=0;n<AliFlowConstants::kSubs;n++)
9777bfcb 3262 {
327288af 3263 for(int k=0;k<AliFlowConstants::kSels;k++)
9777bfcb 3264 {
327288af 3265 for(int j=0;j<AliFlowConstants::kHars;j++) { Float_t aaa = fMultSub[n][k][j] ; cout << aaa << " , " ; }
9777bfcb 3266 cout << endl ;
3267 cout << " " ;
3268 }
3269 }
3270 cout << endl ; cout << "fPsiSub[n][k][j] = " ;
327288af 3271 for(int n=0;n<AliFlowConstants::kSubs;n++)
9777bfcb 3272 {
327288af 3273 for(int k=0;k<AliFlowConstants::kSels;k++)
9777bfcb 3274 {
327288af 3275 for(int j=0;j<AliFlowConstants::kHars;j++) { Float_t aaa = fPsiSub[n][k][j] ; cout << aaa << " , " ; }
9777bfcb 3276 cout << endl ;
3277 cout << " " ;
3278 }
3279 }
3280 cout << endl ; cout << "Delta_PsiSub[k][j] = " ;
327288af 3281 for(int k=0;k<AliFlowConstants::kSels;k++)
9777bfcb 3282 {
327288af 3283 for(int j=0;j<AliFlowConstants::kHars;j++) { Float_t aaa = fPsiSub[0][k][j]-fPsiSub[1][k][j] ; cout << aaa << " , " ; }
9777bfcb 3284 cout << endl ;
3285 cout << " " ;
3286 }
3287 cout << endl ;
92016a03 3288
3289 cout << " N. of selected tracks - SelPart(AliFlowTrack*) = " << fSelParts << endl ; //! n. of tracks selected for correlation analysis
3290 cout << " N. of selected V0s - SelPart(AliFlowV0*) = " << fSelV0s << endl ; //! n. of v0s selected for correlation analysis
3291
3292 cout << endl ;
9777bfcb 3293}
3294//----------------------------------------------------------------------
3295// ###
3296//----------------------------------------------------------------------
3297Bool_t AliFlowAnalyser::Resolution()
3298{
3299 // Calculates the resolution from cos(dPsi) and corrects the observed flow.
3300 // New histograms are then created and filled with the corrected vn
3301 // (see fVnResHistList), and saved in the otput file.
3302
3303 cout << " AliFlowAnalyser::Resolution()" << endl ;
3304
3305 // VnRes histogram collection
327288af 3306 fVnResHistList = new TOrdCollection(AliFlowConstants::kSels*AliFlowConstants::kHars);
9777bfcb 3307
3308 // Calculate resolution from sqrt(fHistCos)
327288af 3309 double cosPair[AliFlowConstants::kSels][AliFlowConstants::kHars];
3310 double cosPairErr[AliFlowConstants::kSels][AliFlowConstants::kHars];
9777bfcb 3311 double content, contentV0;
3312 double error, errorV0;
3313 double totalError;
3314 TString* histTitle;
3315
327288af 3316 for(int k = 0; k < AliFlowConstants::kSels; k++)
9777bfcb 3317 {
3318 // v for tracks (corrected for resolution)
3319 histTitle = new TString("Flow_v_Sel");
3320 *histTitle += k+1;
3321 fHistFull[k].fHistv = fHistFull[k].fHistvObs->ProjectionX(histTitle->Data());
3322 fHistFull[k].fHistv->SetTitle(histTitle->Data());
3323 fHistFull[k].fHistv->SetXTitle("Harmonic");
3324 fHistFull[k].fHistv->SetYTitle("v (%)");
3325 delete histTitle;
3326 fVnResHistList->AddLast(fHistFull[k].fHistv);
3327
3328 // v for v0s (corrected for resolution)
3329 histTitle = new TString("FlowV0_v_Sel");
3330 *histTitle += k+1;
3331 fHistFull[k].fHistV0v = fHistFull[k].fHistV0vObs->ProjectionX(histTitle->Data());
3332 fHistFull[k].fHistV0v->SetTitle(histTitle->Data());
3333 fHistFull[k].fHistV0v->SetXTitle("Harmonic");
3334 fHistFull[k].fHistV0v->SetYTitle("v (%)");
3335 delete histTitle;
3336 fVnResHistList->AddLast(fHistFull[k].fHistV0v);
3337
327288af 3338 for(int j = 0; j < AliFlowConstants::kHars; j++)
9777bfcb 3339 {
3340 double order = (double)(j+1);
3341 cosPair[k][j] = fHistFull[k].fHistCos->GetBinContent(j+1);
3342 cosPairErr[k][j] = fHistFull[k].fHistCos->GetBinError(j+1);
3343 if(cosPair[k][j] > 0.)
3344 {
3345 double resSub = 0. ;
3346 double resSubErr = 0. ;
3347 if(fV1Ep1Ep2 == kTRUE && order == 1) // calculate resolution of second order event plane first
3348 {
3349 double res2 = 0. ;
3350 double res2error = 0. ;
3351 if(fHistFull[k].fHistCos->GetBinContent(2) > 0.)
3352 {
3353 if(fEtaSub) // sub res only
3354 {
3355 res2 = TMath::Sqrt(fHistFull[k].fHistCos->GetBinContent(2));
3356 res2error = fHistFull[k].fHistCos->GetBinError(2) / (2. * res2);
3357 }
3358 else
3359 {
3360 if (fHistFull[k].fHistCos->GetBinContent(2) > 0.92) // resolution saturates
3361 {
3362 res2 = 0.99;
3363 res2error = 0.007;
3364 }
3365 else
3366 {
3367 double deltaRes2Sub = 0.005; // differential for the error propagation
3368 double res2Sub = TMath::Sqrt(fHistFull[k].fHistCos->GetBinContent(2));
3369 double res2SubErr = fHistFull[k].fHistCos->GetBinError(2) / (2. * res2Sub);
3370 double chiSub2 = Chi(res2Sub);
3371 double chiSub2Delta = Chi(res2Sub + deltaRes2Sub);
3372 res2 = ResEventPlane(TMath::Sqrt(2.) * chiSub2); // full event plane res.
3373 double mRes2Delta = ResEventPlane(TMath::Sqrt(2.) * chiSub2Delta);
3374 res2error = res2SubErr * fabs((double)res2 - mRes2Delta) / deltaRes2Sub;
3375 }
3376 }
3377 }
3378 else
3379 {
3380 res2 = 0.;
3381 res2error = 0.;
3382 }
3383 // now put everything together with first order event plane
3384 fRes[k][j] = TMath::Sqrt(cosPair[k][0]*res2);
3385 fResErr[k][j] = 1./(2.*fRes[k][j]) * TMath::Sqrt(cosPairErr[k][0]*cosPairErr[k][0] + res2error*res2error); // Gaussian error propagation
3386 }
3387 else if(fEtaSub) // sub res only
3388 {
3389 resSub = TMath::Sqrt(cosPair[k][j]);
3390 resSubErr = cosPairErr[k][j] / (2. * resSub);
3391 fRes[k][j] = resSub;
3392 fResErr[k][j] = resSubErr;
3393 }
3394 else if(order==4. || order==6.|| order==8.) // 2nd harmonic event plane
3395 {
3396 double deltaResSub = 0.005; // differential for the error propagation
3397 double resSub = TMath::Sqrt(cosPair[k][1]);
3398 double resSubErr = cosPairErr[k][1] / (2. * resSub);
3399 double chiSub = Chi(resSub);
3400 double chiSubDelta = Chi(resSub + deltaResSub);
3401 double mResDelta;
3402 if (order==4.)
3403 {
3404 fRes[k][j] = ResEventPlaneK2(TMath::Sqrt(2.) * chiSub); // full event plane res.
3405 mResDelta = ResEventPlaneK2(TMath::Sqrt(2.) * chiSubDelta);
3406 }
3407 else if(order==6.)
3408 {
3409 fRes[k][j] = ResEventPlaneK3(TMath::Sqrt(2.) * chiSub); // full event plane res.
3410 mResDelta = ResEventPlaneK3(TMath::Sqrt(2.) * chiSubDelta);
3411 }
3412 else
3413 {
3414 fRes[k][j] = ResEventPlaneK4(TMath::Sqrt(2.) * chiSub); // full event plane res.
3415 mResDelta = ResEventPlaneK4(TMath::Sqrt(2.) * chiSubDelta);
3416 }
3417 fResErr[k][j] = resSubErr * fabs((double)fRes[k][j] - mResDelta) / deltaResSub;
3418 }
3419 else
3420 {
3421 if(cosPair[k][j] > 0.92) // resolution saturates
3422 {
3423 fRes[k][j] = 0.99;
3424 fResErr[k][j] = 0.007;
3425 }
3426 else
3427 {
3428 double deltaResSub = 0.005; // differential for the error propagation
3429 double resSub = TMath::Sqrt(cosPair[k][j]);
3430 double resSubErr = cosPairErr[k][j] / (2. * resSub);
3431 double chiSub = Chi(resSub);
3432 double chiSubDelta = Chi(resSub + deltaResSub);
3433 fRes[k][j] = ResEventPlane(TMath::Sqrt(2.) * chiSub); // full event plane res.
3434 double mResDelta = ResEventPlane(TMath::Sqrt(2.) * chiSubDelta);
3435 fResErr[k][j] = resSubErr * TMath::Abs((double)fRes[k][j] - mResDelta) / deltaResSub;
3436 }
3437 }
3438 }
3439 else // subevent correlation must be positive
3440 {
3441 fRes[k][j] = 0.;
3442 fResErr[k][j] = 0.;
3443 }
3444 fHistFull[k].fHistRes->SetBinContent(j+1, fRes[k][j]);
3445 fHistFull[k].fHistRes->SetBinError(j+1, fResErr[k][j]);
3446
3447 // Create the v 2D histogram (Flow corrected for res.) - v,Pt,Eta
3448 histTitle = new TString("Flow_v2D_Sel");
3449 *histTitle += k+1;
3450 histTitle->Append("_Har");
3451 *histTitle += j+1;
3452 fHistFull[k].fHistFullHar[j].fHistv2D = fHistFull[k].fHistFullHar[j].fHistvObs2D->ProjectionXY(histTitle->Data());
3453 fHistFull[k].fHistFullHar[j].fHistv2D->SetTitle(histTitle->Data());
3454 fHistFull[k].fHistFullHar[j].fHistv2D->SetXTitle((char*)fLabel.Data());
3455 fHistFull[k].fHistFullHar[j].fHistv2D->SetYTitle("Pt (GeV/c)");
3456 fHistFull[k].fHistFullHar[j].fHistv2D->SetZTitle("v (%)");
3457 delete histTitle;
3458 fVnResHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistv2D);
3459
3460 // Create the 1D v histograms - v,Eta
3461 histTitle = new TString("Flow_vEta_Sel");
3462 *histTitle += k+1;
3463 histTitle->Append("_Har");
3464 *histTitle += j+1;
3465 fHistFull[k].fHistFullHar[j].fHistvEta = fHistFull[k].fHistFullHar[j].fHistvObsEta->ProjectionX(histTitle->Data());
3466 fHistFull[k].fHistFullHar[j].fHistvEta->SetTitle(histTitle->Data());
3467 fHistFull[k].fHistFullHar[j].fHistvEta->SetXTitle((char*)fLabel.Data());
3468 fHistFull[k].fHistFullHar[j].fHistvEta->SetYTitle("v (%)");
3469 delete histTitle;
3470 fVnResHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistvEta);
3471
3472 // v,Pt
3473 histTitle = new TString("Flow_vPt_Sel");
3474 *histTitle += k+1;
3475 histTitle->Append("_Har");
3476 *histTitle += j+1;
3477 fHistFull[k].fHistFullHar[j].fHistvPt = fHistFull[k].fHistFullHar[j].fHistvObsPt->ProjectionX(histTitle->Data());
3478 fHistFull[k].fHistFullHar[j].fHistvPt->SetTitle(histTitle->Data());
3479 fHistFull[k].fHistFullHar[j].fHistvPt->SetXTitle("Pt (GeV/c)");
3480 fHistFull[k].fHistFullHar[j].fHistvPt->SetYTitle("v (%)");
3481 delete histTitle;
3482 fVnResHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistvPt);
3483
3484 // Create the v 2D histogram (V0s Flow corrected for res.) - v,Pt,Eta
3485 histTitle = new TString("FlowV0_v2D_Sel");
3486 *histTitle += k+1;
3487 histTitle->Append("_Har");
3488 *histTitle += j+1;
3489 fHistFull[k].fHistFullHar[j].fHistV0v2D = fHistFull[k].fHistFullHar[j].fHistV0vObs2D->ProjectionXY(histTitle->Data());
3490 fHistFull[k].fHistFullHar[j].fHistV0v2D->SetTitle(histTitle->Data());
3491 fHistFull[k].fHistFullHar[j].fHistV0v2D->SetXTitle((char*)fLabel.Data());
3492 fHistFull[k].fHistFullHar[j].fHistV0v2D->SetYTitle("Pt (GeV/c)");
3493 fHistFull[k].fHistFullHar[j].fHistV0v2D->SetZTitle("v (%)");
3494 delete histTitle;
3495 fVnResHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistV0v2D);
3496
3497 // Create the 1D v histograms (V0s) - v,Eta
3498 histTitle = new TString("FlowV0_vEta_Sel");
3499 *histTitle += k+1;
3500 histTitle->Append("_Har");
3501 *histTitle += j+1;
3502 fHistFull[k].fHistFullHar[j].fHistV0vEta = fHistFull[k].fHistFullHar[j].fHistV0vObsEta->ProjectionX(histTitle->Data());
3503 fHistFull[k].fHistFullHar[j].fHistV0vEta->SetTitle(histTitle->Data());
3504 fHistFull[k].fHistFullHar[j].fHistV0vEta->SetXTitle((char*)fLabel.Data());
3505 fHistFull[k].fHistFullHar[j].fHistV0vEta->SetYTitle("v (%)");
3506 delete histTitle;
3507 fVnResHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistV0vEta);
3508
3509 // (V0s) v,Pt
3510 histTitle = new TString("FlowV0_vPt_Sel");
3511 *histTitle += k+1;
3512 histTitle->Append("_Har");
3513 *histTitle += j+1;
3514 fHistFull[k].fHistFullHar[j].fHistV0vPt = fHistFull[k].fHistFullHar[j].fHistV0vObsPt->ProjectionX(histTitle->Data());
3515 fHistFull[k].fHistFullHar[j].fHistV0vPt->SetTitle(histTitle->Data());
3516 fHistFull[k].fHistFullHar[j].fHistV0vPt->SetXTitle("Pt (GeV/c)");
3517 fHistFull[k].fHistFullHar[j].fHistV0vPt->SetYTitle("v (%)");
3518 delete histTitle;
3519 fVnResHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistV0vPt);
3520
3521 // Create the v 2D histogram (V0s sidebands Flow corrected for res.) - v,Pt,Eta
3522 histTitle = new TString("FlowV0sb_v2D_Sel");
3523 *histTitle += k+1;
3524 histTitle->Append("_Har");
3525 *histTitle += j+1;
3526 fHistFull[k].fHistFullHar[j].fHistV0sbv2D = fHistFull[k].fHistFullHar[j].fHistV0sbvObs2D->ProjectionXY(histTitle->Data());
3527 fHistFull[k].fHistFullHar[j].fHistV0sbv2D->SetTitle(histTitle->Data());
3528 fHistFull[k].fHistFullHar[j].fHistV0sbv2D->SetXTitle((char*)fLabel.Data());
3529 fHistFull[k].fHistFullHar[j].fHistV0sbv2D->SetYTitle("Pt (GeV/c)");
3530 fHistFull[k].fHistFullHar[j].fHistV0sbv2D->SetZTitle("v (%)");
3531 delete histTitle;
3532 fVnResHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistV0sbv2D);
3533
3534 // Create the 1D v histograms (V0s sidebands) - v,Eta
3535 histTitle = new TString("FlowV0sb_vEta_Sel");
3536 *histTitle += k+1;
3537 histTitle->Append("_Har");
3538 *histTitle += j+1;
3539 fHistFull[k].fHistFullHar[j].fHistV0sbvEta = fHistFull[k].fHistFullHar[j].fHistV0sbvObsEta->ProjectionX(histTitle->Data());
3540 fHistFull[k].fHistFullHar[j].fHistV0sbvEta->SetTitle(histTitle->Data());
3541 fHistFull[k].fHistFullHar[j].fHistV0sbvEta->SetXTitle((char*)fLabel.Data());
3542 fHistFull[k].fHistFullHar[j].fHistV0sbvEta->SetYTitle("v (%)");
3543 delete histTitle;
3544 fVnResHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistV0sbvEta);
3545
3546 // (V0s sidebands) v,Pt
3547 histTitle = new TString("FlowV0sb_vPt_Sel");
3548 *histTitle += k+1;
3549 histTitle->Append("_Har");
3550 *histTitle += j+1;
3551 fHistFull[k].fHistFullHar[j].fHistV0sbvPt = fHistFull[k].fHistFullHar[j].fHistV0sbvObsPt->ProjectionX(histTitle->Data());
3552 fHistFull[k].fHistFullHar[j].fHistV0sbvPt->SetTitle(histTitle->Data());
3553 fHistFull[k].fHistFullHar[j].fHistV0sbvPt->SetXTitle("Pt (GeV/c)");
3554 fHistFull[k].fHistFullHar[j].fHistV0sbvPt->SetYTitle("v (%)");
3555 delete histTitle;
3556 fVnResHistList->AddLast(fHistFull[k].fHistFullHar[j].fHistV0sbvPt);
3557
3558 // Calulate v = vObs / Resolution
3559 if(fRes[k][j])
3560 {
3561 cout << "***** Resolution of the " << j+1 << "th harmonic = " << fRes[k][j] << " +/- " << fResErr[k][j] << endl;
3562 // selected tracks
3563 fHistFull[k].fHistFullHar[j].fHistv2D->Scale(1. / fRes[k][j]);
3564 fHistFull[k].fHistFullHar[j].fHistvEta->Scale(1. / fRes[k][j]);
3565 fHistFull[k].fHistFullHar[j].fHistvPt->Scale(1. / fRes[k][j]);
3566 content = fHistFull[k].fHistv->GetBinContent(j+1);
3567 content /= fRes[k][j];
3568 fHistFull[k].fHistv->SetBinContent(j+1, content);
3569 // selected V0s
3570 fHistFull[k].fHistFullHar[j].fHistV0v2D->Scale(1. / fRes[k][j]);
3571 fHistFull[k].fHistFullHar[j].fHistV0vEta->Scale(1. / fRes[k][j]);
3572 fHistFull[k].fHistFullHar[j].fHistV0vPt->Scale(1. / fRes[k][j]);
3573 contentV0 = fHistFull[k].fHistV0v->GetBinContent(j+1);
3574 contentV0 /= fRes[k][j];
3575 fHistFull[k].fHistV0v->SetBinContent(j+1, contentV0);
3576 // V0s sidebands
3577 fHistFull[k].fHistFullHar[j].fHistV0sbv2D->Scale(1. / fRes[k][j]);
3578 fHistFull[k].fHistFullHar[j].fHistV0sbvEta->Scale(1. / fRes[k][j]);
3579 fHistFull[k].fHistFullHar[j].fHistV0sbvPt->Scale(1. / fRes[k][j]);
3580 // The systematic error of the resolution is folded in.
3581 // tracks
3582 error = fHistFull[k].fHistv->GetBinError(j+1) ;
3583 error /= fRes[k][j];
3584 totalError = fabs(content) * TMath::Sqrt((error/content)*(error/content) + (fResErr[k][j]/fRes[k][j])*(fResErr[k][j]/fRes[k][j]));
3585 fHistFull[k].fHistv->SetBinError(j+1, totalError);
3586 cout << "***** v" << j+1 << "= (" << content << " +/- " << error << " +/- " << totalError << "(with syst.)) %" << endl;
3587 // V0s
3588 errorV0 = fHistFull[k].fHistV0v->GetBinError(j+1) ;
3589 errorV0 /= fRes[k][j];
3590 if (contentV0>0)
3591 totalError = fabs(contentV0) * TMath::Sqrt((errorV0/contentV0)*(errorV0/contentV0) + (fResErr[k][j]/fRes[k][j])*(fResErr[k][j]/fRes[k][j]));
3592 else
3593 totalError = 0;
3594 fHistFull[k].fHistV0v->SetBinError(j+1, totalError);
3595 }
3596 else
3597 {
3598 cout << "##### Resolution of the " << j+1 << "th harmonic was zero." << endl;
3599 // tracks hist
3600 fHistFull[k].fHistFullHar[j].fHistv2D->Reset();
3601 fHistFull[k].fHistFullHar[j].fHistvEta->Reset();
3602 fHistFull[k].fHistFullHar[j].fHistvPt->Reset();
3603 fHistFull[k].fHistv->SetBinContent(j+1, 0.);
3604 fHistFull[k].fHistv->SetBinError(j+1, 0.);
3605 // v0s hist
3606 fHistFull[k].fHistFullHar[j].fHistV0v2D->Reset();
3607 fHistFull[k].fHistFullHar[j].fHistV0vEta->Reset();
3608 fHistFull[k].fHistFullHar[j].fHistV0vPt->Reset();
3609 fHistFull[k].fHistFullHar[j].fHistV0sbv2D->Reset();
3610 fHistFull[k].fHistFullHar[j].fHistV0sbvEta->Reset();
3611 fHistFull[k].fHistFullHar[j].fHistV0sbvPt->Reset();
3612 fHistFull[k].fHistV0v->SetBinContent(j+1, 0.);
3613 fHistFull[k].fHistV0v->SetBinError(j+1, 0.);
3614 }
3615 }
3616 }
3617
3618 return kTRUE ;
3619}
3620//-----------------------------------------------------------------------
3621Double_t AliFlowAnalyser::Chi(Double_t res)
3622{
3623 // chi from the event plane resolution
3624
3625 Double_t chi = 2.0;
3626 Double_t delta = 1.0;
3627 for(int i = 0; i < 15; i++)
3628 {
3629 if(ResEventPlane(chi) < res) { chi = chi + delta ; }
3630 else { chi = chi - delta ; }
3631 delta = delta / 2.;
3632 }
3633
3634 return chi ;
3635}
3636//-----------------------------------------------------------------------
3637Double_t AliFlowAnalyser::ResEventPlane(Double_t chi)
3638{
3639 // plane resolution as function of chi
3640
3641 Double_t con = TMath::Sqrt(TMath::Pi()/2)/2 ; // ~ 0.626657
3642 Double_t arg = chi * chi / 4.;
3643 Double_t res = con * chi * exp(-arg) * (TMath::BesselI0(arg) + TMath::BesselI1(arg));
3644
3645 return res ;
3646}
3647//-----------------------------------------------------------------------
3648Double_t AliFlowAnalyser::ResEventPlaneK2(Double_t chi)
3649{
3650 // ... for the case k=2.
3651
3652 Double_t con = TMath::Sqrt(TMath::Pi()/2)/2 ; // ~ 0.626657
3653 Double_t arg = chi * chi / 4.;
3654
3655 Double_t besselOneHalf = TMath::Sqrt(arg/(TMath::Pi()/2)) * sinh(arg)/arg;
3656 Double_t besselThreeHalfs = TMath::Sqrt(arg/(TMath::Pi()/2)) * (cosh(arg)/arg - sinh(arg)/(arg*arg));
3657 Double_t res = con * chi * exp(-arg) * (besselOneHalf + besselThreeHalfs);
3658
3659 return res ;
3660}
3661//-----------------------------------------------------------------------
3662Double_t AliFlowAnalyser::ResEventPlaneK3(Double_t chi)
3663{
3664 // ...for the case k=3.
3665
3666 Double_t con = TMath::Sqrt(TMath::Pi()/2)/2 ; // ~ 0.626657
3667 Double_t arg = chi * chi / 4.;
3668 Double_t res = con * chi * exp(-arg) * (TMath::BesselI1(arg) + TMath::BesselI(2, arg));
3669
3670 return res ;
3671}
3672//-----------------------------------------------------------------------
3673Double_t AliFlowAnalyser::ResEventPlaneK4(Double_t chi)
3674{
3675 // ...for the case k=4.
3676
3677 Double_t con = TMath::Sqrt(TMath::Pi()/2)/2 ; // ~ 0.626657
3678 Double_t arg = chi * chi / 4.;
3679
3680 Double_t besselOneHalf = TMath::Sqrt(arg/(TMath::Pi()/2)) * sinh(arg)/arg;
3681 Double_t besselThreeHalfs = TMath::Sqrt(arg/(TMath::Pi()/2)) * (cosh(arg)/arg - sinh(arg)/(arg*arg));
3682 Double_t besselFiveHalfs = besselOneHalf - 3*besselThreeHalfs/arg;
3683 Double_t res = con * chi * exp(-arg) * (besselThreeHalfs + besselFiveHalfs);
3684
3685 return res ;
3686}
92016a03 3687//-----------------------------------------------------------------------
9777bfcb 3688// ###
3689//-----------------------------------------------------------------------
3690
3691
3692#endif