]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWGCF/EBYE/macros/drawBalanceFunctionPsi.C
correct normalization for correlation functions (divide by bin width, take into accou...
[u/mrichter/AliRoot.git] / PWGCF / EBYE / macros / drawBalanceFunctionPsi.C
index 3eef2fb15057c3b757b44e4631bc323d4bac5408..8d033a38328f3326be245468b19b65c69e3d908d 100644 (file)
@@ -1,14 +1,20 @@
-const Int_t numberOfCentralityBins = 10;
-TString centralityArray[numberOfCentralityBins] = {"0-10","10-20","20-30","30-40","40-50","50-60","60-70","70-80","0-1","1-2"};
+const Int_t numberOfCentralityBins = 11;
+TString centralityArray[numberOfCentralityBins] = {"0-10","10-20","20-30","30-40","40-50","50-60","60-70","70-80","0-1","1-2","0-100"};
 
 void drawBalanceFunctionPsi(const char* filename = "AnalysisResultsPsi.root", 
                            Int_t gCentrality = 1,
                            Int_t gDeltaEtaDeltaPhi = 2,
+                           Int_t gBit = -1,
+                           const char* gCentralityEstimator = 0x0,
                            Double_t psiMin = -0.5, Double_t psiMax = 0.5,
                            Double_t ptTriggerMin = -1.,
                            Double_t ptTriggerMax = -1.,
                            Double_t ptAssociatedMin = -1.,
-                           Double_t ptAssociatedMax = -1.) {
+                           Double_t ptAssociatedMax = -1.,
+                           Bool_t k2pMethod = kFALSE,
+                           Bool_t k2pMethod2D = kFALSE,
+                           TString eventClass = "EventPlane") //Can be "EventPlane", "Centrality", "Multiplicity"
+{
   //Macro that draws the BF distributions for each centrality bin
   //for reaction plane dependent analysis
   //Author: Panos.Christakoglou@nikhef.nl
@@ -20,10 +26,16 @@ void drawBalanceFunctionPsi(const char* filename = "AnalysisResultsPsi.root",
   gSystem->Load("libPWGTools.so");
   gSystem->Load("libPWGCFebye.so");
 
+  //correction method check
+  if(k2pMethod2D&&!k2pMethod){
+    Printf("Chosen 2D 2particle correction method w/o 2particle correction --> not possible");
+    return;
+  }
+
   //Prepare the objects and return them
-  TList *listBF = GetListOfObjects(filename,gCentrality,0);
-  TList *listBFShuffled = GetListOfObjects(filename,gCentrality,1);
-  TList *listBFMixed = GetListOfObjects(filename,gCentrality,2);
+  TList *listBF = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,0);
+  TList *listBFShuffled = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,1);
+  TList *listBFMixed = GetListOfObjects(filename,gCentrality,gBit,gCentralityEstimator,2);
   if(!listBF) {
     Printf("The TList object was not created");
     return;
@@ -32,33 +44,33 @@ void drawBalanceFunctionPsi(const char* filename = "AnalysisResultsPsi.root",
     draw(listBF,listBFShuffled,listBFMixed,
         gCentrality,gDeltaEtaDeltaPhi,
         psiMin,psiMax,
-        ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);  
+        ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,
+        k2pMethod,k2pMethod2D,eventClass);  
 }
 
 //______________________________________________________//
-TList *GetListOfObjects(const char* filename, 
-                       Int_t gCentrality, Int_t kData = 1) {
+TList *GetListOfObjects(const char* filename,
+                       Int_t gCentrality,
+                       Int_t gBit,
+                       const char *gCentralityEstimator,
+                       Int_t kData = 1) {
   //Get the TList objects (QA, bf, bf shuffled)
-  //kData == 0: data
-  //kData == 1: shuffling
-  //kData == 2: mixing
-  TList *listQA = 0x0;
   TList *listBF = 0x0;
   
   //Open the file
-  TFile *f = TFile::Open(filename);
+  TFile *f = TFile::Open(filename,"UPDATE");
   if((!f)||(!f->IsOpen())) {
     Printf("The file %s is not found. Aborting...",filename);
     return listBF;
   }
-  f->ls();
+  //f->ls();
   
   TDirectoryFile *dir = dynamic_cast<TDirectoryFile *>(f->Get("PWGCFEbyE.outputBalanceFunctionPsiAnalysis"));
   if(!dir) {   
     Printf("The TDirectoryFile is not found. Aborting...",filename);
     return listBF;
   }
-  dir->ls();
+  //dir->ls();
   
   TString listBFName;
   if(kData == 0) {
@@ -74,91 +86,111 @@ TList *GetListOfObjects(const char* filename,
     listBFName = "listBFPsiMixed_";
   }
   listBFName += centralityArray[gCentrality-1];
-  listBFName += "_Bit128_V0M";
-  listBF = dynamic_cast<TList *>(dir->Get(listBFName.Data()));
-  cout<<"======================================================="<<endl;
-  cout<<"List name: "<<listBF->GetName()<<endl;
-  listBF->ls();
-
-  //Get the histograms
-  TString histoName;
-  if(kData == 0)
-    histoName = "fHistPV0M";
-  else if(kData == 1)
-    histoName = "fHistP_shuffleV0M";
-  else if(kData == 2)
-    histoName = "fHistPV0M";
-  AliTHn *fHistP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));  
-  if(!fHistP) {
-    Printf("fHistP %s not found!!!",histoName.Data());
-    break;
-  }
-  fHistP->FillParent(); fHistP->DeleteContainers();
-
-  if(kData == 0)
-    histoName = "fHistNV0M";
-  if(kData == 1)
-    histoName = "fHistN_shuffleV0M";
-  if(kData == 2)
-    histoName = "fHistNV0M";
-  AliTHn *fHistN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
-  if(!fHistN) {
-    Printf("fHistN %s not found!!!",histoName.Data());
-    break;
-  }
-  fHistN->FillParent(); fHistN->DeleteContainers();
+  if(gBit > -1) {
+    listBFName += "_Bit"; listBFName += gBit; }
+  if(gCentralityEstimator) {
+    listBFName += "_"; listBFName += gCentralityEstimator;}
+
+  // histograms were already retrieved (in first iteration)
+  if(dir->Get(Form("%s_histograms",listBFName.Data()))){
+    listBF = dynamic_cast<TList *>(dir->Get(Form("%s_histograms",listBFName.Data())));
+  }
+
+  // histograms were not yet retrieved (this is the first iteration)
+  else{
+
+    listBF = dynamic_cast<TList *>(dir->Get(listBFName.Data()));
+    cout<<"======================================================="<<endl;
+    cout<<"List name (check): "<<listBFName.Data()<<endl;
+    cout<<"List name: "<<listBF->GetName()<<endl;
+    //listBF->ls();
+    
+    //Get the histograms
+    TString histoName;
+    if(kData == 0)
+      histoName = "fHistPV0M";
+    else if(kData == 1)
+      histoName = "fHistP_shuffleV0M";
+    else if(kData == 2)
+      histoName = "fHistPV0M";
+    AliTHn *fHistP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));  
+    if(!fHistP) {
+      Printf("fHistP %s not found!!!",histoName.Data());
+      break;
+    }
+    fHistP->FillParent(); fHistP->DeleteContainers();
+    
+    if(kData == 0)
+      histoName = "fHistNV0M";
+    if(kData == 1)
+      histoName = "fHistN_shuffleV0M";
+    if(kData == 2)
+      histoName = "fHistNV0M";
+    AliTHn *fHistN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
+    if(!fHistN) {
+      Printf("fHistN %s not found!!!",histoName.Data());
+      break;
+    }
+    fHistN->FillParent(); fHistN->DeleteContainers();
     
-  if(kData == 0)
-    histoName = "fHistPNV0M";
-  if(kData == 1)
-    histoName = "fHistPN_shuffleV0M";
-  if(kData == 2)
-    histoName = "fHistPNV0M";
-  AliTHn *fHistPN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
-  if(!fHistPN) {
-    Printf("fHistPN %s not found!!!",histoName.Data());
-    break;
-  }
-  fHistPN->FillParent(); fHistPN->DeleteContainers();
+    if(kData == 0)
+      histoName = "fHistPNV0M";
+    if(kData == 1)
+      histoName = "fHistPN_shuffleV0M";
+    if(kData == 2)
+      histoName = "fHistPNV0M";
+    AliTHn *fHistPN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
+    if(!fHistPN) {
+      Printf("fHistPN %s not found!!!",histoName.Data());
+      break;
+    }
+    fHistPN->FillParent(); fHistPN->DeleteContainers();
+    
+    if(kData == 0)
+      histoName = "fHistNPV0M";
+    if(kData == 1)
+      histoName = "fHistNP_shuffleV0M";
+    if(kData == 2)
+      histoName = "fHistNPV0M";
+    AliTHn *fHistNP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
+    if(!fHistNP) {
+      Printf("fHistNP %s not found!!!",histoName.Data());
+      break;
+    }
+    fHistNP->FillParent(); fHistNP->DeleteContainers();
+    
+    if(kData == 0)
+      histoName = "fHistPPV0M";
+    if(kData == 1)
+      histoName = "fHistPP_shuffleV0M";
+    if(kData == 2)
+      histoName = "fHistPPV0M";
+    AliTHn *fHistPP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
+    if(!fHistPP) {
+      Printf("fHistPP %s not found!!!",histoName.Data());
+      break;
+    }
+    fHistPP->FillParent(); fHistPP->DeleteContainers();
+    
+    if(kData == 0)
+      histoName = "fHistNNV0M";
+    if(kData == 1)
+      histoName = "fHistNN_shuffleV0M";
+    if(kData == 2)
+      histoName = "fHistNNV0M";
+    AliTHn *fHistNN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
+    if(!fHistNN) {
+      Printf("fHistNN %s not found!!!",histoName.Data());
+      break;
+    }
+    fHistNN->FillParent(); fHistNN->DeleteContainers();
+    
+    dir->cd();
+    listBF->Write(Form("%s_histograms",listBFName.Data()), TObject::kSingleKey);
+    
+  }// first iteration
   
