Added Histogram Results/f3pCorrelatorHist holds the value of cos[n(2phi1-phi2-phi3...
authorsnelling <snelling@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 30 Apr 2010 14:50:58 +0000 (14:50 +0000)
committersnelling <snelling@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 30 Apr 2010 14:50:58 +0000 (14:50 +0000)
PWG2/FLOW/AliFlowCommon/AliFlowAnalysisWithMixedHarmonics.cxx
PWG2/FLOW/AliFlowCommon/AliFlowAnalysisWithMixedHarmonics.h
PWG2/FLOW/AliFlowCommon/AliFlowAnalysisWithNestedLoops.cxx
PWG2/FLOW/AliFlowCommon/AliFlowAnalysisWithNestedLoops.h
PWG2/FLOW/AliFlowTasks/AliAnalysisTaskMixedHarmonics.h
PWG2/FLOW/AliFlowTasks/AliAnalysisTaskNestedLoops.h

index b4b6d9f..fee4c22 100644 (file)
@@ -22,7 +22,8 @@
  * can be used to:                                        *
  *                                                        *  
  *  a) Extract subdominant harmonics (like v1 and v4);    *
- *  b) Study strong parity violation.                     * 
+ *  b) Study flow of two-particle resonances;             *
+ *  c) Study strong parity violation.                     * 
  *                                                        * 
  * Author: Ante Bilandzic (abilandzic@gmail.com)          *
  *********************************************************/ 
@@ -55,11 +56,12 @@ fHistList(NULL),
 fHistListName(NULL),
 fAnalysisLabel(NULL),
 fAnalysisSettings(NULL),
-fCorrelatorInteger(2),
+fCorrelatorInteger(1),
 fNoOfMultipicityBins(10),
 fMultipicityBinWidth(2),
 fMinMultiplicity(1),
 fCorrectForDetectorEffects(kTRUE),
+fPrintOnTheScreen(kTRUE),
 fCommonHists(NULL),
 fnBinsPhi(0),
 fPhiMin(0),
@@ -94,13 +96,11 @@ fNonIsotropicTermsVsMPro(NULL),
 fResultsList(NULL),
 f3pCorrelatorHist(NULL),
 fDetectorBiasHist(NULL),
-fNonIsotropicTermsHist(NULL),
-f3pCorrelatorVsMHist(NULL),
-fDetectorBiasVsMHist(NULL),
-fNonIsotropicTermsVsMHist(NULL)
+fDetectorBiasVsMHist(NULL)
 {
  // Constructor. 
-  // Base list to hold all output objects:
+ // Base list to hold all output objects:
  fHistList = new TList();
  fHistListName = new TString("cobjMH");
  fHistList->SetName(fHistListName->Data());
@@ -114,6 +114,9 @@ fNonIsotropicTermsVsMHist(NULL)
  
  // List to hold objects with final results:      
  fResultsList = new TList();
+
+ // Initialize all arrays:  
+ this->InitializeArrays();
  
 } // AliFlowAnalysisWithMixedHarmonics::AliFlowAnalysisWithMixedHarmonics()
  
@@ -139,14 +142,16 @@ void AliFlowAnalysisWithMixedHarmonics::Init()
  // d) Book common control histograms;
  // e) Book all event-by-event quantities;
  // f) Book all all-event quantities;
- // g) Book and fill histograms to hold phi, pt and eta weights.
+ // g) Book and fill histograms to hold phi, pt and eta weights;
+  
  //save old value and prevent histograms from being added to directory
  //to avoid name clashes in case multiple analaysis objects are used
  //in an analysis
  Bool_t oldHistAddStatus = TH1::AddDirectoryStatus();
  TH1::AddDirectory(kFALSE);
  
+ TH1::SetDefaultSumw2();
  this->CrossCheckSettings();
  this->AccessConstants();
  this->BookAndNestAllLists();
@@ -157,6 +162,7 @@ void AliFlowAnalysisWithMixedHarmonics::Init()
  this->BookAndFillWeightsHistograms();
 
  TH1::AddDirectory(oldHistAddStatus);
 } // end of void AliFlowAnalysisWithMixedHarmonics::Init()
 
 //================================================================================================================
@@ -165,15 +171,18 @@ void AliFlowAnalysisWithMixedHarmonics::Make(AliFlowEventSimple* anEvent)
 {
  // Running over data only in this method.
  
- // a) Define local variables;
- // b) Fill common control histograms;
- // c) Loop over data and calculate e-b-e quantities Q_{n,k} and S_{p,k};
- // d) Calculate 3-p azimuthal correlator and non-isotropic terms in terms of Q_{n,k} and S_{p,k};
- // e) Reset all event-by-event quantities.
+ // a) Check all pointers used in this method;
+ // b) Define local variables;
+ // c) Fill common control histograms;
+ // d) Loop over data and calculate e-b-e quantities Q_{n,k} and S_{p,k};
+ // e) Calculate 3-p azimuthal correlator and non-isotropic terms in terms of Q_{n,k} and S_{p,k};
+ // f) Calculate differential 3-p azimuthal correlator cos(2phi1-psi2-psi3) in terms of Q_{2n} and p_{n}:
+ // g) Reset all event-by-event quantities.
  
+ // a) Check all pointers used in this method:
  this->CheckPointersUsedInMake();
  
- // a) Define local variables:
+ // b) Define local variables:
  Double_t dPhi = 0.; // azimuthal angle in the laboratory frame
  Double_t dPt  = 0.; // transverse momentum
  Double_t dEta = 0.; // pseudorapidity
@@ -182,10 +191,10 @@ void AliFlowAnalysisWithMixedHarmonics::Make(AliFlowEventSimple* anEvent)
  Double_t wEta = 1.; // eta weight
  AliFlowTrackSimple *aftsTrack = NULL; // simple track
  
- // b) Fill common control histograms.
+ // c) Fill common control histograms:
  fCommonHists->FillControlHistograms(anEvent);  
  
- // c) Loop over data and calculate e-b-e quantities:
+ // d) Loop over data and calculate e-b-e quantities:
  Int_t nPrim = anEvent->NumberOfTracks();  // nPrim = total number of primary tracks, i.e. nPrim = nRP + nPOI + rest, where:
                                            // nRP   = # of particles used to determine the reaction plane ("Reference Particles");
                                            // nPOI  = # of particles of interest for a detailed flow analysis ("Particles of Interest");
@@ -232,7 +241,27 @@ void AliFlowAnalysisWithMixedHarmonics::Make(AliFlowEventSimple* anEvent)
       (*fSpk)(p,k)+=pow(wPhi*wPt*wEta,k);
      }
     }    
-   } // end of if(aftsTrack->InRPSelection())  
+   } // end of if(aftsTrack->InRPSelection())
+   // POIs:
+   if(aftsTrack->InPOISelection()) // 1st POI
+   {
+    Double_t dPsi1 = aftsTrack->Phi();
+    Double_t dPt1  = aftsTrack->Pt();
+    for(Int_t j=i+1;j<nPrim;j++)
+    {
+     aftsTrack=anEvent->GetTrack(j);
+     if(aftsTrack->InPOISelection()) // 2nd POI
+     {
+      Double_t dPsi2 = aftsTrack->Phi();
+      Double_t dPt2  = aftsTrack->Pt(); 
+      // Fill:
+      fRePEBE[0]->Fill((dPt1+dPt2)/2.,TMath::Cos(n*(dPsi1+dPsi2)),1.);
+      fImPEBE[0]->Fill((dPt1+dPt2)/2.,TMath::Sin(n*(dPsi1+dPsi2)),1.);
+      fRePEBE[1]->Fill(TMath::Abs(dPt1-dPt2),TMath::Cos(n*(dPsi1+dPsi2)),1.);
+      fImPEBE[1]->Fill(TMath::Abs(dPt1-dPt2),TMath::Sin(n*(dPsi1+dPsi2)),1.);
+     }
+    }
+   }  
   } else // to if(aftsTrack)
     {
      cout<<endl;
@@ -250,7 +279,7 @@ void AliFlowAnalysisWithMixedHarmonics::Make(AliFlowEventSimple* anEvent)
   }  
  } 
  
- // d) Calculate 3-p azimuthal correlator in terms of Q_{n,k} and S_{p,k}:
+ // e) Calculate 3-p azimuthal correlator in terms of Q_{n,k} and S_{p,k}:
  if(anEvent->GetEventNSelTracksRP() >= 3) 
  {
   this->Calculate3pCorrelator();
@@ -260,10 +289,14 @@ void AliFlowAnalysisWithMixedHarmonics::Make(AliFlowEventSimple* anEvent)
   this->CalculateNonIsotropicTerms();                          
  }
                  
- // e) Reset all event-by-event quantities: 
+ // f) Calculate differential 3-p azimuthal correlator cos(2phi1-psi2-psi3) in terms of Q_{2n} and p_{n}:
+ this->CalculateDifferential3pCorrelator(); // to be improved - add relevant if statements for the min No of RPs and POIs
+ // g) Reset all event-by-event quantities: 
  this->ResetEventByEventQuantities();
-  
+   
 } // end of AliFlowAnalysisWithMixedHarmonics::Make(AliFlowEventSimple* anEvent)
+
 //================================================================================================================
 
 void AliFlowAnalysisWithMixedHarmonics::Finish()
@@ -272,32 +305,13 @@ void AliFlowAnalysisWithMixedHarmonics::Finish()
  
  // a) Check all pointers used in this method;
  // b) Access settings for analysis with mixed harmonics;
- // c) Calculate errors for non-isotropic terms;
- // d) Correct for detector effects;
- // e) Quantify bias to 3-p correlator coming from detector effects.
+ // c) Correct for detector effects;
+ // d) Print on the screen the final results.
  
  this->CheckPointersUsedInFinish();
  this->AccessSettings();
- this->FinalizeNonIsotropicTerms();
- this->CorrectForDetectorEffects();
- this->QuantifyBiasFromDetectorEffects();
- cout<<endl;
- cout<<"*************************************"<<endl;
- cout<<"*************************************"<<endl;
- cout<<" flow estimates from mixed harmonics "<<endl; 
- cout<<endl;
- if(fCorrelatorInteger!=1)
- {
-  cout<< "  cos["<<fCorrelatorInteger<<"(2phi1-phi2-phi3)] = "<<f3pCorrelatorPro->GetBinContent(1)<<endl;
- } else
-   {
-    cout<< "  cos(2phi1-phi2-phi3) = "<<f3pCorrelatorPro->GetBinContent(1)<<endl;
-   } 
- cout<<endl;
- cout<<"*************************************"<<endl;
- cout<<"*************************************"<<endl;
+ if(fCorrectForDetectorEffects) this->CorrectForDetectorEffects();
+ if(fPrintOnTheScreen) this->PrintOnTheScreen();
                                                                                                                                                                                                                                                                                                                
 } // end of AliFlowAnalysisWithMixedHarmonics::Finish()
 
@@ -378,7 +392,7 @@ void AliFlowAnalysisWithMixedHarmonics::GetPointersForCommonHistograms()
 
 void AliFlowAnalysisWithMixedHarmonics::GetPointersForAllEventProfiles() 
 {
- // Get pointers to histograms holding final results.
+ // Get pointers to profiles holding final results.
  
  TList *profileList = NULL;
  profileList = dynamic_cast<TList*>(fHistList->FindObject("Profiles"));
@@ -389,16 +403,16 @@ void AliFlowAnalysisWithMixedHarmonics::GetPointersForAllEventProfiles()
  }  
  
  TString s3pCorrelatorProName = "f3pCorrelatorPro";
- TProfile *h3pCorrelatorPro = dynamic_cast<TProfile*>(profileList->FindObject(s3pCorrelatorProName.Data()));
- if(h3pCorrelatorPro)
+ TProfile *p3pCorrelatorPro = dynamic_cast<TProfile*>(profileList->FindObject(s3pCorrelatorProName.Data()));
+ if(p3pCorrelatorPro)
  {
-  this->Set3pCorrelatorPro(h3pCorrelatorPro);  
+  this->Set3pCorrelatorPro(p3pCorrelatorPro);  
  }
  TString s3pCorrelatorVsMProName = "f3pCorrelatorVsMPro";
- TProfile *h3pCorrelatorVsMPro = dynamic_cast<TProfile*>(profileList->FindObject(s3pCorrelatorVsMProName.Data()));
- if(h3pCorrelatorVsMPro)
+ TProfile *p3pCorrelatorVsMPro = dynamic_cast<TProfile*>(profileList->FindObject(s3pCorrelatorVsMProName.Data()));
+ if(p3pCorrelatorVsMPro)
  {
-  this->Set3pCorrelatorVsMPro(h3pCorrelatorVsMPro);  
+  this->Set3pCorrelatorVsMPro(p3pCorrelatorVsMPro);  
  }
  TString nonIsotropicTermsProName = "fNonIsotropicTermsPro";
  TProfile *nonIsotropicTermsPro = dynamic_cast<TProfile*>(profileList->FindObject(nonIsotropicTermsProName.Data()));
@@ -411,7 +425,17 @@ void AliFlowAnalysisWithMixedHarmonics::GetPointersForAllEventProfiles()
  if(nonIsotropicTermsVsMPro)
  {
   this->SetNonIsotropicTermsVsMPro(nonIsotropicTermsVsMPro);  
- }
+ } 
+ TString psdFlag[2] = {"PtSum","PtDiff"}; 
+ for(Int_t sd=0;sd<2;sd++)
+ {
+  TProfile *p3pCorrelatorVsPtSumDiffPro = dynamic_cast<TProfile*>(profileList->FindObject(Form("f3pCorrelatorVs%sPro",psdFlag[sd].Data())));
+  if(p3pCorrelatorVsPtSumDiffPro)
+  {
+   this->Set3pCorrelatorVsPtSumDiffPro(p3pCorrelatorVsPtSumDiffPro,sd);  
+  }
+ }  
+  
 } // end of void AliFlowAnalysisWithMixedHarmonics::GetPointersForAllEventProfiles()
 
 //================================================================================================================
