AliPhysics  63e47e1 (63e47e1)
AliAnalysisTaskStudentsML.cxx
Go to the documentation of this file.
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 * template class for student projects *
18 * Marcel Lesch (marcel.lesch@cern.ch) *
19 **********************************************/
20 
21 #include "Riostream.h"
23 #include "AliLog.h"
24 #include "AliAODEvent.h"
25 #include "AliAODInputHandler.h"
26 #include "AliAnalysisManager.h"
27 #include "AliMultSelection.h"
28 #include "TRandom.h"
29 #include "TComplex.h"
30 #include "TProfile.h"
31 #include "TStyle.h"
32 #include "TFile.h"
33 #include "TH1F.h"
34 #include <TArrayF.h>
35 #include <vector>
36 #include "TMath.h"
37 #include "TF1.h"
38 
39 
40 
41 using std::cout;
42 using std::endl;
43 
45 
46 //================================================================================================================
47 
48 AliAnalysisTaskStudentsML::AliAnalysisTaskStudentsML(const char *name, Bool_t useParticleWeights):
49  AliAnalysisTaskSE(name),
50  fHistList(NULL),
51  // Control histograms:
52  fControlHistogramsList(NULL),
53  fPtHist(NULL),
54  fNbins(1000),
55  fMinBin(0.),
56  fMaxBin(10.),
57  fPhiHist(NULL),
58  fEtaHist(NULL),
59  fMultiHisto(NULL),
60  fMaxCorrelator(8),
61  bUseWeights(kFALSE),
62  kNumber(2), //number of correlation
63  kh1(2), kh2(-2), kh3(0), kh4(0), kh5(0), kh6(0), kh7(0), kh8(0), //harmonics
64  kSum((kh1<0?-1*kh1:kh1)+(kh2<0?-1*kh2:kh2)+(kh3<0?-1*kh3:kh3)+(kh4<0?-1*kh4:kh4)
65  + (kh5<0?-1*kh5:kh5)+(kh6<0?-1*kh6:kh6)+(kh7<0?-1*kh7:kh7)+(kh8<0?-1*kh8:kh8)), // We will not go beyond 8-p correlations
66  kMaxHarmonic(kSum+1),
67  kMaxPower(fMaxCorrelator+1),
68  fParticles(0),
69  fCentral(0),
70  fMinCentrality(0.),
71  fMaxCentrality(100.),
72  fAngles(NULL),
73  fWeights(NULL),
74  fBin(NULL),
75  func1(NULL),
76  fCentrality(NULL),
77  fCounterHistogram(NULL),
78  // Final results:
79  fFinalResultsList(NULL)
80  {
81  // Constructor.
82 
83  AliDebug(2,"AliAnalysisTaskStudentsML::AliAnalysisTaskStudentsML(const char *name, Bool_t useParticleWeights)");
84 
85  // Base list:
86  fHistList = new TList();
87  fHistList->SetName("outputStudentAnalysis");
88  fHistList->SetOwner(kTRUE);
89 
90  // Initialize all arrays:
91  this->InitializeArrays();
92 
93  // Define input and output slots here
94  // Input slot #0 works with an AliFlowEventSimple
95  //DefineInput(0, AliFlowEventSimple::Class());
96  // Input slot #1 is needed for the weights input file:
97  //if(useParticleWeights)
98  //{
99  // DefineInput(1, TList::Class());
100  //}
101  // Output slot #0 is reserved
102  // Output slot #1 writes into a TList container
103 
104  DefineOutput(1, TList::Class());
105 
106  if(useParticleWeights)
107  {
108  // not needed for the time being
109  }
110 
111 } // AliAnalysisTaskStudentsML::AliAnalysisTaskStudentsML(const char *name, Bool_t useParticleWeights):
112 
113 //================================================================================================================
114 
117  fHistList(NULL),
118  // Control histograms:
119  fControlHistogramsList(NULL),
120  fPtHist(NULL),
121  fNbins(1000),
122  fMinBin(0.),
123  fMaxBin(10.),
124  fPhiHist(NULL),
125  fEtaHist(NULL),
126  fMultiHisto(NULL),
127  fMaxCorrelator(8),
128  bUseWeights(kFALSE),
129  kNumber(2), //number of correlation
130  kh1(2), kh2(-2), kh3(0), kh4(0), kh5(0), kh6(0), kh7(0), kh8(0), //harmonics
131  kSum((kh1<0?-1*kh1:kh1)+(kh2<0?-1*kh2:kh2)+(kh3<0?-1*kh3:kh3)+(kh4<0?-1*kh4:kh4)
132  + (kh5<0?-1*kh5:kh5)+(kh6<0?-1*kh6:kh6)+(kh7<0?-1*kh7:kh7)+(kh8<0?-1*kh8:kh8)), // We will not go beyond 8-p correlations
133  kMaxHarmonic(kSum+1),
134  kMaxPower(fMaxCorrelator+1),
135  fParticles(0),
136  fCentral(0.),
137  fMinCentrality(0.),
138  fMaxCentrality(100.),
139  fAngles(NULL),
140  fWeights(NULL),
141  fBin(NULL),
142  func1(NULL),
143  // Final results:
144  fCentrality(NULL),
145  fCounterHistogram(NULL),
146  fFinalResultsList(NULL)
147 {
148  // Dummy constructor.
149 
150  AliDebug(2,"AliAnalysisTaskStudentsML::AliAnalysisTaskStudentsML()");
151 
152 } // AliAnalysisTaskStudentsML::AliAnalysisTaskStudentsML():
153 
154 //================================================================================================================
155 
157 {
158  // Destructor.
159 
160  if(fHistList) delete fHistList;
161 
162 } // AliAnalysisTaskStudentsML::~AliAnalysisTaskStudentsML()
163 
164 //================================================================================================================
165 
167 {
168  // Called at every worker node to initialize.
169 
170  // a) Trick to avoid name clashes, part 1;
171  // b) Book and nest all lists;
172  // c) Book all objects;
173  // *) Trick to avoid name clashes, part 2.
174  // d)
175 
176  // a) Trick to avoid name clashes, part 1:
177  Bool_t oldHistAddStatus = TH1::AddDirectoryStatus();
178  TH1::AddDirectory(kFALSE);
179 
180  // b) Book and nest all lists:
181  this->BookAndNestAllLists();
182 
183  // c) Book all objects:
184  this->BookControlHistograms();
186 
187  // *) Trick to avoid name clashes, part 2:
188  TH1::AddDirectory(oldHistAddStatus);
189 
190  PostData(1,fHistList);
191 
192 } // void AliAnalysisTaskStudentsML::UserCreateOutputObjects()
193 
194 //================================================================================================================
195 
197 {
198  // Main loop (called for each event).
199 
200  // a.0) Get pointer to AOD event;
201  // a.1) Centrality;
202  // b) Start analysis over AODs;
203  // c) Reset event-by-event objects;
204  // d) PostData.
205 
206  fCounterHistogram->Fill(0.5); // counter hist 1st bin
207 
208  // a) Get pointer to AOD event:
209  AliAODEvent *aAOD = dynamic_cast<AliAODEvent*>(InputEvent()); // from TaskSE
210  if(!aAOD){return;}
211 
212  fCounterHistogram->Fill(1.5); // counter hist 2nd bin
213 
214  //a.1) Centrality;
215 
216  AliMultSelection *ams = (AliMultSelection*)aAOD->FindListObject("MultSelection");
217  if(!ams){return;}
218 
219  fCounterHistogram->Fill(2.5); // counter hist 3rd bin
220 
221  //func1->SetParameter(2,gRandom->Uniform(0.,TMath::TwoPi())); //*** for testing. for each event psi_2 is different
222 
223 
224  // b) Start analysis over AODs:
225  Int_t nTracks = aAOD->GetNumberOfTracks(); // number of all tracks in current event
226 
227  //fParticles=gRandom->Uniform(50.0,500.0); /*** for testing.
228 
229  Int_t nCounter=0;
230 
231 
232  //first loop: =======================
233 
234  for(Int_t iTrack=0;iTrack<nTracks;iTrack++) // starting a loop over all tracks
235  {
236  AliAODTrack *aTrack = dynamic_cast<AliAODTrack*>(aAOD->GetTrack(iTrack)); // getting a pointer to a track
237  if(!aTrack){continue;} // protection against NULL pointers
238  if(!aTrack->TestFilterBit(128)){continue;} // filter bit 128 denotes TPC-only tracks, use only them for the analysis
239 
240  // example variables for each track:
241  Double_t px = aTrack->Px(); // x-component of momenta
242  Double_t py = aTrack->Py(); // y-component of momenta
243  Double_t pz = aTrack->Pz(); // z-component of momenta
244  Double_t e = aTrack->E(); // energy
245  Double_t phi = aTrack->Phi(); // azimuthal angle
246  Double_t eta = aTrack->Eta(); // pseudorapidity
247  Double_t charge = aTrack->Charge(); // charge
248  Double_t pt = aTrack->Pt(); // Pt (transverse momentum)
249 
250 
251  // apply some cuts: e.g. take for the analysis only particles in -0.8 < eta < 0.8, and 0.2 < pT < 5.0
252  // ... implementation of particle cuts ...
253  if( (-0.8 < eta) && (eta < 0.8) && (0.2 < pt) && (pt < 5.0) )
254  {
255  // Fill some control histograms with the particles which passed cuts:
256  fPtHist->Fill(pt);
257  fPhiHist->Fill(phi);
258  fEtaHist->Fill(eta);
259 
260  fParticles += 1; //one more particle passed the cuts
261 
262 
263  // Do some analysis only with the particles which passed the cuts
264  // ... your analysis code ...
265 
266 
267 
268  } // if( (-0.8 < eta) && (eta < 0.8) && (0.2 < pt) && (pt < 5.0) )
269  } // for(Int_t iTrack=0;iTrack<nTracks;iTrack++) // starting a loop over all tracks
270 
271  fAngles = new TArrayD(fParticles);
272 
273  //second loop: =======================
274 
275  for(Int_t iTrack=0;iTrack<nTracks;iTrack++) // starting a loop over all tracks
276  {
277  AliAODTrack *aTrack = dynamic_cast<AliAODTrack*>(aAOD->GetTrack(iTrack)); // getting a pointer to a track
278  if(!aTrack){continue;} // protection against NULL pointers
279  if(!aTrack->TestFilterBit(128)){continue;} // filter bit 128 denotes TPC-only tracks, use only them for the analysis
280 
281  // example variables for each track:
282  Double_t px = aTrack->Px(); // x-component of momenta
283  Double_t py = aTrack->Py(); // y-component of momenta
284  Double_t pz = aTrack->Pz(); // z-component of momenta
285  Double_t e = aTrack->E(); // energy
286  Double_t phi = aTrack->Phi(); // azimuthal angle
287  Double_t eta = aTrack->Eta(); // pseudorapidity
288  Double_t charge = aTrack->Charge(); // charge
289  Double_t pt = aTrack->Pt(); // Pt (transverse momentum)
290 
291 
292  // apply some cuts: e.g. take for the analysis only particles in -0.8 < eta < 0.8, and 0.2 < pT < 5.0
293  // ... implementation of particle cuts ...
294  if( (-0.8 < eta) && (eta < 0.8) && (0.2 < pt) && (pt < 5.0) )
295  {
296 
297  fAngles->AddAt(phi,nCounter);
298  nCounter += 1;
299 
300  } // if( (-0.8 < eta) && (eta < 0.8) && (0.2 < pt) && (pt < 5.0) )
301  } // for(Int_t iTrack=0;iTrack<nTracks;iTrack++) // starting a loop over all tracks
302 
303 
304  //b.1) analysis
305  if(fParticles>0){fMultiHisto->Fill(fParticles);}
306 
307  //*** for testing.
308  /*fAngles = new TArrayD(fParticles);
309 
310  for(Int_t k=0;k<fParticles;k++)
311  {
312  fAngles->AddAt(func1->GetRandom(),nCounter);
313  nCounter += 1;
314  } */
315 
316  if(fParticles>=kNumber)
317  {
318 
319  if(ams->GetMultiplicityPercentile("V0M") >= fMinCentrality && ams->GetMultiplicityPercentile("V0M") < fMaxCentrality)
320  {
321  if(0. <= (ams->GetMultiplicityPercentile("V0M")) && (ams->GetMultiplicityPercentile("V0M")) < 5. ){fCentral=2.5;}
322  if(5. <= (ams->GetMultiplicityPercentile("V0M")) && (ams->GetMultiplicityPercentile("V0M")) < 10. ){fCentral=7.5;}
323  if(10. <= (ams->GetMultiplicityPercentile("V0M")) && (ams->GetMultiplicityPercentile("V0M")) < 20. ){fCentral=15.5;}
324  if(20. <= (ams->GetMultiplicityPercentile("V0M")) && (ams->GetMultiplicityPercentile("V0M")) < 30. ){fCentral=25.5;}
325  if(30. <= (ams->GetMultiplicityPercentile("V0M")) && (ams->GetMultiplicityPercentile("V0M")) < 40. ){fCentral=35.5;}
326  if(40. <= (ams->GetMultiplicityPercentile("V0M")) && (ams->GetMultiplicityPercentile("V0M")) < 50. ){fCentral=45.5;}
327  if(50. <= (ams->GetMultiplicityPercentile("V0M")) && (ams->GetMultiplicityPercentile("V0M")) < 60. ){fCentral=55.5;}
328  if(60. <= (ams->GetMultiplicityPercentile("V0M")) && (ams->GetMultiplicityPercentile("V0M")) < 70. ){fCentral=65.5;}
329  if(70. <= (ams->GetMultiplicityPercentile("V0M")) && (ams->GetMultiplicityPercentile("V0M")) < 80. ){fCentral=75.5;}
330  }
331  else
332  {
333  return; // this event do not belong to the centrality class specified for this particular analysis
334  }
335 
336 
337  Correlation(); //do the correlation
338  fCentrality->Fill(fCentral,fRecursion[0][kNumber-2]->GetBinContent(1),fRecursion[0][kNumber-2]->GetBinContent(2)); //safe output
339 
340 
341  fRecursion[0][kNumber-2]->Reset(); //Reset
342  fRecursion[1][kNumber-2]->Reset(); //Reset
343 
344 
345  //Testing with nested loops
346 
347  /* for(Int_t m=0;m<fParticles;m++) //making all possible pairs for 500 particles
348  {
349  for(Int_t k=0;k<fParticles;k++)
350  {
351  if(m==k){continue;} //no autocorellation
352  else
353  {
354  fCentrality->Fill(55.5,TMath::Cos(2.0*( (fAngles->GetAt(m)) - (fAngles->GetAt(k)) ) )); //fill TProfile
355  }
356  } //for(Int_t k=0;k<fParticles;k++)
357  } //for(Int_t m=0;m<fParticles;m++)
358  */
359 
360 
361  } //if(fParticles>=kNumber)
362 
363 
364 
365 
366 
367  // c) Reset event-by-event objects:
368  fParticles=0;
369  fCentral=0.;
370  nCounter=0;
371  delete fAngles;
372  fAngles=NULL;
373 
374 
375 
376  // d) PostData:
377  PostData(1,fHistList);
378 
379 } // void AliAnalysisTaskStudentsML::UserExec(Option_t *)
380 
381 //================================================================================================================
382 
384 {
385  // Accessing the merged output list.
386 
387  fHistList = (TList*)GetOutputData(1);
388  if(!fHistList){exit(1);}
389 
390  // Do some calculation in offline mode here:
391  // ...
392 
393  TFile *f = new TFile("AnalysisResults.root","RECREATE");
394  fHistList->Write(fHistList->GetName(),TObject::kSingleKey);
395 
396  delete f;
397 
398 } // end of void AliAnalysisTaskStudentsML::Terminate(Option_t *)
399 
400 //================================================================================================================
401 
403 {
404  // Initialize all data members which are arrays in this method.
405 
406  for(Int_t cs=0;cs<2;cs++)
407  {
408  for(Int_t c=0;c<fMaxCorrelator;c++)
409  {
410 
411  fRecursion[cs][c] = NULL;
412 
413  }
414  } //for(Int_t cs=0;cs<2;cs++)
415 
416  for(Int_t js=0;js<17;js++)
417  {
418  for(Int_t j=0;j<9;j++)
419  {
420 
421  Qvector[js][j] = NULL;
422 
423  }
424  }
425 
426 
427 
428 } // void AliAnalysisTaskStudentsML::InitializeArrays()
429 
430 //=======================================================================================================================
431 
433 {
434  // Book and nest all lists nested in the base list fHistList.
435 
436  // a) Book and nest lists for control histograms;
437  // b) Book and nest lists for final results.
438 
439  TString sMethodName = "void AliAnalysisTaskStudentsML::BookAndNestAllLists()";
440  if(!fHistList){Fatal(sMethodName.Data(),"fHistList is NULL");}
441 
442  // a) Book and nest lists for control histograms:
444  fControlHistogramsList->SetName("ControlHistograms");
445  fControlHistogramsList->SetOwner(kTRUE);
446  fHistList->Add(fControlHistogramsList);
447 
448  // b) Book and nest lists for final results:
449  fFinalResultsList = new TList();
450  fFinalResultsList->SetName("FinalResults");
451  fFinalResultsList->SetOwner(kTRUE);
453 
454 } // void AliAnalysisTaskStudentsML::BookAndNestAllLists()
455 
456 //=======================================================================================================================
457 
459 {
460  // Book all control histograms.
461 
462  // a) Book histogram to hold pt distribution;
463  // b) Book histogram to hold phi distribution;
464  // c) Book histogram to hold eta distribution;
465  // d) Book histogram to hold multiplicty distribution;
466  // e) Book histogram to debug
467 
468 
469 
470 
471  // a) Book histogram to hold pt spectra:
472  fPtHist = new TH1F("fPtHist","atrack->Pt()",fNbins,fMinBin,fMaxBin);
473  fPtHist->SetStats(kFALSE);
474  fPtHist->SetFillColor(kBlue-10);
475  fPtHist->GetXaxis()->SetTitle("p_{t}");
477 
478  // b) Book histogram to hold phi distribution:
479  fPhiHist = new TH1F("fPhiHist","Phi Distribution",1000,0.,6.3);
480  fPhiHist->GetXaxis()->SetTitle("Phi");
481  fPhiHist->SetLineColor(4);
483 
484  // c) Book histogram to hold eta distribution:
485  fEtaHist = new TH1F("fEtaHist","Eta Distribution",1000,-1.,1.);
486  fEtaHist->GetXaxis()->SetTitle("Eta");
487  fEtaHist->SetLineColor(4);
489 
490  // d) Book histogam to hold multiplicty distribution:
491  fMultiHisto = new TH1F("fMultiHisto","Multiplicity",500,0.,500.); //histogram for multiplicity
492  fMultiHisto->GetXaxis()->SetTitle("Multiplicity M");
494 
495  // e) Book histogram to debug
496  fCounterHistogram = new TH1F("fCounterHistogram","Histogram for some checks",3,0.,3.); //histogram for multiplicity
498 
499 } //void AliAnalysisTaskStudentsML::BookControlHistograms()
500 
501 //=======================================================================================================================
502 
504 {
505  // Book all histograms to hold the final results.
506 
507  fCentrality = new TProfile("fCentrality","Result Analysis Centrality Dependence",80,0.,80.); //histogram for multiplicity
508  fCentrality->GetXaxis()->SetTitle("Centrality");
509  fCentrality->GetYaxis()->SetTitle("flow");
511  fCentrality->Sumw2();
512 
513 
514 
515  Cosmetics();
516 
517  //*** for testing
518  /*func1=new TF1("func1","(1/(TMath::TwoPi()))*(1 + 2*[1]*TMath::Cos(2*(x-[2])))",0.,TMath::TwoPi()); //defining function func1 with parameter [1] for v_2 and [2] as pis_2
519  func1->SetParameter(1,0.05); //setting v_2 (aka [1]) fix to 0.05*/
520 
521 } // void AliAnalysisTaskStudentsML::BookFinalResultsHistograms()
522 
523 //=======================================================================================================================
524 
526 {
527  // Book everything here.
528 
529  for(Int_t cs=0;cs<2;cs++)
530  {
531  for(Int_t c=0;c<fMaxCorrelator;c++)
532  {
533 
534  fRecursion[cs][c] = new TProfile("","",2,0.,2.);
535  fRecursion[cs][c]->Sumw2();
536 
537  //NOTE for fRecursion: 1.) [cs] will say if its the real (0) or imaginary part (1)
538  // 2.) [c] gives gives the kind of correlation. [n] is the (n+2)-particle correlation
539  //3.) first bin holds value of correlation (single event). Second bin holds the weight!
540 
541 
542  } // end of for(Int_t c=0;c<fMaxCorrelator;c++)
543  } // end of for(Int_t cs=0;cs<2;cs++)
544 
545 } // void Cosmetics()
546 
547 //=======================================================================================================================
548 
550 {
551  // Calculate Q-vectors.
552 
553  // a) Make sure all Q-vectors are initially zero;
554  // b) Calculate Q-vectors for available angles and weights.
555 
556  // a) Make sure all Q-vectors are initially zero:
557  for(Int_t h=0;h<kMaxHarmonic;h++)
558  {
559  for(Int_t p=0;p<kMaxPower;p++)
560  {
561  Qvector[h][p] = TComplex(0.,0.);
562  } // for(Int_t p=0;p<kMaxPower;p++)
563  } // for(Int_t h=0;h<kMaxHarmonic;h++)
564 
565  // b) Calculate Q-vectors for available angles and weights:
566  Double_t dPhi2 = 0.; // particle angle
567  Double_t wPhi = 1.; // particle weight
568  Double_t wPhiToPowerP = 1.; // particle weight raised to power p
569  for(Int_t i=0;i<fParticles;i++) // loop over particles
570  {
571  dPhi2 = fAngles->GetAt(i);
572  if(bUseWeights){wPhi = fWeights->GetAt(i);}
573  for(Int_t h=0;h<kMaxHarmonic;h++)
574  {
575  for(Int_t p=0;p<kMaxPower;p++)
576  {
577  if(bUseWeights){wPhiToPowerP = pow(wPhi,p);}
578  Qvector[h][p] += TComplex(wPhiToPowerP*TMath::Cos(h*dPhi2),wPhiToPowerP*TMath::Sin(h*dPhi2));
579  } // for(Int_t p=0;p<kMaxPower;p++)
580  } // for(Int_t h=0;h<kMaxHarmonic;h++)
581  } // for(Int_t i=0;i<fParticles;i++) // loop over particles
582 
583 } // void CalculateQvectors()
584 
585 //=======================================================================================================================
586 
587 
589 {
590  // Using the fact that Q{-n,p} = Q{n,p}^*.
591 
592  if(n>=0){return Qvector[n][p];}
593  return TComplex::Conjugate(Qvector[-n][p]);
594 
595 } // TComplex AliAnalysisTaskStudentsML::Q(Int_t n, Int_t p)
596 
597 
598 
599 
600 //========================================================================================================================
601 
602 
603 TComplex AliAnalysisTaskStudentsML::Recursion(Int_t n, Int_t* harmonic, Int_t mult = 1, Int_t skip = 0)
604 {
605  // Calculate multi-particle correlators by using recursion (an improved faster version) originally developed by
606  // Kristjan Gulbrandsen (gulbrand@nbi.dk).
607 
608  Int_t nm1 = n-1;
609  TComplex c(Q(harmonic[nm1], mult));
610  if (nm1 == 0) return c;
611  c *= Recursion(nm1, harmonic);
612  if (nm1 == skip) return c;
613 
614  Int_t multp1 = mult+1;
615  Int_t nm2 = n-2;
616  Int_t counter1 = 0;
617  Int_t hhold = harmonic[counter1];
618  harmonic[counter1] = harmonic[nm2];
619  harmonic[nm2] = hhold + harmonic[nm1];
620  TComplex c2(Recursion(nm1, harmonic, multp1, nm2));
621  Int_t counter2 = n-3;
622  while (counter2 >= skip) {
623  harmonic[nm2] = harmonic[counter1];
624  harmonic[counter1] = hhold;
625  ++counter1;
626  hhold = harmonic[counter1];
627  harmonic[counter1] = harmonic[nm2];
628  harmonic[nm2] = hhold + harmonic[nm1];
629  c2 += Recursion(nm1, harmonic, multp1, counter2);
630  --counter2;
631  }
632  harmonic[nm2] = harmonic[counter1];
633  harmonic[counter1] = hhold;
634 
635  if (mult == 1) return c-c2;
636  return c-Double_t(mult)*c2;
637 
638 } // TComplex AliFlowAnalysisWithMultiparticleCorrelations::Recursion(Int_t n, Int_t* harmonic, Int_t mult = 1, Int_t skip = 0)
639 
640 
641 
642 //========================================================================================================================
643 
644 
646 {
647 
648 
649  // Calculate Q-vectors for available angles and weights;
651 
652  // Calculate n-particle correlations from Q-vectors (using recursion):
653 
654 
655  if(kNumber==2)
656  {
657  Int_t harmonicsTwoNum[2] = {kh1,kh2};
658  Int_t harmonicsTwoDen[2] = {0,0};
659  TComplex twoRecursion = Recursion(2,harmonicsTwoNum)/Recursion(2,harmonicsTwoDen).Re();
660  Double_t wTwoRecursion = Recursion(2,harmonicsTwoDen).Re();
661  fRecursion[0][0]->Fill(0.5,twoRecursion.Re(),wTwoRecursion); // <<cos(h1*phi1+h2*phi2)>>
662  fRecursion[0][0]->Fill(1.5,wTwoRecursion,1.); //weight
663  fRecursion[1][0]->Fill(0.5,twoRecursion.Im(),wTwoRecursion); // <<sin(h1*phi1+h2*phi2)>>
664  fRecursion[1][0]->Fill(1.5,wTwoRecursion,1.); //weight
665  }// 2-p correlation
666 
667  if(kNumber==3)
668  {
669  Int_t harmonicsThreeNum[3] = {kh1,kh2,kh3};
670  Int_t harmonicsThreeDen[3] = {0,0,0};
671  TComplex threeRecursion = Recursion(3,harmonicsThreeNum)/Recursion(3,harmonicsThreeDen).Re();
672  Double_t wThreeRecursion = Recursion(3,harmonicsThreeDen).Re();
673  fRecursion[0][1]->Fill(0.5,threeRecursion.Re(),wThreeRecursion); // <<cos(h1*phi1+h2*phi2+h3*phi3)>>
674  fRecursion[0][1]->Fill(1.5,wThreeRecursion,1.); //weight
675  fRecursion[1][1]->Fill(0.5,threeRecursion.Im(),wThreeRecursion); // <<sin(h1*phi1+h2*phi2+h3*phi3)>>
676  fRecursion[1][1]->Fill(1.5,wThreeRecursion,1.); //weight
677  } // 3-p correlation
678 
679  if(kNumber==4)
680  {
681  Int_t harmonicsFourNum[4] = {kh1,kh2,kh3,kh4};
682  Int_t harmonicsFourDen[4] = {0,0,0,0};
683  TComplex fourRecursion = Recursion(4,harmonicsFourNum)/Recursion(4,harmonicsFourDen).Re();
684  Double_t wFourRecursion = Recursion(4,harmonicsFourDen).Re();
685  fRecursion[0][2]->Fill(0.5,fourRecursion.Re(),wFourRecursion); // <<cos(h1*phi1+h2*phi2+h3*phi3+h4*phi4)>>
686  fRecursion[0][2]->Fill(1.5,wFourRecursion,1.); //weigh
687  fRecursion[1][2]->Fill(0.5,fourRecursion.Im(),wFourRecursion); // <<<sin(h1*phi1+h2*phi2+h3*phi3+h4*phi4)>>
688  fRecursion[1][2]->Fill(1.5,wFourRecursion,1.); //weigh
689  }// 4-p correlation
690 
691  if(kNumber==5)
692  {
693  Int_t harmonicsFiveNum[5] = {kh1,kh2,kh3,kh4,kh5};
694  Int_t harmonicsFiveDen[5] = {0,0,0,0,0};
695  TComplex fiveRecursion = Recursion(5,harmonicsFiveNum)/Recursion(5,harmonicsFiveDen).Re();
696  Double_t wFiveRecursion = Recursion(5,harmonicsFiveDen).Re();
697  fRecursion[0][3]->Fill(0.5,fiveRecursion.Re(),wFiveRecursion); // <<cos(h1*phi1+h2*phi2+h3*phi3+h4*phi4+h5*phi5)>>
698  fRecursion[0][3]->Fill(1.5,wFiveRecursion,1.);
699  fRecursion[1][3]->Fill(0.5,fiveRecursion.Im(),wFiveRecursion); // <<<sin(h1*phi1+h2*phi2+h3*phi3+h4*phi4+h5*phi5)>>
700  fRecursion[1][3]->Fill(1.5,wFiveRecursion,1.);
701  }// 5-p correlation
702 
703  if(kNumber==6)
704  {
705  Int_t harmonicsSixNum[6] = {kh1,kh2,kh3,kh4,kh5,kh6};
706  Int_t harmonicsSixDen[6] = {0,0,0,0,0,0};
707  TComplex sixRecursion = Recursion(6,harmonicsSixNum)/Recursion(6,harmonicsSixDen).Re();
708  Double_t wSixRecursion = Recursion(6,harmonicsSixDen).Re();
709  fRecursion[0][4]->Fill(0.5,sixRecursion.Re(),wSixRecursion); // <<cos(h1*phi1+h2*phi2+h3*phi3+h4*phi4+h5*phi5+h6*phi6)>>
710  fRecursion[0][4]->Fill(1.5,wSixRecursion,1.);
711  fRecursion[1][4]->Fill(0.5,sixRecursion.Im(),wSixRecursion); // <<<sin(h1*phi1+h2*phi2+h3*phi3+h4*phi4+h5*phi5+h6*phi6)>>
712  fRecursion[1][4]->Fill(1.5,wSixRecursion,1.);
713  }// 6-p correlation
714 
715 
716  if(kNumber==7)
717  {
718  Int_t harmonicsSevenNum[7] = {kh1,kh2,kh3,kh4,kh5,kh6,kh7};
719  Int_t harmonicsSevenDen[7] = {0,0,0,0,0,0,0};
720  TComplex sevenRecursion = Recursion(7,harmonicsSevenNum)/Recursion(7,harmonicsSevenDen).Re();
721  Double_t wSevenRecursion = Recursion(7,harmonicsSevenDen).Re();
722  fRecursion[0][5]->Fill(0.5,sevenRecursion.Re(),wSevenRecursion); // <<cos(h1*phi1+h2*phi2+h3*phi3+h4*phi4+h5*phi5+h6*phi6+h7*phi7)>>
723  fRecursion[0][5]->Fill(1.5,wSevenRecursion,1.);
724  fRecursion[1][5]->Fill(0.5,sevenRecursion.Im(),wSevenRecursion); // <<<sin(h1*phi1+h2*phi2+h3*phi3+h4*phi4+h5*phi5+h6*phi6+h7*phi7)>>
725  fRecursion[1][5]->Fill(1.5,wSevenRecursion,1.);
726  }// 7-p correlation
727 
728 
729  if(kNumber==8)
730  {
731  Int_t harmonicsEightNum[8] = {kh1,kh2,kh3,kh4,kh5,kh6,kh7,kh8};
732  Int_t harmonicsEightDen[8] = {0,0,0,0,0,0,0,0};
733  TComplex eightRecursion = Recursion(8,harmonicsEightNum)/Recursion(8,harmonicsEightDen).Re();
734  Double_t wEightRecursion = Recursion(8,harmonicsEightDen).Re();
735  fRecursion[0][6]->Fill(0.5,eightRecursion.Re(),wEightRecursion); // <<cos(h1*phi1+h2*phi2+h3*phi3+h4*phi4+h5*phi5+h6*phi6+h7*phi7+h8*phi8)>>
736  fRecursion[0][6]->Fill(1.5,wEightRecursion,1.);
737  fRecursion[1][6]->Fill(0.5,eightRecursion.Im(),wEightRecursion); // <<<sin(h1*phi1+h2*phi2+h3*phi3+h4*phi4+h5*phi5+h6*phi6+h7*phi7+h8*phi8)>>
738  fRecursion[1][6]->Fill(1.5,wEightRecursion,1.);
739  }// 8-p correlation
740 
741  if(kNumber!=2 && kNumber!=3 && kNumber!=4 && kNumber!=5 && kNumber!=6 && kNumber!=7 && kNumber!=8)
742  {
743  return;
744  }
745 
746 
747 
748  }//void Correlation()
749 
750 
751 //========================================================================================================================
752 
753 
Int_t charge
double Double_t
Definition: External.C:58
Bool_t bUseWeights
//how can i set the 8 as fMaxCorrelator??????
virtual void Terminate(Option_t *)
TProfile * fCentrality
//[fMaxHarmonic*fMaxCorrelator+1][fMaxCorrelator+1]
TCanvas * c
Definition: TestFitELoss.C:172
TArrayD * fWeights
Azimuthal angles.
int Int_t
Definition: External.C:63
const char Option_t
Definition: External.C:48
bool Bool_t
Definition: External.C:53
TComplex Recursion(Int_t n, Int_t *harmonic, Int_t mult, Int_t skip)