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 * In this class azimuthal correlators in mixed harmonics *
20 * are implemented in terms of Q-vectors. This approach *
21 * doesn't require evaluation of nested loops. This class *
24 * a) Extract subdominant harmonics (like v1 and v4); *
25 * b) Study flow of two-particle resonances; *
26 * c) Study strong parity violation. *
28 * Author: Ante Bilandzic (abilandzic@gmail.com) *
29 *********************************************************/
31 #define AliFlowAnalysisWithMixedHarmonics_cxx
33 #include "Riostream.h"
34 #include "AliFlowCommonConstants.h"
35 #include "AliFlowCommonHist.h"
36 #include "AliFlowCommonHistResults.h"
42 #include "TProfile2D.h"
44 #include "AliFlowEventSimple.h"
45 #include "AliFlowTrackSimple.h"
46 #include "AliFlowAnalysisWithMixedHarmonics.h"
50 ClassImp(AliFlowAnalysisWithMixedHarmonics)
52 //================================================================================================================
54 AliFlowAnalysisWithMixedHarmonics::AliFlowAnalysisWithMixedHarmonics():
58 fAnalysisSettings(NULL),
59 fCorrelatorInteger(1),
60 fNoOfMultipicityBins(10),
61 fMultipicityBinWidth(2),
63 fCorrectForDetectorEffects(kTRUE),
64 fPrintOnTheScreen(kTRUE),
79 fUsePhiWeights(kFALSE),
80 fUsePtWeights(kFALSE),
81 fUseEtaWeights(kFALSE),
82 fUseParticleWeights(NULL),
90 f3pCorrelatorPro(NULL),
91 fNonIsotropicTermsPro(NULL),
92 f3pCorrelatorVsMPro(NULL),
93 fNonIsotropicTermsVsMPro(NULL),
95 f3pCorrelatorHist(NULL),
96 fDetectorBiasHist(NULL),
97 f3pCorrelatorVsMHist(NULL),
98 fDetectorBiasVsMHist(NULL)
102 // Base list to hold all output objects:
103 fHistList = new TList();
104 fHistListName = new TString("cobjMH");
105 fHistList->SetName(fHistListName->Data());
106 fHistList->SetOwner(kTRUE);
108 // List to hold histograms with phi, pt and eta weights:
109 fWeightsList = new TList();
111 // List to hold all all-event profiles:
112 fProfileList = new TList();
114 // List to hold objects with final results:
115 fResultsList = new TList();
117 // Initialize all arrays:
118 this->InitializeArrays();
120 } // AliFlowAnalysisWithMixedHarmonics::AliFlowAnalysisWithMixedHarmonics()
122 //================================================================================================================
124 AliFlowAnalysisWithMixedHarmonics::~AliFlowAnalysisWithMixedHarmonics()
130 } // end of AliFlowAnalysisWithMixedHarmonics::~AliFlowAnalysisWithMixedHarmonics()
132 //================================================================================================================
134 void AliFlowAnalysisWithMixedHarmonics::Init()
136 // Initialize and book all objects.
138 // a) Cross check if the user settings make sense before starting;
139 // b) Access all common constants;
140 // c) Book and nest all lists in the base list fHistList;
141 // d) Book common control histograms;
142 // e) Book all event-by-event quantities;
143 // f) Book all all-event quantities;
144 // g) Book and fill histograms to hold phi, pt and eta weights;
146 //save old value and prevent histograms from being added to directory
147 //to avoid name clashes in case multiple analaysis objects are used
149 Bool_t oldHistAddStatus = TH1::AddDirectoryStatus();
150 TH1::AddDirectory(kFALSE);
152 TH1::SetDefaultSumw2();
154 this->CrossCheckSettings();
155 this->AccessConstants();
156 this->BookAndNestAllLists();
157 this->BookProfileHoldingSettings();
158 this->BookCommonHistograms();
159 this->BookAllEventByEventQuantities();
160 this->BookAllAllEventQuantities();
161 this->BookAndFillWeightsHistograms();
163 TH1::AddDirectory(oldHistAddStatus);
165 } // end of void AliFlowAnalysisWithMixedHarmonics::Init()
167 //================================================================================================================
169 void AliFlowAnalysisWithMixedHarmonics::Make(AliFlowEventSimple* anEvent)
171 // Running over data only in this method.
173 // a) Check all pointers used in this method;
174 // b) Define local variables;
175 // c) Fill common control histograms;
176 // d) Loop over data and calculate e-b-e quantities Q_{n,k} and S_{p,k};
177 // e) Calculate 3-p azimuthal correlator and non-isotropic terms in terms of Q_{n,k} and S_{p,k};
178 // f) Calculate differential 3-p azimuthal correlator cos(2phi1-psi2-psi3) in terms of Q_{2n} and p_{n}:
179 // g) Reset all event-by-event quantities.
181 // a) Check all pointers used in this method:
182 this->CheckPointersUsedInMake();
184 // b) Define local variables:
185 Double_t dPhi = 0.; // azimuthal angle in the laboratory frame
186 Double_t dPt = 0.; // transverse momentum
187 Double_t dEta = 0.; // pseudorapidity
188 Double_t wPhi = 1.; // phi weight
189 Double_t wPt = 1.; // pt weight
190 Double_t wEta = 1.; // eta weight
191 AliFlowTrackSimple *aftsTrack = NULL; // simple track
193 // c) Fill common control histograms:
194 fCommonHists->FillControlHistograms(anEvent);
196 // d) Loop over data and calculate e-b-e quantities:
197 Int_t nPrim = anEvent->NumberOfTracks(); // nPrim = total number of primary tracks, i.e. nPrim = nRP + nPOI + rest, where:
198 // nRP = # of particles used to determine the reaction plane ("Reference Particles");
199 // nPOI = # of particles of interest for a detailed flow analysis ("Particles of Interest");
200 // rest = # of particles which are not niether RPs nor POIs.
201 // Start loop over data:
202 for(Int_t i=0;i<nPrim;i++)
204 aftsTrack=anEvent->GetTrack(i);
207 if(!(aftsTrack->InRPSelection() || aftsTrack->InPOISelection())) continue; // consider only tracks which are either RPs or POIs
208 Int_t n = fCorrelatorInteger; // integer n in the correlator cos[n(2phi1-phi2-phi3)]
209 if(aftsTrack->InRPSelection()) // checking RP condition:
211 dPhi = aftsTrack->Phi();
212 dPt = aftsTrack->Pt();
213 dEta = aftsTrack->Eta();
214 if(fUsePhiWeights && fPhiWeights && fnBinsPhi) // determine phi-weight for this particle:
216 wPhi = fPhiWeights->GetBinContent(1+(Int_t)(TMath::Floor(dPhi*fnBinsPhi/TMath::TwoPi())));
218 if(fUsePtWeights && fPtWeights && fnBinsPt) // determine pt-weight for this particle:
220 wPt = fPtWeights->GetBinContent(1+(Int_t)(TMath::Floor((dPt-fPtMin)/fPtBinWidth)));
222 if(fUseEtaWeights && fEtaWeights && fEtaBinWidth) // determine eta-weight for this particle:
224 wEta = fEtaWeights->GetBinContent(1+(Int_t)(TMath::Floor((dEta-fEtaMin)/fEtaBinWidth)));
226 // Calculate Re[Q_{m*n,k}] and Im[Q_{m*n,k}], (m = 1,2 and k = 0,1,2,3) for this event:
227 for(Int_t m=0;m<2;m++)
229 for(Int_t k=0;k<4;k++) // to be improved (what is maximum k that I need?)
231 (*fReQnk)(m,k)+=pow(wPhi*wPt*wEta,k)*TMath::Cos((m+1)*n*dPhi);
232 (*fImQnk)(m,k)+=pow(wPhi*wPt*wEta,k)*TMath::Sin((m+1)*n*dPhi);
235 // Calculate partially S_{p,k} for this event (final calculation of S_{p,k} follows after the loop over data bellow):
236 for(Int_t p=0;p<4;p++) // to be improved (what is maximum p that I need?)
238 for(Int_t k=0;k<4;k++) // to be improved (what is maximum k that I need?)
240 (*fSpk)(p,k)+=pow(wPhi*wPt*wEta,k);
243 } // end of if(aftsTrack->InRPSelection())
245 if(aftsTrack->InPOISelection()) // 1st POI
247 Double_t dPsi1 = aftsTrack->Phi();
248 Double_t dPt1 = aftsTrack->Pt();
249 for(Int_t j=i+1;j<nPrim;j++)
251 aftsTrack=anEvent->GetTrack(j);
252 if(aftsTrack->InPOISelection()) // 2nd POI
254 Double_t dPsi2 = aftsTrack->Phi();
255 Double_t dPt2 = aftsTrack->Pt();
257 fRePEBE[0]->Fill((dPt1+dPt2)/2.,TMath::Cos(n*(dPsi1+dPsi2)),1.);
258 fImPEBE[0]->Fill((dPt1+dPt2)/2.,TMath::Sin(n*(dPsi1+dPsi2)),1.);
259 fRePEBE[1]->Fill(TMath::Abs(dPt1-dPt2),TMath::Cos(n*(dPsi1+dPsi2)),1.);
260 fImPEBE[1]->Fill(TMath::Abs(dPt1-dPt2),TMath::Sin(n*(dPsi1+dPsi2)),1.);
264 } else // to if(aftsTrack)
267 cout<<" WARNING (MH): No particle! (i.e. aftsTrack is a NULL pointer in AFAWMH::Make().)"<<endl;
270 } // end of for(Int_t i=0;i<nPrim;i++)
272 // Calculate the final expressions for S_{p,k}:
273 for(Int_t p=0;p<4;p++) // to be improved (what is maximum p that I need?)
275 for(Int_t k=0;k<4;k++) // to be improved (what is maximum that I need?)
277 (*fSpk)(p,k)=pow((*fSpk)(p,k),p+1);
281 // e) Calculate 3-p azimuthal correlator in terms of Q_{n,k} and S_{p,k}:
282 if(anEvent->GetEventNSelTracksRP() >= 3)
284 this->Calculate3pCorrelator();
286 if(anEvent->GetEventNSelTracksRP() >= 0) // to be improved (is this correct if statement?)
288 this->CalculateNonIsotropicTerms();
291 // f) Calculate differential 3-p azimuthal correlator cos(2phi1-psi2-psi3) in terms of Q_{2n} and p_{n}:
292 this->CalculateDifferential3pCorrelator(); // to be improved - add relevant if statements for the min No of RPs and POIs
294 // g) Reset all event-by-event quantities:
295 this->ResetEventByEventQuantities();
297 } // end of AliFlowAnalysisWithMixedHarmonics::Make(AliFlowEventSimple* anEvent)
299 //================================================================================================================
301 void AliFlowAnalysisWithMixedHarmonics::Finish()
303 // Calculate the final results.
305 // a) Check all pointers used in this method;
306 // b) Access settings for analysis with mixed harmonics;
307 // c) Correct for detector effects;
308 // d) Print on the screen the final results.
310 this->CheckPointersUsedInFinish();
311 this->AccessSettings();
312 if(fCorrectForDetectorEffects)
314 this->CorrectForDetectorEffects();
315 this->CorrectForDetectorEffectsVsM();
317 if(fPrintOnTheScreen) this->PrintOnTheScreen();
319 } // end of AliFlowAnalysisWithMixedHarmonics::Finish()
321 //================================================================================================================
323 void AliFlowAnalysisWithMixedHarmonics::GetOutputHistograms(TList *outputListHistos)
325 // Get pointers to all objects saved in the output file.
327 // a) Get pointers for common control histograms.
330 this->SetHistList(outputListHistos);
334 cout<<" WARNING (MH): fHistList is NULL in AFAWMH::GOH() !!!!"<<endl;
338 this->GetPointersForBaseHistograms();
339 this->GetPointersForCommonHistograms();
340 this->GetPointersForAllEventProfiles();
341 this->GetPointersForResultsHistograms();
345 cout<<" WARNING (MH): outputListHistos is NULL in AFAWMH::GOH() !!!!"<<endl;
350 } // end of void AliFlowAnalysisWithMixedHarmonics::GetOutputHistograms(TList *outputListHistos)
352 //================================================================================================================
354 void AliFlowAnalysisWithMixedHarmonics::GetPointersForBaseHistograms()
356 // Get pointers to base histograms.
358 TString analysisSettingsName = "fAnalysisSettings";
359 TProfile *analysisSettings = dynamic_cast<TProfile*>(fHistList->FindObject(analysisSettingsName.Data()));
362 this->SetAnalysisSettings(analysisSettings);
366 cout<<" WARNING (MH): analysisSettings is NULL in AFAWMH::GPFBH() !!!!"<<endl;
371 } // end of void AliFlowAnalysisWithMixedHarmonics::GetPointersForBaseHistograms()
373 //================================================================================================================
375 void AliFlowAnalysisWithMixedHarmonics::GetPointersForCommonHistograms()
377 // Get pointers to common control histograms.
379 TString commonHistsName = "AliFlowCommonHistMH";
380 AliFlowCommonHist *commonHist = dynamic_cast<AliFlowCommonHist*>(fHistList->FindObject(commonHistsName.Data()));
383 this->SetCommonHists(commonHist);
387 cout<<" WARNING (MH): commonHist is NULL in AFAWMH::GPFCH() !!!!"<<endl;
392 } // end of void AliFlowAnalysisWithMixedHarmonics::GetPointersForCommonHistograms()
394 //================================================================================================================
396 void AliFlowAnalysisWithMixedHarmonics::GetPointersForAllEventProfiles()
398 // Get pointers to profiles holding final results.
400 TList *profileList = NULL;
401 profileList = dynamic_cast<TList*>(fHistList->FindObject("Profiles"));
404 cout<<"WARNING: profileList is NULL in AFAWMH::GPFAEP() !!!!"<<endl;
408 TString s3pCorrelatorProName = "f3pCorrelatorPro";
409 TProfile *p3pCorrelatorPro = dynamic_cast<TProfile*>(profileList->FindObject(s3pCorrelatorProName.Data()));
412 this->Set3pCorrelatorPro(p3pCorrelatorPro);
414 TString s3pCorrelatorVsMProName = "f3pCorrelatorVsMPro";
415 TProfile *p3pCorrelatorVsMPro = dynamic_cast<TProfile*>(profileList->FindObject(s3pCorrelatorVsMProName.Data()));
416 if(p3pCorrelatorVsMPro)
418 this->Set3pCorrelatorVsMPro(p3pCorrelatorVsMPro);
420 TString nonIsotropicTermsProName = "fNonIsotropicTermsPro";
421 TProfile *nonIsotropicTermsPro = dynamic_cast<TProfile*>(profileList->FindObject(nonIsotropicTermsProName.Data()));
422 if(nonIsotropicTermsPro)
424 this->SetNonIsotropicTermsPro(nonIsotropicTermsPro);
426 TString nonIsotropicTermsVsMProName = "fNonIsotropicTermsVsMPro";
427 TProfile2D *nonIsotropicTermsVsMPro = dynamic_cast<TProfile2D*>(profileList->FindObject(nonIsotropicTermsVsMProName.Data()));
428 if(nonIsotropicTermsVsMPro)
430 this->SetNonIsotropicTermsVsMPro(nonIsotropicTermsVsMPro);
432 TString psdFlag[2] = {"PtSum","PtDiff"};
433 for(Int_t sd=0;sd<2;sd++)
435 TProfile *p3pCorrelatorVsPtSumDiffPro = dynamic_cast<TProfile*>(profileList->FindObject(Form("f3pCorrelatorVs%sPro",psdFlag[sd].Data())));
436 if(p3pCorrelatorVsPtSumDiffPro)
438 this->Set3pCorrelatorVsPtSumDiffPro(p3pCorrelatorVsPtSumDiffPro,sd);
442 } // end of void AliFlowAnalysisWithMixedHarmonics::GetPointersForAllEventProfiles()
444 //================================================================================================================
446 void AliFlowAnalysisWithMixedHarmonics::GetPointersForResultsHistograms()
448 // Get pointers to histograms holding final results.
450 TList *resultsList = NULL;
451 resultsList = dynamic_cast<TList*>(fHistList->FindObject("Results"));
454 cout<<"WARNING: resultsList is NULL in AFAWMH::GPFRH() !!!!"<<endl;
457 TString s3pCorrelatorHistName = "f3pCorrelatorHist";
458 TH1D *h3pCorrelatorHist = dynamic_cast<TH1D*>(resultsList->FindObject(s3pCorrelatorHistName.Data()));
459 if(h3pCorrelatorHist)
461 this->Set3pCorrelatorHist(h3pCorrelatorHist);
463 TString s3pCorrelatorVsMHistName = "f3pCorrelatorVsMHist";
464 TH1D *h3pCorrelatorVsMHist = dynamic_cast<TH1D*>(resultsList->FindObject(s3pCorrelatorVsMHistName.Data()));
465 if(h3pCorrelatorVsMHist)
467 this->Set3pCorrelatorVsMHist(h3pCorrelatorVsMHist);
469 TString detectorBiasHistName = "fDetectorBiasHist";
470 TH1D *detectorBiasHist = dynamic_cast<TH1D*>(resultsList->FindObject(detectorBiasHistName.Data()));
473 this->SetDetectorBiasHist(detectorBiasHist);
475 TString detectorBiasVsMHistName = "fDetectorBiasVsMHist";
476 TH1D *detectorBiasVsMHist = dynamic_cast<TH1D*>(resultsList->FindObject(detectorBiasVsMHistName.Data()));
477 if(detectorBiasVsMHist)
479 this->SetDetectorBiasVsMHist(detectorBiasVsMHist);
482 } // end of void AliFlowAnalysisWithMixedHarmonics::GetPointersForResultsHistograms()
484 //================================================================================================================
486 void AliFlowAnalysisWithMixedHarmonics::WriteHistograms(TString outputFileName)
488 // Store the final results in output .root file.
489 TFile *output = new TFile(outputFileName.Data(),"RECREATE");
490 fHistList->Write(fHistList->GetName(),TObject::kSingleKey);
494 //================================================================================================================
496 void AliFlowAnalysisWithMixedHarmonics::WriteHistograms(TDirectoryFile *outputFileName)
498 // Store the final results in output .root file.
499 fHistList->SetName("cobjMH");
500 fHistList->SetOwner(kTRUE);
501 outputFileName->Add(fHistList);
502 outputFileName->Write(outputFileName->GetName(),TObject::kSingleKey);
505 //================================================================================================================
507 void AliFlowAnalysisWithMixedHarmonics::InitializeArrays()
509 // Initialize arrays.
511 for(Int_t sd=0;sd<2;sd++)
515 f3pCorrelatorVsPtSumDiffPro[sd] = NULL;
518 } // end of AliFlowAnalysisWithMixedHarmonics::InitializeArrays()
520 //================================================================================================================
522 void AliFlowAnalysisWithMixedHarmonics::BookAndNestAllLists()
524 // Book and nest all list in base list fHistList.
527 fWeightsList->SetName("Weights");
528 fWeightsList->SetOwner(kTRUE);
529 fHistList->Add(fWeightsList);
531 fProfileList->SetName("Profiles");
532 fProfileList->SetOwner(kTRUE);
533 fHistList->Add(fProfileList);
535 fResultsList->SetName("Results");
536 fResultsList->SetOwner(kTRUE);
537 fHistList->Add(fResultsList);
539 } // end of void AliFlowAnalysisWithMixedHarmonics::BookAndNestAllLists()
541 //================================================================================================================
543 void AliFlowAnalysisWithMixedHarmonics::BookProfileHoldingSettings()
545 // Book profile to hold all analysis settings.
547 TString analysisSettingsName = "fAnalysisSettings";
548 fAnalysisSettings = new TProfile(analysisSettingsName.Data(),"Settings for analysis with mixed harmonics",6,0,6);
549 fAnalysisSettings->GetXaxis()->SetLabelSize(0.025);
550 fAnalysisSettings->GetXaxis()->SetBinLabel(1,"Integer n in cos(n(2#phi_{1}-#phi_{2}-#phi_{3}))");
551 fAnalysisSettings->Fill(0.5,fCorrelatorInteger);
552 fAnalysisSettings->GetXaxis()->SetBinLabel(2,"Corrected for detector effects?");
553 fAnalysisSettings->Fill(1.5,(Int_t)fCorrectForDetectorEffects);
554 fAnalysisSettings->GetXaxis()->SetBinLabel(3,"# of multiplicity bins");
555 fAnalysisSettings->Fill(2.5,fNoOfMultipicityBins);
556 fAnalysisSettings->GetXaxis()->SetBinLabel(4,"Width of multiplicity bins");
557 fAnalysisSettings->Fill(3.5,fMultipicityBinWidth);
558 fAnalysisSettings->GetXaxis()->SetBinLabel(5,"Minimal multiplicity");
559 fAnalysisSettings->Fill(4.5,fMinMultiplicity);
560 fAnalysisSettings->GetXaxis()->SetBinLabel(6,"Print on the screen?");
561 fAnalysisSettings->Fill(5.5,(Int_t)fPrintOnTheScreen);
562 fHistList->Add(fAnalysisSettings);
564 } // end of void AliFlowAnalysisWithMixedHarmonics::BookProfileHoldingSettings()
566 //================================================================================================================
568 void AliFlowAnalysisWithMixedHarmonics::BookCommonHistograms()
570 // Book common control histograms and common histograms for final results.
572 TString commonHistsName = "AliFlowCommonHistMH";
573 fCommonHists = new AliFlowCommonHist(commonHistsName.Data());
574 fHistList->Add(fCommonHists);
576 } // end of void AliFlowAnalysisWithMixedHarmonics::BookCommonHistograms()
578 //================================================================================================================
580 void AliFlowAnalysisWithMixedHarmonics::BookAllEventByEventQuantities()
582 // Book all event-by-event quantitites.
584 // Q_{n,k} and S{p,k}:
585 fReQnk = new TMatrixD(2,4); // to be improved(bound on k)
586 fImQnk = new TMatrixD(2,4); // to be improved(bound on k)
587 fSpk = new TMatrixD(4,4); // to be improved(bound on p and k)
589 // p_n vs [(p1+p2)/2,|p1-p2|]
590 TString psdFlag[2] = {"PtSum","PtDiff"};
591 for(Int_t sd=0;sd<2;sd++)
593 // to be improved: hardwired ,fnBinsPt,0.,fPtMax):
594 fRePEBE[sd] = new TProfile(Form("fRePEBE%s",psdFlag[sd].Data()),"",fnBinsPt,0.,fPtMax);
595 fImPEBE[sd] = new TProfile(Form("fImPEBE%s",psdFlag[sd].Data()),"",fnBinsPt,0.,fPtMax);
598 } // end fo void AliFlowAnalysisWithMixedHarmonics::BookAllEventByEventQuantities()
600 //================================================================================================================
602 void AliFlowAnalysisWithMixedHarmonics::BookAllAllEventQuantities()
604 // Book all all-event quantitites.
606 // a) Book quantites without multiplicity binning;
607 // b) Book quantites with multiplicity binning;
608 // c) Book quantities with binning in (p1+p2)/2, |p1-p2|, (eta1+eta2)/2 or |eta1-eta2|.
610 // a) Book quantites without multiplicity binning:
611 // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> for all events:
612 TString s3pCorrelatorProName = "f3pCorrelatorPro";
613 f3pCorrelatorPro = new TProfile(s3pCorrelatorProName.Data(),"",1,0,1);
614 f3pCorrelatorPro->GetXaxis()->SetLabelOffset(0.01);
615 f3pCorrelatorPro->GetXaxis()->SetLabelSize(0.05);
616 if(fCorrelatorInteger == 1)
618 f3pCorrelatorPro->GetXaxis()->SetBinLabel(1,"cos(2#phi_{1}-#phi_{2}-#phi_{3})");
621 f3pCorrelatorPro->GetXaxis()->SetBinLabel(1,Form("cos[%d(2#phi_{1}-#phi_{2}-#phi_{3})]",fCorrelatorInteger));
623 fProfileList->Add(f3pCorrelatorPro);
624 // Non-isotropic terms in the decomposition of <cos[n(2phi1-phi2-phi3)]:
625 TString nonIsotropicTermsProName = "fNonIsotropicTermsPro";
626 fNonIsotropicTermsPro = new TProfile(nonIsotropicTermsProName.Data(),"",10,0,10);
627 if(fCorrelatorInteger == 1)
629 fNonIsotropicTermsPro->SetTitle("Non-isotropic terms in decomposition of cos(2#phi_{1}-#phi_{2}-#phi_{3})");
630 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(1,"cos(#phi_{1})");
631 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(2,"sin(#phi_{1})");
632 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(3,"cos(2#phi_{1})");
633 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(4,"sin(2#phi_{1})");
634 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(5,"cos(#phi_{1}+#phi_{2})");
635 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(6,"sin(#phi_{1}+#phi_{2})");
636 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(7,"cos(2#phi_{1}-#phi_{2})");
637 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(8,"sin(2#phi_{1}-#phi_{2})");
638 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(9,"cos(#phi_{1}-#phi_{2}-#phi_{3})");
639 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(10,"sin(#phi_{1}-#phi_{2}-#phi_{3})");
642 fNonIsotropicTermsPro->SetTitle(Form("Non-isotropic terms in decomposition of cos[%d(2#phi_{1}-#phi_{2}-#phi_{3})]",fCorrelatorInteger));
643 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(1,Form("cos(%d#phi_{1})",fCorrelatorInteger));
644 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(2,Form("sin(%d#phi_{1})",fCorrelatorInteger));
645 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(3,Form("cos(%d#phi_{1})",2*fCorrelatorInteger));
646 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(4,Form("sin(%d#phi_{1})",2*fCorrelatorInteger));
647 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(5,Form("cos(%d(#phi_{1}+#phi_{2}))",fCorrelatorInteger));
648 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(6,Form("sin(%d(#phi_{1}+#phi_{2}))",fCorrelatorInteger));
649 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(7,Form("cos(%d(2#phi_{1}-#phi_{2}))",fCorrelatorInteger));
650 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(8,Form("sin(%d(2#phi_{1}-#phi_{2}))",fCorrelatorInteger));
651 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(9,Form("cos(%d(#phi_{1}-#phi_{2}-#phi_{3}))",fCorrelatorInteger));
652 fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(10,Form("sin(%d(#phi_{1}-#phi_{2}-#phi_{3}))",fCorrelatorInteger));
654 fProfileList->Add(fNonIsotropicTermsPro);
655 // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> corrected for detector effects:
656 TString s3pCorrelatorHistName = "f3pCorrelatorHist";
657 f3pCorrelatorHist = new TH1D(s3pCorrelatorHistName.Data(),"",1,0,1);
658 f3pCorrelatorHist->GetXaxis()->SetLabelOffset(0.01);
659 f3pCorrelatorHist->GetXaxis()->SetLabelSize(0.05);
660 if(fCorrelatorInteger == 1)
662 f3pCorrelatorHist->GetXaxis()->SetBinLabel(1,"cos(2#phi_{1}-#phi_{2}-#phi_{3})");
665 f3pCorrelatorHist->GetXaxis()->SetBinLabel(1,Form("cos[%d(2#phi_{1}-#phi_{2}-#phi_{3})]",fCorrelatorInteger));
667 fResultsList->Add(f3pCorrelatorHist);
668 // Quantified bias comming from detector inefficiencies to 3-p correlator <<cos[n(2phi1-phi2-phi3)]>>:
669 TString detectorBiasHistName = "fDetectorBiasHist";
670 fDetectorBiasHist = new TH1D(detectorBiasHistName.Data(),"Bias coming from detector inefficiences",1,0,1);
671 fDetectorBiasHist->GetXaxis()->SetBinLabel(1,"#frac{corrected}{measured}");
672 fResultsList->Add(fDetectorBiasHist);
674 // b) Book quantites with multiplicity binning.
675 // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> versus multiplicity:
676 TString s3pCorrelatorVsMProName = "f3pCorrelatorVsMPro";
677 f3pCorrelatorVsMPro = new TProfile(s3pCorrelatorVsMProName.Data(),"",fNoOfMultipicityBins+2,0,fNoOfMultipicityBins+2);
678 if(fCorrelatorInteger == 1)
680 f3pCorrelatorVsMPro->SetTitle("cos(2#phi_{1}-#phi_{2}-#phi_{3}) #font[72]{vs} M");
683 f3pCorrelatorVsMPro->SetTitle(Form("cos[%d(2#phi_{1}-#phi_{2}-#phi_{3})] #font[72]{vs} M",fCorrelatorInteger));
685 f3pCorrelatorVsMPro->GetXaxis()->SetBinLabel(1,Form("M < %d",(Int_t)fMinMultiplicity));
686 for(Int_t b=2;b<=fNoOfMultipicityBins+1;b++)
688 f3pCorrelatorVsMPro->GetXaxis()->SetBinLabel(b,Form("%d #leq M < %d",(Int_t)(fMinMultiplicity+(b-2)*fMultipicityBinWidth),(Int_t)(fMinMultiplicity+(b-1)*fMultipicityBinWidth)));
690 f3pCorrelatorVsMPro->GetXaxis()->SetBinLabel(fNoOfMultipicityBins+2,Form(" M #geq %d",(Int_t)(fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)));
691 fProfileList->Add(f3pCorrelatorVsMPro);
692 // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> versus multiplicity corrected for detector effects:
693 TString s3pCorrelatorVsMHistName = "f3pCorrelatorVsMHist";
694 f3pCorrelatorVsMHist = new TH1D(s3pCorrelatorVsMHistName.Data(),"",fNoOfMultipicityBins+2,0,fNoOfMultipicityBins+2);
695 if(fCorrelatorInteger == 1)
697 f3pCorrelatorVsMHist->SetTitle("cos(2#phi_{1}-#phi_{2}-#phi_{3}) #font[72]{vs} M");
700 f3pCorrelatorVsMHist->SetTitle(Form("cos[%d(2#phi_{1}-#phi_{2}-#phi_{3})] #font[72]{vs} M",fCorrelatorInteger));
702 f3pCorrelatorVsMHist->GetXaxis()->SetBinLabel(1,Form("M < %d",(Int_t)fMinMultiplicity));
703 for(Int_t b=2;b<=fNoOfMultipicityBins+1;b++)
705 f3pCorrelatorVsMHist->GetXaxis()->SetBinLabel(b,Form("%d #leq M < %d",(Int_t)(fMinMultiplicity+(b-2)*fMultipicityBinWidth),(Int_t)(fMinMultiplicity+(b-1)*fMultipicityBinWidth)));
707 f3pCorrelatorVsMHist->GetXaxis()->SetBinLabel(fNoOfMultipicityBins+2,Form(" M #geq %d",(Int_t)(fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)));
708 fResultsList->Add(f3pCorrelatorVsMHist);
709 // Non-isotropic terms in the decomposition of <cos[n(2phi1-phi2-phi3)] vs multiplicity:
710 TString nonIsotropicTermsVsMProName = "fNonIsotropicTermsVsMPro";
711 fNonIsotropicTermsVsMPro = new TProfile2D(nonIsotropicTermsVsMProName.Data(),"Non-isotropic terms in the decomposition of cos[n(2phi1-phi2-phi3)] #font[72]{vs} M",10,0,10,fNoOfMultipicityBins+2,0,fNoOfMultipicityBins+2);
712 if(fCorrelatorInteger == 1)
714 fNonIsotropicTermsVsMPro->SetTitle("Non-isotropic terms in decomposition of cos(2#phi_{1}-#phi_{2}-#phi_{3}) #font[72]{vs} M");
715 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(1,"cos(#phi_{1})");
716 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(2,"sin(#phi_{1})");
717 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(3,"cos(2#phi_{1})");
718 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(4,"sin(2#phi_{1})");
719 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(5,"cos(#phi_{1}+#phi_{2})");
720 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(6,"sin(#phi_{1}+#phi_{2})");
721 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(7,"cos(2#phi_{1}-#phi_{2})");
722 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(8,"sin(2#phi_{1}-#phi_{2})");
723 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(9,"cos(#phi_{1}-#phi_{2}-#phi_{3})");
724 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(10,"sin(#phi_{1}-#phi_{2}-#phi_{3})");
727 fNonIsotropicTermsVsMPro->SetTitle(Form("Non-isotropic terms in decomposition of cos[%d(2#phi_{1}-#phi_{2}-#phi_{3})]",fCorrelatorInteger));
728 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(1,Form("cos(%d#phi_{1})",fCorrelatorInteger));
729 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(2,Form("sin(%d#phi_{1})",fCorrelatorInteger));
730 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(3,Form("cos(%d#phi_{1})",2*fCorrelatorInteger));
731 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(4,Form("sin(%d#phi_{1})",2*fCorrelatorInteger));
732 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(5,Form("cos(%d(#phi_{1}+#phi_{2}))",fCorrelatorInteger));
733 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(6,Form("sin(%d(#phi_{1}+#phi_{2}))",fCorrelatorInteger));
734 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(7,Form("cos(%d(2#phi_{1}-#phi_{2}))",fCorrelatorInteger));
735 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(8,Form("sin(%d(2#phi_{1}-#phi_{2}))",fCorrelatorInteger));
736 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(9,Form("cos(%d(#phi_{1}-#phi_{2}-#phi_{3}))",fCorrelatorInteger));
737 fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(10,Form("sin(%d(#phi_{1}-#phi_{2}-#phi_{3}))",fCorrelatorInteger));
739 fNonIsotropicTermsVsMPro->GetYaxis()->SetBinLabel(1,Form("M < %d",(Int_t)fMinMultiplicity));
740 for(Int_t b=2;b<=fNoOfMultipicityBins+1;b++)
742 fNonIsotropicTermsVsMPro->GetYaxis()->SetBinLabel(b,Form("%d #leq M < %d",(Int_t)(fMinMultiplicity+(b-2)*fMultipicityBinWidth),(Int_t)(fMinMultiplicity+(b-1)*fMultipicityBinWidth)));
744 fNonIsotropicTermsVsMPro->GetYaxis()->SetBinLabel(fNoOfMultipicityBins+2,Form(" M #geq %d",(Int_t)(fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)));
745 fProfileList->Add(fNonIsotropicTermsVsMPro);
746 // Quantified bias comming from detector inefficiencies to 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> vs multiplicity
747 TString detectorBiasVsMHistName = "fDetectorBiasVsMHist";
748 fDetectorBiasVsMHist = new TH1D(detectorBiasVsMHistName.Data(),"",fNoOfMultipicityBins+2,0,fNoOfMultipicityBins+2);
749 if(fCorrelatorInteger == 1)
751 fDetectorBiasVsMHist->SetTitle("#frac{corrected}{measured} cos(2#phi_{1}-#phi_{2}-#phi_{3}) #font[72]{vs} M"); // to be improved (title)
754 fDetectorBiasVsMHist->SetTitle(Form("#frac{corrected}{measured} cos[%d(2#phi_{1}-#phi_{2}-#phi_{3})] #font[72]{vs} M",fCorrelatorInteger)); // to be improved (title)
756 fDetectorBiasVsMHist->GetXaxis()->SetBinLabel(1,Form("M < %d",(Int_t)fMinMultiplicity));
757 for(Int_t b=2;b<=fNoOfMultipicityBins+1;b++)
759 fDetectorBiasVsMHist->GetXaxis()->SetBinLabel(b,Form("%d #leq M < %d",(Int_t)(fMinMultiplicity+(b-2)*fMultipicityBinWidth),(Int_t)(fMinMultiplicity+(b-1)*fMultipicityBinWidth)));
761 fDetectorBiasVsMHist->GetXaxis()->SetBinLabel(fNoOfMultipicityBins+2,Form(" M #geq %d",(Int_t)(fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)));
762 fResultsList->Add(fDetectorBiasVsMHist);
764 // c) Book quantities with binning in (p1+p2)/2, |p1-p2|, (eta1+eta2)/2 or |eta1-eta2|:
765 TString psdFlag[2] = {"PtSum","PtDiff"};
766 TString psdTitleFlag[2] = {"(p_{t,1}+p_{t,2})/2","#left|p_{t,1}-p_{t,2}#right|"};
767 //TString s3pCorrelatorVsPtSumDiffProName = "f3pCorrelatorVsPtSumDiffPro";
768 for(Int_t sd=0;sd<2;sd++)
770 // to be improved: hardwired ,fnBinsPt,0.,fPtMax):
771 f3pCorrelatorVsPtSumDiffPro[sd] = new TProfile(Form("f3pCorrelatorVs%sPro",psdFlag[sd].Data()),"",fnBinsPt,0.,fPtMax);
772 //f3pCorrelatorVsPtSumDiffPro[sd]->SetLabelSize(0.05);
773 //f3pCorrelatorVsPtSumDiffPro[sd]->SetMarkerStyle(25);
774 f3pCorrelatorVsPtSumDiffPro[sd]->GetXaxis()->SetTitle(psdTitleFlag[sd].Data());
775 fProfileList->Add(f3pCorrelatorVsPtSumDiffPro[sd]);
778 } // end of void AliFlowAnalysisWithMixedHarmonics::BookAllAllEventQuantities()
780 //================================================================================================================
782 void AliFlowAnalysisWithMixedHarmonics::AccessConstants()
784 // Access needed common constants from AliFlowCommonConstants.
786 fnBinsPhi = AliFlowCommonConstants::GetMaster()->GetNbinsPhi();
787 fPhiMin = AliFlowCommonConstants::GetMaster()->GetPhiMin();
788 fPhiMax = AliFlowCommonConstants::GetMaster()->GetPhiMax();
789 if(fnBinsPhi) fPhiBinWidth = (fPhiMax-fPhiMin)/fnBinsPhi;
790 fnBinsPt = AliFlowCommonConstants::GetMaster()->GetNbinsPt();
791 fPtMin = AliFlowCommonConstants::GetMaster()->GetPtMin();
792 fPtMax = AliFlowCommonConstants::GetMaster()->GetPtMax();
793 if(fnBinsPt) fPtBinWidth = (fPtMax-fPtMin)/fnBinsPt;
794 fnBinsEta = AliFlowCommonConstants::GetMaster()->GetNbinsEta();
795 fEtaMin = AliFlowCommonConstants::GetMaster()->GetEtaMin();
796 fEtaMax = AliFlowCommonConstants::GetMaster()->GetEtaMax();
797 if(fnBinsEta) fEtaBinWidth = (fEtaMax-fEtaMin)/fnBinsEta;
799 } // end of void AliFlowAnalysisWithMixedHarmonics::AccessConstants()
801 //================================================================================================================
803 void AliFlowAnalysisWithMixedHarmonics::CrossCheckSettings()
805 // Cross-check if the user settings make sense.
807 } // end of void AliFlowAnalysisWithMixedHarmonics::CrossCheckSettings()
809 //================================================================================================================
811 void AliFlowAnalysisWithMixedHarmonics::BookAndFillWeightsHistograms()
813 // Book and fill (by accessing file "weights.root") histograms which hold phi, pt and eta weights.
817 cout<<"WARNING: fWeightsList is NULL in AFAWMH::BAFWH() !!!!"<<endl;
820 // Profile to hold flags for weights:
821 TString fUseParticleWeightsName = "fUseParticleWeightsMH";
822 fUseParticleWeights = new TProfile(fUseParticleWeightsName.Data(),"0 = particle weight not used, 1 = particle weight used ",3,0,3);
823 fUseParticleWeights->SetLabelSize(0.06);
824 (fUseParticleWeights->GetXaxis())->SetBinLabel(1,"w_{#phi}");
825 (fUseParticleWeights->GetXaxis())->SetBinLabel(2,"w_{p_{T}}");
826 (fUseParticleWeights->GetXaxis())->SetBinLabel(3,"w_{#eta}");
827 fUseParticleWeights->Fill(0.5,(Int_t)fUsePhiWeights);
828 fUseParticleWeights->Fill(1.5,(Int_t)fUsePtWeights);
829 fUseParticleWeights->Fill(2.5,(Int_t)fUseEtaWeights);
830 fWeightsList->Add(fUseParticleWeights);
834 if(fWeightsList->FindObject("phi_weights"))
836 fPhiWeights = dynamic_cast<TH1F*>(fWeightsList->FindObject("phi_weights"));
837 if(TMath::Abs(fPhiWeights->GetBinWidth(1)-fPhiBinWidth)>pow(10.,-6.))
840 cout<<"WARNING (MH): Inconsistent binning in histograms for phi-weights throughout the code."<<endl;
846 cout<<"WARNING (MH): fWeightsList->FindObject(\"phi_weights\") is NULL in AFAWMH::BAFWH() !!!!"<<endl;
849 } // end of if(fUsePhiWeights)
853 if(fWeightsList->FindObject("pt_weights"))
855 fPtWeights = dynamic_cast<TH1D*>(fWeightsList->FindObject("pt_weights"));
856 if(TMath::Abs(fPtWeights->GetBinWidth(1)-fPtBinWidth)>pow(10.,-6.))
859 cout<<"WARNING (MH): Inconsistent binning in histograms for pt-weights throughout the code."<<endl;
865 cout<<"WARNING (MH): fWeightsList->FindObject(\"pt_weights\") is NULL in AFAWMH::BAFWH() !!!!"<<endl;
868 } // end of if(fUsePtWeights)
872 if(fWeightsList->FindObject("eta_weights"))
874 fEtaWeights = dynamic_cast<TH1D*>(fWeightsList->FindObject("eta_weights"));
875 if(TMath::Abs(fEtaWeights->GetBinWidth(1)-fEtaBinWidth)>pow(10.,-6.))
878 cout<<"WARNING (MH): Inconsistent binning in histograms for eta-weights throughout the code."<<endl;
884 cout<<"WARNING: fUseEtaWeights && fWeightsList->FindObject(\"eta_weights\") is NULL in AFAWMH::BAFWH() !!!!"<<endl;
887 } // end of if(fUseEtaWeights)
889 } // end of AliFlowAnalysisWithMixedHarmonics::BookAndFillWeightsHistograms()
891 //================================================================================================================
893 void AliFlowAnalysisWithMixedHarmonics::CheckPointersUsedInMake()
895 // Check pointers used in method Make().
897 if(!fReQnk || !fImQnk || !fSpk )
900 cout<<" WARNING (MH::Make()): (!fReQnk || !fImQnk || !fSpk) is NULL !!!!"<<endl;
904 if(!f3pCorrelatorPro)
907 cout<<" WARNING (MH::Make()): (!f3pCorrelatorPro) is NULL !!!!"<<endl;
911 if(!fNonIsotropicTermsPro)
914 cout<<" WARNING (MH::Make()): !fNonIsotropicTermsPro is NULL !!!!"<<endl;
918 if(!f3pCorrelatorVsMPro)
921 cout<<" WARNING (MH::Make()): !f3pCorrelatorVsMPro is NULL !!!!"<<endl;
925 if(!fNonIsotropicTermsVsMPro)
928 cout<<" WARNING (MH::Make()): !fNonIsotropicTermsVsMPro is NULL !!!!"<<endl;
932 for(Int_t sd=0;sd<2;sd++)
934 if(!(f3pCorrelatorVsPtSumDiffPro[sd]))
937 cout<<" WARNING (MH::Make()): !"<<Form("f3pCorrelatorVsPtSumDiffPro[%d]",sd)<<" is NULL !!!!"<<endl;
943 } // end of AliFlowAnalysisWithMixedHarmonics::CheckPointersUsedInMake()
945 //================================================================================================================
947 void AliFlowAnalysisWithMixedHarmonics::CheckPointersUsedInFinish()
949 // Check pointers used in method Finish().
951 if(!fAnalysisSettings)
954 cout<<" WARNING (MH::Finish()): !fAnalysisSettings is NULL !!!!"<<endl;
958 if(!f3pCorrelatorPro)
961 cout<<" WARNING (MH::Finish()): !f3pCorrelatorPro is NULL !!!!"<<endl;
965 if(!fNonIsotropicTermsPro)
968 cout<<" WARNING (MH::Finish()): !fNonIsotropicTermsPro is NULL !!!!"<<endl;
972 if(!f3pCorrelatorVsMPro)
975 cout<<" WARNING (MH::Finish()): !f3pCorrelatorVsMPro is NULL !!!!"<<endl;
979 if(!f3pCorrelatorVsMHist)
982 cout<<" WARNING (MH::Finish()): !f3pCorrelatorVsMHist is NULL !!!!"<<endl;
986 if(!fNonIsotropicTermsVsMPro)
989 cout<<" WARNING (MH::Finish()): !fNonIsotropicTermsVsMPro is NULL !!!!"<<endl;
993 if(!f3pCorrelatorHist)
996 cout<<" WARNING (MH::Finish()): !f3pCorrelatorHist is NULL !!!!"<<endl;
1000 if(!fDetectorBiasHist)
1003 cout<<" WARNING (MH::Finish()): !fDetectorBiasHist is NULL !!!!"<<endl;
1007 /* to be improved - enabled eventually
1008 if(!fDetectorBiasVsMHist)
1011 cout<<" WARNING (MH::Finish()): !fDetectorBiasVsMHist is NULL !!!!"<<endl;
1016 for(Int_t sd=0;sd<2;sd++)
1018 if(!(f3pCorrelatorVsPtSumDiffPro[sd]))
1021 cout<<" WARNING (MH::Finish()): !"<<Form("f3pCorrelatorVsPtSumDiffPro[%d]",sd)<<" is NULL !!!!"<<endl;
1027 } // end of AliFlowAnalysisWithMixedHarmonics::CheckPointersUsedInFinish()
1029 //================================================================================================================
1031 void AliFlowAnalysisWithMixedHarmonics::PrintOnTheScreen()
1033 // Print the final results on the screen.
1036 cout<<"*******************************************************"<<endl;
1037 cout<<"*******************************************************"<<endl;
1038 cout<<" Mixed Harmonics "<<endl;
1040 if(fCorrelatorInteger!=1)
1042 cout<<" cos["<<fCorrelatorInteger<<"(2phi1-phi2-phi3)] = "<<f3pCorrelatorHist->GetBinContent(1)<<
1043 " +/- "<<f3pCorrelatorHist->GetBinError(1)<<endl;
1046 cout<<" cos(2phi1-phi2-phi3) = "<<f3pCorrelatorHist->GetBinContent(1)<<" +/- "<<f3pCorrelatorHist->GetBinError(1)<<endl;
1048 cout<<" Detector Bias = "<<fDetectorBiasHist->GetBinContent(1)<<endl;
1050 cout<<"*******************************************************"<<endl;
1051 cout<<"*******************************************************"<<endl;
1053 } // end of void AliFlowAnalysisWithMixedHarmonics::PrintOnTheScreen()
1055 //================================================================================================================
1057 void AliFlowAnalysisWithMixedHarmonics::AccessSettings()
1059 // Access the settings for analysis with mixed harmonics.
1061 fCorrelatorInteger = (Int_t)fAnalysisSettings->GetBinContent(1);
1062 fCorrectForDetectorEffects = (Bool_t)fAnalysisSettings->GetBinContent(2);
1063 fNoOfMultipicityBins = (Int_t)fAnalysisSettings->GetBinContent(3);
1064 fMultipicityBinWidth = (Double_t)fAnalysisSettings->GetBinContent(4);
1065 fMinMultiplicity = (Double_t)fAnalysisSettings->GetBinContent(5);
1066 fPrintOnTheScreen = (Bool_t)fAnalysisSettings->GetBinContent(6);
1068 } // end of AliFlowAnalysisWithMixedHarmonics::AccessSettings()
1070 //================================================================================================================
1072 void AliFlowAnalysisWithMixedHarmonics::CorrectForDetectorEffects()
1074 // Correct measured 3-p correlator cos[n(2phi1-phi2-phi3)] for detector effects.
1076 Double_t measured3pCorrelator = f3pCorrelatorPro->GetBinContent(1); // biased by detector effects
1077 Double_t corrected3pCorrelator = 0.; // corrected for detector effects
1078 Double_t nonIsotropicTerms[10] = {0.}; // there are 10 distinct non-isotropic terms
1079 for(Int_t nit=0;nit<10;nit++)
1081 nonIsotropicTerms[nit] = fNonIsotropicTermsPro->GetBinContent(nit+1);
1083 // Calculate corrected 3-p correlator:
1084 corrected3pCorrelator = measured3pCorrelator
1085 - nonIsotropicTerms[2]*nonIsotropicTerms[4]
1086 - nonIsotropicTerms[3]*nonIsotropicTerms[5]
1087 - 2.*nonIsotropicTerms[0]*nonIsotropicTerms[6]
1088 - 2.*nonIsotropicTerms[1]*nonIsotropicTerms[7]
1089 + 2.*nonIsotropicTerms[2]*(pow(nonIsotropicTerms[0],2.)-pow(nonIsotropicTerms[1],2.))
1090 + 4.*nonIsotropicTerms[3]*nonIsotropicTerms[0]*nonIsotropicTerms[1];
1091 // Store corrected correlator:
1092 f3pCorrelatorHist->SetBinContent(1,corrected3pCorrelator);
1093 f3pCorrelatorHist->SetBinError(1,f3pCorrelatorPro->GetBinError(1)); // to be improved (propagate error for non-isotropic terms)
1094 // Quantify bias from detector inefficiences to 3-p correlator. Remark: Bias is quantified as a
1095 // ratio between corrected and measured 3-p correlator:
1096 // bias = corrected/measured
1097 // This bias is stored in histogram fDetectorBias.
1099 if(measured3pCorrelator)
1101 bias = corrected3pCorrelator/measured3pCorrelator;
1102 fDetectorBiasHist->SetBinContent(1,bias);
1105 } // end of AliFlowAnalysisWithMixedHarmonics::CorrectForDetectorEffects()
1107 //================================================================================================================
1109 void AliFlowAnalysisWithMixedHarmonics::CorrectForDetectorEffectsVsM()
1111 // Correct measured 3-p correlator cos[n(2phi1-phi2-phi3)] vs M for detector effects.
1113 for(Int_t b=1;b<=fNoOfMultipicityBins+2;b++)
1115 Double_t measured3pCorrelator = f3pCorrelatorVsMPro->GetBinContent(b); // biased by detector effects
1116 Double_t corrected3pCorrelator = 0.; // corrected for detector effects
1117 Double_t nonIsotropicTerms[10] = {0.}; // there are 10 distinct non-isotropic terms
1118 for(Int_t nit=0;nit<10;nit++)
1120 nonIsotropicTerms[nit] = fNonIsotropicTermsVsMPro->GetBinContent(fNonIsotropicTermsVsMPro->GetBin(nit+1,b));
1122 // Calculate corrected 3-p correlator:
1123 corrected3pCorrelator = measured3pCorrelator
1124 - nonIsotropicTerms[2]*nonIsotropicTerms[4]
1125 - nonIsotropicTerms[3]*nonIsotropicTerms[5]
1126 - 2.*nonIsotropicTerms[0]*nonIsotropicTerms[6]
1127 - 2.*nonIsotropicTerms[1]*nonIsotropicTerms[7]
1128 + 2.*nonIsotropicTerms[2]*(pow(nonIsotropicTerms[0],2.)-pow(nonIsotropicTerms[1],2.))
1129 + 4.*nonIsotropicTerms[3]*nonIsotropicTerms[0]*nonIsotropicTerms[1];
1130 // Store corrected correlator:
1131 f3pCorrelatorVsMHist->SetBinContent(b,corrected3pCorrelator);
1132 f3pCorrelatorVsMHist->SetBinError(b,f3pCorrelatorVsMPro->GetBinError(b)); // to be improved (propagate error for non-isotropic terms)
1133 // Quantify bias from detector inefficiences to 3-p correlator. Remark: Bias is quantified as a
1134 // ratio between corrected and measured 3-p correlator:
1135 // bias = corrected/measured
1136 // This bias is stored in histogram fDetectorBias.
1138 if(measured3pCorrelator)
1140 bias = corrected3pCorrelator/measured3pCorrelator;
1141 fDetectorBiasVsMHist->SetBinContent(b,bias);
1143 } // end of for(Int_t b=1;b<=fNoOfMultipicityBins;b++)
1145 } // end of AliFlowAnalysisWithMixedHarmonics::CorrectForDetectorEffectsVsM()
1147 //================================================================================================================
1149 void AliFlowAnalysisWithMixedHarmonics::ResetEventByEventQuantities()
1151 // Reset all event by event quantities.
1157 for(Int_t sd=0;sd<2;sd++)
1159 fRePEBE[sd]->Reset();
1160 fImPEBE[sd]->Reset();
1163 } // end of void AliFlowAnalysisWithMixedHarmonics::ResetEventByEventQuantities()
1165 //================================================================================================================
1167 void AliFlowAnalysisWithMixedHarmonics::Calculate3pCorrelator()
1169 // Calculate 3-p azimuthal correlator cos[n(2phi1-phi2-phi3)] in terms of Q_{n,k} and S_{p,k}.
1171 // a) Calculate 3-p correlator without using particle weights;
1172 // b) Calculate 3-p correlator with using particle weights.
1174 // a) Calculate 3-p correlator without using particle weights:
1175 if(!(fUsePhiWeights || fUsePtWeights || fUseEtaWeights))
1177 // Multiplicity (number of RPs):
1178 Double_t dMult = (*fSpk)(0,0);
1179 // Real and imaginary parts of non-weighted Q-vectors (Q_{n,0}) evaluated in harmonics n and 2n:
1180 Double_t dReQ1n = (*fReQnk)(0,0);
1181 Double_t dReQ2n = (*fReQnk)(1,0);
1182 Double_t dImQ1n = (*fImQnk)(0,0);
1183 Double_t dImQ2n = (*fImQnk)(1,0);
1184 // 3-particle azimuthal correlator <cos(n*(2.*phi1-phi2-phi3))>:
1185 Double_t three2n1n1n = (pow(dReQ1n,2.)*dReQ2n + 2.*dReQ1n*dImQ1n*dImQ2n - pow(dImQ1n,2.)*dReQ2n
1186 - 2.*(pow(dReQ1n,2.)+pow(dImQ1n,2.))
1187 - (pow(dReQ2n,2.)+pow(dImQ2n,2.))+2.*dMult)
1188 / (dMult*(dMult-1.)*(dMult-2.));
1189 // Fill all-events profile:
1190 f3pCorrelatorPro->Fill(0.5,three2n1n1n,dMult*(dMult-1.)*(dMult-2.));
1192 // 3-particle azimuthal correlator <cos(n*(2.*phi1-phi2-phi3))> vs multiplicity:
1193 if(dMult<fMinMultiplicity)
1195 f3pCorrelatorVsMPro->Fill(0.5,three2n1n1n,dMult*(dMult-1.)*(dMult-2.));
1196 } else if(dMult>=fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)
1198 f3pCorrelatorVsMPro->Fill(0.5+fNoOfMultipicityBins+1,three2n1n1n,dMult*(dMult-1.)*(dMult-2.));
1201 f3pCorrelatorVsMPro->Fill(1.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),three2n1n1n,dMult*(dMult-1.)*(dMult-2.));
1204 } // end of if(!(fUsePhiWeights || fUsePtWeights || fUseEtaWeights))
1206 // b) Calculate 3-p correlator without using particle weights:
1207 if(fUsePhiWeights || fUsePtWeights || fUseEtaWeights)
1210 } // end of if(fUsePhiWeights || fUsePtWeights || fUseEtaWeights)
1212 } // end of void AliFlowAnalysisWithMixedHarmonics::Calculate3pCorrelator()
1214 //================================================================================================================
1216 void AliFlowAnalysisWithMixedHarmonics::CalculateNonIsotropicTerms()
1218 // Calculate non-isotropic terms which appear in the decomposition of 3-p correlator <cos[n(2phi1-phi2-phi3)]>.
1220 // For detector with uniform acceptance all these terms vanish. These non-isotropic terms are stored in fNonIsotropicTermsPro.
1221 // Binning of fNonIsotropicTermsPro is organized as follows:
1222 // 1st bin: <<cos(n*phi1)>>
1223 // 2nd bin: <<sin(n*phi1)>>
1224 // 3rd bin: <<cos(2n*phi1)>>
1225 // 4th bin: <<sin(2n*phi1)>>
1226 // 5th bin: <<cos(n*(phi1+phi2)>>
1227 // 6th bin: <<sin(n*(phi1+phi2)>>
1228 // 7th bin: <<cos(n*(2phi1-phi2)>>
1229 // 8th bin: <<sin(n*(2phi1-phi2)>>
1230 // 9th bin: <<cos(n*(phi1-phi2-phi3)>>
1231 // 10th bin: <<sin(n*(phi1-phi2-phi3)>>
1233 // a) Calculate using particle weights;
1234 // b) Calculate without using particle weights.
1236 // a) Calculate using particle weights:
1237 if(!(fUsePhiWeights || fUsePtWeights || fUseEtaWeights))
1239 // Multiplicity (number of RPs):
1240 Double_t dMult = (*fSpk)(0,0);
1241 // Real and imaginary parts of non-weighted Q-vectors (Q_{n,0}) evaluated in harmonics n and 2n:
1242 Double_t dReQ1n = (*fReQnk)(0,0);
1243 Double_t dReQ2n = (*fReQnk)(1,0);
1244 Double_t dImQ1n = (*fImQnk)(0,0);
1245 Double_t dImQ2n = (*fImQnk)(1,0);
1246 // 1-particle terms:
1247 Double_t cosP1n = 0.; // <cos(n*(phi1))>
1248 Double_t sinP1n = 0.; // <sin(n*(phi1))>
1249 Double_t cosP2n = 0.; // <cos(2n*(phi1))>
1250 Double_t sinP2n = 0.; // <sin(2n*(phi1))>
1253 cosP1n = dReQ1n/dMult;
1254 sinP1n = dImQ1n/dMult;
1255 cosP2n = dReQ2n/dMult;
1256 sinP2n = dImQ2n/dMult;
1257 // All-events avarages:
1258 fNonIsotropicTermsPro->Fill(0.5,cosP1n,dMult); // <<cos(n*(phi1))>>
1259 fNonIsotropicTermsPro->Fill(1.5,sinP1n,dMult); // <<sin(n*(phi1))>>
1260 fNonIsotropicTermsPro->Fill(2.5,cosP2n,dMult); // <<cos(2n*(phi1))>>
1261 fNonIsotropicTermsPro->Fill(3.5,sinP2n,dMult); // <<sin(2n*(phi1))>>
1262 // All-events avarages vs M:
1263 if(dMult<fMinMultiplicity)
1265 fNonIsotropicTermsVsMPro->Fill(0.5,0.5,cosP1n,dMult); // <<cos(n*(phi1))>>
1266 fNonIsotropicTermsVsMPro->Fill(1.5,0.5,sinP1n,dMult); // <<sin(n*(phi1))>>
1267 fNonIsotropicTermsVsMPro->Fill(2.5,0.5,cosP2n,dMult); // <<cos(2n*(phi1))>>
1268 fNonIsotropicTermsVsMPro->Fill(3.5,0.5,sinP2n,dMult); // <<sin(2n*(phi1))>>
1269 } else if(dMult>=fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)
1271 fNonIsotropicTermsVsMPro->Fill(0.5,0.5+fNoOfMultipicityBins+1,cosP1n,dMult); // <<cos(n*(phi1))>>
1272 fNonIsotropicTermsVsMPro->Fill(1.5,0.5+fNoOfMultipicityBins+1,sinP1n,dMult); // <<sin(n*(phi1))>>
1273 fNonIsotropicTermsVsMPro->Fill(2.5,0.5+fNoOfMultipicityBins+1,cosP2n,dMult); // <<cos(2n*(phi1))>>
1274 fNonIsotropicTermsVsMPro->Fill(3.5,0.5+fNoOfMultipicityBins+1,sinP2n,dMult); // <<sin(2n*(phi1))>>
1277 fNonIsotropicTermsVsMPro->Fill(0.5,1.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),cosP1n,dMult); // <<cos(n*(phi1))>>
1278 fNonIsotropicTermsVsMPro->Fill(1.5,1.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),sinP1n,dMult); // <<sin(n*(phi1))>>
1279 fNonIsotropicTermsVsMPro->Fill(2.5,1.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),cosP2n,dMult); // <<cos(2n*(phi1))>>
1280 fNonIsotropicTermsVsMPro->Fill(3.5,1.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),sinP2n,dMult); // <<sin(2n*(phi1))>>
1282 } // end of if(dMult>0)
1283 // 2-particle terms:
1284 Double_t cosP1nP1n = 0.; // <cos(n*(phi1+phi2))>
1285 Double_t sinP1nP1n = 0.; // <sin(n*(phi1+phi2))>
1286 Double_t cosP2nM1n = 0.; // <cos(n*(2phi1-phi2))>
1287 Double_t sinP2nM1n = 0.; // <sin(n*(2phi1-phi2))>
1290 cosP1nP1n = (pow(dReQ1n,2)-pow(dImQ1n,2)-dReQ2n)/(dMult*(dMult-1));
1291 sinP1nP1n = (2.*dReQ1n*dImQ1n-dImQ2n)/(dMult*(dMult-1));
1292 cosP2nM1n = (dReQ2n*dReQ1n+dImQ2n*dImQ1n-dReQ1n)/(dMult*(dMult-1));
1293 sinP2nM1n = (dImQ2n*dReQ1n-dReQ2n*dImQ1n-dImQ1n)/(dMult*(dMult-1));
1294 // All-events avarages:
1295 fNonIsotropicTermsPro->Fill(4.5,cosP1nP1n,dMult*(dMult-1.)); // <<cos(n*(phi1+phi2))>>
1296 fNonIsotropicTermsPro->Fill(5.5,sinP1nP1n,dMult*(dMult-1.)); // <<sin(n*(phi1+phi2))>>
1297 fNonIsotropicTermsPro->Fill(6.5,cosP2nM1n,dMult*(dMult-1.)); // <<cos(n*(2phi1-phi2))>>
1298 fNonIsotropicTermsPro->Fill(7.5,sinP2nM1n,dMult*(dMult-1.)); // <<sin(n*(2phi1-phi2))>>
1299 // All-events avarages vs M:
1300 if(dMult<fMinMultiplicity)
1302 fNonIsotropicTermsVsMPro->Fill(4.5,0.5,cosP1nP1n,dMult*(dMult-1.)); // <<cos(n*(phi1+phi2))>>
1303 fNonIsotropicTermsVsMPro->Fill(5.5,0.5,sinP1nP1n,dMult*(dMult-1.)); // <<sin(n*(phi1+phi2))>>
1304 fNonIsotropicTermsVsMPro->Fill(6.5,0.5,cosP2nM1n,dMult*(dMult-1.)); // <<cos(n*(2phi1-phi2))>>
1305 fNonIsotropicTermsVsMPro->Fill(7.5,0.5,sinP2nM1n,dMult*(dMult-1.)); // <<sin(n*(2phi1-phi2))>>
1306 } else if(dMult>=fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)
1308 fNonIsotropicTermsVsMPro->Fill(4.5,0.5+fNoOfMultipicityBins+1,cosP1nP1n,dMult*(dMult-1.)); // <<cos(n*(phi1+phi2))>>
1309 fNonIsotropicTermsVsMPro->Fill(5.5,0.5+fNoOfMultipicityBins+1,sinP1nP1n,dMult*(dMult-1.)); // <<sin(n*(phi1+phi2))>>
1310 fNonIsotropicTermsVsMPro->Fill(6.5,0.5+fNoOfMultipicityBins+1,cosP2nM1n,dMult*(dMult-1.)); // <<cos(n*(2phi1-phi2))>>
1311 fNonIsotropicTermsVsMPro->Fill(7.5,0.5+fNoOfMultipicityBins+1,sinP2nM1n,dMult*(dMult-1.)); // <<sin(n*(2phi1-phi2))>>
1314 fNonIsotropicTermsVsMPro->Fill(4.5,1.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),cosP1nP1n,dMult*(dMult-1.)); // <<cos(n*(phi1+phi2))>>
1315 fNonIsotropicTermsVsMPro->Fill(5.5,1.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),sinP1nP1n,dMult*(dMult-1.)); // <<sin(n*(phi1+phi2))>>
1316 fNonIsotropicTermsVsMPro->Fill(6.5,1.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),cosP2nM1n,dMult*(dMult-1.)); // <<cos(n*(2phi1-phi2))>>
1317 fNonIsotropicTermsVsMPro->Fill(7.5,1.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),sinP2nM1n,dMult*(dMult-1.)); // <<sin(n*(2phi1-phi2))>>
1319 } // end of if(dMult>1)
1321 Double_t cosP1nM1nM1n = 0.; // <cos(n*(phi1-phi2-phi3))>
1322 Double_t sinP1nM1nM1n = 0.; // <sin(n*(phi1-phi2-phi3))>
1325 cosP1nM1nM1n = (dReQ1n*(pow(dReQ1n,2)+pow(dImQ1n,2))-dReQ1n*dReQ2n-dImQ1n*dImQ2n-2.*(dMult-1)*dReQ1n)
1326 / (dMult*(dMult-1)*(dMult-2));
1327 sinP1nM1nM1n = (-dImQ1n*(pow(dReQ1n,2)+pow(dImQ1n,2))+dReQ1n*dImQ2n-dImQ1n*dReQ2n+2.*(dMult-1)*dImQ1n)
1328 / (dMult*(dMult-1)*(dMult-2));
1329 // All-events avarages:
1330 fNonIsotropicTermsPro->Fill(8.5,cosP1nM1nM1n,dMult*(dMult-1.)*(dMult-2.)); // <<cos(n*(phi1-phi2-phi3))>>
1331 fNonIsotropicTermsPro->Fill(9.5,sinP1nM1nM1n,dMult*(dMult-1.)*(dMult-2.)); // <<sin(n*(phi1-phi2-phi3))>>
1332 // All-events avarages vs M:
1333 if(dMult<fMinMultiplicity)
1335 fNonIsotropicTermsVsMPro->Fill(8.5,0.5,cosP1nM1nM1n,dMult*(dMult-1.)*(dMult-2.)); // <<cos(n*(phi1-phi2-phi3))>>
1336 fNonIsotropicTermsVsMPro->Fill(9.5,0.5,sinP1nM1nM1n,dMult*(dMult-1.)*(dMult-2.)); // <<sin(n*(phi1-phi2-phi3))>>
1337 } else if(dMult>=fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)
1339 fNonIsotropicTermsVsMPro->Fill(8.5,0.5+fNoOfMultipicityBins+1,cosP1nM1nM1n,dMult*(dMult-1.)*(dMult-2.)); // <<cos(n*(phi1-phi2-phi3))>>
1340 fNonIsotropicTermsVsMPro->Fill(9.5,0.5+fNoOfMultipicityBins+1,sinP1nM1nM1n,dMult*(dMult-1.)*(dMult-2.)); // <<sin(n*(phi1-phi2-phi3))>>
1343 fNonIsotropicTermsVsMPro->Fill(8.5,1.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),
1344 cosP1nM1nM1n,dMult*(dMult-1.)*(dMult-2.)); // <<cos(n*(phi1-phi2-phi3))>>
1345 fNonIsotropicTermsVsMPro->Fill(9.5,1.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),
1346 sinP1nM1nM1n,dMult*(dMult-1.)*(dMult-2.)); // <<sin(n*(phi1-phi2-phi3))>>
1348 } // end of if(dMult>2)
1349 } // end of if(!(fUsePhiWeights || fUsePtWeights || fUseEtaWeights))
1351 // b) Calculate without using particle weights:
1352 if(fUsePhiWeights || fUsePtWeights || fUseEtaWeights)
1355 } // end of if(fUsePhiWeights || fUsePtWeights || fUseEtaWeights)
1357 } // end of void AliFlowAnalysisWithMixedHarmonics::CalculateNonIsotropicTerms()
1359 //================================================================================================================
1361 void AliFlowAnalysisWithMixedHarmonics::CalculateDifferential3pCorrelator()
1363 // Calculate differential 3-p azimuthal correlator cos[n(2phi1-psi2-psi3)] in terms of Q_{2n} and p_{n}.
1365 // a) Calculate differential 3-p correlator without using particle weights;
1366 // b) Calculate differential 3-p correlator with using particle weights.
1368 // a) Calculate differential 3-p correlator without using particle weights:
1369 if(!(fUsePhiWeights || fUsePtWeights || fUseEtaWeights))
1371 // Multiplicity (number of RPs):
1372 Double_t dMult = (*fSpk)(0,0);
1373 // Real and imaginary parts of non-weighted Q-vectors (Q_{n,0}) evaluated in harmonic 2n:
1374 Double_t dReQ2n = (*fReQnk)(1,0);
1375 Double_t dImQ2n = (*fImQnk)(1,0);
1376 for(Int_t sd=0;sd<2;sd++) // [(p1+p2)/2,|p1-p2|]
1378 // looping over all bins and calculating reduced correlations:
1379 for(Int_t b=1;b<=fnBinsPt;b++)
1381 // real and imaginary parts of p_{n}:
1382 Double_t p1nRe = fRePEBE[sd]->GetBinContent(b)*fRePEBE[sd]->GetBinEntries(b);
1383 Double_t p1nIm = fImPEBE[sd]->GetBinContent(b)*fImPEBE[sd]->GetBinEntries(b);
1384 // number of pairs of POIs in particular (p1+p2)/2 or |p1-p2| bin:
1385 Double_t mp = fRePEBE[sd]->GetBinEntries(b);
1386 Double_t cosP2nphi1M1npsi2M1npsi2 = 0; // cos[n(2phi1-psi2-psi3)]
1389 cosP2nphi1M1npsi2M1npsi2 = (p1nRe*dReQ2n+p1nIm*dImQ2n)/(mp*dMult);
1391 f3pCorrelatorVsPtSumDiffPro[sd]->Fill(fPtMin+(b-1)*fPtBinWidth,cosP2nphi1M1npsi2M1npsi2,mp*dMult);
1392 } // end of for(Int_t b=1;b<=fnBinsPt;b++)
1393 } // end of for(Int_t sd=0;sd<2;sd++)
1394 } // end of if(!(fUsePhiWeights || fUsePtWeights || fUseEtaWeights))
1396 // b) Calculate differential 3-p correlator by using particle weights:
1397 if(fUsePhiWeights || fUsePtWeights || fUseEtaWeights)
1400 } // end of if(fUsePhiWeights || fUsePtWeights || fUseEtaWeights)
1402 } // end of void AliFlowAnalysisWithMixedHarmonics::CalculateDifferential3pCorrelator()
1404 //================================================================================================================