@@ -427,19 +451,12 @@ void AliFlowAnalysisWithMixedHarmonics::GetPointersForResultsHistograms()
   cout<<"WARNING: resultsList is NULL in AFAWMH::GPFRH() !!!!"<<endl;
   exit(0); 
  }  
  TString s3pCorrelatorHistName = "f3pCorrelatorHist";
  TH1D *h3pCorrelatorHist = dynamic_cast<TH1D*>(resultsList->FindObject(s3pCorrelatorHistName.Data()));
  if(h3pCorrelatorHist)
  {
   this->Set3pCorrelatorHist(h3pCorrelatorHist);  
  }
- TString s3pCorrelatorVsMHistName = "f3pCorrelatorVsMHist";
- TH1D *h3pCorrelatorVsMHist = dynamic_cast<TH1D*>(resultsList->FindObject(s3pCorrelatorVsMHistName.Data()));
- if(h3pCorrelatorVsMHist)
- {
-  this->Set3pCorrelatorVsMHist(h3pCorrelatorVsMHist);  
- }
  TString detectorBiasHistName = "fDetectorBiasHist";
  TH1D *detectorBiasHist = dynamic_cast<TH1D*>(resultsList->FindObject(detectorBiasHistName.Data()));
  if(detectorBiasHist)
@@ -452,18 +469,6 @@ void AliFlowAnalysisWithMixedHarmonics::GetPointersForResultsHistograms()
  {
   this->SetDetectorBiasVsMHist(detectorBiasVsMHist);  
  }
- TString nonIsotropicTermsHistName = "fNonIsotropicTermsHist";
- TH1D *nonIsotropicTermsHist = dynamic_cast<TH1D*>(resultsList->FindObject(nonIsotropicTermsHistName.Data()));
- if(nonIsotropicTermsHist)
- {
-  this->SetNonIsotropicTermsHist(nonIsotropicTermsHist);  
- }  
- TString nonIsotropicTermsVsMHistName = "fNonIsotropicTermsVsMHist";
- TH2D *nonIsotropicTermsVsMHist = dynamic_cast<TH2D*>(resultsList->FindObject(nonIsotropicTermsVsMHistName.Data()));
- if(nonIsotropicTermsVsMHist)
- {
-  this->SetNonIsotropicTermsVsMHist(nonIsotropicTermsVsMHist);  
- }
 
 } // end of void AliFlowAnalysisWithMixedHarmonics::GetPointersForResultsHistograms()
 
@@ -490,6 +495,21 @@ void AliFlowAnalysisWithMixedHarmonics::WriteHistograms(TDirectoryFile *outputFi
 
 //================================================================================================================
 
+void AliFlowAnalysisWithMixedHarmonics::InitializeArrays()
+{
+ // Initialize arrays.
+ for(Int_t sd=0;sd<2;sd++)
+ {
+  fRePEBE[sd] = NULL;
+  fImPEBE[sd] = NULL;
+  f3pCorrelatorVsPtSumDiffPro[sd] = NULL;
+ }
+  
+} // end of AliFlowAnalysisWithMixedHarmonics::InitializeArrays()
+
+//================================================================================================================  
+
 void AliFlowAnalysisWithMixedHarmonics::BookAndNestAllLists()
 {
  // Book and nest all list in base list fHistList.
@@ -516,7 +536,8 @@ void AliFlowAnalysisWithMixedHarmonics::BookProfileHoldingSettings()
  // Book profile to hold all analysis settings.
 
  TString analysisSettingsName = "fAnalysisSettings";
- fAnalysisSettings = new TProfile(analysisSettingsName.Data(),"Settings for analysis with mixed harmonics",5,0,5);
+ fAnalysisSettings = new TProfile(analysisSettingsName.Data(),"Settings for analysis with mixed harmonics",6,0,6);
+ fAnalysisSettings->GetXaxis()->SetLabelSize(0.025);
  fAnalysisSettings->GetXaxis()->SetBinLabel(1,"Integer n in cos(n(2#phi_{1}-#phi_{2}-#phi_{3}))");
  fAnalysisSettings->Fill(0.5,fCorrelatorInteger);
  fAnalysisSettings->GetXaxis()->SetBinLabel(2,"Corrected for detector effects?");
@@ -527,6 +548,8 @@ void AliFlowAnalysisWithMixedHarmonics::BookProfileHoldingSettings()
  fAnalysisSettings->Fill(3.5,fMultipicityBinWidth);  
  fAnalysisSettings->GetXaxis()->SetBinLabel(5,"Minimal multiplicity");
  fAnalysisSettings->Fill(4.5,fMinMultiplicity);
+ fAnalysisSettings->GetXaxis()->SetBinLabel(6,"Print on the screen?");
+ fAnalysisSettings->Fill(5.5,(Int_t)fPrintOnTheScreen);
  fHistList->Add(fAnalysisSettings);
  
 } // end of void AliFlowAnalysisWithMixedHarmonics::BookProfileHoldingSettings()
@@ -559,8 +582,18 @@ void AliFlowAnalysisWithMixedHarmonics::BookAllEventByEventQuantities()
  // Correction terms for non-uniform acceptance to <cos[n(2phi1-phi2-phi3)]> for single event:
  TString nonIsotropicTermsEBEName = "fNonIsotropicTermsEBE";
  fNonIsotropicTermsEBE = new TH1D(nonIsotropicTermsEBEName.Data(),"Non-isotropic terms for single event",10,0,10);
-
+ // p_n vs [(p1+p2)/2,|p1-p2|]
+ TString psdFlag[2] = {"PtSum","PtDiff"};
+ for(Int_t sd=0;sd<2;sd++)
+ {
+  // to be improved: hardwired ,fnBinsPt,0.,fPtMax):
+  fRePEBE[sd] = new TProfile(Form("fRePEBE%s",psdFlag[sd].Data()),"",fnBinsPt,0.,fPtMax);
+  fImPEBE[sd] = new TProfile(Form("fImPEBE%s",psdFlag[sd].Data()),"",fnBinsPt,0.,fPtMax);
+ }  
 } // end fo void AliFlowAnalysisWithMixedHarmonics::BookAllEventByEventQuantities()
+
 //================================================================================================================
 
 void AliFlowAnalysisWithMixedHarmonics::BookAllAllEventQuantities()
@@ -568,118 +601,152 @@ void AliFlowAnalysisWithMixedHarmonics::BookAllAllEventQuantities()
  // Book all all-event quantitites.
  
  // a) Book quantites without multiplicity binning;
- // b) Book quantites with multiplicity binning.
+ // b) Book quantites with multiplicity binning;
+ // c) Book quantities with binning in (p1+p2)/2, |p1-p2|, (eta1+eta2)/2 or |eta1-eta2|.
+  
  // a) Book quantites without multiplicity binning:
- // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> for all events (with wrong error):
+ // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> for all events:
  TString s3pCorrelatorProName = "f3pCorrelatorPro";
- f3pCorrelatorPro = new TProfile(s3pCorrelatorProName.Data(),"<<cos[n(2#phi_{1}-#phi_{2}-#phi_{3})]>> for all events (errors are wrong here!)",1,0,1,"s");
+ f3pCorrelatorPro = new TProfile(s3pCorrelatorProName.Data(),"",1,0,1);
+ f3pCorrelatorPro->GetXaxis()->SetLabelOffset(0.01);
+ f3pCorrelatorPro->GetXaxis()->SetLabelSize(0.05);
+ if(fCorrelatorInteger == 1)
+ {
+  f3pCorrelatorPro->GetXaxis()->SetBinLabel(1,"cos(2#phi_{1}-#phi_{2}-#phi_{3})");
+ } else
+   {
+    f3pCorrelatorPro->GetXaxis()->SetBinLabel(1,Form("cos[%d(2#phi_{1}-#phi_{2}-#phi_{3})]",fCorrelatorInteger)); 
+   }
  fProfileList->Add(f3pCorrelatorPro);
- // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> for all events (with correct error):
- TString s3pCorrelatorHistName = "f3pCorrelatorHist";
- f3pCorrelatorHist = new TH1D(s3pCorrelatorHistName.Data(),"<<cos[n(2#phi_{1}-#phi_{2}-#phi_{3})]>> for all events",1,0,1);
- fResultsList->Add(f3pCorrelatorHist);
- // Correction terms for non-uniform acceptance to <cos[n(2phi1-phi2-phi3)]> for all events (with wrong errors):
+ // Non-isotropic terms in the decomposition of <cos[n(2phi1-phi2-phi3)]:
  TString nonIsotropicTermsProName = "fNonIsotropicTermsPro";
- fNonIsotropicTermsPro = new TProfile(nonIsotropicTermsProName.Data(),"Non-isotropic terms for all events (errors are wrong here!)",10,0,10,"s");
- fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(1,"cos(n(#phi_{1}))");
- fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(2,"sin(n(#phi_{1}))");
- fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(3,"cos(2n(#phi_{1}))");
- fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(4,"sin(2n(#phi_{1}))");
- fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(5,"cos(n(#phi_{1}+#phi_{2}))");
- fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(6,"sin(n(#phi_{1}+#phi_{2}))");
- fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(7,"cos(n(2#phi_{1}-#phi_{2}))");
- fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(8,"sin(n(2#phi_{1}-#phi_{2}))");
- fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(9,"cos(n(#phi_{1}-#phi_{2}-#phi_{3}))");
- fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(10,"sin(n(#phi_{1}-#phi_{2}-#phi_{3})");  
+ fNonIsotropicTermsPro = new TProfile(nonIsotropicTermsProName.Data(),"",10,0,10);
+ if(fCorrelatorInteger == 1)
+ {
+  fNonIsotropicTermsPro->SetTitle("Non-isotropic terms in decomposition of cos(2#phi_{1}-#phi_{2}-#phi_{3})");
+  fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(1,"cos(#phi_{1})");
+  fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(2,"sin(#phi_{1})");
+  fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(3,"cos(2#phi_{1})");
+  fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(4,"sin(2#phi_{1})");
+  fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(5,"cos(#phi_{1}+#phi_{2})");
+  fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(6,"sin(#phi_{1}+#phi_{2})");
+  fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(7,"cos(2#phi_{1}-#phi_{2})");
+  fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(8,"sin(2#phi_{1}-#phi_{2})");
+  fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(9,"cos(#phi_{1}-#phi_{2}-#phi_{3})");
+  fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(10,"sin(#phi_{1}-#phi_{2}-#phi_{3})");  
+ } else
+   {
+    fNonIsotropicTermsPro->SetTitle(Form("Non-isotropic terms in decomposition of cos[%d(2#phi_{1}-#phi_{2}-#phi_{3})]",fCorrelatorInteger));
+    fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(1,Form("cos(%d#phi_{1})",fCorrelatorInteger));
+    fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(2,Form("sin(%d#phi_{1})",fCorrelatorInteger));
+    fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(3,Form("cos(%d#phi_{1})",2*fCorrelatorInteger));
+    fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(4,Form("sin(%d#phi_{1})",2*fCorrelatorInteger));
+    fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(5,Form("cos(%d(#phi_{1}+#phi_{2}))",fCorrelatorInteger));
+    fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(6,Form("sin(%d(#phi_{1}+#phi_{2}))",fCorrelatorInteger));
+    fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(7,Form("cos(%d(2#phi_{1}-#phi_{2}))",fCorrelatorInteger));
+    fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(8,Form("sin(%d(2#phi_{1}-#phi_{2}))",fCorrelatorInteger));
+    fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(9,Form("cos(%d(#phi_{1}-#phi_{2}-#phi_{3}))",fCorrelatorInteger));
+    fNonIsotropicTermsPro->GetXaxis()->SetBinLabel(10,Form("sin(%d(#phi_{1}-#phi_{2}-#phi_{3}))",fCorrelatorInteger));  
+   } 
  fProfileList->Add(fNonIsotropicTermsPro);
- // Correction terms for non-uniform acceptance to <cos[n(2phi1-phi2-phi3)]> for all events (with correct errors):
- TString nonIsotropicTermsHistName = "fNonIsotropicTermsHist";
- fNonIsotropicTermsHist = new TH1D(nonIsotropicTermsHistName.Data(),"Non-isotropic terms for all events",10,0,10);
- fNonIsotropicTermsHist->GetXaxis()->SetBinLabel(1,"cos(n(#phi_{1}))");
- fNonIsotropicTermsHist->GetXaxis()->SetBinLabel(2,"sin(n(#phi_{1}))");
- fNonIsotropicTermsHist->GetXaxis()->SetBinLabel(3,"cos(2n(#phi_{1}))");
- fNonIsotropicTermsHist->GetXaxis()->SetBinLabel(4,"sin(2n(#phi_{1}))");
- fNonIsotropicTermsHist->GetXaxis()->SetBinLabel(5,"cos(n(#phi_{1}+#phi_{2}))");
- fNonIsotropicTermsHist->GetXaxis()->SetBinLabel(6,"sin(n(#phi_{1}+#phi_{2}))");
- fNonIsotropicTermsHist->GetXaxis()->SetBinLabel(7,"cos(n(2#phi_{1}-#phi_{2}))");
- fNonIsotropicTermsHist->GetXaxis()->SetBinLabel(8,"sin(n(2#phi_{1}-#phi_{2}))");
- fNonIsotropicTermsHist->GetXaxis()->SetBinLabel(9,"cos(n(#phi_{1}-#phi_{2}-#phi_{3}))");
- fNonIsotropicTermsHist->GetXaxis()->SetBinLabel(10,"sin(n(#phi_{1}-#phi_{2}-#phi_{3})");  
- fResultsList->Add(fNonIsotropicTermsHist);
- // Quantified bias comming from detector inefficiencies to 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> (in %)
+ // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> corrected for detector effects:
+ TString s3pCorrelatorHistName = "f3pCorrelatorHist";
+ f3pCorrelatorHist = new TH1D(s3pCorrelatorHistName.Data(),"",1,0,1);
+ f3pCorrelatorHist->GetXaxis()->SetLabelOffset(0.01);
+ f3pCorrelatorHist->GetXaxis()->SetLabelSize(0.05);
+ if(fCorrelatorInteger == 1)
+ {
+  f3pCorrelatorHist->GetXaxis()->SetBinLabel(1,"cos(2#phi_{1}-#phi_{2}-#phi_{3})");
+ } else
+   {
+    f3pCorrelatorHist->GetXaxis()->SetBinLabel(1,Form("cos[%d(2#phi_{1}-#phi_{2}-#phi_{3})]",fCorrelatorInteger)); 
+   }
+ fResultsList->Add(f3pCorrelatorHist);
+ // Quantified bias comming from detector inefficiencies to 3-p correlator <<cos[n(2phi1-phi2-phi3)]>>:
  TString detectorBiasHistName = "fDetectorBiasHist";
  fDetectorBiasHist = new TH1D(detectorBiasHistName.Data(),"Bias coming from detector inefficiences",1,0,1);