-  if(kData == 0)
-    histoName = "fHistNPV0M";
-  if(kData == 1)
-    histoName = "fHistNP_shuffleV0M";
-  if(kData == 2)
-    histoName = "fHistNPV0M";
-  AliTHn *fHistNP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
-  if(!fHistNP) {
-    Printf("fHistNP %s not found!!!",histoName.Data());
-    break;
-  }
-  fHistNP->FillParent(); fHistNP->DeleteContainers();
-
-  if(kData == 0)
-    histoName = "fHistPPV0M";
-  if(kData == 1)
-    histoName = "fHistPP_shuffleV0M";
-  if(kData == 2)
-    histoName = "fHistPPV0M";
-  AliTHn *fHistPP = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
-  if(!fHistPP) {
-    Printf("fHistPP %s not found!!!",histoName.Data());
-    break;
-  }
-  fHistPP->FillParent(); fHistPP->DeleteContainers();
-
-  if(kData == 0)
-    histoName = "fHistNNV0M";
-  if(kData == 1)
-    histoName = "fHistNN_shuffleV0M";
-  if(kData == 2)
-    histoName = "fHistNNV0M";
-  AliTHn *fHistNN = dynamic_cast<AliTHn *>(listBF->FindObject(histoName.Data()));
-  if(!fHistNN) {
-    Printf("fHistNN %s not found!!!",histoName.Data());
-    break;
-  }
-  fHistNN->FillParent(); fHistNN->DeleteContainers();
+  f->Close();
   
   return listBF;
 }
