1 /*************************************************************************
2 * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 /***************************************************************
19 * Only in this class nested loops are used for flow analysis. *
20 * Nested loops are used to evaluate: *
22 * a) Distribution of relative angle difference (phi1-phi2); *
23 * b) Cross-check the results for mixed harmonics. *
25 * Author: Ante Bilandzic (abilandzic@gmail.com) *
26 ***************************************************************/
28 #define AliFlowAnalysisWithNestedLoops_cxx
30 #include "Riostream.h"
31 #include "AliFlowCommonConstants.h"
32 #include "AliFlowCommonHist.h"
33 #include "AliFlowCommonHistResults.h"
40 #include "AliFlowEventSimple.h"
41 #include "AliFlowTrackSimple.h"
42 #include "AliFlowAnalysisWithNestedLoops.h"
49 ClassImp(AliFlowAnalysisWithNestedLoops)
51 //================================================================================================================
53 AliFlowAnalysisWithNestedLoops::AliFlowAnalysisWithNestedLoops():
58 fAnalysisSettings(NULL),
59 fOppositeChargesPOI(kFALSE),
60 fEvaluateDifferential3pCorrelator(kFALSE),
61 fPrintOnTheScreen(kTRUE),
76 fUsePhiWeights(kFALSE),
77 fUsePtWeights(kFALSE),
78 fUseEtaWeights(kFALSE),
79 fUseParticleWeights(NULL),
84 fEvaluateNestedLoopsForRAD(kTRUE),
85 fRelativeAngleDistribution(NULL),
88 fEvaluateNestedLoopsForQC(kFALSE),
90 fEvaluateNestedLoopsForMH(kFALSE),
91 f3pCorrelatorPro(NULL),
92 f5pCorrelatorPro(NULL)
96 // Base list to hold all output objects:
97 fHistList = new TList();
98 fHistListName = new TString("cobjNL");
99 fHistList->SetName(fHistListName->Data());
100 fHistList->SetOwner(kTRUE);
102 // List to hold histograms with phi, pt and eta weights:
103 fWeightsList = new TList();
105 // List to hold objects relevant for relative angle distributions:
106 fListRAD = new TList();
108 // List holding objects relevant for debugging and cross-checking of Q-cumulants class:
109 fListQC = new TList();
111 // List holding objects relevant for debugging and cross-checking of MH class:
112 fListMH = new TList();
114 // Initialize all arrays:
115 this->InitializeArraysForMH();
117 } // AliFlowAnalysisWithNestedLoops::AliFlowAnalysisWithNestedLoops()
119 //================================================================================================================
121 AliFlowAnalysisWithNestedLoops::~AliFlowAnalysisWithNestedLoops()
127 } // end of AliFlowAnalysisWithNestedLoops::~AliFlowAnalysisWithNestedLoops()
129 //================================================================================================================
131 void AliFlowAnalysisWithNestedLoops::Init()
133 // Initialize and book all objects.
135 // a) Cross check if the user settings make sense before starting;
136 // b) Access all common constants;
137 // c) Book and nest all lists in the base list fHistList;
138 // d) Book profile holding seetings for analysis with nested loops;
139 // e) Book common control histograms;
140 // f) Book all objects relevant for distributions;
141 // g) Book and fill histograms to hold phi, pt and eta weights;
142 // h) Store harmonic n.
145 //save old value and prevent histograms from being added to directory
146 //to avoid name clashes in case multiple analaysis objects are used
148 Bool_t oldHistAddStatus = TH1::AddDirectoryStatus();
149 TH1::AddDirectory(kFALSE);
151 TH1::SetDefaultSumw2();
153 this->CrossCheckSettings();
154 this->AccessConstants();
155 this->BookAndNestAllLists();
156 this->BookAndFillProfileHoldingSettings();
157 this->BookCommonHistograms();
158 this->BookEverythingForRAD();
159 this->BookEverythingForMH();
160 this->BookAndFillWeightsHistograms();
161 this->StoreHarmonic();
164 TH1::AddDirectory(oldHistAddStatus);
165 } // end of void AliFlowAnalysisWithNestedLoops::Init()
167 //================================================================================================================
169 void AliFlowAnalysisWithNestedLoops::Make(AliFlowEventSimple* anEvent)
171 // Running over data only in this method.
173 // a) Check all pointers used in this method;
174 // b) Fill common control histograms;
175 // c) Evaluate nested loops for relative angle distribution;
176 // d) Evaluate nested loops for mixed harmonics.
178 this->CheckPointersUsedInMake();
179 fCommonHists->FillControlHistograms(anEvent);
180 if(fEvaluateNestedLoopsForRAD) this->EvaluateNestedLoopsForRAD(anEvent);
181 if(fEvaluateNestedLoopsForMH) this->EvaluateNestedLoopsForMH(anEvent);
183 } // end of AliFlowAnalysisWithNestedLoops::Make(AliFlowEventSimple* anEvent)
185 //================================================================================================================
187 void AliFlowAnalysisWithNestedLoops::Finish()
189 // Calculate the final results.
191 // a) Access settings for analysis with mixed harmonics;
192 // b) Check all pointers used in this method;
193 // c) Print on the screen.
195 this->AccessSettings();
196 this->CheckPointersUsedInFinish();
197 if(fPrintOnTheScreen) this->PrintOnTheScreen();
199 } // end of AliFlowAnalysisWithNestedLoops::Finish()
201 //================================================================================================================
203 void AliFlowAnalysisWithNestedLoops::GetOutputHistograms(TList *outputListHistos)
205 // Get pointers to all objects saved in the output file.
207 // a) Get pointers for common control histograms.
210 this->SetHistList(outputListHistos);
214 cout<<" WARNING (NL): fHistList is NULL in AFAWNL::GOH() !!!!"<<endl;
218 this->GetPointersForBaseHistograms();
219 this->GetPointersForCommonHistograms();
220 this->AccessSettings();
221 if(fEvaluateNestedLoopsForRAD) this->GetPointersForRAD();
222 if(fEvaluateNestedLoopsForMH) this->GetPointersForMH();
226 cout<<" WARNING (NL): outputListHistos is NULL in AFAWNL::GOH() !!!!"<<endl;
231 } // end of void AliFlowAnalysisWithNestedLoops::GetOutputHistograms(TList *outputListHistos)
233 //================================================================================================================
235 void AliFlowAnalysisWithNestedLoops::GetPointersForBaseHistograms()
237 // Get pointers to base histograms.
239 TString analysisSettingsName = "fAnalysisSettings";
240 TProfile *analysisSettings = dynamic_cast<TProfile*>(fHistList->FindObject(analysisSettingsName.Data()));
243 this->SetAnalysisSettings(analysisSettings);
247 cout<<" WARNING (NL): analysisSettings is NULL in AFAWNL::GPFBH() !!!!"<<endl;
252 } // end of void AliFlowAnalysisWithNestedLoops::GetPointersForBaseHistograms()
254 //================================================================================================================
256 void AliFlowAnalysisWithNestedLoops::GetPointersForCommonHistograms()
258 // Get pointers to common control histograms.
260 TString commonHistsName = "AliFlowCommonHistNL";
261 AliFlowCommonHist *commonHist = dynamic_cast<AliFlowCommonHist*>(fHistList->FindObject(commonHistsName.Data()));
264 this->SetCommonHists(commonHist);
268 cout<<" WARNING (NL): commonHist is NULL in AFAWNL::GPFCH() !!!!"<<endl;
273 } // end of void AliFlowAnalysisWithNestedLoops::GetPointersForCommonHistograms()
275 //================================================================================================================
277 void AliFlowAnalysisWithNestedLoops::GetPointersForRAD()
279 // Get pointers to objects relevant for relative angle distributions.
281 TList *listRAD = NULL;
282 listRAD = dynamic_cast<TList*>(fHistList->FindObject("Relative Angle Distribution"));
285 cout<<"WARNING: listRAD is NULL in AFAWNL::GPFRAD() !!!!"<<endl;
289 TString relativeAngleDistributionName = "fRelativeAngleDistribution";
290 TH1D *relativeAngleDistribution = dynamic_cast<TH1D*>(listRAD->FindObject(relativeAngleDistributionName.Data()));
291 if(relativeAngleDistribution)
293 this->SetRelativeAngleDistribution(relativeAngleDistribution);
296 TString chargeName = "fCharge";
297 TH1D *charge = dynamic_cast<TH1D*>(listRAD->FindObject(chargeName.Data()));
300 this->SetCharge(charge);
303 } // end of void AliFlowAnalysisWithNestedLoops::GetPointersForRAD()
305 //================================================================================================================
307 void AliFlowAnalysisWithNestedLoops::GetPointersForMH()
309 // Get pointers to objects evaluated with nested loops.
311 TList *listMH = NULL;
312 listMH = dynamic_cast<TList*>(fHistList->FindObject("Mixed Harmonics"));
315 cout<<"WARNING: listMH is NULL in AFAWNL::GPFMH() !!!!"<<endl;
319 TString s3pCorrelatorProName = "f3pCorrelatorPro";
320 TProfile *p3pCorrelatorPro = dynamic_cast<TProfile*>(listMH->FindObject(s3pCorrelatorProName.Data()));
323 this->Set3pCorrelatorPro(p3pCorrelatorPro);
325 TString s5pCorrelatorProName = "f5pCorrelatorPro";
326 TProfile *p5pCorrelatorPro = dynamic_cast<TProfile*>(listMH->FindObject(s5pCorrelatorProName.Data()));
329 this->Set5pCorrelatorPro(p5pCorrelatorPro);
331 if(!fEvaluateDifferential3pCorrelator){return;}
332 TString psdFlag[2] = {"PtSum","PtDiff"};
333 for(Int_t sd=0;sd<2;sd++)
335 TProfile *p3pCorrelatorVsPtSumDiffDirectPro = dynamic_cast<TProfile*>(listMH->FindObject(Form("f3pCorrelatorDirectVs%s",psdFlag[sd].Data())));
336 if(p3pCorrelatorVsPtSumDiffDirectPro)
338 this->Set3pCorrelatorVsPtSumDiffDirectPro(p3pCorrelatorVsPtSumDiffDirectPro,sd);
340 } // end of for(Int_t sd=0;sd<2;sd++)
342 } // end of void AliFlowAnalysisWithNestedLoops::GetPointersForMH()
344 //================================================================================================================
346 void AliFlowAnalysisWithNestedLoops::WriteHistograms(TString outputFileName)
348 // Store the final results in output .root file.
349 TFile *output = new TFile(outputFileName.Data(),"RECREATE");
350 fHistList->Write(fHistList->GetName(),TObject::kSingleKey);
354 //================================================================================================================
356 void AliFlowAnalysisWithNestedLoops::WriteHistograms(TDirectoryFile *outputFileName)
358 // Store the final results in output .root file.
359 fHistList->SetName("cobjNL");
360 fHistList->SetOwner(kTRUE);
361 outputFileName->Add(fHistList);
362 outputFileName->Write(outputFileName->GetName(),TObject::kSingleKey);
365 //================================================================================================================
367 void AliFlowAnalysisWithNestedLoops::StoreHarmonic()
369 // Store harmonic n used in cos[n*(phi1+phi2-2phi3)] and cos[n*(psi1+psi2-2phi3)].
371 (fCommonHists->GetHarmonic())->Fill(0.5,fHarmonic);
373 } // end of void AliFlowAnalysisWithNestedLoops::StoreHarmonic()
375 //================================================================================================================
377 void AliFlowAnalysisWithNestedLoops::BookAndNestAllLists()
379 // Book and nest all list in base list fHistList.
382 fWeightsList->SetName("Weights");
383 fWeightsList->SetOwner(kTRUE);
384 fHistList->Add(fWeightsList);
385 // List for Relative Angle Distribution:
386 fListRAD->SetName("Relative Angle Distribution");
387 fListRAD->SetOwner(kTRUE);
388 if(fEvaluateNestedLoopsForRAD) fHistList->Add(fListRAD);
389 // List for Q-cumulants:
390 fListQC->SetName("Q-cumulants");
391 fListQC->SetOwner(kTRUE);
392 if(fEvaluateNestedLoopsForQC) fHistList->Add(fListQC);
393 // List for Mixed Harmonics:
394 fListMH->SetName("Mixed Harmonics");
395 fListMH->SetOwner(kTRUE);
396 if(fEvaluateNestedLoopsForMH) fHistList->Add(fListMH);
398 } // end of void AliFlowAnalysisWithNestedLoops::BookAndNestAllLists()
400 //================================================================================================================
402 void AliFlowAnalysisWithNestedLoops::BookAndFillProfileHoldingSettings()
404 // Book profile to hold all analysis settings.
406 TString analysisSettingsName = "fAnalysisSettings";
407 fAnalysisSettings = new TProfile(analysisSettingsName.Data(),"Settings for analysis with nested loops",7,0,7);
408 fAnalysisSettings->GetXaxis()->SetLabelSize(0.035);
409 fAnalysisSettings->GetXaxis()->SetBinLabel(1,"Nested loops for RAD?");
410 fAnalysisSettings->Fill(0.5,(Int_t)fEvaluateNestedLoopsForRAD);
411 fAnalysisSettings->GetXaxis()->SetBinLabel(2,"Nested loops for QC?");
412 fAnalysisSettings->Fill(1.5,(Int_t)fEvaluateNestedLoopsForQC);
413 fAnalysisSettings->GetXaxis()->SetBinLabel(3,"Nested loops for MH?");
414 fAnalysisSettings->Fill(2.5,(Int_t)fEvaluateNestedLoopsForMH);
415 fAnalysisSettings->GetXaxis()->SetBinLabel(4,"fHarmonic");
416 fAnalysisSettings->Fill(3.5,(Int_t)fHarmonic);
417 fAnalysisSettings->GetXaxis()->SetBinLabel(5,"Print on the screen?");
418 fAnalysisSettings->Fill(4.5,(Int_t)fPrintOnTheScreen);
419 fAnalysisSettings->GetXaxis()->SetBinLabel(6,"fOppositeChargesPOI");
420 fAnalysisSettings->Fill(5.5,fOppositeChargesPOI);
421 fAnalysisSettings->GetXaxis()->SetBinLabel(7,"fEvaluateDifferential3pCorrelator");
422 fAnalysisSettings->Fill(6.5,fEvaluateDifferential3pCorrelator);
423 fHistList->Add(fAnalysisSettings);
425 } // end of void AliFlowAnalysisWithNestedLoops::BookAndFillProfileHoldingSettings()
427 //================================================================================================================
429 void AliFlowAnalysisWithNestedLoops::BookCommonHistograms()
431 // Book common control histograms and common histograms for final results.
433 TString commonHistsName = "AliFlowCommonHistNL";
434 fCommonHists = new AliFlowCommonHist(commonHistsName.Data());
435 fHistList->Add(fCommonHists);
437 } // end of void AliFlowAnalysisWithNestedLoops::BookCommonHistograms()
439 //================================================================================================================
441 void AliFlowAnalysisWithNestedLoops::BookEverythingForRAD()
443 // Book all objects relevant calculation of relative angle distribution.
445 TString relativeAngleDistributionName = "fRelativeAngleDistribution";
446 fRelativeAngleDistribution = new TH1D(relativeAngleDistributionName.Data(),"Relative angle distribution",720,-TMath::TwoPi(),TMath::TwoPi());
447 fRelativeAngleDistribution->GetYaxis()->SetTitle("#frac{dN}{#Delta #phi}");
448 fRelativeAngleDistribution->GetXaxis()->SetTitle("#Delta #phi");
449 fListRAD->Add(fRelativeAngleDistribution);
451 TString chargeName = "fCharge";
452 fCharge = new TH1D(chargeName.Data(),"Charges",3,0,3);
454 fCharge->GetXaxis()->SetBinLabel(1,"+");
455 fCharge->GetXaxis()->SetBinLabel(2,"-");
456 fCharge->GetXaxis()->SetBinLabel(3,"rest");
457 fListRAD->Add(fCharge);
459 } // end fo void AliFlowAnalysisWithNestedLoops::BookEverythingForRAD()
461 //================================================================================================================
463 void AliFlowAnalysisWithNestedLoops::AccessConstants()
465 // Access needed common constants from AliFlowCommonConstants.
467 fnBinsPhi = AliFlowCommonConstants::GetMaster()->GetNbinsPhi();
468 fPhiMin = AliFlowCommonConstants::GetMaster()->GetPhiMin();
469 fPhiMax = AliFlowCommonConstants::GetMaster()->GetPhiMax();
470 if(fnBinsPhi) fPhiBinWidth = (fPhiMax-fPhiMin)/fnBinsPhi;
471 fnBinsPt = AliFlowCommonConstants::GetMaster()->GetNbinsPt();
472 fPtMin = AliFlowCommonConstants::GetMaster()->GetPtMin();
473 fPtMax = AliFlowCommonConstants::GetMaster()->GetPtMax();
474 if(fnBinsPt) fPtBinWidth = (fPtMax-fPtMin)/fnBinsPt;
475 fnBinsEta = AliFlowCommonConstants::GetMaster()->GetNbinsEta();
476 fEtaMin = AliFlowCommonConstants::GetMaster()->GetEtaMin();
477 fEtaMax = AliFlowCommonConstants::GetMaster()->GetEtaMax();
478 if(fnBinsEta) fEtaBinWidth = (fEtaMax-fEtaMin)/fnBinsEta;
480 } // end of void AliFlowAnalysisWithNestedLoops::AccessConstants()
482 //================================================================================================================
484 void AliFlowAnalysisWithNestedLoops::CrossCheckSettings()
486 // Cross-check if the user settings make sense.
490 } // end of void AliFlowAnalysisWithNestedLoops::CrossCheckSettings()
492 //================================================================================================================
494 void AliFlowAnalysisWithNestedLoops::BookAndFillWeightsHistograms()
496 // Book and fill (by accessing file "weights.root") histograms which hold phi, pt and eta weights.
500 cout<<"WARNING: fWeightsList is NULL in AFAWNL::BAFWH() !!!!"<<endl;
503 // Profile to hold flags for weights:
504 TString fUseParticleWeightsName = "fUseParticleWeightsNL";
505 fUseParticleWeights = new TProfile(fUseParticleWeightsName.Data(),"0 = particle weight not used, 1 = particle weight used ",3,0,3);
506 fUseParticleWeights->SetLabelSize(0.06);
507 (fUseParticleWeights->GetXaxis())->SetBinLabel(1,"w_{#phi}");
508 (fUseParticleWeights->GetXaxis())->SetBinLabel(2,"w_{p_{T}}");
509 (fUseParticleWeights->GetXaxis())->SetBinLabel(3,"w_{#eta}");
510 fUseParticleWeights->Fill(0.5,(Int_t)fUsePhiWeights);
511 fUseParticleWeights->Fill(1.5,(Int_t)fUsePtWeights);
512 fUseParticleWeights->Fill(2.5,(Int_t)fUseEtaWeights);
513 fWeightsList->Add(fUseParticleWeights);
517 if(fWeightsList->FindObject("phi_weights"))
519 fPhiWeights = dynamic_cast<TH1F*>(fWeightsList->FindObject("phi_weights"));
522 printf("WARNING: no phi weights object, bye!\n");
525 if(TMath::Abs(fPhiWeights->GetBinWidth(1)-fPhiBinWidth)>pow(10.,-6.))
528 cout<<"WARNING (NL): Inconsistent binning in histograms for phi-weights throughout the code."<<endl;
534 cout<<"WARNING (NL): fWeightsList->FindObject(\"phi_weights\") is NULL in AFAWNL::BAFWH() !!!!"<<endl;
537 } // end of if(fUsePhiWeights)
541 if(fWeightsList->FindObject("pt_weights"))
543 fPtWeights = dynamic_cast<TH1D*>(fWeightsList->FindObject("pt_weights"));
546 printf("WARNING: no pt weights object, bye!\n");
549 if(TMath::Abs(fPtWeights->GetBinWidth(1)-fPtBinWidth)>pow(10.,-6.))
552 cout<<"WARNING (NL): Inconsistent binning in histograms for pt-weights throughout the code."<<endl;
558 cout<<"WARNING (NL): fWeightsList->FindObject(\"pt_weights\") is NULL in AFAWNL::BAFWH() !!!!"<<endl;
561 } // end of if(fUsePtWeights)
565 if(fWeightsList->FindObject("eta_weights"))
567 fEtaWeights = dynamic_cast<TH1D*>(fWeightsList->FindObject("eta_weights"));
570 printf("WARNING: no eta weights object, bye!\n");
573 if(TMath::Abs(fEtaWeights->GetBinWidth(1)-fEtaBinWidth)>pow(10.,-6.))
576 cout<<"WARNING (NL): Inconsistent binning in histograms for eta-weights throughout the code."<<endl;
582 cout<<"WARNING: fUseEtaWeights && fWeightsList->FindObject(\"eta_weights\") is NULL in AFAWNL::BAFWH() !!!!"<<endl;
585 } // end of if(fUseEtaWeights)
587 } // end of AliFlowAnalysisWithNestedLoops::BookAndFillWeightsHistograms()
589 //================================================================================================================
591 void AliFlowAnalysisWithNestedLoops::CheckPointersUsedInMake()
593 // Check pointers used in method Make().
595 if(fEvaluateNestedLoopsForRAD) CheckPointersForRAD("Make");
596 if(fEvaluateNestedLoopsForMH) CheckPointersForMH("Make");
598 } // end of AliFlowAnalysisWithNestedLoops::CheckPointersUsedInMake()
600 //================================================================================================================
602 void AliFlowAnalysisWithNestedLoops::CheckPointersUsedInFinish()
604 // Check pointers used in method Finish().
606 if(fEvaluateNestedLoopsForRAD) CheckPointersForRAD("Finish");
607 if(fEvaluateNestedLoopsForMH) CheckPointersForMH("Finish");
609 } // end of AliFlowAnalysisWithNestedLoops::CheckPointersUsedInFinish()
611 //================================================================================================================
613 void AliFlowAnalysisWithNestedLoops::CheckPointersForRAD(TString where)
615 // Check pointers relevant for calculation of relative angle distribution.
617 if(!fRelativeAngleDistribution)
620 cout<<" WARNING (NL): fRelativeAngleDistribution is NULL in "<<where.Data()<<"() !!!!"<<endl;
625 if(strcmp(where.Data(),"Make") == 0)
627 // Check pointers used only in method Make():
630 else if(strcmp(where.Data(),"Finish") == 0)
632 // Check pointers used only in method Finish():
636 } // end of void AliFlowAnalysisWithNestedLoops::CheckPointersForRAD(TString where)
638 //================================================================================================================
640 void AliFlowAnalysisWithNestedLoops::CheckPointersForMH(TString where)
642 // Check pointers relevant for calculation of mixed harmonics.
644 if(strcmp(where.Data(),"Make") == 0)
646 // Check pointers used only in method Make():
647 if(!f3pCorrelatorPro)
650 cout<<" WARNING (NL): f3pCorrelatorPro is NULL in Make() !!!!"<<endl;
654 if(!f5pCorrelatorPro)
657 cout<<" WARNING (NL): f5pCorrelatorPro is NULL in Make() !!!!"<<endl;
661 if(!fEvaluateDifferential3pCorrelator){return;}
662 for(Int_t sd=0;sd<2;sd++)
664 if(!(f3pCorrelatorVsPtSumDiffDirectPro[sd]))
667 cout<<" WARNING (NL): "<<Form("f3pCorrelatorVsPtSumDiffDirectPro[%d]",sd)<<" is NULL in Make() !!!!"<<endl;
671 } // end of for(Int_t sd=0;sd<2;sd++)
672 } // if(strcmp(where.Data(),"Make") == 0)
673 else if(strcmp(where.Data(),"Finish") == 0)
675 // Check pointers used only in method Finish():
676 if(!f3pCorrelatorPro)
679 cout<<" WARNING (NL): f3pCorrelatorPro is NULL in Finish() !!!!"<<endl;
683 if(!f5pCorrelatorPro)
686 cout<<" WARNING (NL): f5pCorrelatorPro is NULL in Finish() !!!!"<<endl;
690 } // else if(strcmp(where.Data(),"Finish") == 0)
692 } // end of void AliFlowAnalysisWithNestedLoops::CheckPointersForMH(TString where)
694 //================================================================================================================
696 void AliFlowAnalysisWithNestedLoops::AccessSettings()
698 // Access the settings for analysis.
700 fEvaluateNestedLoopsForRAD = (Bool_t)fAnalysisSettings->GetBinContent(1);
701 fEvaluateNestedLoopsForQC = (Bool_t)fAnalysisSettings->GetBinContent(2);
702 fEvaluateNestedLoopsForMH = (Bool_t)fAnalysisSettings->GetBinContent(3);
703 fHarmonic = (Int_t)fAnalysisSettings->GetBinContent(4);
704 fPrintOnTheScreen = (Bool_t)fAnalysisSettings->GetBinContent(5);
705 fOppositeChargesPOI = (Bool_t)fAnalysisSettings->GetBinContent(6);
706 fEvaluateDifferential3pCorrelator = (Bool_t)fAnalysisSettings->GetBinContent(7);
708 } // end of AliFlowAnalysisWithNestedLoops::AccessSettings()
710 //================================================================================================================
712 void AliFlowAnalysisWithNestedLoops::InitializeArraysForMH()
714 // Initialize arrays mixed harmonics calculations.
716 for(Int_t sd=0;sd<2;sd++) // sum or difference
718 f3pCorrelatorVsPtSumDiffDirectPro[sd] = NULL;
721 } // end of AliFlowAnalysisWithNestedLoops::InitializeArraysForMH()
723 //================================================================================================================
725 void AliFlowAnalysisWithNestedLoops::BookEverythingForMH()
727 // Book all objects relevant for mixed harmonics.
729 if(fEvaluateNestedLoopsForMH)
732 TString s3pCorrelatorProName = "f3pCorrelatorPro";
733 f3pCorrelatorPro = new TProfile(s3pCorrelatorProName.Data(),"",1,0,1);
734 f3pCorrelatorPro->SetStats(kFALSE);
735 f3pCorrelatorPro->GetXaxis()->SetLabelOffset(0.01);
736 f3pCorrelatorPro->GetXaxis()->SetLabelSize(0.05);
739 f3pCorrelatorPro->GetXaxis()->SetBinLabel(1,"#LT#LTcos(#phi_{1}+#phi_{2}-2#phi_{3})#GT#GT");
742 f3pCorrelatorPro->GetXaxis()->SetBinLabel(1,Form("#LT#LTcos[%i(#phi_{1}+#phi_{2}-2#phi_{3})]#GT#GT",fHarmonic));
744 fListMH->Add(f3pCorrelatorPro);
746 TString s5pCorrelatorProName = "f5pCorrelatorPro";
747 f5pCorrelatorPro = new TProfile(s5pCorrelatorProName.Data(),"",1,0,1);
748 f5pCorrelatorPro->SetStats(kFALSE);
749 f5pCorrelatorPro->GetXaxis()->SetLabelOffset(0.01);
750 f5pCorrelatorPro->GetXaxis()->SetLabelSize(0.05);
753 f5pCorrelatorPro->GetXaxis()->SetBinLabel(1,"#LT#LTcos(2#phi_{1}+2#phi_{2}+2#phi_{3}-3#phi_{4}-3#phi_{5})#GT#GT");
756 f5pCorrelatorPro->GetXaxis()->SetBinLabel(1,Form("#LT#LTcos[%i(2#phi_{1}+2#phi_{2}+2#phi_{3}-3#phi_{4}-3#phi_{5})]#GT#GT",fHarmonic));
758 fListMH->Add(f5pCorrelatorPro);
759 if(fEvaluateDifferential3pCorrelator)
761 TString psdFlag[2] = {"PtSum","PtDiff"};
762 TString psdTitleFlag[2] = {"(p_{t,1}+p_{t,2})/2","#left|p_{t,1}-p_{t,2}#right|"};
763 //TString s3pCorrelatorVsPtSumDiffDirectProName = "f3pCorrelatorVsPtSumDiffDirectPro";
764 for(Int_t sd=0;sd<2;sd++)
766 // to be improved: hardwired ,fnBinsPt,0.,fPtMax):
767 f3pCorrelatorVsPtSumDiffDirectPro[sd] = new TProfile(Form("f3pCorrelatorDirectVs%s",psdFlag[sd].Data()),"",fnBinsPt,0.,fPtMax);
768 //f3pCorrelatorVsPtSumDiffDirectPro[sd]->SetLabelSize(0.05);
769 //f3pCorrelatorVsPtSumDiffDirectPro[sd]->SetMarkerStyle(25);
770 f3pCorrelatorVsPtSumDiffDirectPro[sd]->GetXaxis()->SetTitle(psdTitleFlag[sd].Data());
771 fListMH->Add(f3pCorrelatorVsPtSumDiffDirectPro[sd]);
773 } // end of if(fEvaluateDifferential3pCorrelator)
774 } // end of if(fEvaluateNestedLoopsForMH)
776 } // end of AliFlowAnalysisWithNestedLoops::BookEverythingForMH()
778 //================================================================================================================
780 void AliFlowAnalysisWithNestedLoops::EvaluateNestedLoopsForRAD(AliFlowEventSimple *anEvent)
782 // Evaluate nested loops needed for calculation of relative angle distribution.
784 Double_t dPhi1=0., dPhi2=0.; // azimuthal angles in the laboratory frame
785 AliFlowTrackSimple *aftsTrack = NULL; // simple track
787 // Loop over data and store for each distinct pair phi1-phi2 in fRelativeAngleDistribution:
788 Int_t nPrim = anEvent->NumberOfTracks(); // nPrim = total number of primary tracks, i.e. nPrim = nRP + nPOI + rest, where:
789 // nRP = # of particles used to determine the reaction plane ("Reference Particles");
790 // nPOI = # of particles of interest for a detailed flow analysis ("Particles of Interest");
791 // rest = # of particles which are not niether RPs nor POIs.
792 // Start nested loops over data:
793 for(Int_t i=0;i<nPrim;i++)
795 aftsTrack=anEvent->GetTrack(i);
798 if(!aftsTrack->InRPSelection()) continue; // consider only tracks which are RPs
799 dPhi1 = aftsTrack->Phi();
800 for(Int_t j=0;j<nPrim;j++)
802 if(j==i) continue; // eliminating trivial contribution from autocorrelation
803 aftsTrack=anEvent->GetTrack(j);
806 if(!aftsTrack->InRPSelection()) continue; // consider only tracks which are RPs
807 dPhi2 = aftsTrack->Phi();
808 // Fill the histogram:
809 fRelativeAngleDistribution->Fill(dPhi1-dPhi2);
811 } // end of for(Int_t j=0;j<nPrim;j++)
812 } else // to if(aftsTrack)
815 cout<<" WARNING (NL): No particle! (i.e. aftsTrack is a NULL pointer in AFAWNL::Make().)"<<endl;
818 } // end of for(Int_t i=0;i<nPrim;i++)
820 } // end of void AliFlowAnalysisWithNestedLoops::EvaluateNestedLoopsForRAD(AliFlowEventSimple *anEvent)
822 //================================================================================================================
824 void AliFlowAnalysisWithNestedLoops::EvaluateNestedLoopsForMH(AliFlowEventSimple *anEvent)
826 // Evaluate nested loops needed for mixed harmonics.
827 // Remark: phi labels the azimuthal angle of RP particle and psi labels azimuthal angle of POI particle.
829 Int_t nPrim = anEvent->NumberOfTracks();
830 Int_t nRP = anEvent->GetEventNSelTracksRP();
831 AliFlowTrackSimple *aftsTrack = NULL;
832 Double_t phi1=0.,phi2=0.,phi3=0.,phi4=0.,phi5=0.; // azimuthal angles of RPs
833 Double_t psi1=0.,psi2=0.; // azimuthal angles of POIs
834 Int_t charge1=0,charge2=0; // charge of POIs
835 Double_t pt1=0.,pt2=0.; // transverse momenta of POIs
838 // Evaluting correlator cos[n(phi1+phi2-2*phi3)] with three nested loops:
841 for(Int_t i1=0;i1<nPrim;i1++)
843 aftsTrack=anEvent->GetTrack(i1);
844 if(!(aftsTrack->InRPSelection())) continue;
845 phi1 = aftsTrack->Phi();
846 for(Int_t i2=0;i2<nPrim;i2++)
849 aftsTrack = anEvent->GetTrack(i2);
850 if(!(aftsTrack->InRPSelection())) continue;
851 phi2 = aftsTrack->Phi();
852 for(Int_t i3=0;i3<nPrim;i3++)
854 if(i3==i1||i3==i2) continue;
855 aftsTrack=anEvent->GetTrack(i3);
856 if(!(aftsTrack->InRPSelection())) continue;
857 phi3 = aftsTrack->Phi();
858 f3pCorrelatorPro->Fill(0.5,cos(n*(phi1+phi2-2.*phi3)),1.);
859 } // end of for(Int_t i3=0;i3<nPrim;i3++)
860 } // end of for(Int_t i2=0;i2<nPrim;i2++)
861 } // end of for(Int_t i1=0;i1<nPrim;i1++)
864 // Evaluting correlator cos[n(2*phi1+2*phi2+2*phi3-3*phi4-3*phi5)] with five nested loops:
867 for(Int_t i1=0;i1<nPrim;i1++)
869 aftsTrack=anEvent->GetTrack(i1);
870 if(!(aftsTrack->InRPSelection())) continue;
871 phi1=aftsTrack->Phi();
872 for(Int_t i2=0;i2<nPrim;i2++)
875 aftsTrack=anEvent->GetTrack(i2);
876 if(!(aftsTrack->InRPSelection())) continue;
877 phi2=aftsTrack->Phi();
878 for(Int_t i3=0;i3<nPrim;i3++)
880 if(i3==i1||i3==i2)continue;
881 aftsTrack=anEvent->GetTrack(i3);
882 if(!(aftsTrack->InRPSelection())) continue;
883 phi3=aftsTrack->Phi();
884 for(Int_t i4=0;i4<nPrim;i4++)
886 if(i4==i1||i4==i2||i4==i3)continue;
887 aftsTrack=anEvent->GetTrack(i4);
888 if(!(aftsTrack->InRPSelection())) continue;
889 phi4=aftsTrack->Phi();
890 for(Int_t i5=0;i5<nPrim;i5++)
892 if(i5==i1||i5==i2||i5==i3||i5==i4)continue;
893 aftsTrack=anEvent->GetTrack(i5);
894 if(!(aftsTrack->InRPSelection())) continue;
895 phi5=aftsTrack->Phi();
896 f5pCorrelatorPro->Fill(0.5,cos(2.*n*phi1+2.*n*phi2+2.*n*phi3-3.*n*phi4-3.*n*phi5),1.);
897 } // end of for(Int_t i5=0;i5<nPrim;i5++)
898 } // end of for(Int_t i4=0;i4<nPrim;i4++)
899 } // end of for(Int_t i3=0;i3<nPrim;i3++)
900 } // end of for(Int_t i2=0;i2<nPrim;i2++)
901 } // end of for(Int_t i1=0;i1<nPrim;i1++)
902 } // end of if(nPrim>=5)
904 // Evaluting correlator cos[n(psi1+psi2-2*phi3)] with three nested loops:
905 if(!fEvaluateDifferential3pCorrelator){return;}
906 for(Int_t i1=0;i1<nPrim;i1++)
908 aftsTrack=anEvent->GetTrack(i1);
909 if(!(aftsTrack->InPOISelection())){continue;}
910 psi1 = aftsTrack->Phi();
911 pt1 = aftsTrack->Pt();
912 charge1 = aftsTrack->Charge();
913 for(Int_t i2=0;i2<nPrim;i2++)
915 if(i2==i1){continue;}
916 aftsTrack = anEvent->GetTrack(i2);
917 if(!(aftsTrack->InPOISelection())){continue;}
918 psi2 = aftsTrack->Phi();
919 pt2 = aftsTrack->Pt();
920 charge2 = aftsTrack->Charge();
921 if(fOppositeChargesPOI && charge1 == charge2){continue;}
922 for(Int_t i3=0;i3<nPrim;i3++)
924 if(i3==i1||i3==i2){continue;}
925 aftsTrack=anEvent->GetTrack(i3);
926 if(!(aftsTrack->InRPSelection())){continue;}
927 phi3 = aftsTrack->Phi();
928 // Evaluate and store differential correlator cos[n(psi1+psi2-2*phi3)]:
929 Double_t ptSum = (pt1+pt2)/2.;
930 Double_t ptDiff = TMath::Abs(pt1-pt2);
931 Double_t diff3pCorrelator = TMath::Cos(n*(psi1+psi2-2.*phi3));
932 f3pCorrelatorVsPtSumDiffDirectPro[0]->Fill(ptSum,diff3pCorrelator,1.);
933 f3pCorrelatorVsPtSumDiffDirectPro[1]->Fill(ptDiff,diff3pCorrelator,1.);
934 } // end of for(Int_t i3=0;i3<nPrim;i3++)
935 } // end of for(Int_t i2=0;i2<nPrim;i2++)
936 } // end of for(Int_t i1=0;i1<nPrim;i1++)
938 } // end of void AliFlowAnalysisWithNestedLoops::EvaluateNestedLoopsForMH(AliFlowEventSimple *anEvent)
940 //================================================================================================================
942 void AliFlowAnalysisWithNestedLoops::PrintOnTheScreen()
944 // Print on the screen.
947 cout<<"********************************************************************"<<endl;
948 cout<<"********************************************************************"<<endl;
949 cout<<" Nested Loops "<<endl;
952 if(fEvaluateNestedLoopsForRAD)
954 cout<<" Evaluated for relative angle distribution."<<endl;
957 if(fEvaluateNestedLoopsForMH)
959 cout<<" Evaluated for mixed harmonics:"<<endl;
962 cout<<" cos["<<fHarmonic<<"(phi1+phi2-2phi3)] = "<<f3pCorrelatorPro->GetBinContent(1)<<" +/- " <<f3pCorrelatorPro->GetBinError(1)<<endl;
963 cout<<" cos["<<fHarmonic<<"(2phi1+2phi2+2phi3-3phi4-3phi5)] = "<<f5pCorrelatorPro->GetBinContent(1)<<" +/- " <<f5pCorrelatorPro->GetBinError(1)<<endl;
966 cout<<" cos(phi1+phi2-2phi3) = "<<f3pCorrelatorPro->GetBinContent(1)<<" +/- " <<f3pCorrelatorPro->GetBinError(1)<<endl;
967 cout<<" cos[(2phi1+2phi2+2phi3-3phi4-3phi5)] = "<<f5pCorrelatorPro->GetBinContent(1)<<" +/- " <<f5pCorrelatorPro->GetBinError(1)<<endl;
969 if(fEvaluateDifferential3pCorrelator)
971 cout<<" Evaluated also for differential 3p correlator "<<endl;
974 cout<<" cos["<<fHarmonic<<"(psi1+psi2-2phi3)]. "<<endl;
977 cout<<" cos(psi1+psi2-2phi3). "<<endl;
979 } // end of if(fEvaluateDifferential3pCorrelator)
980 } // end of if(fEvaluateNestedLoopsForMH)
982 if(!fEvaluateNestedLoopsForRAD && !fEvaluateNestedLoopsForMH)
984 cout<<" Not evaluated."<<endl;
987 cout<<"********************************************************************"<<endl;
988 cout<<"********************************************************************"<<endl;
991 } // end of void AliFlowAnalysisWithNestedLoops::PrintOnTheScreen()