- fDetectorBiasHist->GetXaxis()->SetBinLabel(1,"#left|1 - #frac{corrected}{not-corrected}#right| (in %)");
+ fDetectorBiasHist->GetXaxis()->SetBinLabel(1,"#frac{corrected}{measured}");
  fResultsList->Add(fDetectorBiasHist);
  
  // b) Book quantites with multiplicity binning.
- // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> for all events versus multiplicity (with wrong error):
+ // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> versus multiplicity:
  TString s3pCorrelatorVsMProName = "f3pCorrelatorVsMPro";
- f3pCorrelatorVsMPro = new TProfile(s3pCorrelatorVsMProName.Data(),"<<cos[n(2#phi_{1}-#phi_{2}-#phi_{3})]>> vs M for all events (errors are wrong here!)",fNoOfMultipicityBins+1,0,fNoOfMultipicityBins+1,"s");
- for(Int_t b=1;b<=fNoOfMultipicityBins;b++)
+ f3pCorrelatorVsMPro = new TProfile(s3pCorrelatorVsMProName.Data(),"",fNoOfMultipicityBins+2,0,fNoOfMultipicityBins+2);
+ if(fCorrelatorInteger == 1)
  {
-  f3pCorrelatorVsMPro->GetXaxis()->SetBinLabel(b,Form("%d #leq M < %d",(Int_t)(fMinMultiplicity+(b-1)*fMultipicityBinWidth),(Int_t)(fMinMultiplicity+b*fMultipicityBinWidth)));
- }
- f3pCorrelatorVsMPro->GetXaxis()->SetBinLabel(fNoOfMultipicityBins+1,Form(" M #geq %d",(Int_t)(fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)));
- fProfileList->Add(f3pCorrelatorVsMPro);
- // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> for all events (with correct error):
- TString s3pCorrelatorVsMHistName = "f3pCorrelatorVsMHist";
- f3pCorrelatorVsMHist = new TH1D(s3pCorrelatorVsMHistName.Data(),"<<cos[n(2#phi_{1}-#phi_{2}-#phi_{3})]>> vs M for all events",fNoOfMultipicityBins+1,0,fNoOfMultipicityBins+1);
- for(Int_t b=1;b<=fNoOfMultipicityBins;b++)
+  f3pCorrelatorVsMPro->SetTitle("cos(2#phi_{1}-#phi_{2}-#phi_{3}) #font[72]{vs} M");
+ } else
+   {
+    f3pCorrelatorVsMPro->SetTitle(Form("cos[%d(2#phi_{1}-#phi_{2}-#phi_{3})] #font[72]{vs} M",fCorrelatorInteger)); 
+   }
+ f3pCorrelatorVsMPro->GetXaxis()->SetBinLabel(1,Form("M < %d",(Int_t)fMinMultiplicity));
+ for(Int_t b=2;b<=fNoOfMultipicityBins+1;b++)
  {
-  f3pCorrelatorVsMHist->GetXaxis()->SetBinLabel(b,Form("%d #leq M < %d",(Int_t)(fMinMultiplicity+(b-1)*fMultipicityBinWidth),(Int_t)(fMinMultiplicity+b*fMultipicityBinWidth)));
+  f3pCorrelatorVsMPro->GetXaxis()->SetBinLabel(b,Form("%d #leq M < %d",(Int_t)(fMinMultiplicity+(b-2)*fMultipicityBinWidth),(Int_t)(fMinMultiplicity+(b-1)*fMultipicityBinWidth)));
  }
- f3pCorrelatorVsMHist->GetXaxis()->SetBinLabel(fNoOfMultipicityBins+1,Form(" M #geq %d",(Int_t)(fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)));
- fResultsList->Add(f3pCorrelatorVsMHist);
- // Correction terms for non-uniform acceptance to <cos[n(2phi1-phi2-phi3)]> for all events vs M (with wrong errors):
+ f3pCorrelatorVsMPro->GetXaxis()->SetBinLabel(fNoOfMultipicityBins+2,Form(" M #geq %d",(Int_t)(fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)));
+ fProfileList->Add(f3pCorrelatorVsMPro);
+ // Non-isotropic terms in the decomposition of <cos[n(2phi1-phi2-phi3)] vs multiplicity:
  TString nonIsotropicTermsVsMProName = "fNonIsotropicTermsVsMPro";
- fNonIsotropicTermsVsMPro = new TProfile2D(nonIsotropicTermsVsMProName.Data(),"Non-isotropic terms for all events vs M (errors are wrong here!)",10,0,10,fNoOfMultipicityBins+1,0,fNoOfMultipicityBins+1,"s");
- fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(1,"cos(n(#phi_{1}))");
- fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(2,"sin(n(#phi_{1}))");
- fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(3,"cos(2n(#phi_{1}))");
- fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(4,"sin(2n(#phi_{1}))");
- fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(5,"cos(n(#phi_{1}+#phi_{2}))");
- fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(6,"sin(n(#phi_{1}+#phi_{2}))");
- fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(7,"cos(n(2#phi_{1}-#phi_{2}))");
- fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(8,"sin(n(2#phi_{1}-#phi_{2}))");
- fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(9,"cos(n(#phi_{1}-#phi_{2}-#phi_{3}))");
- fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(10,"sin(n(#phi_{1}-#phi_{2}-#phi_{3})");  
- for(Int_t b=1;b<=fNoOfMultipicityBins;b++)
+ 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);
+ if(fCorrelatorInteger == 1)
  {
-  fNonIsotropicTermsVsMPro->GetYaxis()->SetBinLabel(b,Form("%d #leq M < %d",(Int_t)(fMinMultiplicity+(b-1)*fMultipicityBinWidth),(Int_t)(fMinMultiplicity+b*fMultipicityBinWidth)));
- }
- fNonIsotropicTermsVsMPro->GetYaxis()->SetBinLabel(fNoOfMultipicityBins+1,Form(" M #geq %d",(Int_t)(fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)));
- fProfileList->Add(fNonIsotropicTermsVsMPro);
- // Correction terms for non-uniform acceptance to <cos[n(2phi1-phi2-phi3)]> for all events (with correct errors):
- TString nonIsotropicTermsVsMHistName = "fNonIsotropicTermsVsMHist";
- fNonIsotropicTermsVsMHist = new TH2D(nonIsotropicTermsVsMHistName.Data(),"Non-isotropic terms for all events vs M ",10,0,10,fNoOfMultipicityBins+1,0,fNoOfMultipicityBins+1);
- fNonIsotropicTermsVsMHist->GetXaxis()->SetBinLabel(1,"cos(n(#phi_{1}))");
- fNonIsotropicTermsVsMHist->GetXaxis()->SetBinLabel(2,"sin(n(#phi_{1}))");
- fNonIsotropicTermsVsMHist->GetXaxis()->SetBinLabel(3,"cos(2n(#phi_{1}))");
- fNonIsotropicTermsVsMHist->GetXaxis()->SetBinLabel(4,"sin(2n(#phi_{1}))");
- fNonIsotropicTermsVsMHist->GetXaxis()->SetBinLabel(5,"cos(n(#phi_{1}+#phi_{2}))");
- fNonIsotropicTermsVsMHist->GetXaxis()->SetBinLabel(6,"sin(n(#phi_{1}+#phi_{2}))");
- fNonIsotropicTermsVsMHist->GetXaxis()->SetBinLabel(7,"cos(n(2#phi_{1}-#phi_{2}))");
- fNonIsotropicTermsVsMHist->GetXaxis()->SetBinLabel(8,"sin(n(2#phi_{1}-#phi_{2}))");
- fNonIsotropicTermsVsMHist->GetXaxis()->SetBinLabel(9,"cos(n(#phi_{1}-#phi_{2}-#phi_{3}))");
- fNonIsotropicTermsVsMHist->GetXaxis()->SetBinLabel(10,"sin(n(#phi_{1}-#phi_{2}-#phi_{3})");  
- for(Int_t b=1;b<=fNoOfMultipicityBins;b++)
+  fNonIsotropicTermsVsMPro->SetTitle("Non-isotropic terms in decomposition of cos(2#phi_{1}-#phi_{2}-#phi_{3}) #font[72]{vs} M");
+  fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(1,"cos(#phi_{1})");
+  fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(2,"sin(#phi_{1})");
+  fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(3,"cos(2#phi_{1})");
+  fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(4,"sin(2#phi_{1})");
+  fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(5,"cos(#phi_{1}+#phi_{2})");
+  fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(6,"sin(#phi_{1}+#phi_{2})");
+  fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(7,"cos(2#phi_{1}-#phi_{2})");
+  fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(8,"sin(2#phi_{1}-#phi_{2})");
+  fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(9,"cos(#phi_{1}-#phi_{2}-#phi_{3})");
+  fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(10,"sin(#phi_{1}-#phi_{2}-#phi_{3})");  
+ } else
+   {
+    fNonIsotropicTermsVsMPro->SetTitle(Form("Non-isotropic terms in decomposition of cos[%d(2#phi_{1}-#phi_{2}-#phi_{3})]",fCorrelatorInteger));
+    fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(1,Form("cos(%d#phi_{1})",fCorrelatorInteger));
+    fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(2,Form("sin(%d#phi_{1})",fCorrelatorInteger));
+    fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(3,Form("cos(%d#phi_{1})",2*fCorrelatorInteger));
+    fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(4,Form("sin(%d#phi_{1})",2*fCorrelatorInteger));
+    fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(5,Form("cos(%d(#phi_{1}+#phi_{2}))",fCorrelatorInteger));
+    fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(6,Form("sin(%d(#phi_{1}+#phi_{2}))",fCorrelatorInteger));
+    fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(7,Form("cos(%d(2#phi_{1}-#phi_{2}))",fCorrelatorInteger));
+    fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(8,Form("sin(%d(2#phi_{1}-#phi_{2}))",fCorrelatorInteger));
+    fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(9,Form("cos(%d(#phi_{1}-#phi_{2}-#phi_{3}))",fCorrelatorInteger));
+    fNonIsotropicTermsVsMPro->GetXaxis()->SetBinLabel(10,Form("sin(%d(#phi_{1}-#phi_{2}-#phi_{3}))",fCorrelatorInteger));  
+   } 
+ fNonIsotropicTermsVsMPro->GetYaxis()->SetBinLabel(1,Form("M < %d",(Int_t)fMinMultiplicity));
+ for(Int_t b=2;b<=fNoOfMultipicityBins+1;b++)
  {
-  fNonIsotropicTermsVsMHist->GetYaxis()->SetBinLabel(b,Form("%d #leq M < %d",(Int_t)(fMinMultiplicity+(b-1)*fMultipicityBinWidth),(Int_t)(fMinMultiplicity+b*fMultipicityBinWidth)));
+ fNonIsotropicTermsVsMPro->GetYaxis()->SetBinLabel(b,Form("%d #leq M < %d",(Int_t)(fMinMultiplicity+(b-2)*fMultipicityBinWidth),(Int_t)(fMinMultiplicity+(b-1)*fMultipicityBinWidth)));
  }
- fNonIsotropicTermsVsMHist->GetYaxis()->SetBinLabel(fNoOfMultipicityBins+1,Form(" M #geq %d",(Int_t)(fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)));
- fResultsList->Add(fNonIsotropicTermsVsMHist);
+ fNonIsotropicTermsVsMPro->GetYaxis()->SetBinLabel(fNoOfMultipicityBins+2,Form(" M #geq %d",(Int_t)(fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)));
+ fProfileList->Add(fNonIsotropicTermsVsMPro); 
  // Quantified bias comming from detector inefficiencies to 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> (in %)
  TString detectorBiasVsMHistName = "fDetectorBiasVsMHist";
- fDetectorBiasVsMHist = new TH1D(detectorBiasVsMHistName.Data(),"#left|1 - #frac{corrected}{not-corrected}#right| (in %)",fNoOfMultipicityBins+1,0,fNoOfMultipicityBins+1);
+ fDetectorBiasVsMHist = new TH1D(detectorBiasVsMHistName.Data(),"",fNoOfMultipicityBins+1,0,fNoOfMultipicityBins+1);
  for(Int_t b=1;b<=fNoOfMultipicityBins;b++)
  {
-  fDetectorBiasVsMHist->GetXaxis()->SetBinLabel(b,Form("%d #leq M < %d",(Int_t)(fMinMultiplicity+(b-1)*fMultipicityBinWidth),(Int_t)(fMinMultiplicity+b*fMultipicityBinWidth)));
+  fDetectorBiasVsMHist->GetXaxis()->SetBinLabel(b,Form("%d #leq M < %d",(Int_t)(fMinMultiplicity+(b-1)*fMultipicityBinWidth),(Int_t)(fMinMultiplicity+b*fMultipicityBinWidth))); // to be imroved - wrong labeling - see other vs M profiles for correct version
  }
  fDetectorBiasVsMHist->GetXaxis()->SetBinLabel(fNoOfMultipicityBins+1,Form(" M #geq %d",(Int_t)(fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)));
- fResultsList->Add(fDetectorBiasVsMHist);
+ //fResultsList->Add(fDetectorBiasVsMHist); // to be improved - calculated and added eventually
  