@@ -168,7 +200,8 @@ void draw(TList *listBF, TList *listBFShuffled, TList *listBFMixed,
          Int_t gCentrality, Int_t gDeltaEtaDeltaPhi, 
          Double_t psiMin, Double_t psiMax,
          Double_t ptTriggerMin, Double_t ptTriggerMax,
-         Double_t ptAssociatedMin, Double_t ptAssociatedMax) {
+         Double_t ptAssociatedMin, Double_t ptAssociatedMax,
+         Bool_t k2pMethod = kFALSE,Bool_t k2pMethod2D = kFALSE, TString eventClass="EventPlane") {
   gROOT->LoadMacro("~/SetPlotStyle.C");
   SetPlotStyle();
   gStyle->SetPalette(1,0);
@@ -192,6 +225,7 @@ void draw(TList *listBF, TList *listBFShuffled, TList *listBFMixed,
   hNN = (AliTHn*) listBF->FindObject("fHistNNV0M");
 
   AliBalancePsi *b = new AliBalancePsi();
+  b->SetEventClass(eventClass);
   b->SetHistNp(hP);
   b->SetHistNn(hN);
   b->SetHistNpn(hPN);
@@ -215,6 +249,7 @@ void draw(TList *listBF, TList *listBFShuffled, TList *listBFMixed,
   hNNShuffled = (AliTHn*) listBFShuffled->FindObject("fHistNN_shuffleV0M");
 
   AliBalancePsi *bShuffled = new AliBalancePsi();
+  bShuffled->SetEventClass(eventClass);
   bShuffled->SetHistNp(hPShuffled);
   bShuffled->SetHistNn(hNShuffled);
   bShuffled->SetHistNpn(hPNShuffled);
@@ -238,6 +273,7 @@ void draw(TList *listBF, TList *listBFShuffled, TList *listBFMixed,
   hNNMixed = (AliTHn*) listBFMixed->FindObject("fHistNNV0M");
 
   AliBalancePsi *bMixed = new AliBalancePsi();
+  bMixed->SetEventClass(eventClass);
   bMixed->SetHistNp(hPMixed);
   bMixed->SetHistNn(hNMixed);
   bMixed->SetHistNpn(hPNMixed);
@@ -252,32 +288,114 @@ void draw(TList *listBF, TList *listBFShuffled, TList *listBFMixed,
   TString histoTitle, pngName;
   TLegend *legend;
   
-  histoTitle = "Centrality: "; 
-  histoTitle += centralityArray[gCentrality-1]; 
-  histoTitle += "%";
-  if((psiMin == -0.5)&&(psiMax == 0.5))
-    histoTitle += " (-7.5^{o} < #phi - #Psi_{2} < 7.5^{o})"; 
-  else if((psiMin == 0.5)&&(psiMax == 1.5))
-    histoTitle += " (37.5^{o} < #phi - #Psi_{2} < 52.5^{o})"; 
-  else if((psiMin == 1.5)&&(psiMax == 2.5))
-    histoTitle += " (82.5^{o} < #phi - #Psi_{2} < 97.5^{o})"; 
-  else 
-    histoTitle += " (0^{o} < #phi - #Psi_{2} < 180^{o})"; 
-  
+  if(eventClass == "Centrality"){
+    histoTitle = "Centrality: ";
+    histoTitle += psiMin;
+    histoTitle += " - ";
+    histoTitle += psiMax;
+    histoTitle += " % ";
+    histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})"; 
+  }
+  else if(eventClass == "Multiplicity"){
+    histoTitle = "Multiplicity: ";
+    histoTitle += psiMin;
+    histoTitle += " - ";
+    histoTitle += psiMax;
+    histoTitle += " tracks";
+    histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})"; 
+  }
+  else{ // "EventPlane" (default)
+    histoTitle = "Centrality: ";
+    histoTitle += centralityArray[gCentrality-1]; 
+    histoTitle += "%";
+    if((psiMin == -0.5)&&(psiMax == 0.5))
+      histoTitle += " (-7.5^{o} < #varphi - #Psi_{2} < 7.5^{o})"; 
+    else if((psiMin == 0.5)&&(psiMax == 1.5))
+      histoTitle += " (37.5^{o} < #varphi - #Psi_{2} < 52.5^{o})"; 
+    else if((psiMin == 1.5)&&(psiMax == 2.5))
+      histoTitle += " (82.5^{o} < #varphi - #Psi_{2} < 97.5^{o})"; 
+    else 
+      histoTitle += " (0^{o} < #varphi - #Psi_{2} < 180^{o})"; 
+  }
+
   //Raw balance function
-  gHistBalanceFunction = b->GetBalanceFunctionHistogram(0,gDeltaEtaDeltaPhi,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
+  if(k2pMethod){ 
+    if(bMixed){
+      gHistBalanceFunction = b->GetBalanceFunctionHistogram2pMethod(0,gDeltaEtaDeltaPhi,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
+    }
+    else{
+      cerr<<"RAW: NO MIXED BF BUT REQUESTED CORRECTING WITH IT! --> FAIL"<<endl;
+      return;
+    }
+  }
+  else if(k2pMethod2D){ 
+    if(bMixed){
+      if(gDeltaEtaDeltaPhi==1) //Delta eta
+       gHistBalanceFunction = b->GetBalanceFunction1DFrom2D2pMethod(0,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
+      else //Delta phi
+       gHistBalanceFunction = b->GetBalanceFunction1DFrom2D2pMethod(1,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
+    }
+    else{
+      cerr<<"RAW: NO MIXED BF BUT REQUESTED CORRECTING WITH IT! --> FAIL"<<endl;
+      return;
+    }
+  }
+  else
+    gHistBalanceFunction = b->GetBalanceFunctionHistogram(0,gDeltaEtaDeltaPhi,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
   gHistBalanceFunction->SetMarkerStyle(20);
   gHistBalanceFunction->SetTitle(histoTitle.Data());
   gHistBalanceFunction->GetYaxis()->SetTitleOffset(1.3);
   gHistBalanceFunction->SetName("gHistBalanceFunction");
   
   //Shuffled balance function
-  gHistBalanceFunctionShuffled = bShuffled->GetBalanceFunctionHistogram(0,gDeltaEtaDeltaPhi,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
+  if(k2pMethod){ 
+    if(bMixed)
+      gHistBalanceFunctionShuffled = bShuffled->GetBalanceFunctionHistogram2pMethod(0,gDeltaEtaDeltaPhi,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
+    else{
+      cerr<<"SHUFFLE: NO MIXED BF BUT REQUESTED CORRECTING WITH IT! --> FAIL"<<endl;
+      return;
+    }
+  }
+  else if(k2pMethod2D){ 
+    if(bMixed){
+      if(gDeltaEtaDeltaPhi==1) //Delta eta
+       gHistBalanceFunctionShuffled = bShuffled->GetBalanceFunction1DFrom2D2pMethod(0,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
+      else //Delta phi
+       gHistBalanceFunctionShuffled = bShuffled->GetBalanceFunction1DFrom2D2pMethod(1,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
+    }
+    else{
+      cerr<<"SHUFFLE: NO MIXED BF BUT REQUESTED CORRECTING WITH IT! --> FAIL"<<endl;
+      return;
+    }
+  }
+  else
+    gHistBalanceFunctionShuffled = bShuffled->GetBalanceFunctionHistogram(0,gDeltaEtaDeltaPhi,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
   gHistBalanceFunctionShuffled->SetMarkerStyle(24);
   gHistBalanceFunctionShuffled->SetName("gHistBalanceFunctionShuffled");
 
   //Mixed balance function
-  gHistBalanceFunctionMixed = bMixed->GetBalanceFunctionHistogram(0,gDeltaEtaDeltaPhi,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
+  if(k2pMethod){ 
+    if(bMixed)
+      gHistBalanceFunctionMixed = bMixed->GetBalanceFunctionHistogram2pMethod(0,gDeltaEtaDeltaPhi,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
+    else{
+      cerr<<"MIXED: NO MIXED BF BUT REQUESTED CORRECTING WITH IT! --> FAIL"<<endl;
+      return;
+    }
+  }
+  else if(k2pMethod2D){ 
+    if(bMixed){
+      if(gDeltaEtaDeltaPhi==1) //Delta eta
+       gHistBalanceFunctionMixed = bMixed->GetBalanceFunction1DFrom2D2pMethod(0,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
+      else //Delta phi
+       gHistBalanceFunctionMixed = bMixed->GetBalanceFunction1DFrom2D2pMethod(1,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax,bMixed);
+    }
+    else{
+      cerr<<"MIXED: NO MIXED BF BUT REQUESTED CORRECTING WITH IT! --> FAIL"<<endl;
+      return;
+    }
+  }
+  else
+    gHistBalanceFunctionMixed = bMixed->GetBalanceFunctionHistogram(0,gDeltaEtaDeltaPhi,psiMin,psiMax,ptTriggerMin,ptTriggerMax,ptAssociatedMin,ptAssociatedMax);
   gHistBalanceFunctionMixed->SetMarkerStyle(25);
   gHistBalanceFunctionMixed->SetName("gHistBalanceFunctionMixed");
   
@@ -285,6 +403,7 @@ void draw(TList *listBF, TList *listBFShuffled, TList *listBFMixed,
   gHistBalanceFunctionSubtracted = dynamic_cast<TH1D *>(gHistBalanceFunction->Clone());
   gHistBalanceFunctionSubtracted->Add(gHistBalanceFunctionMixed,-1);
   gHistBalanceFunctionSubtracted->Rebin(gRebin);
+  gHistBalanceFunctionSubtracted->Scale(1./(Double_t)(gRebin));    
   gHistBalanceFunctionSubtracted->SetMarkerStyle(20);
   gHistBalanceFunctionSubtracted->SetTitle(histoTitle.Data());
   gHistBalanceFunctionSubtracted->GetYaxis()->SetTitleOffset(1.3);
@@ -311,10 +430,36 @@ void draw(TList *listBF, TList *listBFShuffled, TList *listBFMixed,
   legend->AddEntry(gHistBalanceFunctionMixed,"Mixed data","lp");
   legend->Draw();
   
-  pngName = "BalanceFunctionDeltaEta.Centrality"; 
-  pngName += centralityArray[gCentrality-1]; 
-  pngName += ".Psi"; //pngName += psiMin; pngName += "To"; pngName += psiMax;
+  pngName = "BalanceFunction."; 
+  if(eventClass == "Centrality"){
+    pngName += Form("Centrality%.1fTo%.1f",psiMin,psiMax);
+    if(gDeltaEtaDeltaPhi == 1) pngName += ".InDeltaEta.PsiAll.PttFrom";
+    else if(gDeltaEtaDeltaPhi == 2) pngName += ".InDeltaPhi.PsiAll.PttFrom";
+  }
+  else if(eventClass == "Multiplicity"){
+    pngName += Form("Multiplicity%.0fTo%.0f",psiMin,psiMax);
+    if(gDeltaEtaDeltaPhi == 1) pngName += ".InDeltaEta.PsiAll.PttFrom";
+    else if(gDeltaEtaDeltaPhi == 2) pngName += ".InDeltaPhi.PsiAll.PttFrom";  
+  }
+  else{ // "EventPlane" (default)
+    pngName += "Centrality";
+    pngName += gCentrality; 
+    if(gDeltaEtaDeltaPhi == 1) pngName += ".InDeltaEta.Psi";
+    else if(gDeltaEtaDeltaPhi == 2) pngName += ".InDeltaPhi.Psi";
+    if((psiMin == -0.5)&&(psiMax == 0.5)) pngName += "InPlane.Ptt";
+    else if((psiMin == 0.5)&&(psiMax == 1.5)) pngName += "Intermediate.Ptt";
+    else if((psiMin == 1.5)&&(psiMax == 2.5)) pngName += "OutOfPlane.Ptt";
+    else if((psiMin == 2.5)&&(psiMax == 3.5)) pngName += "Rest.PttFrom";
+    else pngName += "All.PttFrom";
+  }  
+  pngName += Form("%.1f",ptTriggerMin); pngName += "To"; 
+  pngName += Form("%.1f",ptTriggerMax); pngName += "PtaFrom";
+  pngName += Form("%.1f",ptAssociatedMin); pngName += "To"; 
+  pngName += Form("%.1f",ptAssociatedMax); 
+  if(k2pMethod2D) pngName += "_2pMethod2D";
+  else if(k2pMethod) pngName += "_2pMethod";
   pngName += ".png";
+
   c1->SaveAs(pngName.Data());
   
   GetWeightedMean(gHistBalanceFunction);
@@ -360,19 +505,34 @@ void draw(TList *listBF, TList *listBFShuffled, TList *listBFMixed,
   latex->DrawLatex(0.64,0.77,skewnessLatex.Data());
   latex->DrawLatex(0.64,0.73,kurtosisLatex.Data());
 
-  TString newFileName = "balanceFunction.Centrality"; 
-  newFileName += gCentrality; newFileName += ".In";
-  if(gDeltaEtaDeltaPhi == 1) newFileName += "DeltaEta.Psi";
-  else if(gDeltaEtaDeltaPhi == 2) newFileName += "DeltaPhi.Psi";
-  if((psiMin == -0.5)&&(psiMax == 0.5)) newFileName += "InPlane.Ptt";
-  else if((psiMin == 0.5)&&(psiMax == 1.5)) newFileName += "Intermediate.Ptt";
-  else if((psiMin == 1.5)&&(psiMax == 2.5)) newFileName += "OutOfPlane.Ptt";
-  else if((psiMin == 2.5)&&(psiMax == 3.5)) newFileName += "Rest.PttFrom";
-  else newFileName += "All.PttFrom";
+  TString newFileName = "balanceFunction."; 
+  if(eventClass == "Centrality"){
+    newFileName += Form("Centrality%.1fTo%.1f",psiMin,psiMax);
+    if(gDeltaEtaDeltaPhi == 1) newFileName += ".InDeltaEta.PsiAll.PttFrom";
+    else if(gDeltaEtaDeltaPhi == 2) newFileName += ".InDeltaPhi.PsiAll.PttFrom";
+  }
+  else if(eventClass == "Multiplicity"){
+    newFileName += Form("Multiplicity%.0fTo%.0f",psiMin,psiMax);
+    if(gDeltaEtaDeltaPhi == 1) newFileName += ".InDeltaEta.PsiAll.PttFrom";
+    else if(gDeltaEtaDeltaPhi == 2) newFileName += ".InDeltaPhi.PsiAll.PttFrom";  
+  }
+  else{ // "EventPlane" (default)
+    newFileName += "Centrality";
+    newFileName += gCentrality; 
+    if(gDeltaEtaDeltaPhi == 1) newFileName += ".InDeltaEta.Psi";
+    else if(gDeltaEtaDeltaPhi == 2) newFileName += ".InDeltaPhi.Psi";
+    if((psiMin == -0.5)&&(psiMax == 0.5)) newFileName += "InPlane.Ptt";
+    else if((psiMin == 0.5)&&(psiMax == 1.5)) newFileName += "Intermediate.Ptt";
+    else if((psiMin == 1.5)&&(psiMax == 2.5)) newFileName += "OutOfPlane.Ptt";
+    else if((psiMin == 2.5)&&(psiMax == 3.5)) newFileName += "Rest.PttFrom";
+    else newFileName += "All.PttFrom";
+  }  
   newFileName += Form("%.1f",ptTriggerMin); newFileName += "To"; 
   newFileName += Form("%.1f",ptTriggerMax); newFileName += "PtaFrom";
   newFileName += Form("%.1f",ptAssociatedMin); newFileName += "To"; 
   newFileName += Form("%.1f",ptAssociatedMax); 
+  if(k2pMethod2D) newFileName += "_2pMethod2D";
+  else if(k2pMethod) newFileName += "_2pMethod";
   newFileName += ".root";
 
   TFile *fOutput = new TFile(newFileName.Data(),"recreate");
@@ -397,15 +557,15 @@ void GetWeightedMean(TH1D *gHistBalance, Int_t fStartBin = 1) {
   Int_t fNumberOfBins = gHistBalance->GetNbinsX();
   Double_t fP2Step    = gHistBalance->GetBinWidth(1); // assume equal binning!
   
-  cout<<"=================================================="<<endl;
-  cout<<"RECALCULATION OF BF WIDTH (StartBin = "<<fStartBin<<")"<<endl;
-  cout<<"HISTOGRAM has "<<fNumberOfBins<<" bins with bin size of "<<fP2Step<<endl;
-  cout<<"=================================================="<<endl;
+  //cout<<"=================================================="<<endl;
+  //cout<<"RECALCULATION OF BF WIDTH (StartBin = "<<fStartBin<<")"<<endl;
+  //cout<<"HISTOGRAM has "<<fNumberOfBins<<" bins with bin size of "<<fP2Step<<endl;
+  //cout<<"=================================================="<<endl;
   for(Int_t i = 1; i <= fNumberOfBins; i++) {
     // this is to simulate |Delta eta| or |Delta phi|
     if(fNumberOfBins/2 - fStartBin + 1 < i && i < fNumberOfBins/2 + fStartBin ) continue;
 
-    cout<<"B: "<<gHistBalance->GetBinContent(i)<<"\t Error: "<<gHistBalance->GetBinError(i)<<"\t bin: "<<TMath::Abs(gHistBalance->GetBinCenter(i))<<endl;
+    //cout<<"B: "<<gHistBalance->GetBinContent(i)<<"\t Error: "<<gHistBalance->GetBinError(i)<<"\t bin: "<<TMath::Abs(gHistBalance->GetBinCenter(i))<<endl;
 
     gSumXi += TMath::Abs(gHistBalance->GetBinCenter(i)); // this is to simulate |Delta eta| or |Delta phi|
     gSumBi += gHistBalance->GetBinContent(i);