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