+ // c) Book quantities with binning in (p1+p2)/2, |p1-p2|, (eta1+eta2)/2 or |eta1-eta2|: 
+ TString psdFlag[2] = {"PtSum","PtDiff"};
+ TString psdTitleFlag[2] = {"(p_{t,1}+p_{t,2})/2","#left|p_{t,1}-p_{t,2}#right|"};
+ //TString s3pCorrelatorVsPtSumDiffProName = "f3pCorrelatorVsPtSumDiffPro";
+ for(Int_t sd=0;sd<2;sd++)
+ {
+  // to be improved: hardwired ,fnBinsPt,0.,fPtMax):
+  f3pCorrelatorVsPtSumDiffPro[sd] = new TProfile(Form("f3pCorrelatorVs%sPro",psdFlag[sd].Data()),"",fnBinsPt,0.,fPtMax);
+  //f3pCorrelatorVsPtSumDiffPro[sd]->SetLabelSize(0.05);
+  //f3pCorrelatorVsPtSumDiffPro[sd]->SetMarkerStyle(25);
+  f3pCorrelatorVsPtSumDiffPro[sd]->GetXaxis()->SetTitle(psdTitleFlag[sd].Data());
+  fProfileList->Add(f3pCorrelatorVsPtSumDiffPro[sd]);
+ }  
+  
 } // end of void AliFlowAnalysisWithMixedHarmonics::BookAllAllEventQuantities()
 
 //================================================================================================================
@@ -802,17 +869,48 @@ void AliFlowAnalysisWithMixedHarmonics::CheckPointersUsedInMake()
  if(!fReQnk || !fImQnk || !fSpk || !f3pCorrelatorEBE)
  {                        
   cout<<endl;
-  cout<<" WARNING (MH): (!fReQnk || !fImQnk || !fSpk || !f3pCorrelatorEBE) is NULL !!!!"<<endl;
+  cout<<" WARNING (MH::Make()): (!fReQnk || !fImQnk || !fSpk || !f3pCorrelatorEBE) is NULL !!!!"<<endl;
   cout<<endl;
   exit(0);
  }
  if(!f3pCorrelatorPro)
  {
   cout<<endl;
-  cout<<" WARNING (MH): (!f3pCorrelatorPro) is NULL !!!!"<<endl;
+  cout<<" WARNING (MH::Make()): (!f3pCorrelatorPro) is NULL !!!!"<<endl;
   cout<<endl;
   exit(0); 
  }
+ if(!fNonIsotropicTermsPro)
+ {                        
+  cout<<endl;
+  cout<<" WARNING (MH::Make()): !fNonIsotropicTermsPro is NULL !!!!"<<endl;
+  cout<<endl;
+  exit(0);
+ } 
+ if(!f3pCorrelatorVsMPro)
+ {                        
+  cout<<endl;
+  cout<<" WARNING (MH::Make()): !f3pCorrelatorVsMPro is NULL !!!!"<<endl;
+  cout<<endl;
+  exit(0);
+ }
+ if(!fNonIsotropicTermsVsMPro)
+ {                        
+  cout<<endl;
+  cout<<" WARNING (MH::Make()): !fNonIsotropicTermsVsMPro is NULL !!!!"<<endl;
+  cout<<endl;
+  exit(0);
+ }
+ for(Int_t sd=0;sd<2;sd++)
+ {
+  if(!(f3pCorrelatorVsPtSumDiffPro[sd]))
+  {
+   cout<<endl;
+   cout<<" WARNING (MH::Make()): !"<<Form("f3pCorrelatorVsPtSumDiffPro[%d]",sd)<<" is NULL !!!!"<<endl;
+   cout<<endl;
+   exit(0);   
+  } 
+ }
                                                                                                                                                                                                                                                                                                                                    
 } // end of AliFlowAnalysisWithMixedHarmonics::CheckPointersUsedInMake()
 
@@ -825,85 +923,102 @@ void AliFlowAnalysisWithMixedHarmonics::CheckPointersUsedInFinish()
  if(!fAnalysisSettings)
  {                        
   cout<<endl;
-  cout<<" WARNING (MH): !fAnalysisSettings is NULL !!!!"<<endl;
+  cout<<" WARNING (MH::Finish()): !fAnalysisSettings is NULL !!!!"<<endl;
   cout<<endl;
   exit(0);
  }
  if(!f3pCorrelatorPro)
  {                        
   cout<<endl;
-  cout<<" WARNING (MH): !f3pCorrelatorPro is NULL !!!!"<<endl;
+  cout<<" WARNING (MH::Finish()): !f3pCorrelatorPro is NULL !!!!"<<endl;
   cout<<endl;
   exit(0);
  }
  if(!fNonIsotropicTermsPro)
  {                        
   cout<<endl;
-  cout<<" WARNING (MH): !fNonIsotropicTermsPro is NULL !!!!"<<endl;
+  cout<<" WARNING (MH::Finish()): !fNonIsotropicTermsPro is NULL !!!!"<<endl;
   cout<<endl;
   exit(0);
  } 
  if(!f3pCorrelatorVsMPro)
  {                        
   cout<<endl;
-  cout<<" WARNING (MH): !f3pCorrelatorVsMPro is NULL !!!!"<<endl;
+  cout<<" WARNING (MH::Finish()): !f3pCorrelatorVsMPro is NULL !!!!"<<endl;
   cout<<endl;
   exit(0);
  }
  if(!fNonIsotropicTermsVsMPro)
  {                        
   cout<<endl;
-  cout<<" WARNING (MH): !fNonIsotropicTermsVsMPro is NULL !!!!"<<endl;
+  cout<<" WARNING (MH::Finish()): !fNonIsotropicTermsVsMPro is NULL !!!!"<<endl;
   cout<<endl;
   exit(0);
  }
  if(!f3pCorrelatorHist)
  {                        
   cout<<endl;
-  cout<<" WARNING (MH): !f3pCorrelatorHist is NULL !!!!"<<endl;
+  cout<<" WARNING (MH::Finish()): !f3pCorrelatorHist is NULL !!!!"<<endl;
   cout<<endl;
   exit(0);
- }                                                                                                                                                                                                                                                                                                                                   
- if(!fNonIsotropicTermsHist)
- {                        
-  cout<<endl;
-  cout<<" WARNING (MH): !fNonIsotropicTermsHist is NULL !!!!"<<endl;
-  cout<<endl;
-  exit(0);
- }
+ }   
  if(!fDetectorBiasHist)
  {                        
   cout<<endl;
-  cout<<" WARNING (MH): !fDetectorBiasHist is NULL !!!!"<<endl;
+  cout<<" WARNING (MH::Finish()): !fDetectorBiasHist is NULL !!!!"<<endl;
   cout<<endl;
   exit(0);
  }   
- if(!f3pCorrelatorVsMHist)
- {                        
-  cout<<endl;
-  cout<<" WARNING (MH): !f3pCorrelatorVsMHist is NULL !!!!"<<endl;
-  cout<<endl;
-  exit(0);
- }                                                                                                                                                                                                                                                                                                                                   
- if(!fNonIsotropicTermsVsMHist)
- {                        
-  cout<<endl;
-  cout<<" WARNING (MH): !fNonIsotropicTermsVsMHist is NULL !!!!"<<endl;
-  cout<<endl;
-  exit(0);
- }
+ /* to be improved - enabled eventually
  if(!fDetectorBiasVsMHist)
  {                        
   cout<<endl;
-  cout<<" WARNING (MH): !fDetectorBiasVsMHist is NULL !!!!"<<endl;
+  cout<<" WARNING (MH::Finish()): !fDetectorBiasVsMHist is NULL !!!!"<<endl;
   cout<<endl;
   exit(0);
- }   
+ } 
+ */  
+ for(Int_t sd=0;sd<2;sd++)
+ {
+  if(!(f3pCorrelatorVsPtSumDiffPro[sd]))
+  {
+   cout<<endl;
+   cout<<" WARNING (MH::Finish()): !"<<Form("f3pCorrelatorVsPtSumDiffPro[%d]",sd)<<" is NULL !!!!"<<endl;
+   cout<<endl;
+   exit(0);   
+  } 
+ }
 
 } // end of AliFlowAnalysisWithMixedHarmonics::CheckPointersUsedInFinish()
 
 //================================================================================================================
 
+void AliFlowAnalysisWithMixedHarmonics::PrintOnTheScreen()
+{
+ // Print the final results on the screen.
+ cout<<endl;
+ cout<<"*******************************************************"<<endl;
+ cout<<"*******************************************************"<<endl;
+ cout<<"                    Mixed Harmonics                      "<<endl; 
+ cout<<endl;
+ if(fCorrelatorInteger!=1)
+ {
+  cout<<"  cos["<<fCorrelatorInteger<<"(2phi1-phi2-phi3)] = "<<f3pCorrelatorHist->GetBinContent(1)<<
+  " +/- "<<f3pCorrelatorHist->GetBinError(1)<<endl;
+ } else
+   {
+    cout<<"  cos(2phi1-phi2-phi3) = "<<f3pCorrelatorHist->GetBinContent(1)<<" +/- "<<f3pCorrelatorHist->GetBinError(1)<<endl;
+   }
+ cout<<"  Detector Bias = "<<fDetectorBiasHist->GetBinContent(1)<<endl;
+ cout<<endl;
+ cout<<"*******************************************************"<<endl;
+ cout<<"*******************************************************"<<endl;
+
+} // end of void AliFlowAnalysisWithMixedHarmonics::PrintOnTheScreen()
+
+//================================================================================================================
+
 void AliFlowAnalysisWithMixedHarmonics::AccessSettings()
 {
  // Access the settings for analysis with mixed harmonics.
@@ -913,6 +1028,7 @@ void AliFlowAnalysisWithMixedHarmonics::AccessSettings()
  fNoOfMultipicityBins = (Int_t)fAnalysisSettings->GetBinContent(3);
  fMultipicityBinWidth = (Double_t)fAnalysisSettings->GetBinContent(4);
  fMinMultiplicity = (Double_t)fAnalysisSettings->GetBinContent(5);
+ fPrintOnTheScreen = (Bool_t)fAnalysisSettings->GetBinContent(6);
                                                                                                                                                                                                                                                                                                                                    
 } // end of AliFlowAnalysisWithMixedHarmonics::AccessSettings()
 
@@ -921,7 +1037,7 @@ void AliFlowAnalysisWithMixedHarmonics::AccessSettings()
 void AliFlowAnalysisWithMixedHarmonics::CorrectForDetectorEffects()
 {
  // Correct measured 3-p correlator cos[n(2phi1-phi2-phi3)] for detector effects.
+  
  Double_t measured3pCorrelator = f3pCorrelatorPro->GetBinContent(1); // biased by detector effects
  Double_t corrected3pCorrelator = 0.; // corrected for detector effects
  Double_t nonIsotropicTerms[10] = {0.}; // there are 10 distinct non-isotropic terms
@@ -936,44 +1052,25 @@ void AliFlowAnalysisWithMixedHarmonics::CorrectForDetectorEffects()
                        - 2.*nonIsotropicTerms[0]*nonIsotropicTerms[6]                                       
                        - 2.*nonIsotropicTerms[1]*nonIsotropicTerms[7]                                       
                        + 2.*nonIsotropicTerms[2]*(pow(nonIsotropicTerms[0],2.)-pow(nonIsotropicTerms[1],2.))                                       
-                       + 4.*nonIsotropicTerms[3]*nonIsotropicTerms[0]*nonIsotropicTerms[1];                                       
- // Store corrected 3-p correlator:
- f3pCorrelatorHist->SetBinContent(1,corrected3pCorrelator);                 
+                       + 4.*nonIsotropicTerms[3]*nonIsotropicTerms[0]*nonIsotropicTerms[1]; 
+ // Store corrected correlator:
+ f3pCorrelatorHist->SetBinContent(1,corrected3pCorrelator);
+ f3pCorrelatorHist->SetBinError(1,f3pCorrelatorPro->GetBinError(1)); // to be improved (propagate error for non-isotropic terms)
+ // Quantify bias from detector inefficiences to 3-p correlator. Remark: Bias is quantified as a 
+ // ratio between corrected and measured 3-p correlator:
+ //              bias = corrected/measured
+ // This bias is stored in histogram fDetectorBias.
+ Double_t bias = 0.;
+ if(measured3pCorrelator)
+ {
+  bias = corrected3pCorrelator/measured3pCorrelator;
+  fDetectorBiasHist->SetBinContent(1,bias);                                                          
+ } 
                                                                                                                                                                                                                                                                                                                                    
 } // end of AliFlowAnalysisWithMixedHarmonics::CorrectForDetectorEffects()
 
 //================================================================================================================
 
-void AliFlowAnalysisWithMixedHarmonics::FinalizeNonIsotropicTerms()
-{
- // Calculate correctly the errors of non-isotropic terms.
- // Transfer mean values from profiles to histograms: // to be improved
- for(Int_t mv=1;mv<=10;mv++)
- {
-  fNonIsotropicTermsHist->SetBinContent(mv,fNonIsotropicTermsPro->GetBinContent(mv));
- }
-} // end of void AliFlowAnalysisWithMixedHarmonics::FinalizeNonIsotropicTerms() 
-
-//================================================================================================================
-
-void AliFlowAnalysisWithMixedHarmonics::QuantifyBiasFromDetectorEffects()
-{
- // Quantify bias from detector inefficiences to 3-p correlator.
- // Remark: Bias is quantified as 1 - |corrected/nonCorrected| (in %) and stored in histogram fDetectorBias
- Double_t bias = 0.;
- if(f3pCorrelatorPro->GetBinContent(1))
- {
-  bias = 100.*TMath::Abs(1.-f3pCorrelatorHist->GetBinContent(1)/f3pCorrelatorPro->GetBinContent(1));
- }
- fDetectorBiasHist->SetBinContent(1,bias);
-} // end of void AliFlowAnalysisWithMixedHarmonics::QuantifyBiasFromDetectorEffects()
-
-//================================================================================================================
 void AliFlowAnalysisWithMixedHarmonics::ResetEventByEventQuantities()
 {
  // Reset all event by event quantities.
@@ -983,7 +1080,13 @@ void AliFlowAnalysisWithMixedHarmonics::ResetEventByEventQuantities()
  fSpk->Zero();
  
  f3pCorrelatorEBE->Reset();
-
+ for(Int_t sd=0;sd<2;sd++)
+ {
+  fRePEBE[sd]->Reset();
+  fImPEBE[sd]->Reset();
+ }
 } // end of void AliFlowAnalysisWithMixedHarmonics::ResetEventByEventQuantities()
 
 //================================================================================================================
