9f95aacf32d47dfcc8d9964541272a9e1c005535
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowAnalysisWithCumulants.cxx
1 /*************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
3 *                                                                        *
4 * Author: The ALICE Off-line Project.                                    *
5 * Contributors are mentioned in the code where appropriate.              *
6 *                                                                        *
7 * Permission to use, copy, modify and distribute this software and its   *
8 * documentation strictly for non-commercial purposes is hereby granted   *
9 * without fee, provided that the above copyright notice appears in all   *
10 * copies and that both the copyright notice and this permission notice   *
11 * appear in the supporting documentation. The authors make no claims     *
12 * about the suitability of this software for any purpose. It is          *
13 * provided "as is" without express or implied warranty.                  * 
14 **************************************************************************/
15
16 /******************************** 
17  * flow analysis with cumulants * 
18  *                              *
19  * author: Ante Bilandzic       * 
20  *          (anteb@nikhef.nl)   *
21  *******************************/ 
22
23 #define AliFlowAnalysisWithCumulants_cxx
24
25 #include "Riostream.h"
26 #include "AliFlowCommonConstants.h"
27 #include "AliFlowCommonHist.h"
28 #include "AliFlowCommonHistResults.h"
29 #include "TChain.h"
30 #include "TFile.h"
31 #include "TList.h" //NEW
32 #include "TParticle.h"
33 #include "TRandom3.h"
34 #include "TProfile.h"
35 #include "TProfile2D.h" 
36 #include "TProfile3D.h"
37 #include "AliFlowEventSimple.h"
38 #include "AliFlowTrackSimple.h"
39 #include "AliFlowAnalysisWithCumulants.h"
40 #include "AliFlowCumuConstants.h"
41 #include "AliCumulantsFunctions.h"
42
43 class TH1;
44 class TGraph;
45 class TPave;
46 class TLatex;
47 class TMarker;
48 class TRandom3;
49 class TObjArray;
50 class TList;
51 class TCanvas;
52 class TSystem;
53 class TROOT;
54 class AliFlowVector;
55 class TVector;
56
57 //================================================================================================================
58
59 ClassImp(AliFlowAnalysisWithCumulants)
60
61 AliFlowAnalysisWithCumulants::AliFlowAnalysisWithCumulants():  
62  fTrack(NULL),
63  fHistList(NULL),
64  fR0(0),
65  fPtMax(0),
66  fPtMin(0),
67  fBinWidthPt(0),
68  fgknBinsPt(0),
69  fEtaMax(0),
70  fEtaMin(0),
71  fBinWidthEta(0),
72  fgknBinsEta(0),
73  fAvQx(0),
74  fAvQy(0),
75  fAvQ2x(0),
76  fAvQ2y(0),
77  fAvMultIntFlowGFC(NULL),
78  fQVectorComponentsGFC(NULL),
79  fIntFlowResultsGFC(NULL),
80  fDiffFlowResults2ndOrderGFC(NULL),
81  fDiffFlowResults4thOrderGFC(NULL),
82  fDiffFlowResults6thOrderGFC(NULL),
83  fDiffFlowResults8thOrderGFC(NULL),
84  fCommonHistsResults2nd(NULL),
85  fCommonHistsResults4th(NULL),
86  fCommonHistsResults6th(NULL),
87  fCommonHistsResults8th(NULL),
88  fIntFlowGenFun(NULL),
89  fIntFlowGenFun4(NULL),//(only for other system of Eq.)
90  fIntFlowGenFun6(NULL),//(only for other system of Eq.)
91  fIntFlowGenFun8(NULL),//(only for other system of Eq.)
92  fIntFlowGenFun16(NULL),//(only for other system of Eq.)
93  fAvMultIntFlow4GFC(NULL),//(only for other system of Eq.)
94  fAvMultIntFlow6GFC(NULL),//(only for other system of Eq.)
95  fAvMultIntFlow8GFC(NULL),//(only for other system of Eq.)
96  fAvMultIntFlow16GFC(NULL),//(only for other system of Eq.)
97  fDiffFlowPtRPGenFunRe(NULL),
98  fDiffFlowPtRPGenFunIm(NULL),
99  fPtBinRPNoOfParticles(NULL),
100  fDiffFlowEtaRPGenFunRe(NULL),
101  fDiffFlowEtaRPGenFunIm(NULL),
102  fEtaBinRPNoOfParticles(NULL),
103  fDiffFlowPtPOIGenFunRe(NULL),
104  fDiffFlowPtPOIGenFunIm(NULL),
105  fPtBinPOINoOfParticles(NULL),
106  fDiffFlowEtaPOIGenFunRe(NULL),
107  fDiffFlowEtaPOIGenFunIm(NULL),
108  fEtaBinPOINoOfParticles(NULL),
109  /*
110  fDiffFlowGenFunRe0(NULL),
111  fDiffFlowGenFunRe1(NULL),
112  fDiffFlowGenFunRe2(NULL),
113  fDiffFlowGenFunRe3(NULL),
114  fDiffFlowGenFunRe4(NULL),
115  fDiffFlowGenFunRe5(NULL),
116  fDiffFlowGenFunRe6(NULL),
117  fDiffFlowGenFunRe7(NULL),
118  fDiffFlowGenFunIm0(NULL),
119  fDiffFlowGenFunIm1(NULL),
120  fDiffFlowGenFunIm2(NULL),
121  fDiffFlowGenFunIm3(NULL),
122  fDiffFlowGenFunIm4(NULL),
123  fDiffFlowGenFunIm5(NULL),
124  fDiffFlowGenFunIm6(NULL),
125  fDiffFlowGenFunIm7(NULL),
126  */
127  fCommonHists(NULL),
128  fOtherEquations(kFALSE)
129 {
130  //constructor 
131  fHistList = new TList();
132  fR0=AliFlowCumuConstants::fgR0;
133  //Pt:
134  fPtMax=AliFlowCommonConstants::GetPtMax(); 
135  fPtMin=AliFlowCommonConstants::GetPtMin();
136  fgknBinsPt=AliFlowCommonConstants::GetNbinsPt();
137  if(fgknBinsPt)
138  {
139   fBinWidthPt=(fPtMax-fPtMin)/fgknBinsPt;  
140  } 
141  //Eta: 
142  fEtaMax=AliFlowCommonConstants::GetEtaMax(); 
143  fEtaMin=AliFlowCommonConstants::GetEtaMin();
144  fgknBinsEta=AliFlowCommonConstants::GetNbinsEta();
145  if(fgknBinsEta)
146  {
147   fBinWidthEta=(fEtaMax-fEtaMin)/fgknBinsEta;  
148  } 
149  
150  fOtherEquations=AliFlowCumuConstants::fgOtherEquations;
151 }
152
153 AliFlowAnalysisWithCumulants::~AliFlowAnalysisWithCumulants()
154 {
155  //desctructor
156  delete fHistList; 
157 }
158
159 //================================================================================================================
160
161 void AliFlowAnalysisWithCumulants::CreateOutputObjects()
162 {
163  //various output histograms
164  
165  //average multiplicity 
166  fAvMultIntFlowGFC = new TProfile("fAvMultIntFlowGFC","Average Multiplicity",1,0,1,"s");
167  fAvMultIntFlowGFC->SetXTitle("");
168  fAvMultIntFlowGFC->SetYTitle("");
169  fAvMultIntFlowGFC->SetLabelSize(0.06);
170  fAvMultIntFlowGFC->SetMarkerStyle(25);
171  fAvMultIntFlowGFC->SetLabelOffset(0.01);
172  (fAvMultIntFlowGFC->GetXaxis())->SetBinLabel(1,"Average Multiplicity");
173  fHistList->Add(fAvMultIntFlowGFC);
174  
175  //averages of Q-vector components (1st bin: <Q_x>, 2nd bin: <Q_y>, 3rd bin: <(Q_x)^2>, 4th bin: <(Q_y)^2>)
176  fQVectorComponentsGFC = new TProfile("fQVectorComponentsGFC","Average of Q-vector components",4,0.,4.);
177  fQVectorComponentsGFC->SetXTitle("");
178  fQVectorComponentsGFC->SetYTitle("");
179  fQVectorComponentsGFC->SetLabelSize(0.06);
180  fQVectorComponentsGFC->SetMarkerStyle(25);
181  (fQVectorComponentsGFC->GetXaxis())->SetBinLabel(1,"<Q_{x}>");
182  (fQVectorComponentsGFC->GetXaxis())->SetBinLabel(2,"<Q_{y}>");
183  (fQVectorComponentsGFC->GetXaxis())->SetBinLabel(3,"<Q_{x}^{2}>");
184  (fQVectorComponentsGFC->GetXaxis())->SetBinLabel(4,"<Q_{y}^{2}>");
185  fHistList->Add(fQVectorComponentsGFC);
186  
187  //final results for integrated flow (v_n{2}, v_n{4},..., v_n{16}) from cumulants (by default n=2) 
188  fIntFlowResultsGFC = new TH1D("fIntFlowResultsGFC","Integrated Flow From Cumulants (Generating Function)",4,0,4);
189  fIntFlowResultsGFC->SetXTitle("");
190  fIntFlowResultsGFC->SetYTitle("");
191  fIntFlowResultsGFC->SetLabelSize(0.06);
192  //fIntFlowResultsGFC->SetTickLength(1);
193  fIntFlowResultsGFC->SetMarkerStyle(25);
194  (fIntFlowResultsGFC->GetXaxis())->SetBinLabel(1,"v_{n}{2}");
195  (fIntFlowResultsGFC->GetXaxis())->SetBinLabel(2,"v_{n}{4}");
196  (fIntFlowResultsGFC->GetXaxis())->SetBinLabel(3,"v_{n}{6}");
197  (fIntFlowResultsGFC->GetXaxis())->SetBinLabel(4,"v_{n}{8}");
198  fHistList->Add(fIntFlowResultsGFC);
199   
200  //final results for differential flow v_p/n{2} (by default p=n=2)
201  fDiffFlowResults2ndOrderGFC = new TH1D("fDiffFlowResults2ndOrderGFC","v'_2/2{2}",fgknBinsPt,fPtMin,fPtMax);
202  fDiffFlowResults2ndOrderGFC->SetXTitle("pt [GeV]");
203  fDiffFlowResults2ndOrderGFC->SetYTitle("");
204  fHistList->Add(fDiffFlowResults2ndOrderGFC);
205
206  //final results for differential flow v_p/n{4} (by default p=n=2) 
207  fDiffFlowResults4thOrderGFC = new TH1D("fDiffFlowResults4thOrderGFC","v'_2/2{4}",fgknBinsPt,fPtMin,fPtMax);
208  fDiffFlowResults4thOrderGFC->SetXTitle("pt [GeV]");
209  fDiffFlowResults4thOrderGFC->SetYTitle("");
210  fHistList->Add(fDiffFlowResults4thOrderGFC);
211  
212  //final results for differential flow v_p/n{6} (by default p=n=2)  
213  fDiffFlowResults6thOrderGFC = new TH1D("fDiffFlowResults6thOrderGFC","v'_2/2{6}",fgknBinsPt,fPtMin,fPtMax);
214  fDiffFlowResults6thOrderGFC->SetXTitle("pt [GeV]");
215  fDiffFlowResults6thOrderGFC->SetYTitle("");
216  fHistList->Add(fDiffFlowResults6thOrderGFC);
217  
218  //final results for differential flow v_p/n{8} (by default p=n=2)
219  fDiffFlowResults8thOrderGFC = new TH1D("fDiffFlowResults8thOrderGFC","v'_2/2{8}",fgknBinsPt,fPtMin,fPtMax);
220  fDiffFlowResults8thOrderGFC->SetXTitle("pt [GeV]");
221  fDiffFlowResults8thOrderGFC->SetYTitle("");
222  fHistList->Add(fDiffFlowResults8thOrderGFC);
223   
224  //avarage of the generating function for integrated flow <G[p][q]>
225  fIntFlowGenFun = new TProfile2D("fIntFlowGenFun","<G[p][q]>",fgkPmax,0.,(Double_t)fgkPmax,fgkQmax,0.,(Double_t)fgkQmax);
226  fIntFlowGenFun->SetXTitle("p");
227  fIntFlowGenFun->SetYTitle("q");
228  fHistList->Add(fIntFlowGenFun);
229
230  if(fOtherEquations)
231  {
232   //avarage of the generating function for integrated flow <G[p][q]> (only for other system of Eq. - up to 4th order)
233   fIntFlowGenFun4 = new TProfile2D("fIntFlowGenFun4","<G4[p4][q4]>",fgkPmax4,0.,(Double_t)fgkPmax4,fgkQmax4,0.,(Double_t)fgkQmax4);
234   fIntFlowGenFun4->SetXTitle("p4");
235   fIntFlowGenFun4->SetYTitle("q4");
236   fHistList->Add(fIntFlowGenFun4);
237
238   //avarage of the generating function for integrated flow <G[p][q]> (only for other system of Eq. - up to 6th order) 
239   fIntFlowGenFun6 = new TProfile2D("fIntFlowGenFun6","<G6[p6][q6]>",fgkPmax6,0.,(Double_t)fgkPmax6,fgkQmax6,0.,(Double_t)fgkQmax6);
240   fIntFlowGenFun6->SetXTitle("p6");
241   fIntFlowGenFun6->SetYTitle("q6");
242   fHistList->Add(fIntFlowGenFun6);
243
244   //avarage of the generating function for integrated flow <G[p][q]> (only for other system of Eq. - up to 8th order)
245   fIntFlowGenFun8 = new TProfile2D("fIntFlowGenFun8","<G8[p8][q8]>",fgkPmax8,0.,(Double_t)fgkPmax8,fgkQmax8,0.,(Double_t)fgkQmax8);
246   fIntFlowGenFun8->SetXTitle("p8");
247   fIntFlowGenFun8->SetYTitle("q8");
248   fHistList->Add(fIntFlowGenFun8);
249
250   //avarage of the generating function for integrated flow <G[p][q]> (only for other system of Eq. - up to 16th order)
251   fIntFlowGenFun16 = new TProfile2D("fIntFlowGenFun16","<G16[p16][q16]>",fgkPmax16,0.,(Double_t)fgkPmax16,fgkQmax16,0.,(Double_t)fgkQmax16);
252   fIntFlowGenFun16->SetXTitle("p16");
253   fIntFlowGenFun16->SetYTitle("q16");
254   fHistList->Add(fIntFlowGenFun16);
255  
256   //average multiplicity (only for other system of Eq. - up to 4th order)
257   fAvMultIntFlow4GFC = new TProfile("fAvMultIntFlow4GFC","Average Multiplicity",1,0,1,"s");
258   fAvMultIntFlow4GFC->SetXTitle("");
259   fAvMultIntFlow4GFC->SetYTitle("");
260   fAvMultIntFlow4GFC->SetLabelSize(0.06);
261   fAvMultIntFlow4GFC->SetMarkerStyle(25);
262   fAvMultIntFlow4GFC->SetLabelOffset(0.01);
263   (fAvMultIntFlow4GFC->GetXaxis())->SetBinLabel(1,"Average Multiplicity");
264   fHistList->Add(fAvMultIntFlow4GFC);
265  
266   //average multiplicity (only for other system of Eq. - up to 6th order)
267   fAvMultIntFlow6GFC = new TProfile("fAvMultIntFlow6GFC","Average Multiplicity",1,0,1,"s");
268   fAvMultIntFlow6GFC->SetXTitle("");
269   fAvMultIntFlow6GFC->SetYTitle("");
270   fAvMultIntFlow6GFC->SetLabelSize(0.06);
271   fAvMultIntFlow6GFC->SetMarkerStyle(25);
272   fAvMultIntFlow6GFC->SetLabelOffset(0.01);
273   (fAvMultIntFlow6GFC->GetXaxis())->SetBinLabel(1,"Average Multiplicity");
274   fHistList->Add(fAvMultIntFlow6GFC);
275  
276   //average multiplicity (only for other system of Eq. - up to 8th order)
277   fAvMultIntFlow8GFC = new TProfile("fAvMultIntFlow8GFC","Average Multiplicity",1,0,1,"s");
278   fAvMultIntFlow8GFC->SetXTitle("");
279   fAvMultIntFlow8GFC->SetYTitle("");
280   fAvMultIntFlow8GFC->SetLabelSize(0.06);
281   fAvMultIntFlow8GFC->SetMarkerStyle(25);
282   fAvMultIntFlow8GFC->SetLabelOffset(0.01);
283   (fAvMultIntFlow8GFC->GetXaxis())->SetBinLabel(1,"Average Multiplicity");
284   fHistList->Add(fAvMultIntFlow8GFC);
285  
286   //average multiplicity (only for other system of Eq. - up to 16th order)
287   fAvMultIntFlow16GFC = new TProfile("fAvMultIntFlow16GFC","Average Multiplicity",1,0,1,"s");
288   fAvMultIntFlow16GFC->SetXTitle("");
289   fAvMultIntFlow16GFC->SetYTitle("");
290   fAvMultIntFlow16GFC->SetLabelSize(0.06);
291   fAvMultIntFlow16GFC->SetMarkerStyle(25);
292   fAvMultIntFlow16GFC->SetLabelOffset(0.01);
293   (fAvMultIntFlow16GFC->GetXaxis())->SetBinLabel(1,"Average Multiplicity");
294   fHistList->Add(fAvMultIntFlow16GFC);
295  }
296  
297  //avarage of the real part of generating function for differential flow in Pt <Re(D[b][p][q])>
298  fDiffFlowPtRPGenFunRe = new TProfile3D("fDiffFlowPtRPGenFunRe","<Re(D[b][p][q])>",fgknBinsPt,(Double_t)(fPtMin/fBinWidthPt),(Double_t)(fPtMax/fBinWidthPt),fgkPmax,0.,(Double_t)fgkPmax,fgkQmax,0.,(Double_t)fgkQmax);
299  fDiffFlowPtRPGenFunRe->SetXTitle("b");
300  fDiffFlowPtRPGenFunRe->SetYTitle("p");
301  fDiffFlowPtRPGenFunRe->SetZTitle("q");
302  fDiffFlowPtRPGenFunRe->SetTitleOffset(1.44,"X");
303  fDiffFlowPtRPGenFunRe->SetTitleOffset(1.44,"Y");
304  fHistList->Add(fDiffFlowPtRPGenFunRe);
305  
306  //avarage of the imaginary part of generating function for differential flow in Pt <Im(D[b][p][q])>
307  fDiffFlowPtRPGenFunIm = new TProfile3D("fDiffFlowPtRPGenFunIm","<Im(D[b][p][q])>",fgknBinsPt,(Double_t)(fPtMin/fBinWidthPt),(Double_t)(fPtMax/fBinWidthPt),fgkPmax,0.,(Double_t)fgkPmax,fgkQmax,0.,(Double_t)fgkQmax);
308  fDiffFlowPtRPGenFunIm->SetXTitle("b");
309  fDiffFlowPtRPGenFunIm->SetYTitle("p");
310  fDiffFlowPtRPGenFunIm->SetZTitle("q");
311  fDiffFlowPtRPGenFunIm->SetTitleOffset(1.44,"X");
312  fDiffFlowPtRPGenFunIm->SetTitleOffset(1.44,"Y");
313  fHistList->Add(fDiffFlowPtRPGenFunIm);
314  
315  //number of particles per pt bin
316  fPtBinRPNoOfParticles = new TProfile("fPtBinRPNoOfParticles","Number of particles per #p_{t} bin",fgknBinsPt,fPtMin,fPtMax);
317  fPtBinRPNoOfParticles->SetXTitle("pt [GeV]");
318  fPtBinRPNoOfParticles->SetYTitle("");
319  fHistList->Add(fPtBinRPNoOfParticles);
320  
321  //avarage of the real part of generating function for differential flow in Eta <Re(D[b][p][q])>
322  fDiffFlowEtaRPGenFunRe = new TProfile3D("fDiffFlowEtaRPGenFunRe","<Re(D[b][p][q])>",fgknBinsEta,(Double_t)(fEtaMin/fBinWidthEta),(Double_t)(fEtaMax/fBinWidthEta),fgkPmax,0.,(Double_t)fgkPmax,fgkQmax,0.,(Double_t)fgkQmax);
323  fDiffFlowEtaRPGenFunRe->SetXTitle("b");
324  fDiffFlowEtaRPGenFunRe->SetYTitle("p");
325  fDiffFlowEtaRPGenFunRe->SetZTitle("q");
326  fDiffFlowEtaRPGenFunRe->SetTitleOffset(1.44,"X");
327  fDiffFlowEtaRPGenFunRe->SetTitleOffset(1.44,"Y");
328  fHistList->Add(fDiffFlowEtaRPGenFunRe);
329  
330  //avarage of the imaginary part of generating function for differential flow in Eta <Im(D[b][p][q])>
331  fDiffFlowEtaRPGenFunIm = new TProfile3D("fDiffFlowEtaRPGenFunIm","<Im(D[b][p][q])>",fgknBinsEta,(Double_t)(fEtaMin/fBinWidthEta),(Double_t)(fEtaMax/fBinWidthEta),fgkPmax,0.,(Double_t)fgkPmax,fgkQmax,0.,(Double_t)fgkQmax);
332  fDiffFlowEtaRPGenFunIm->SetXTitle("b");
333  fDiffFlowEtaRPGenFunIm->SetYTitle("p");
334  fDiffFlowEtaRPGenFunIm->SetZTitle("q");
335  fDiffFlowEtaRPGenFunIm->SetTitleOffset(1.44,"X");
336  fDiffFlowEtaRPGenFunIm->SetTitleOffset(1.44,"Y");
337  fHistList->Add(fDiffFlowEtaRPGenFunIm);
338  
339  //number of particles per eta bin
340  fEtaBinRPNoOfParticles = new TProfile("fEtaBinRPNoOfParticles","Number of particles per #eta bin",fgknBinsEta,fEtaMin,fEtaMax);
341  fEtaBinRPNoOfParticles->SetXTitle("#eta");
342  fEtaBinRPNoOfParticles->SetYTitle("");
343  fHistList->Add(fEtaBinRPNoOfParticles);
344  
345  //avarage of the real part of generating function for differential flow in Pt <Re(D[b][p][q])>
346  fDiffFlowPtPOIGenFunRe = new TProfile3D("fDiffFlowPtPOIGenFunRe","<Re(D[b][p][q])>",fgknBinsPt,(Double_t)(fPtMin/fBinWidthPt),(Double_t)(fPtMax/fBinWidthPt),fgkPmax,0.,(Double_t)fgkPmax,fgkQmax,0.,(Double_t)fgkQmax);
347  fDiffFlowPtPOIGenFunRe->SetXTitle("b");
348  fDiffFlowPtPOIGenFunRe->SetYTitle("p");
349  fDiffFlowPtPOIGenFunRe->SetZTitle("q");
350  fDiffFlowPtPOIGenFunRe->SetTitleOffset(1.44,"X");
351  fDiffFlowPtPOIGenFunRe->SetTitleOffset(1.44,"Y");
352  fHistList->Add(fDiffFlowPtPOIGenFunRe);
353  
354  //avarage of the imaginary part of generating function for differential flow in Pt <Im(D[b][p][q])>
355  fDiffFlowPtPOIGenFunIm = new TProfile3D("fDiffFlowPtPOIGenFunIm","<Im(D[b][p][q])>",fgknBinsPt,(Double_t)(fPtMin/fBinWidthPt),(Double_t)(fPtMax/fBinWidthPt),fgkPmax,0.,(Double_t)fgkPmax,fgkQmax,0.,(Double_t)fgkQmax);
356  fDiffFlowPtPOIGenFunIm->SetXTitle("b");
357  fDiffFlowPtPOIGenFunIm->SetYTitle("p");
358  fDiffFlowPtPOIGenFunIm->SetZTitle("q");
359  fDiffFlowPtPOIGenFunIm->SetTitleOffset(1.44,"X");
360  fDiffFlowPtPOIGenFunIm->SetTitleOffset(1.44,"Y");
361  fHistList->Add(fDiffFlowPtPOIGenFunIm);
362  
363  //number of particles per pt bin
364  fPtBinPOINoOfParticles = new TProfile("fPtBinPOINoOfParticles","Number of particles per #p_{t} bin",fgknBinsPt,fPtMin,fPtMax);
365  fPtBinPOINoOfParticles->SetXTitle("pt [GeV]");
366  fPtBinPOINoOfParticles->SetYTitle("");
367  fHistList->Add(fPtBinPOINoOfParticles);
368  
369  //avarage of the real part of generating function for differential flow in Eta <Re(D[b][p][q])>
370  fDiffFlowEtaPOIGenFunRe = new TProfile3D("fDiffFlowEtaPOIGenFunRe","<Re(D[b][p][q])>",fgknBinsEta,(Double_t)(fEtaMin/fBinWidthEta),(Double_t)(fEtaMax/fBinWidthEta),fgkPmax,0.,(Double_t)fgkPmax,fgkQmax,0.,(Double_t)fgkQmax);
371  fDiffFlowEtaPOIGenFunRe->SetXTitle("b");
372  fDiffFlowEtaPOIGenFunRe->SetYTitle("p");
373  fDiffFlowEtaPOIGenFunRe->SetZTitle("q");
374  fDiffFlowEtaPOIGenFunRe->SetTitleOffset(1.44,"X");
375  fDiffFlowEtaPOIGenFunRe->SetTitleOffset(1.44,"Y");
376  fHistList->Add(fDiffFlowEtaPOIGenFunRe);
377  
378  //avarage of the imaginary part of generating function for differential flow in Eta <Im(D[b][p][q])>
379  fDiffFlowEtaPOIGenFunIm = new TProfile3D("fDiffFlowEtaPOIGenFunIm","<Im(D[b][p][q])>",fgknBinsEta,(Double_t)(fEtaMin/fBinWidthEta),(Double_t)(fEtaMax/fBinWidthEta),fgkPmax,0.,(Double_t)fgkPmax,fgkQmax,0.,(Double_t)fgkQmax);
380  fDiffFlowEtaPOIGenFunIm->SetXTitle("b");
381  fDiffFlowEtaPOIGenFunIm->SetYTitle("p");
382  fDiffFlowEtaPOIGenFunIm->SetZTitle("q");
383  fDiffFlowEtaPOIGenFunIm->SetTitleOffset(1.44,"X");
384  fDiffFlowEtaPOIGenFunIm->SetTitleOffset(1.44,"Y");
385  fHistList->Add(fDiffFlowEtaPOIGenFunIm);
386  
387  //number of particles per eta bin
388  fEtaBinPOINoOfParticles = new TProfile("fEtaBinPOINoOfParticles","Number of particles per #eta bin",fgknBinsEta,fEtaMin,fEtaMax);
389  fEtaBinPOINoOfParticles->SetXTitle("#eta");
390  fEtaBinPOINoOfParticles->SetYTitle("");
391  fHistList->Add(fEtaBinPOINoOfParticles);
392  
393  /*
394  fDiffFlowGenFunRe0 = new TProfile2D("fDiffFlowGenFunRe0","Re(<D[b][0][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
395  fDiffFlowGenFunRe0->SetXTitle("b");
396  fDiffFlowGenFunRe0->SetYTitle("q");
397  fHistList->Add(fDiffFlowGenFunRe0);
398
399  fDiffFlowGenFunIm0 = new TProfile2D("fDiffFlowGcout<<"HEY M1"<<endl;enFunIm0","Im(<D[b][0][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
400  fDiffFlowGenFunIm0->SetXTitle("b");
401  fDiffFlowGenFunIm0->SetYTitle("q");
402  fHistList->Add(fDiffFlowGenFunIm0);
403  
404  fDiffFlowGenFunRe1 = new TProfile2D("fDiffFlowGenFunRe1","Re(<D[b][1][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
405  fDiffFlowGenFunRe1->SetXTitle("b");
406  fDiffFlowGenFunRe1->SetYTitle("q");
407  fHistList->Add(fDiffFlowGenFunRe1);
408  
409  fDiffFlowGenFunIm1 = new TProfile2D("fDiffFlowGenFunIm1","Im(<D[b][1][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
410  fDiffFlowGenFunIm1->SetXTitle("b");
411  fDiffFlowGenFunIm1->SetYTitle("q");
412  fHistList->Add(fDiffFlowGenFunIm1);
413  
414  fDiffFlowGenFunRe2 = new TProfile2D("fDiffFlowGenFunRe2","Re(<D[b][2][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
415  fDiffFlowGenFunRe2->SetXTitle("b");
416  fDiffFlowGenFunRe2->SetYTitle("q");
417  fHistList->Add(fDiffFlowGenFunRe2);
418  
419  fDiffFlowGenFunIm2 = new TProfile2D("fDiffFlowGenFunIm2","Im(<D[b][2][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
420  fDiffFlowGenFunIm2->SetXTitle("b");
421  fDiffFlowGenFunIm2->SetYTitle("q");
422  fHistList->Add(fDiffFlowGenFunIm2);
423  
424  fDiffFlowGenFunRe3 = new TProfile2D("fDiffFlowGenFunRe3","Re(<D[b][3][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
425  fDiffFlowGenFunRe3->SetXTitle("b");
426  fDiffFlowGenFunRe3->SetYTitle("q");
427  fHistList->Add(fDiffFlowGenFunRe3);
428  
429  fDiffFlowGenFunIm3 = new TProfile2D("fDiffFlowGenFunIm3","Im(<D[b][3][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
430  fDiffFlowGenFunIm3->SetXTitle("b");
431  fDiffFlowGenFunIm3->SetYTitle("q");
432  fHistList->Add(fDiffFlowGenFunIm3);
433  
434  fDiffFlowGenFunRe4 = new TProfile2D("fDiffFlowGenFunRe4","Re(<D[b][4][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
435  fDiffFlowGenFunRe4->SetXTitle("b");
436  fDiffFlowGenFunRe4->SetYTitle("q");
437  fHistList->Add(fDiffFlowGenFunRe4);
438  
439  fDiffFlowGenFunIm4 = new TProfile2D("fDiffFlowGenFunIm4","Im(<D[b][4][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
440  fDiffFlowGenFunIm4->SetXTitle("b");
441  fDiffFlowGenFunIm4->SetYTitle("q");
442  fHistList->Add(fDiffFlowGenFunIm4);
443  
444  fDiffFlowGenFunRe5 = new TProfile2D("fDiffFlowGenFunRe5","Re(<D[b][5][q]>)",fgkQmax,0.,(Double_t)fgkQmax,fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth));
445  fDiffFlowGenFunRe5->SetXTitle("b");
446  fDiffFlowGenFunRe5->SetYTitle("q");
447  fHistList->Add(fDiffFlowGenFunRe5);
448  
449  fDiffFlowGenFunIm5 = new TProfile2D("fDiffFlowGenFunIm5","Im(<D[b][5][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
450  fDiffFlowGenFunIm5->SetXTitle("b");
451  fDiffFlowGenFunIm5->SetYTitle("q");
452  fHistList->Add(fDiffFlowGenFunIm5);
453  
454  fDiffFlowGenFunRe6 = new TProfile2D("fDiffFlowGenFunRe6","Re(<D[b][6][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
455  fDiffFlowGenFunRe6->SetXTitle("b");
456  fDiffFlowGenFunRe6->SetYTitle("q");
457  fHistList->Add(fDiffFlowGenFunRe6);
458  
459  fDiffFlowGenFunIm6 = new TProfile2D("fDiffFlowGenFunIm6","Im(<D[b][6][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
460  fDiffFlowGenFunIm6->SetXTitle("b");
461  fDiffFlowGenFunIm6->SetYTitle("q");
462  fHistList->Add(fDiffFlowGenFunIm6);
463  
464  fDiffFlowGenFunRe7 = new TProfile2D("fDiffFlowGenFunRe7","Re(<D[b][7][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
465  fDiffFlowGenFunRe7->SetXTitle("b");
466  fDiffFlowGenFunRe7->SetYTitle("q");
467  fHistList->Add(fDiffFlowGenFunRe7);
468  
469  fDiffFlowGenFunIm7 = new TProfile2D("fDiffFlowGenFunIm7","Im(<D[b][7][q]>)",fgknBins,(Double_t)(fPtMin/fBinWidth),(Double_t)(fPtMax/fBinWidth),fgkQmax,0.,(Double_t)fgkQmax);
470  fDiffFlowGenFunIm7->SetXTitle("b");
471  fDiffFlowGenFunIm7->SetYTitle("q");
472  fHistList->Add(fDiffFlowGenFunIm7);
473  */
474  
475  //common control histograms
476  fCommonHists = new AliFlowCommonHist("AliFlowCommonHistGFC");
477  fHistList->Add(fCommonHists);  
478  
479  //common histograms for final results (2nd order)
480  fCommonHistsResults2nd = new AliFlowCommonHistResults("AliFlowCommonHistResults2ndOrderGFC");
481  fHistList->Add(fCommonHistsResults2nd); 
482  
483  //common histograms for final results (4th order)
484  fCommonHistsResults4th = new AliFlowCommonHistResults("AliFlowCommonHistResults4thOrderGFC");
485  fHistList->Add(fCommonHistsResults4th);
486  
487  //common histograms for final results (6th order)
488  fCommonHistsResults6th = new AliFlowCommonHistResults("AliFlowCommonHistResults6thOrderGFC");
489  fHistList->Add(fCommonHistsResults6th); 
490  
491  //common histograms for final results (8th order)
492  fCommonHistsResults8th = new AliFlowCommonHistResults("AliFlowCommonHistResults8thOrderGFC");
493  fHistList->Add(fCommonHistsResults8th);
494   
495 }//end of CreateOutputObjects()
496
497 //================================================================================================================
498
499 void AliFlowAnalysisWithCumulants::Make(AliFlowEventSimple* anEvent)
500 {
501  //running over data:
502  Int_t nPrim = anEvent->NumberOfTracks(); //total multiplicity
503  
504  Int_t nEventNSelTracksIntFlow = anEvent->GetEventNSelTracksIntFlow(); //selected multiplicity (particles used for int. flow)
505   
506  if(nEventNSelTracksIntFlow>9) //generating function formalism applied here make sense only for selected multiplicity >= 10 
507  { 
508  //fill the common control histograms
509  fCommonHists->FillControlHistograms(anEvent);   
510   
511  //initializing the generating function G[p][q] for integrated flow 
512  Double_t G[fgkPmax][fgkQmax];
513  
514  for(Int_t p=0;p<fgkPmax;p++)
515   {
516    for(Int_t q=0;q<fgkQmax;q++)
517    {
518     G[p][q]=1.;
519    }   
520   }
521
522  Int_t nSelTracksIntFlow = 0; //cross-checking the selected multiplicity
523  
524  for(Int_t i=0;i<nPrim;i++)
525  {
526   fTrack=anEvent->GetTrack(i);
527   if(fTrack && fTrack->UseForIntegratedFlow())
528   {
529    nSelTracksIntFlow++;
530    for(Int_t p=0;p<fgkPmax;p++)
531    {
532     for(Int_t q=0;q<fgkQmax;q++)
533     {
534      G[p][q]*=(1.+(2.*fR0*sqrt(p+1.)/nEventNSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax));
535     }
536    } 
537   }
538  }  
539  
540  //storing the value of G[p][q] in 2D profile in order to get automatically the avarage <G[p][q]>
541  for(Int_t p=0;p<fgkPmax;p++)
542  {
543   for(Int_t q=0;q<fgkQmax;q++)
544   {
545    fIntFlowGenFun->Fill((Double_t)p,(Double_t)q,G[p][q],1);
546   }
547  } 
548  
549  //storing the selected multiplicity (if fAvMultIntFlow is not filled here then you had wrongly selected the particles used for integrated flow)
550  if(nSelTracksIntFlow==nEventNSelTracksIntFlow)
551  {
552   fAvMultIntFlowGFC->Fill(0.,nSelTracksIntFlow,1);
553  }
554  
555  //calculating Q-vector of event (needed for errors)
556  AliFlowVector fQVector;
557  fQVector.Set(0.,0.);
558  fQVector.SetMult(0);
559  fQVector=anEvent->GetQ();                               //get the Q vector for this event
560  fQVectorComponentsGFC->Fill(0.,fQVector.X(),1);         //in the 1st bin fill Q_x
561  fQVectorComponentsGFC->Fill(1.,fQVector.Y(),1);         //in the 2nd bin fill Q_y
562  fQVectorComponentsGFC->Fill(2.,pow(fQVector.X(),2.),1); //in the 3rd bin fill (Q_x)^2
563  fQVectorComponentsGFC->Fill(3.,pow(fQVector.Y(),2.),1); //in the 4th bin fill (Q_y)^2
564  
565  //3D profiles for differential flow in pt and eta
566  //second loop over data: evaluating the generating function D[b][p][q] for differential flow 
567  //remark 0: D[b][p][q] is a complex number => real and imaginary part are calculated separately
568  //remark 1: note that bellow G[p][q] is needed, the value of generating function for integrated flow for the CURRENT event
569  //remark 2: results are stored in 3D profiles in order to automatically get <Re(D[b][p][q])> and <Im(D[b][p][q])>
570  for(Int_t i=0;i<nPrim;i++)
571  {
572   fTrack=anEvent->GetTrack(i);
573   if(fTrack)
574   {
575    if(fTrack->UseForIntegratedFlow() && fTrack->UseForDifferentialFlow())
576    {
577     fPtBinPOINoOfParticles->Fill(fTrack->Pt(),fTrack->Pt(),1.);
578     fEtaBinPOINoOfParticles->Fill(fTrack->Eta(),fTrack->Eta(),1.);
579     for(Int_t p=0;p<fgkPmax;p++)
580     {
581      for(Int_t q=0;q<fgkQmax;q++)
582      {
583       //real part (Pt)
584       fDiffFlowPtPOIGenFunRe->Fill(fTrack->Pt()/fBinWidthPt,(Double_t)p,(Double_t)q,G[p][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(p+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
585       //imaginary part (Pt)
586       fDiffFlowPtPOIGenFunIm->Fill(fTrack->Pt()/fBinWidthPt,(Double_t)p,(Double_t)q,G[p][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(p+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.); 
587       //real part (Eta)
588       fDiffFlowEtaPOIGenFunRe->Fill(fTrack->Eta()/fBinWidthEta,(Double_t)p,(Double_t)q,G[p][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(p+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
589       //imaginary part (Eta)
590       fDiffFlowEtaPOIGenFunIm->Fill(fTrack->Eta()/fBinWidthEta,(Double_t)p,(Double_t)q,G[p][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(p+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.); 
591      }
592     }
593    } 
594    else if(fTrack->UseForDifferentialFlow() && !(fTrack->UseForIntegratedFlow()))
595    {
596     fPtBinPOINoOfParticles->Fill(fTrack->Pt(),fTrack->Pt(),1.);
597     fEtaBinPOINoOfParticles->Fill(fTrack->Eta(),fTrack->Eta(),1.);
598     for(Int_t p=0;p<fgkPmax;p++)
599     {
600      for(Int_t q=0;q<fgkQmax;q++)
601      {
602       //real part (Pt)
603       fDiffFlowPtPOIGenFunRe->Fill(fTrack->Pt()/fBinWidthPt,(Double_t)p,(Double_t)q,G[p][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi()),1.);
604       //imaginary part (Pt)
605       fDiffFlowPtPOIGenFunIm->Fill(fTrack->Pt()/fBinWidthPt,(Double_t)p,(Double_t)q,G[p][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi()),1.); 
606       //real part (Eta)
607       fDiffFlowEtaPOIGenFunRe->Fill(fTrack->Eta()/fBinWidthEta,(Double_t)p,(Double_t)q,G[p][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi()),1.);
608       //imaginary part (Eta)
609       fDiffFlowEtaPOIGenFunIm->Fill(fTrack->Eta()/fBinWidthEta,(Double_t)p,(Double_t)q,G[p][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi()),1.);                     
610      } 
611     }
612    }
613    //RP:
614    if(fTrack->UseForIntegratedFlow())       
615    {
616     fPtBinRPNoOfParticles->Fill(fTrack->Pt(),fTrack->Pt(),1.);
617     fEtaBinRPNoOfParticles->Fill(fTrack->Eta(),fTrack->Eta(),1.);
618     for(Int_t p=0;p<fgkPmax;p++)
619     {
620      for(Int_t q=0;q<fgkQmax;q++)
621      {
622       //real part (Pt)
623       fDiffFlowPtRPGenFunRe->Fill(fTrack->Pt()/fBinWidthPt,(Double_t)p,(Double_t)q,G[p][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(p+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
624       //imaginary part (Pt)
625       fDiffFlowPtRPGenFunIm->Fill(fTrack->Pt()/fBinWidthPt,(Double_t)p,(Double_t)q,G[p][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(p+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.); 
626       //real part (Eta)
627       fDiffFlowEtaRPGenFunRe->Fill(fTrack->Eta()/fBinWidthEta,(Double_t)p,(Double_t)q,G[p][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(p+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
628       //imaginary part (Eta)
629       fDiffFlowEtaRPGenFunIm->Fill(fTrack->Eta()/fBinWidthEta,(Double_t)p,(Double_t)q,G[p][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(p+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.); 
630      }
631     }
632    }//end of if(fTrack->UseForIntegratedFlow())                   
633   }//end of if(fTrack)  
634  }//ending the second loop over data    
635  
636  
637  
638  /*
639  //sixteen 2D profiles for differential flow          
640  for(Int_t i=0;i<nPrim;i++){
641   fTrack=anEvent->GetTrack(i);
642   if (fTrack && fTrack->UseForDifferentialFlow()){
643    //for(Int_t p=0;p<fgkPmax;p++){
644     for(Int_t q=0;q<fgkQmax;q++){
645      //real part
646      fDiffFlowGenFunRe0->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[0][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(0.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
647      //imaginary part
648      fDiffFlowGenFunIm0->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[0][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(0.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.); 
649      //-----------------------------------------------------------------------
650      //real part
651      fDiffFlowGenFunRe1->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[1][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(1.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
652      //imaginary part
653      fDiffFlowGenFunIm1->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[1][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(1.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);             
654     //-----------------------------------------------------------------------
655      //real part
656      fDiffFlowGenFunRe2->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[2][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(2.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
657      //imaginary part
658      fDiffFlowGenFunIm2->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[2][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(2.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);  
659      //-----------------------------------------------------------------------
660      //real part
661      fDiffFlowGenFunRe3->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[3][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(3.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
662      //imaginary part
663      fDiffFlowGenFunIm3->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[3][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(3.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);  
664      //-----------------------------------------------------------------------
665      //real part
666      fDiffFlowGenFunRe4->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[4][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(4.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
667      //imaginary part
668      fDiffFlowGenFunIm4->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[4][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(4.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);  
669      //-----------------------------------------------------------------------
670      //real part
671      fDiffFlowGenFunRe5->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[5][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(5.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
672      //imaginary part
673      fDiffFlowGenFunIm5->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[5][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(5.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);  
674      //-----------------------------------------------------------------------
675      //real part
676      fDiffFlowGenFunRe6->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[6][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(6.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
677      //imaginary part
678      fDiffFlowGenFunIm6->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[6][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(6.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
679      //-----------------------------------------------------------------------
680      //real part
681      fDiffFlowGenFunRe7->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[7][q]*cos(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(7.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);
682      //imaginary part
683      fDiffFlowGenFunIm7->Fill(fTrack->Pt()/fBinWidth,(Double_t)q,G[7][q]*sin(fgkMltpl*fgkFlow*fTrack->Phi())/(1.+(2.*fR0*sqrt(7.+1.)/nSelTracksIntFlow)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax)),1.);    
684    }
685    //}
686   }  
687  }//ending the second loop over data            
688  */
689       
690  }//end of if(nEventNSelTracksIntFlow>9)                   
691  
692
693  
694   
695    
696     
697      
698       
699        
700         
701          
702           
703            
704  //off the record: numerical equations for cumulants solved up to different highest order  
705  if(fOtherEquations)
706  {
707   //running over data
708   Int_t nPrimOE = anEvent->NumberOfTracks();//total multiplicity 
709   
710   Int_t nEventNSelTracksIntFlowOE = anEvent->GetEventNSelTracksIntFlow();
711   
712   Double_t G4[fgkPmax4][fgkQmax4];
713   Double_t G6[fgkPmax6][fgkQmax6];
714   Double_t G8[fgkPmax8][fgkQmax8];
715   Double_t G16[fgkPmax16][fgkQmax16];
716   for(Int_t p=0;p<fgkPmax16;p++)
717   {
718    for(Int_t q=0;q<fgkQmax16;q++)
719    {
720     G16[p][q]=1.;
721     if(p<fgkPmax8 && q<fgkQmax8)
722     {
723      G8[p][q]=1.;
724      if(p<fgkPmax6 && q<fgkQmax6)
725      {
726       G6[p][q]=1.;
727       if(p<fgkPmax4 && q<fgkQmax4)
728       {
729        G4[p][q]=1.;
730       }
731      }
732     } 
733    }
734   } 
735    
736   //multiplicities: 
737   if(nEventNSelTracksIntFlowOE>15) fAvMultIntFlow16GFC->Fill(0.,nEventNSelTracksIntFlowOE,1);
738   if(nEventNSelTracksIntFlowOE>7) fAvMultIntFlow8GFC->Fill(0.,nEventNSelTracksIntFlowOE,1);
739   if(nEventNSelTracksIntFlowOE>5) fAvMultIntFlow6GFC->Fill(0.,nEventNSelTracksIntFlowOE,1);
740   if(nEventNSelTracksIntFlowOE>3) fAvMultIntFlow4GFC->Fill(0.,nEventNSelTracksIntFlowOE,1);  
741   
742   //first loop over data: evaluating the generating function G[p][q] for integrated flow 
743   for(Int_t i=0;i<nPrimOE;i++)
744   {
745    fTrack=anEvent->GetTrack(i);
746    if(fTrack && fTrack->UseForIntegratedFlow())
747    {
748     for(Int_t p=0;p<fgkPmax16;p++)
749     {
750      for(Int_t q=0;q<fgkQmax16;q++)
751      {
752       if(nEventNSelTracksIntFlowOE>15)
753       {
754        G16[p][q]*=(1.+(2.*fR0*sqrt(p+1.)/nEventNSelTracksIntFlowOE)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax16));
755       }       
756       if(p<fgkPmax8 && q<fgkQmax8)
757       {
758        if(nEventNSelTracksIntFlowOE>7)
759        { 
760         G8[p][q]*=(1.+(2.*fR0*sqrt(p+1.)/nEventNSelTracksIntFlowOE)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax8));
761        }
762        if(p<fgkPmax6 && q<fgkQmax6)
763        {
764         if(nEventNSelTracksIntFlowOE>5) 
765         {
766          G6[p][q]*=(1.+(2.*fR0*sqrt(p+1.)/nEventNSelTracksIntFlowOE)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax6));
767         }
768         if(p<fgkPmax4 && q<fgkQmax4)
769         {
770          if(nEventNSelTracksIntFlowOE>3)
771          {
772           G4[p][q]*=(1.+(2.*fR0*sqrt(p+1.)/nEventNSelTracksIntFlowOE)*cos(fgkFlow*fTrack->Phi()-2.*q*TMath::Pi()/fgkQmax4));
773          } 
774         }
775        }
776       }
777      } 
778     }  
779    }//end of if(fTrack && fTrack->UseForIntegratedFlow())
780   }//ending the loop over data
781  
782  //storing the value of G[p][q] in 2D profile in order to get automatically the avarage <G[p][q]>
783  for(Int_t p=0;p<fgkPmax16;p++)
784  {
785   for(Int_t q=0;q<fgkQmax16;q++)
786   {
787    if(nEventNSelTracksIntFlowOE>15) fIntFlowGenFun16->Fill((Double_t)p,(Double_t)q,G16[p][q],1);
788    if(p<fgkPmax8 && q<fgkQmax8)
789    {
790     if(nEventNSelTracksIntFlowOE>7) fIntFlowGenFun8->Fill((Double_t)p,(Double_t)q,G8[p][q],1);
791     if(p<fgkPmax6 && q<fgkQmax6)
792     {
793      if(nEventNSelTracksIntFlowOE>5) fIntFlowGenFun6->Fill((Double_t)p,(Double_t)q,G6[p][q],1);
794      if(p<fgkPmax4 && q<fgkQmax4)
795      {
796       if(nEventNSelTracksIntFlowOE>3) fIntFlowGenFun4->Fill((Double_t)p,(Double_t)q,G4[p][q],1);
797      }
798     }
799    } 
800   }
801  }
802 }//end of if(fOtherEquations)                  
803                                                                                                                                                                                                                                                                                                                                                                                                                                                                        
804 }//end of Make()
805
806 //================================================================================================================
807
808 void AliFlowAnalysisWithCumulants::Finish()
809 {
810  //calculate the final results
811  //AliCumulantsFunctions finalResults(fIntFlowGenFun,NULL,NULL, fIntFlowResults,fDiffFlowResults2,fDiffFlowResults4,fDiffFlowResults6,fDiffFlowResults8,fAvMultIntFlow,fQVectorComponents,  fQDist,fDiffFlowGenFunRe0,fDiffFlowGenFunRe1,fDiffFlowGenFunRe2, fDiffFlowGenFunRe3,fDiffFlowGenFunRe4,fDiffFlowGenFunRe5,fDiffFlowGenFunRe6,fDiffFlowGenFunRe7,fDiffFlowGenFunIm0,fDiffFlowGenFunIm1, fDiffFlowGenFunIm2,fDiffFlowGenFunIm3,fDiffFlowGenFunIm4,fDiffFlowGenFunIm5,fDiffFlowGenFunIm6,fDiffFlowGenFunIm7);
812
813  AliCumulantsFunctions finalResults(fIntFlowGenFun,fIntFlowGenFun4,fIntFlowGenFun6,fIntFlowGenFun8,fIntFlowGenFun16,fAvMultIntFlow4GFC, fAvMultIntFlow6GFC,fAvMultIntFlow8GFC,fAvMultIntFlow16GFC,fDiffFlowPtRPGenFunRe,fDiffFlowPtRPGenFunIm,fPtBinRPNoOfParticles, fDiffFlowEtaRPGenFunRe,fDiffFlowEtaRPGenFunIm,fEtaBinRPNoOfParticles,fDiffFlowPtPOIGenFunRe,fDiffFlowPtPOIGenFunIm,fPtBinPOINoOfParticles, fDiffFlowEtaPOIGenFunRe,fDiffFlowEtaPOIGenFunIm,fEtaBinPOINoOfParticles,fIntFlowResultsGFC,fDiffFlowResults2ndOrderGFC,fDiffFlowResults4thOrderGFC,fDiffFlowResults6thOrderGFC,fDiffFlowResults8thOrderGFC, fAvMultIntFlowGFC,fQVectorComponentsGFC,fCommonHistsResults2nd, fCommonHistsResults4th,fCommonHistsResults6th,fCommonHistsResults8th,fCommonHists);
814                            
815  finalResults.Calculate();  
816 }
817
818 //================================================================================================================
819
820 void AliFlowAnalysisWithCumulants::WriteHistograms(TString* outputFileName)
821 {
822  //store the final results in output .root file
823  TFile *output = new TFile(outputFileName->Data(),"RECREATE");
824  output->WriteObject(fHistList, "cobjGFC","SingleKey");
825  delete output;
826 }
827
828 //================================================================================================================
829
830
831
832