@@ -1013,9 +1116,19 @@ void AliFlowAnalysisWithMixedHarmonics::Calculate3pCorrelator()
   // Fill event-by-event histogram and all-events profile:                     
   f3pCorrelatorEBE->SetBinContent(1,three2n1n1n);
   f3pCorrelatorPro->Fill(0.5,three2n1n1n,dMult*(dMult-1.)*(dMult-2.));
-
-  // Correlator vs multiplicity: // to be improved
-  f3pCorrelatorVsMPro->Fill(0.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),dMult*(dMult-1.)*(dMult-2.)); // to be improved (cross-checked)
+  
+  // 3-particle azimuthal correlator <cos(n*(2.*phi1-phi2-phi3))> vs multiplicity:
+  if(dMult<fMinMultiplicity) 
+  {
+   f3pCorrelatorVsMPro->Fill(0.5,three2n1n1n,dMult*(dMult-1.)*(dMult-2.));
+  } else if(dMult>=fMinMultiplicity+fNoOfMultipicityBins*fMultipicityBinWidth)
+    {
+     f3pCorrelatorVsMPro->Fill(0.5+fNoOfMultipicityBins+1,three2n1n1n,dMult*(dMult-1.)*(dMult-2.));  
+    } else
+      {
+       f3pCorrelatorVsMPro->Fill(1.5+(Int_t)((dMult-fMinMultiplicity)/fMultipicityBinWidth),three2n1n1n,dMult*(dMult-1.)*(dMult-2.));
+      }
+      
  } // end of if(!(fUsePhiWeights || fUsePtWeights || fUseEtaWeights)) 
 
  // b) Calculate 3-p correlator without using particle weights: 
@@ -1030,14 +1143,10 @@ void AliFlowAnalysisWithMixedHarmonics::Calculate3pCorrelator()
 
 void AliFlowAnalysisWithMixedHarmonics::CalculateNonIsotropicTerms()
 {
- // Calculate correction terms for non-uniform acceptance to 3-p correlator <cos[n(2phi1-phi2-phi3)]>.
- // a) Calculate using particle weights; 
- // b) Calculate without using particle weights. 
+ // Calculate non-isotropic terms which appear in the decomposition of 3-p correlator <cos[n(2phi1-phi2-phi3)]>.
  
- // These non-isotropic terms are stored in fNonIsotropicTermsPro. The same results but with
- // correct errors are stored in fNonIsotropicTermsHist. Binning of fNonIsotropicTermsPro and 
- // fNonIsotropicTermsHist is organized as follows:
+ // For detector with uniform acceptance all these terms vanish. These non-isotropic terms are stored in fNonIsotropicTermsPro.
+ // Binning of fNonIsotropicTermsPro is organized as follows:
  //  1st bin: <<cos(n*phi1)>>
  //  2nd bin: <<sin(n*phi1)>>
  //  3rd bin: <<cos(2n*phi1)>>
@@ -1049,6 +1158,9 @@ void AliFlowAnalysisWithMixedHarmonics::CalculateNonIsotropicTerms()
  //  9th bin: <<cos(n*(phi1-phi2-phi3)>>
  // 10th bin: <<sin(n*(phi1-phi2-phi3)>> 
  
+ // a) Calculate using particle weights; 
+ // b) Calculate without using particle weights. 
  // a) Calculate using particle weights: 
  if(!(fUsePhiWeights || fUsePtWeights || fUseEtaWeights))
  {
@@ -1130,3 +1242,50 @@ void AliFlowAnalysisWithMixedHarmonics::CalculateNonIsotropicTerms()
  
 } // end of void AliFlowAnalysisWithMixedHarmonics::CalculateNonIsotropicTerms()
 
+//================================================================================================================
+
+void AliFlowAnalysisWithMixedHarmonics::CalculateDifferential3pCorrelator()
+{
+ // Calculate differential 3-p azimuthal correlator cos[n(2phi1-psi2-psi3)] in terms of Q_{2n} and p_{n}.
+ // a) Calculate differential 3-p correlator without using particle weights;
+ // b) Calculate differential 3-p correlator with using particle weights.
+
+ // a) Calculate differential 3-p correlator without using particle weights: 
+ if(!(fUsePhiWeights || fUsePtWeights || fUseEtaWeights))
+ {
+  // Multiplicity (number of RPs):
+  Double_t dMult = (*fSpk)(0,0);
+  // Real and imaginary parts of non-weighted Q-vectors (Q_{n,0}) evaluated in harmonic 2n: 
+  Double_t dReQ2n = (*fReQnk)(1,0);
+  Double_t dImQ2n = (*fImQnk)(1,0);
+  for(Int_t sd=0;sd<2;sd++) // [(p1+p2)/2,|p1-p2|]
+  {
+   // looping over all bins and calculating reduced correlations: 
+   for(Int_t b=1;b<=fnBinsPt;b++)
+   {
+    // real and imaginary parts of p_{n}: 
+    Double_t p1nRe = fRePEBE[sd]->GetBinContent(b)*fRePEBE[sd]->GetBinEntries(b);
+    Double_t p1nIm = fImPEBE[sd]->GetBinContent(b)*fImPEBE[sd]->GetBinEntries(b);
+    // number of pairs of POIs in particular (p1+p2)/2 or |p1-p2| bin:
+    Double_t mp = fRePEBE[sd]->GetBinEntries(b);
+    Double_t cosP2nphi1M1npsi2M1npsi2 = 0; // cos[n(2phi1-psi2-psi3)]
+    if(mp*dMult>0.)
+    {
+     cosP2nphi1M1npsi2M1npsi2 = (p1nRe*dReQ2n+p1nIm*dImQ2n)/(mp*dMult);
+    }
+    f3pCorrelatorVsPtSumDiffPro[sd]->Fill(fPtMin+(b-1)*fPtBinWidth,cosP2nphi1M1npsi2M1npsi2,mp*dMult);
+   } // end of for(Int_t b=1;b<=fnBinsPt;b++)
+  } // end of for(Int_t sd=0;sd<2;sd++)      
+ } // end of if(!(fUsePhiWeights || fUsePtWeights || fUseEtaWeights)) 
+
+ // b) Calculate differential 3-p correlator by using particle weights: 
+ if(fUsePhiWeights || fUsePtWeights || fUseEtaWeights)
+ {
+ } // end of if(fUsePhiWeights || fUsePtWeights || fUseEtaWeights)
+} // end of void AliFlowAnalysisWithMixedHarmonics::CalculateDifferential3pCorrelator() 
+
+//================================================================================================================
+
index 5b652af..c350c3e 100644 (file)
@@ -11,7 +11,8 @@
  * can be used to:                                        *
  *                                                        *  
  *  a) Extract subdominant harmonics (like v1 and v4);    *
- *  b) Study strong parity violation.                     * 
+ *  b) Study flow of two-particle resonances;             *
+ *  c) Study strong parity violation.                     * 
  *                                                        * 
  * Author: Ante Bilandzic (abilandzic@gmail.com)          *
  *********************************************************/ 
@@ -44,6 +45,7 @@ class AliFlowAnalysisWithMixedHarmonics
   AliFlowAnalysisWithMixedHarmonics();
   virtual ~AliFlowAnalysisWithMixedHarmonics(); 
   // 0.) Methods called in the constructor:
+    virtual void InitializeArrays();
   // 1.) Method Init() and methods called within Init():
   virtual void Init();
     virtual void CrossCheckSettings();
@@ -54,19 +56,19 @@ class AliFlowAnalysisWithMixedHarmonics
     virtual void BookAllEventByEventQuantities();
     virtual void BookAllAllEventQuantities();
     virtual void BookAndFillWeightsHistograms();
-  // 2.) Method Make() and methods called within Make():
+    // 2.) Method Make() and methods called within Make():
   virtual void Make(AliFlowEventSimple *anEvent);
     virtual void CheckPointersUsedInMake();
     virtual void Calculate3pCorrelator();
     virtual void CalculateNonIsotropicTerms();
+    virtual void CalculateDifferential3pCorrelator();
     virtual void ResetEventByEventQuantities();
   // 3.) Method Finish() and methods called within Finish():
   virtual void Finish();  
     virtual void CheckPointersUsedInFinish(); 
     virtual void AccessSettings();       
-    virtual void FinalizeNonIsotropicTerms();       
     virtual void CorrectForDetectorEffects();
-    virtual void QuantifyBiasFromDetectorEffects();
+    virtual void PrintOnTheScreen();  
   // 4.) Method GetOutputHistograms and method called within it:
   virtual void GetOutputHistograms(TList *outputListHistos);
     virtual void GetPointersForBaseHistograms();
@@ -95,6 +97,8 @@ class AliFlowAnalysisWithMixedHarmonics
   Double_t GetMinMultiplicity() const {return this->fMinMultiplicity;}; 
   void SetCorrectForDetectorEffects(Bool_t const cfde) {this->fCorrectForDetectorEffects = cfde;};
   Bool_t GetCorrectForDetectorEffects() const {return this->fCorrectForDetectorEffects;}; 
+  void SetPrintOnTheScreen(Bool_t const pots) {this->fPrintOnTheScreen = pots;};
+  Bool_t GetPrintOnTheScreen() const {return this->fPrintOnTheScreen;};   
   void SetCommonHists(AliFlowCommonHist* const ch) {this->fCommonHists = ch;};
   AliFlowCommonHist* GetCommonHists() const {return this->fCommonHists;};
   void SetWeightsList(TList* const wl) {this->fWeightsList = (TList*)wl->Clone();}
@@ -124,20 +128,16 @@ class AliFlowAnalysisWithMixedHarmonics
   void SetNonIsotropicTermsVsMPro(TProfile2D* const nitVsMPro) {this->fNonIsotropicTermsVsMPro = nitVsMPro;};
   TProfile2D* GetNonIsotropicTermsVsMPro() const {return this->fNonIsotropicTermsVsMPro;};
   void SetResultsList(TList* const rlist) {this->fResultsList = rlist;}
-  TList* GetResultsList() const {return this->fResultsList;}  
-  void Set3pCorrelatorHist(TH1D* const s3pcHist) {this->f3pCorrelatorHist = s3pcHist;};
-  TH1D* Get3pCorrelatorHist() const {return this->f3pCorrelatorHist;};
+  TList* GetResultsList() const {return this->fResultsList;}    
+  void Set3pCorrelatorHist(TH1D* const s3pHist) {this->f3pCorrelatorHist = s3pHist;};
+  TH1D* Get3pCorrelatorHist() const {return this->f3pCorrelatorHist;};  
   void SetDetectorBiasHist(TH1D* const dbHist) {this->fDetectorBiasHist = dbHist;};
   TH1D* GetDetectorBiasHist() const {return this->fDetectorBiasHist;};  
-  void SetNonIsotropicTermsHist(TH1D* const nitHist) {this->fNonIsotropicTermsHist = nitHist;};
-  TH1D* GetNonIsotropicTermsHist() const {return this->fNonIsotropicTermsHist;};  
-  void Set3pCorrelatorVsMHist(TH1D* const s3pcVsMHist) {this->f3pCorrelatorVsMHist = s3pcVsMHist;};
-  TH1D* Get3pCorrelatorVsMHist() const {return this->f3pCorrelatorVsMHist;};
   void SetDetectorBiasVsMHist(TH1D* const dbVsMHist) {this->fDetectorBiasVsMHist = dbVsMHist;};
   TH1D* GetDetectorBiasVsMHist() const {return this->fDetectorBiasVsMHist;};  
-  void SetNonIsotropicTermsVsMHist(TH2D* const nitVsMHist) {this->fNonIsotropicTermsVsMHist = nitVsMHist;};
-  TH2D* GetNonIsotropicTermsVsMHist() const {return this->fNonIsotropicTermsVsMHist;};
-  
+  void Set3pCorrelatorVsPtSumDiffPro(TProfile* const s3pcvpsd, Int_t const sd) {this->f3pCorrelatorVsPtSumDiffPro[sd] = s3pcvpsd;};
+  TProfile* Get3pCorrelatorVsPtSumDiffPro(Int_t sd) const {return this->f3pCorrelatorVsPtSumDiffPro[sd];};
+    
  private:
   AliFlowAnalysisWithMixedHarmonics(const AliFlowAnalysisWithMixedHarmonics& afawQc);
   AliFlowAnalysisWithMixedHarmonics& operator=(const AliFlowAnalysisWithMixedHarmonics& afawQc); 
@@ -151,6 +151,7 @@ class AliFlowAnalysisWithMixedHarmonics
   Double_t fMultipicityBinWidth; // width of multiplicity bin
   Double_t fMinMultiplicity; // minimal multiplicity
   Bool_t fCorrectForDetectorEffects; // correct 3-p correlator for detector effects
+  Bool_t fPrintOnTheScreen; // print or not the final results on the screen
   // 1.) Common:
   AliFlowCommonHist *fCommonHists; // common control histograms (filled only with events with 3 or more tracks for 3-p correlators) 
   Int_t fnBinsPhi; // number of phi bins
@@ -180,21 +181,21 @@ class AliFlowAnalysisWithMixedHarmonics
   TMatrixD *fSpk; // fS[p][k] = S_{p,k} = (sum_{i=1}^{M} w_{i}^{k})^{p+1} // note p+1 in the power to use 0th index in p in non-trivial way
   TH1D *f3pCorrelatorEBE; // 3-p correlator <cos[n(2phi1-phi2-phi3)]> for single event
   TH1D *fNonIsotropicTermsEBE; // correction terms to 3-p correlator <cos[n(2phi1-phi2-phi3)]> for single event
+  TProfile *fRePEBE[2]; // real part of p_n vs [(p1+p2)/2,|p1-p2|]
+  TProfile *fImPEBE[2]; // imaginary part of p_n vs [(p1+p2)/2,|p1-p2|]
   // 4.) Profiles:
   TList *fProfileList; // list holding all all-event profiles 
-  TProfile *f3pCorrelatorPro; // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> for all events (with wrong errors)
-  TProfile *fNonIsotropicTermsPro; // correction terms to 3-p correlator <cos[n(2phi1-phi2-phi3)]> for all events (with wrong errors)
-  TProfile *f3pCorrelatorVsMPro; // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> for all events (with wrong errors) versus multiplicity
-  TProfile2D *fNonIsotropicTermsVsMPro; // correction terms to <cos[n(2phi1-phi2-phi3)]> for all events (with wrong errors) versus multiplicity
+  TProfile *f3pCorrelatorPro; // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> not corrected for detector effects
+  TProfile *fNonIsotropicTermsPro; // non-isotropic terms in the decomposition of 3-p correlator <<cos[n(2phi1-phi2-phi3)]>>
+  TProfile *f3pCorrelatorVsMPro; // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> vs multiplicity
+  TProfile2D *fNonIsotropicTermsVsMPro; // non-isotropic terms in the decomposition of <cos[n(2phi1-phi2-phi3)]> vs multiplicity
+  TProfile *f3pCorrelatorVsPtSumDiffPro[2]; // differential 3-p correlator <<cos[n(2phi1-psi2-psi3)]>> vs [(p1+p2)/2,|p1-p2|]
   // 5.) Final results:
   TList *fResultsList; // list holding objects with final results 
-  TH1D *f3pCorrelatorHist; // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> for all events (with correct errors)
+  TH1D *f3pCorrelatorHist; // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> corrected for detector effects
   TH1D *fDetectorBiasHist; // bias comming from detector inefficiencies to 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> (in %)
-  TH1D *fNonIsotropicTermsHist; // correction terms to 3-p correlator <cos[n(2phi1-phi2-phi3)]> for all events (with correct errors)
-  TH1D *f3pCorrelatorVsMHist; // 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> for all events (with correct errors) versus multiplicity
   TH1D *fDetectorBiasVsMHist; // bias comming from detector inefficiencies to 3-p correlator <<cos[n(2phi1-phi2-phi3)]>> (in %) versus multiplicity
-  TH2D *fNonIsotropicTermsVsMHist; // correction terms to <cos[n(2phi1-phi2-phi3)]> for all events (with correct errors) versus multiplicity
-  
+
   ClassDef(AliFlowAnalysisWithMixedHarmonics, 0);
 
 };
index 175928a..24f17e0 100644 (file)
@@ -19,7 +19,8 @@
  * Only in this class nested loops are used for flow analysis. *
  * Nested loops are used to evaluate:                          *
  *                                                             *  
- *  a) Distribution of relative angle difference (phi1-phi2).  *
+ *  a) Distribution of relative angle difference (phi1-phi2);  *
+ *  b) Cross-check the results for mixed harmonics.            *
  *                                                             *
  *       Author: Ante Bilandzic (abilandzic@gmail.com)         *
  ***************************************************************/ 
@@ -51,6 +52,7 @@ fHistList(NULL),
 fHistListName(NULL),
 fAnalysisLabel(NULL),
 fAnalysisSettings(NULL),
+fPrintOnTheScreen(kTRUE),
 fCommonHists(NULL),
 fnBinsPhi(0),
 fPhiMin(0),
@@ -72,11 +74,18 @@ fUseParticleWeights(NULL),
 fPhiWeights(NULL),
 fPtWeights(NULL),
 fEtaWeights(NULL),
-fResultsList(NULL),
-fRelativeAngleDistribution(NULL)
+fListRAD(NULL),
+fEvaluateNestedLoopsForRAD(kTRUE),
+fRelativeAngleDistribution(NULL),
+fListMH(NULL),
+fEvaluateNestedLoopsForMH(kFALSE),
+fCorrelatorIntegerMH(1),
+fCrossCheckInPtSumBinNo(4),  
+fCrossCheckInPtDiffBinNo(5) 
 {
  // Constructor. 
-  // Base list to hold all output objects:
+ // Base list to hold all output objects:
  fHistList = new TList();
  fHistListName = new TString("cobjNL");
  fHistList->SetName(fHistListName->Data());
@@ -85,8 +94,14 @@ fRelativeAngleDistribution(NULL)
  // List to hold histograms with phi, pt and eta weights:      
  fWeightsList = new TList();
  
- // List to hold objects with final results:      
- fResultsList = new TList();
+ // List to hold objects relevant for relative angle distributions:      
+ fListRAD = new TList();
+ // List holding objects relevant for debugging and cross-checking of MH class: 
+ fListMH = new TList();
+ // Initialize all arrays: 
+ this->InitializeArraysForMH();
  
 } // AliFlowAnalysisWithNestedLoops::AliFlowAnalysisWithNestedLoops()
  
@@ -119,13 +134,16 @@ void AliFlowAnalysisWithNestedLoops::Init()
  //in an analysis
  Bool_t oldHistAddStatus = TH1::AddDirectoryStatus();
  TH1::AddDirectory(kFALSE);
+ TH1::SetDefaultSumw2();
   
  this->CrossCheckSettings();
  this->AccessConstants();
  this->BookAndNestAllLists();
- this->BookProfileHoldingSettings();
+ this->BookAndFillProfileHoldingSettings();
  this->BookCommonHistograms();
- this->BookEverythingForDistributions();
+ this->BookEverythingForRAD();
+ this->BookEverythingForMH();
  this->BookAndFillWeightsHistograms();
 
  //restore old status
@@ -138,56 +156,18 @@ void AliFlowAnalysisWithNestedLoops::Make(AliFlowEventSimple* anEvent)
 {
  // Running over data only in this method.
  
- // a) Check all pointers used in method Make();
- // b) Define local variables;
- // c) Fill common control histograms;
- // d) Loop over data and store for each distinct pair phi1-phi2 in fRelativeAngleDistribution.
+ // a) Check all pointers used in this method;
+ // b) Fill common control histograms;
+ // c) Evaluate nested loops for relative angle distribution;
+ // d) Evaluate nested loops for mixed harmonics.
  
- // a) Check all pointers used in method Make():
  this->CheckPointersUsedInMake();
- // b) Define local variables:
- Double_t dPhi1=0., dPhi2=0.; // azimuthal angles in the laboratory frame
- AliFlowTrackSimple *aftsTrack = NULL; // simple track
- // c) Fill common control histograms:
  fCommonHists->FillControlHistograms(anEvent);  
- // d) Loop over data and store for each distinct pair phi1-phi2 in fRelativeAngleDistribution:
- Int_t nPrim = anEvent->NumberOfTracks();  // nPrim = total number of primary tracks, i.e. nPrim = nRP + nPOI + rest, where:
-                                           // nRP   = # of particles used to determine the reaction plane ("Reference Particles");
-                                           // nPOI  = # of particles of interest for a detailed flow analysis ("Particles of Interest");
-                                           // rest  = # of particles which are not niether RPs nor POIs.  
- // Start nested loops over data:
- for(Int_t i=0;i<nPrim;i++) 
- { 
-  aftsTrack=anEvent->GetTrack(i);
-  if(aftsTrack)
-  {
-   if(!aftsTrack->InRPSelection()) continue; // consider only tracks which are RPs 
-   dPhi1 = aftsTrack->Phi();
-   for(Int_t j=0;j<nPrim;j++) 
-   { 
-    if(j==i) continue; // eliminating trivial contribution from autocorrelation
-    aftsTrack=anEvent->GetTrack(j);
-    if(aftsTrack)
-    {
-     if(!aftsTrack->InRPSelection()) continue; // consider only tracks which are RPs 
-     dPhi2 = aftsTrack->Phi();
-     // Fill the histogram:
-     fRelativeAngleDistribution->Fill(dPhi1-dPhi2);
-    }
-   } // end of for(Int_t j=0;j<nPrim;j++)
-  } else // to if(aftsTrack)
-    {
-     cout<<endl;
-     cout<<" WARNING (NL): No particle! (i.e. aftsTrack is a NULL pointer in AFAWNL::Make().)"<<endl;
-     cout<<endl;       
-    }
-    
- } // end of for(Int_t i=0;i<nPrim;i++) 
+ if(fEvaluateNestedLoopsForRAD) this->EvaluateNestedLoopsForRAD(anEvent);
+ if(fEvaluateNestedLoopsForMH) this->EvaluateNestedLoopsForMH(anEvent);
   
 } // end of AliFlowAnalysisWithNestedLoops::Make(AliFlowEventSimple* anEvent)
+
 //================================================================================================================
 
 void AliFlowAnalysisWithNestedLoops::Finish()
@@ -196,9 +176,11 @@ void AliFlowAnalysisWithNestedLoops::Finish()
  
  // a) Check all pointers used in this method;
  // b) Access settings for analysis with mixed harmonics;
+ // c) Print on the screen.
  
  this->CheckPointersUsedInFinish();
  this->AccessSettings();
+ if(fPrintOnTheScreen) this->PrintOnTheScreen();
                                                                                                                                                                                                                                                                                                                
 } // end of AliFlowAnalysisWithNestedLoops::Finish()
 
@@ -219,8 +201,12 @@ void AliFlowAnalysisWithNestedLoops::GetOutputHistograms(TList *outputListHistos
    cout<<endl;
    exit(0);
   }
+  this->GetPointersForBaseHistograms();
   this->GetPointersForCommonHistograms();
-  this->GetPointersForResultsHistograms();
+  Bool_t bEvaluateNestedLoopsForRAD = (Bool_t) fAnalysisSettings->GetBinContent(1); // to be improved (not needed here?)
+  Bool_t bEvaluateNestedLoopsForMH = (Bool_t) fAnalysisSettings->GetBinContent(2); // to be improved (not needed here?)
+  if(bEvaluateNestedLoopsForRAD) this->GetPointersForRAD();
+  if(bEvaluateNestedLoopsForMH) this->GetPointersForMH();
  } else 
    {
     cout<<endl;
@@ -233,6 +219,27 @@ void AliFlowAnalysisWithNestedLoops::GetOutputHistograms(TList *outputListHistos
 
 //================================================================================================================
 
+void AliFlowAnalysisWithNestedLoops::GetPointersForBaseHistograms() 
+{
+ // Get pointers to base histograms.
+ TString analysisSettingsName = "fAnalysisSettings";
+ TProfile *analysisSettings = dynamic_cast<TProfile*>(fHistList->FindObject(analysisSettingsName.Data()));
+ if(analysisSettings) 
+ {
+  this->SetAnalysisSettings(analysisSettings); 
+ } else
+   {
+    cout<<endl;
+    cout<<" WARNING (NL): analysisSettings is NULL in AFAWNL::GPFBH() !!!!"<<endl;
+    cout<<endl;
+    exit(0);  
+   }
+} // end of void AliFlowAnalysisWithNestedLoops::GetPointersForBaseHistograms()
+
+//================================================================================================================
+
 void AliFlowAnalysisWithNestedLoops::GetPointersForCommonHistograms() 
 {
  // Get pointers to common control histograms.
@@ -254,26 +261,52 @@ void AliFlowAnalysisWithNestedLoops::GetPointersForCommonHistograms()
 
 //================================================================================================================
 
-void AliFlowAnalysisWithNestedLoops::GetPointersForResultsHistograms() 
+void AliFlowAnalysisWithNestedLoops::GetPointersForRAD() 
 {
- // Get pointers to histograms holding final results.
+ // Get pointers to objects relevant for relative angle distributions.
  
- TList *resultsList = NULL;
- resultsList = dynamic_cast<TList*>(fHistList->FindObject("Results"));
- if(!resultsList) 
+ TList *listRAD = NULL;
+ listRAD = dynamic_cast<TList*>(fHistList->FindObject("Relative Angle Distribution"));
+ if(!listRAD) 
  {
-  cout<<"WARNING: resultsList is NULL in AFAWNL::GPFRH() !!!!"<<endl;
+  cout<<"WARNING: listRAD is NULL in AFAWNL::GPFRAD() !!!!"<<endl;
   exit(0); 
  }  
+
  TString relativeAngleDistributionName = "fRelativeAngleDistribution";
- TH1D *relativeAngleDistribution = dynamic_cast<TH1D*>(resultsList->FindObject(relativeAngleDistributionName.Data()));
+ TH1D *relativeAngleDistribution = dynamic_cast<TH1D*>(listRAD->FindObject(relativeAngleDistributionName.Data()));
  if(relativeAngleDistribution)
  {
   this->SetRelativeAngleDistribution(relativeAngleDistribution);  
  }
   
-} // end of void AliFlowAnalysisWithNestedLoops::GetPointersForResultsHistograms()
+} // end of void AliFlowAnalysisWithNestedLoops::GetPointersForRAD()
+
+//================================================================================================================
+
+void AliFlowAnalysisWithNestedLoops::GetPointersForMH() 
+{
+ // Get pointers to objects evaluated with nested loops.
+ TList *listMH = NULL;
+ listMH = dynamic_cast<TList*>(fHistList->FindObject("Mixed Harmonics"));
+ if(!listMH) 
+ {
+  cout<<"WARNING: listMH is NULL in AFAWNL::GPFMH() !!!!"<<endl;
+  exit(0); 
+ }  
+  
+ TString psdFlag[2] = {"PtSum","PtDiff"};
+ for(Int_t sd=0;sd<2;sd++)
+ {
+  TProfile *p3pCorrelatorVsPtSumDiffDirectPro = dynamic_cast<TProfile*>(listMH->FindObject(Form("f3pCorrelatorDirectVs%s",psdFlag[sd].Data())));
+  if(p3pCorrelatorVsPtSumDiffDirectPro)
+  {
+   this->Set3pCorrelatorVsPtSumDiffDirectPro(p3pCorrelatorVsPtSumDiffDirectPro,sd);  
+  }  
+ } // end of for(Int_t sd=0;sd<2;sd++)
+} // end of void AliFlowAnalysisWithNestedLoops::GetPointersForMH() 
 
 //================================================================================================================
 
@@ -306,26 +339,41 @@ void AliFlowAnalysisWithNestedLoops::BookAndNestAllLists()
  fWeightsList->SetName("Weights");
  fWeightsList->SetOwner(kTRUE);   
  fHistList->Add(fWeightsList); 
- // Results:
- fResultsList->SetName("Results");
- fResultsList->SetOwner(kTRUE);   
- fHistList->Add(fResultsList); 
+ // List for Relative Angle Distribution:
+ fListRAD->SetName("Relative Angle Distribution");
+ fListRAD->SetOwner(kTRUE);   
+ if(fEvaluateNestedLoopsForRAD) fHistList->Add(fListRAD); 
+ // List for Mixed Harmonics:
+ fListMH->SetName("Mixed Harmonics");
+ fListMH->SetOwner(kTRUE);   
+ if(fEvaluateNestedLoopsForMH) fHistList->Add(fListMH); 
 
 } // end of void AliFlowAnalysisWithNestedLoops::BookAndNestAllLists()
 
 //================================================================================================================
 
-void AliFlowAnalysisWithNestedLoops::BookProfileHoldingSettings()
+void AliFlowAnalysisWithNestedLoops::BookAndFillProfileHoldingSettings()
 {
  // Book profile to hold all analysis settings.
 
  TString analysisSettingsName = "fAnalysisSettings";
- fAnalysisSettings = new TProfile(analysisSettingsName.Data(),"Settings for analysis with nested loops",1,0,1);
- //fAnalysisSettings->GetXaxis()->SetBinLabel(1," ... ");
- //fAnalysisSettings->Fill(0.5, ... );
+ fAnalysisSettings = new TProfile(analysisSettingsName.Data(),"Settings for analysis with nested loops",6,0,6);
+ fAnalysisSettings->GetXaxis()->SetLabelSize(0.035);
+ fAnalysisSettings->GetXaxis()->SetBinLabel(1,"Nested loops for RAD?");
+ fAnalysisSettings->Fill(0.5,(Int_t)fEvaluateNestedLoopsForRAD);
+ fAnalysisSettings->GetXaxis()->SetBinLabel(2,"Nested loops for MH?");
+ fAnalysisSettings->Fill(1.5,(Int_t)fEvaluateNestedLoopsForMH);
+ fAnalysisSettings->GetXaxis()->SetBinLabel(3,"Integer n in cos(n(2#phi_{1}-#psi_{2}-#psi_{3}))");
+ fAnalysisSettings->Fill(2.5,(Int_t)fCorrelatorIntegerMH);
+ fAnalysisSettings->GetXaxis()->SetBinLabel(4,"Print on the screen?");
+ fAnalysisSettings->Fill(3.5,(Int_t)fPrintOnTheScreen); 
+ fAnalysisSettings->GetXaxis()->SetBinLabel(5,"fCrossCheckInPtSumBinNo");
+ fAnalysisSettings->Fill(4.5,fCrossCheckInPtSumBinNo);
+ fAnalysisSettings->GetXaxis()->SetBinLabel(6,"fCrossCheckInPtDiffBinNo");
+ fAnalysisSettings->Fill(5.5,fCrossCheckInPtDiffBinNo);
  fHistList->Add(fAnalysisSettings);
  
-} // end of void AliFlowAnalysisWithNestedLoops::BookProfileHoldingSettings()
+} // end of void AliFlowAnalysisWithNestedLoops::BookAndFillProfileHoldingSettings()
 
 //================================================================================================================
 
@@ -341,18 +389,17 @@ void AliFlowAnalysisWithNestedLoops::BookCommonHistograms()
 
 //================================================================================================================
 
-void AliFlowAnalysisWithNestedLoops::BookEverythingForDistributions()
+void AliFlowAnalysisWithNestedLoops::BookEverythingForRAD()
 {
- // Book all objects relevant for distributions.
+ // Book all objects relevant calculation of relative angle distribution.
  
- // Histogram to hold distribution of phi1-phi2:
  TString relativeAngleDistributionName = "fRelativeAngleDistribution";
  fRelativeAngleDistribution = new TH1D(relativeAngleDistributionName.Data(),"Relative angle distribution",720,-TMath::TwoPi(),TMath::TwoPi());
  fRelativeAngleDistribution->GetYaxis()->SetTitle("#frac{dN}{#Delta #phi}"); 
  fRelativeAngleDistribution->GetXaxis()->SetTitle("#Delta #phi");
- fResultsList->Add(fRelativeAngleDistribution);
+ fListRAD->Add(fRelativeAngleDistribution);
 
-} // end fo void AliFlowAnalysisWithNestedLoops::BookEverythingForDistributions()
+} // end fo void AliFlowAnalysisWithNestedLoops::BookEverythingForRAD()
 
 //================================================================================================================
 
@@ -471,13 +518,8 @@ void AliFlowAnalysisWithNestedLoops::CheckPointersUsedInMake()
 {
  // Check pointers used in method Make().
                         
- if(!fRelativeAngleDistribution)
- {
-  cout<<endl;
-  cout<<" WARNING (NL): !fRelativeAngleDistribution is NULL !!!!"<<endl;
-  cout<<endl;
-  exit(0); 
- }
+ if(fEvaluateNestedLoopsForRAD) CheckPointersForRAD("Make");
+ if(fEvaluateNestedLoopsForMH) CheckPointersForMH("Make"); 
                                                                                                                                                                                                                                                                                                                                    
 } // end of AliFlowAnalysisWithNestedLoops::CheckPointersUsedInMake()
 
@@ -486,15 +528,256 @@ void AliFlowAnalysisWithNestedLoops::CheckPointersUsedInMake()
 void AliFlowAnalysisWithNestedLoops::CheckPointersUsedInFinish()
 {
  // Check pointers used in method Finish().
+ if(fEvaluateNestedLoopsForRAD) CheckPointersForRAD("Finish");
+ if(fEvaluateNestedLoopsForMH) CheckPointersForMH("Finish"); 
                                                                                                                                                                                                                                                                                                                                    
 } // end of AliFlowAnalysisWithNestedLoops::CheckPointersUsedInFinish()
 
 //================================================================================================================
 
+void AliFlowAnalysisWithNestedLoops::CheckPointersForRAD(TString where)
+{
+ // Check pointers relevant for calculation of relative angle distribution.
+ if(!fRelativeAngleDistribution)
+ {
+  cout<<endl;
+  cout<<" WARNING (NL): !fRelativeAngleDistribution is NULL in "<<where.Data()<<"() !!!!"<<endl;
+  cout<<endl;
+  exit(0); 
+ }
+ if(strcmp(where.Data(),"Make") == 0)
+ {
+  // Check pointers used only in method Make():
+  // ...
+ }
+ else if(strcmp(where.Data(),"Finish") == 0)
+ {
+  // Check pointers used only in method Finish():
+  // ...
+ }
+
+} // end of void AliFlowAnalysisWithNestedLoops::CheckPointersForRAD(TString where)
+
+//================================================================================================================
+
+void AliFlowAnalysisWithNestedLoops::CheckPointersForMH(TString where)
+{
+ // Check pointers relevant for calculation of mixed harmonics.
+ for(Int_t sd=0;sd<2;sd++)
+ {
+  if(!(f3pCorrelatorVsPtSumDiffDirectPro[sd]))
+  {
+   cout<<endl;
+   cout<<" WARNING (NL): !"<<Form("f3pCorrelatorVsPtSumDiffDirectPro[%d]",sd)<<" is NULL in "<<where.Data()<<"() !!!!"<<endl;
+   cout<<endl;
+   exit(0);   
+  } 
+ }
+
+ if(strcmp(where.Data(),"Make") == 0)
+ {
+  // Check pointers used only in method Make():
+  // ...
+ }
+ else if(strcmp(where.Data(),"Finish") == 0)
+ {
+  // Check pointers used only in method Finish():
+  // ...
+ }
+
+} // end of void AliFlowAnalysisWithNestedLoops::CheckPointersForMH(TString where)
+
+//================================================================================================================
+
 void AliFlowAnalysisWithNestedLoops::AccessSettings()
 {
- // Access the settings for analysis with mixed harmonics.
+ // Access the settings for analysis.
+ fEvaluateNestedLoopsForRAD = (Bool_t)fAnalysisSettings->GetBinContent(1);
+ fEvaluateNestedLoopsForMH = (Bool_t)fAnalysisSettings->GetBinContent(2);
+ fCorrelatorIntegerMH = (Int_t)fAnalysisSettings->GetBinContent(3);
+ fPrintOnTheScreen = (Bool_t)fAnalysisSettings->GetBinContent(4);
+ fCrossCheckInPtSumBinNo = (Int_t)fAnalysisSettings->GetBinContent(5);
+ fCrossCheckInPtDiffBinNo = (Int_t)fAnalysisSettings->GetBinContent(6);
                                                                                                                                                                                                                                                                                                                                    
 } // end of AliFlowAnalysisWithNestedLoops::AccessSettings()
 
 //================================================================================================================
+
+void AliFlowAnalysisWithNestedLoops::InitializeArraysForMH()
+{
+ // Initialize arrays mixed harmonics calculations.
+ for(Int_t sd=0;sd<2;sd++) // sum or difference
+ {
+  f3pCorrelatorVsPtSumDiffDirectPro[sd] = NULL;
+ }
+  
+} // end of AliFlowAnalysisWithNestedLoops::InitializeArraysForMH()
+
+//================================================================================================================  
+
+void AliFlowAnalysisWithNestedLoops::BookEverythingForMH()
+{
+ // Book all objects relevant for mixed harmonics.
+ if(fEvaluateNestedLoopsForMH)
+ {
+  TString psdFlag[2] = {"PtSum","PtDiff"};
+  TString psdTitleFlag[2] = {"(p_{t,1}+p_{t,2})/2","#left|p_{t,1}-p_{t,2}#right|"};
+  //TString s3pCorrelatorVsPtSumDiffDirectProName = "f3pCorrelatorVsPtSumDiffDirectPro";
+  for(Int_t sd=0;sd<2;sd++)
+  {
+   // to be improved: hardwired ,fnBinsPt,0.,fPtMax):
+   f3pCorrelatorVsPtSumDiffDirectPro[sd] = new TProfile(Form("f3pCorrelatorDirectVs%s",psdFlag[sd].Data()),"",fnBinsPt,0.,fPtMax);
+   //f3pCorrelatorVsPtSumDiffDirectPro[sd]->SetLabelSize(0.05);
+   //f3pCorrelatorVsPtSumDiffDirectPro[sd]->SetMarkerStyle(25);
+   f3pCorrelatorVsPtSumDiffDirectPro[sd]->GetXaxis()->SetTitle(psdTitleFlag[sd].Data());
+   fListMH->Add(f3pCorrelatorVsPtSumDiffDirectPro[sd]);
+  }  
+ } // end of if(fEvaluateNestedLoopsForMH)
+
+} // end of AliFlowAnalysisWithNestedLoops::BookEverythingForMH()
+
+//================================================================================================================
+
+void AliFlowAnalysisWithNestedLoops::EvaluateNestedLoopsForRAD(AliFlowEventSimple *anEvent)
+{
+ // Evaluate nested loops needed for calculation of relative angle distribution.
+ Double_t dPhi1=0., dPhi2=0.; // azimuthal angles in the laboratory frame
+ AliFlowTrackSimple *aftsTrack = NULL; // simple track
+  
+ // Loop over data and store for each distinct pair phi1-phi2 in fRelativeAngleDistribution:
+ Int_t nPrim = anEvent->NumberOfTracks();  // nPrim = total number of primary tracks, i.e. nPrim = nRP + nPOI + rest, where:
+                                           // nRP   = # of particles used to determine the reaction plane ("Reference Particles");
+                                           // nPOI  = # of particles of interest for a detailed flow analysis ("Particles of Interest");
+                                           // rest  = # of particles which are not niether RPs nor POIs.  
+ // Start nested loops over data:
+ for(Int_t i=0;i<nPrim;i++) 
+ { 
+  aftsTrack=anEvent->GetTrack(i);
+  if(aftsTrack)
+  {
+   if(!aftsTrack->InRPSelection()) continue; // consider only tracks which are RPs 
+   dPhi1 = aftsTrack->Phi();
+   for(Int_t j=0;j<nPrim;j++) 
+   { 
+    if(j==i) continue; // eliminating trivial contribution from autocorrelation
+    aftsTrack=anEvent->GetTrack(j);
+    if(aftsTrack)
+    {
+     if(!aftsTrack->InRPSelection()) continue; // consider only tracks which are RPs 
+     dPhi2 = aftsTrack->Phi();
+     // Fill the histogram:
+     fRelativeAngleDistribution->Fill(dPhi1-dPhi2);
+    }
+   } // end of for(Int_t j=0;j<nPrim;j++)
+  } else // to if(aftsTrack)
+    {
+     cout<<endl;
+     cout<<" WARNING (NL): No particle! (i.e. aftsTrack is a NULL pointer in AFAWNL::Make().)"<<endl;
+     cout<<endl;       
+    }
+ } // end of for(Int_t i=0;i<nPrim;i++) 
+} // end of void AliFlowAnalysisWithNestedLoops::EvaluateNestedLoopsForRAD(AliFlowEventSimple *anEvent)
+
+//================================================================================================================
+
+void AliFlowAnalysisWithNestedLoops::EvaluateNestedLoopsForMH(AliFlowEventSimple *anEvent)
+{
+ // Evaluate nested loops needed for calculation of mixed harmonics.
+ // Remark: phi label azimuthal angle of RP particle and psi label azimuthal angle of POI particle.
+ Int_t nPrim = anEvent->NumberOfTracks(); 
+ AliFlowTrackSimple *aftsTrack = NULL;
+ Double_t phi1=0.,psi2=0.,psi3=0.; // angles in the correlator cos[n(2phi1-psi2-psi3)] 
+ Double_t pt2=0.,pt3=0.; // transverse momenta of psi2 and psi3
+ Int_t n = fCorrelatorIntegerMH; 
+ // Evaluting differential correlator cos[n(2phi1-psi2-psi3)] with three nested loops:
+ for(Int_t i1=0;i1<nPrim;i1++)
+ {
+  aftsTrack=anEvent->GetTrack(i1);
+  // RP condition (first particle in the correlator must be RP): 
+  if(!(aftsTrack->InRPSelection())) continue;
+  phi1 = aftsTrack->Phi();  
+  for(Int_t i2=0;i2<nPrim;i2++)
+  {
+   if(i2==i1) continue;
+   aftsTrack = anEvent->GetTrack(i2);
+   // POI condition (second particle in the correlator must be POI):
+   if(!(aftsTrack->InPOISelection())) continue;
+   psi2 = aftsTrack->Phi();
+   pt2 = aftsTrack->Pt();
+   for(Int_t i3=0;i3<nPrim;i3++)
+   {
+    if(i3==i1||i3==i2) continue;
+    aftsTrack=anEvent->GetTrack(i3);
+    // POI condition (third particle in the correlator must be POI):
+    if(!(aftsTrack->InPOISelection())) continue;
+    psi3 = aftsTrack->Phi();
+    pt3 = aftsTrack->Pt();   
+    // Evaluate and store differential correlator cos[n(2phi1-psi2-psi3)]:
+    Double_t ptSum = (pt2+pt3)/2.;
+    Double_t ptDiff = TMath::Abs(pt2-pt3);
+    Double_t diff3pCorrelator = TMath::Cos(n*(2.*phi1-psi2-psi3));
+    f3pCorrelatorVsPtSumDiffDirectPro[0]->Fill(ptSum,diff3pCorrelator,1.);
+    f3pCorrelatorVsPtSumDiffDirectPro[1]->Fill(ptDiff,diff3pCorrelator,1.);
+   } // end of for(Int_t i3=0;i3<nPrim;i3++)  
+  } // end of for(Int_t i2=0;i2<nPrim;i2++)  
+ } // end of for(Int_t i1=0;i1<nPrim;i1++)
+
+} // end of void AliFlowAnalysisWithNestedLoops::EvaluateNestedLoopsForMH(AliFlowEventSimple *anEvent)
+
+//================================================================================================================
+
+void AliFlowAnalysisWithNestedLoops::PrintOnTheScreen()
+{
+ // Print on the screen.
+ cout<<endl;
+ cout<<"****************************************************"<<endl;
+ cout<<"****************************************************"<<endl;
+ cout<<"                  Nested Loops                 "<<endl; 
+ cout<<endl;
+ if(fEvaluateNestedLoopsForRAD) 
+ {
+  cout<<"  Evaluated for relative angle distribution."<<endl;
+ }
+ if(fEvaluateNestedLoopsForMH) 
+ {
+  cout<<"  Evaluated for mixed harmonics."<<endl;
+  if(fCorrelatorIntegerMH!=1)
+  {
+   cout<< "  cos["<<fCorrelatorIntegerMH<<"(2phi1-psi2-psi3)] = "<<endl;
+  } else
+    {
+     cout<< "  cos(2phi1-psi2-psi3) = "<<endl;
+    } 
+  cout<< "  a) in pt sum bin "<<fCrossCheckInPtSumBinNo<<": "<<
+  f3pCorrelatorVsPtSumDiffDirectPro[0]->GetBinContent(fCrossCheckInPtSumBinNo)<<
+  " +/- "<<f3pCorrelatorVsPtSumDiffDirectPro[0]->GetBinError(fCrossCheckInPtSumBinNo)<<endl;
+  cout<< "  b) in pt diff bin "<<fCrossCheckInPtDiffBinNo<<": "<<
+  f3pCorrelatorVsPtSumDiffDirectPro[1]->GetBinContent(fCrossCheckInPtDiffBinNo)<<
+  " +/- "<<f3pCorrelatorVsPtSumDiffDirectPro[1]->GetBinError(fCrossCheckInPtDiffBinNo)<<endl;
+ }
+ if(!fEvaluateNestedLoopsForRAD && !fEvaluateNestedLoopsForMH)
+ {
+  cout<<"  Not evaluated."<<endl;
+ }
+ cout<<endl;
+ cout<<"****************************************************"<<endl;
+ cout<<"****************************************************"<<endl;
+ cout<<endl;
+} // end of void AliFlowAnalysisWithNestedLoops::PrintOnTheScreen()
+
+
index 36e897d..026c9cd 100644 (file)
@@ -8,7 +8,8 @@
  * Only in this class nested loops are used for flow analysis. *
  * Nested loops are used to evaluate:                          *
  *                                                             *  
- *  a) Distribution of relative angle difference (phi1-phi2).  *
+ *  a) Distribution of relative angle difference (phi1-phi2);  *
+ *  b) Cross-check the results for mixed harmonics.            *
  *                                                             *
  *       Author: Ante Bilandzic (abilandzic@gmail.com)         *
  ***************************************************************/ 
@@ -36,29 +37,38 @@ class AliFlowAnalysisWithNestedLoops
   AliFlowAnalysisWithNestedLoops();
   virtual ~AliFlowAnalysisWithNestedLoops(); 
   // 0.) Methods called in the constructor:
+    virtual void InitializeArraysForMH();
   // 1.) Method Init() and methods called within Init():
   virtual void Init();
     virtual void CrossCheckSettings();
     virtual void AccessConstants();
     virtual void BookAndNestAllLists();
-    virtual void BookProfileHoldingSettings();
+    virtual void BookAndFillProfileHoldingSettings();
     virtual void BookCommonHistograms();
-    virtual void BookEverythingForDistributions();
+    virtual void BookEverythingForRAD(); // RAD = relative angle distribution phi1-phi2
+    virtual void BookEverythingForMH(); // MH = Mixed Harmonics
     virtual void BookAndFillWeightsHistograms();
   // 2.) Method Make() and methods called within Make():
   virtual void Make(AliFlowEventSimple *anEvent);
     virtual void CheckPointersUsedInMake();
+    virtual void EvaluateNestedLoopsForRAD(AliFlowEventSimple *anEvent);
+    virtual void EvaluateNestedLoopsForMH(AliFlowEventSimple *anEvent);
   // 3.) Method Finish() and methods called within Finish():
   virtual void Finish();  
     virtual void CheckPointersUsedInFinish(); 
-    virtual void AccessSettings();       
+    virtual void AccessSettings();  
+    virtual void PrintOnTheScreen();     
   // 4.) Method GetOutputHistograms and method called within it:
   virtual void GetOutputHistograms(TList *outputListHistos);
+    virtual void GetPointersForBaseHistograms();
     virtual void GetPointersForCommonHistograms();
-    virtual void GetPointersForResultsHistograms();
+    virtual void GetPointersForRAD();
+    virtual void GetPointersForMH();
   // 5.) Other methods:   
   virtual void WriteHistograms(TString outputFileName);
   virtual void WriteHistograms(TDirectoryFile *outputFileName);  
+  virtual void CheckPointersForRAD(TString where);
+  virtual void CheckPointersForMH(TString where);
   // 6.) Setters and getters:
   void SetHistList(TList* const hl) {this->fHistList = hl;}
   TList* GetHistList() const {return this->fHistList;}  
@@ -68,6 +78,8 @@ class AliFlowAnalysisWithNestedLoops
   TString *GetAnalysisLabel() const {return this->fAnalysisLabel;};
   void SetAnalysisSettings(TProfile* const as) {this->fAnalysisSettings = as;};
   TProfile* GetAnalysisSettings() const {return this->fAnalysisSettings;};
+  void SetPrintOnTheScreen(Bool_t const pots) {this->fPrintOnTheScreen = pots;};
+  Bool_t GetPrintOnTheScreen() const {return this->fPrintOnTheScreen;};   
   void SetCommonHists(AliFlowCommonHist* const ch) {this->fCommonHists = ch;};
   AliFlowCommonHist* GetCommonHists() const {return this->fCommonHists;};
   void SetWeightsList(TList* const wl) {this->fWeightsList = (TList*)wl->Clone();}
@@ -86,10 +98,27 @@ class AliFlowAnalysisWithNestedLoops
   TH1D* GetPtWeights() const {return this->fPtWeights;};
   void SetEtaWeights(TH1D* const histEtaWeights) {this->fEtaWeights = histEtaWeights;};
   TH1D* GetEtaWeights() const {return this->fEtaWeights;};
-  void SetResultsList(TList* const rlist) {this->fResultsList = rlist;}
-  TList* GetResultsList() const {return this->fResultsList;}  
+  void SetListRAD(TList* const lRAD) {this->fListRAD = lRAD;}
+  TList* GetListRAD() const {return this->fListRAD;}  
+  void SetEvaluateNestedLoopsForRAD(Bool_t const enlfRAD) {this->fEvaluateNestedLoopsForRAD = enlfRAD;};
+  Bool_t GetEvaluateNestedLoopsForRAD() const {return this->fEvaluateNestedLoopsForRAD;};
   void SetRelativeAngleDistribution(TH1D* const rad) {this->fRelativeAngleDistribution = rad;};
-  TH1D* GetRelativeAngleDistribution() const {return this->fRelativeAngleDistribution;};
+  TH1D* GetRelativeAngleDistribution() const {return this->fRelativeAngleDistribution;}; 
+  // QC:
+  // ...
+  // MH:
+  void SetListMH(TList* const lMH) {this->fListMH = lMH;}
+  TList* GetListMH() const {return this->fListMH;}  
+  void SetEvaluateNestedLoopsForMH(Bool_t const enlfMH) {this->fEvaluateNestedLoopsForMH = enlfMH;};
+  Bool_t GetEvaluateNestedLoopsForMH() const {return this->fEvaluateNestedLoopsForMH;};
+  void SetCorrelatorIntegerMH(Int_t const ciMH) {this->fCorrelatorIntegerMH = ciMH;};
+  Int_t GetCorrelatorIntegerMH() const {return this->fCorrelatorIntegerMH;}; 
+  void Set3pCorrelatorVsPtSumDiffDirectPro(TProfile* const s3pcvpsdd, Int_t const sd) {this->f3pCorrelatorVsPtSumDiffDirectPro[sd] = s3pcvpsdd;};
+  TProfile* Get3pCorrelatorVsPtSumDiffDirectPro(Int_t sd) const {return this->f3pCorrelatorVsPtSumDiffDirectPro[sd];};
+  void SetCrossCheckInPtSumBinNo(Int_t const ccipsbn) {this->fCrossCheckInPtSumBinNo = ccipsbn;};
+  Int_t GetCrossCheckInPtSumBinNo() const {return this->fCrossCheckInPtSumBinNo;}; 
+  void SetCrossCheckInPtDiffBinNo(Int_t const ccipdbn) {this->fCrossCheckInPtDiffBinNo = ccipdbn;};
+  Int_t GetCrossCheckInPtDiffBinNo() const {return this->fCrossCheckInPtDiffBinNo;};  
   
  private:
   AliFlowAnalysisWithNestedLoops(const AliFlowAnalysisWithNestedLoops& afawQc);
@@ -99,6 +128,7 @@ class AliFlowAnalysisWithNestedLoops
   TString *fHistListName; // name of base list
   TString *fAnalysisLabel; // analysis label 
   TProfile *fAnalysisSettings; // profile to hold analysis settings
+  Bool_t fPrintOnTheScreen; // print or not on the screen
   // 1.) Common:
   AliFlowCommonHist *fCommonHists; // common control histograms (filled only with events with 3 or more tracks for 3-p correlators) 
   Int_t fnBinsPhi; // number of phi bins
@@ -122,12 +152,21 @@ class AliFlowAnalysisWithNestedLoops
   TH1F *fPhiWeights; // histogram holding phi weights
   TH1D *fPtWeights; // histogram holding phi weights
   TH1D *fEtaWeights; // histogram holding phi weights 
-  // 3.) Final results:
-  TList *fResultsList; // list holding objects with final results 
+  // 3.) Relative angle distribution (RAD):
+  TList *fListRAD; // list holding objects for calculation of relative angle distribution phi1-phi2 
+  Bool_t fEvaluateNestedLoopsForRAD; // evaluate nested loops for relative angle distribution
   TH1D *fRelativeAngleDistribution; // distribution of phi1-phi2 for all distinct pairs of particles
+  // 4.) Debugging and cross-checking QC:
+  // ...
+  // 5.) Debugging and cross-checking MH:
+  TList *fListMH; // list holding objects relevant for debugging and cross-checking of MH class
+  Bool_t fEvaluateNestedLoopsForMH; // evaluate nested loops for mixed harmonics
+  Int_t fCorrelatorIntegerMH; // integer n in cos[n(2phi1-psi2-psi3)]
+  TProfile *f3pCorrelatorVsPtSumDiffDirectPro[2]; // differential 3-p correlator cos[n(2phi1-psi2-psi3)] vs [(p1+p2)/2,|p1-p2|]
+  Int_t fCrossCheckInPtSumBinNo; // print on the screen result of cos[n(2phi1-psi2-psi3)] vs (p1+p2)/2 in this bin number 
+  Int_t fCrossCheckInPtDiffBinNo; // print on the screen result of cos[n(2phi1-psi2-psi3)] vs |p1-p2| in this bin number 
   
   ClassDef(AliFlowAnalysisWithNestedLoops, 0);
-
 };
 
 //================================================================================================================
index 44967b3..b9a43a7 100644 (file)
@@ -76,7 +76,7 @@ class AliAnalysisTaskMixedHarmonics : public AliAnalysisTaskSE{
   Bool_t fUseEtaWeights; // use eta weights  
   TList *fWeightsList; // list with weights
   
-  ClassDef(AliAnalysisTaskMixedHarmonics, 0); 
+  ClassDef(AliAnalysisTaskMixedHarmonics, 1); 
 };
 
 //================================================================================================================
index 23fc41b..eb784c4 100644 (file)
@@ -59,7 +59,7 @@ class AliAnalysisTaskNestedLoops : public AliAnalysisTaskSE{
   Bool_t fUseEtaWeights; // use eta weights  
   TList *fWeightsList; // list with weights
   
-  ClassDef(AliAnalysisTaskNestedLoops, 0); 
+  ClassDef(AliAnalysisTaskNestedLoops, 1); 
 };
 
 //================================================================================================================