PWGCF/EBYE/macros/readBalanceFunction.C

   1 const TString gBFAnalysisType[7] = {"y","eta","qlong","qout","qside","qinv","phi"};
   2
   3 const Int_t nrOfCentralities = 9;
   4 const Double_t centralityArray[nrOfCentralities+1] = {0.,5.,10.,20.,30.,40.,50.,60.,70.,80.};  // in centrality percentile (0-5,5-10,10-20,20-30,...,70-80)
   5 //const Double_t centralityArray[nrOfCentralities+1] = {0.,1.,2.,3.,4.,6.,10.,20.,30.,80.};  // in centrality percentile (0-5,5-10,10-20,20-30,...,70-80)
   6 const Double_t cent[nrOfCentralities]  = {382.8,329.7,260.5,186.4,128.9,85.,52.8,30.,15.8};   // hard coded at the moment for centrality percentiles
   7 const Double_t centE[nrOfCentralities] = {3.1,4.6,4.4,3.9,3.3,2.6,2.0,1.3,0.6};               // (0-5,5-10,10-20,20-30,...,70-80)
   8
   9 void readBalanceFunction(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResults.root",Int_t fStartBinBFWidth = 3, Int_t fRebin = 2,Int_t fStartBinBFWidthPhi = 2, Int_t fRebinPhi = 2,TString centEst = "V0M",Double_t etaWindow = -1) {
  10   // Macro to read the output of the BF analysis:  MW: CHANGE THIS!!!!
  11   //i) Prints and draws the final BF output
  12   //ii) Plots the QA part of the analysis
  13   //iii) store BF in output file
  14   //Author: Panos.Christakoglou@cern.ch, m.weber@cern.ch
  15   //Loading the needed libraries
  16   gSystem->Load("libProofPlayer.so");
  17   gSystem->Load("libANALYSIS.so");
  18   gSystem->Load("libANALYSISalice.so");
  19   gSystem->Load("libPWGCFebye.so");
  20
  21   //Draw BF
  22   drawBF(bHistos,inFile, fStartBinBFWidth, fRebin,fStartBinBFWidthPhi, fRebinPhi,centEst, "",  etaWindow);
  23
  24
  25   //Merge the output
  26   //mergeOutput("/alice/cern.ch/user/p/pchrist/Balance/pp/7TeV/LHC10b/output/");
  27 }
  28
  29 //___________________________________________________________//
  30 void drawBF(Bool_t bHistos = kFALSE, TString inFile = "AnalysisResults.root", Int_t fStartBinBFWidth = 1, Int_t fRebin = 1, Int_t fStartBinBFWidthPhi = 1, Int_t fRebinPhi = 1, TString centEst = "V0M",TString extraString = "", Double_t etaWindow = -1) {
  31   //Function to draw the BF objects and write them into the output file
  32
  33   Int_t maximumCanvases = 13;
  34   Int_t iCanvas         = 0;
  35   Int_t iList           = -1;
  36   TCanvas *cQA[13][10];
  37   TCanvas *cQAV0M = new TCanvas("cQAV0M","V0M multiplicities");
  38   cQAV0M->Divide(4,3);
  39   TCanvas *cQARef = new TCanvas("cQARef","reference track multiplicities");
  40   cQARef->Divide(4,3);
  41   TCanvas *cBF[13][10];
  42   TCanvas *cBFS[13][10];
  43
  44   // get the file
  45   TFile *f = TFile::Open(inFile.Data());
  46   if(!f) {
  47     Printf("File not found!!!");
  48     break;
  49   }
  50
  51   // get the BF output directory
  52   TDirectoryFile *dir = dynamic_cast<TDirectoryFile *>(f->Get("PWGCFEbyE.outputBalanceFunctionAnalysis"));
  53   if(!dir) {
  54     Printf("Output directory not found!!!");
  55     break;
  56   }
  57
  58   // loop over all lists and plot the BF and QA
  59   TList *list = NULL;
  60   TString listName;
  61   TIter nextkey( dir->GetListOfKeys() );
  62   TKey *key;
  63
  64   AliBalance *bf[13][7];
  65   AliBalance *bfs[13][7];
  66   TH1D *gbf[13][10][7];
  67   TH1D *gbfs[13][10][7];
  68
  69   for(Int_t i = 0; i < 13; i++){
  70     for(Int_t j = 0; j < 10; j++){
  71       for(Int_t k = 0; k < 7; k++){
  72         gbf[i][j][k]  = NULL;
  73         gbfs[i][j][k] = NULL;
  74       }
  75     }
  76   }
  77
  78   TH2D *fHistP[13][7]; //N+
  79   TH2D *fHistN[13][7]; //N-
  80   TH2D *fHistPN[13][7]; //N+-
  81   TH2D *fHistNP[13][7]; //N-+
  82   TH2D *fHistPP[13][7]; //N++
  83   TH2D *fHistNN[13][7]; //N--
  84
  85   TH2D *fHistPS[13][7]; //N+
  86   TH2D *fHistNS[13][7]; //N-
  87   TH2D *fHistPNS[13][7]; //N+-
  88   TH2D *fHistNPS[13][7]; //N-+
  89   TH2D *fHistPPS[13][7]; //N++
  90   TH2D *fHistNNS[13][7]; //N--
  91
  92   Double_t WM[13][10];     // weighted mean for eta (recalculated from fStartBin)
  93   Double_t WME[13][10];    // error
  94   Double_t WMS[13][10];     // weighted mean for eta (recalculated from fStartBin) (shuffled)
  95   Double_t WMSE[13][10];    // error (shuffled)
  96
  97   Double_t WMP[13][10];     // weighted mean for phi (recalculated from fStartBinPhi)
  98   Double_t WMPE[13][10];    // error
  99   Double_t WMPS[13][10];     // weighted mean for phi (recalculated from fStartBin) (shuffled)
 100   Double_t WMPSE[13][10];    // error (shuffled)
 101
 102   Double_t integ[13][10];     // integral for eta (calculated from bin 1)
 103   Double_t integE[13][10];    // error
 104   Double_t integS[13][10];     // integral for eta (calculated from bin 1) (shuffled)
 105   Double_t integSE[13][10];    // error (shuffled)
 106
 107   Double_t integP[13][10];     // integral for phi (calculated from bin 1)
 108   Double_t integPE[13][10];    // error
 109   Double_t integPS[13][10];     // integral for phi (calculated from bin 1) (shuffled)
 110   Double_t integPSE[13][10];    // error (shuffled)
 111
 112
 113   while ( (key = (TKey*)nextkey())) {
 114
 115     list = (TList*)key->ReadObj();
 116     listName = TString(list->GetName());
 117
 118     cout<<"Processing list "<<listName<<endl;
 119
 120
 121     // ----------------------------------------------------
 122     // plot QA histograms
 123     if(listName.Contains("QA")){
 124
 125       if(iList<13) iList++;
 126       else{
 127         cerr<<"TOO MANY LISTS!!!"<<endl;
 128         return;
 129       }
 130
 131       cQA[iList][iCanvas] = new TCanvas(listName,listName);
 132       cQA[iList][iCanvas]->Divide(4,3);
 133
 134      cQA[iList][iCanvas]->cd(1);
 135       TH1F* histEventStats = (TH1F*)list->FindObject("fHistEventStats");
 136       if(histEventStats){
 137         histEventStats->SetFillColor(9);
 138         histEventStats->Draw();
 139       }
 140
 141       cQA[iList][iCanvas]->cd(2);
 142       TH1F* histTriggerStats = (TH1F*)list->FindObject("fHistTriggerStats");
 143       if(histTriggerStats){
 144         histTriggerStats->SetFillColor(9);
 145         histTriggerStats->Draw();
 146       }
 147
 148      cQA[iList][iCanvas]->cd(3);
 149       TH1F* histTrackStats = (TH1F*)list->FindObject("fHistTrackStats");
 150       if(histTrackStats){
 151         histTrackStats->SetFillColor(9);
 152         histTrackStats->Draw();
 153       }
 154
 155       cQA[iList][iCanvas]->cd(4);
 156       TH1F* histCentStats = (TH1F*)list->FindObject("fHistCentStats");
 157       if(histCentStats){
 158         histCentStats->SetFillColor(9);
 159         histCentStats->Draw("colz");
 160       }
 161
 162
 163
 164       cQA[iList][iCanvas]->cd(5);
 165       TH1F* histVx = (TH1F*)list->FindObject("fHistVx");
 166       if(histVx){
 167         histVx->SetFillColor(9);
 168         histVx->Draw();
 169       }
 170       cQA[iList][iCanvas]->cd(6);
 171       TH1F* histVy = (TH1F*)list->FindObject("fHistVy");
 172       if(histVy){
 173         histVy->SetFillColor(9);
 174         histVy->Draw();
 175       }
 176
 177       cQA[iList][iCanvas]->cd(7);
 178       TH1F* histVz = (TH1F*)list->FindObject("fHistVz");
 179       if(histVz){
 180         histVz->SetFillColor(9);
 181         histVz->Draw();
 182       }
 183
 184       cQA[iList][iCanvas]->cd(8);
 185       cQA[iList][iCanvas]->cd(8)->SetLogz();
 186       TH2F* histDCA = (TH2F*)list->FindObject("fHistDCA");
 187       if(histDCA) histDCA->Draw("colz");
 188
 189
 190       cQA[iList][iCanvas]->cd(9);
 191       cQA[iList][iCanvas]->cd(9)->SetLogz();
 192       TH2F* histClus = (TH2F*)list->FindObject("fHistClus");
 193       if(histClus) histClus->Draw("colz");
 194
 195       cQA[iList][iCanvas]->cd(10);
 196       TH1F* histPt = (TH1F*)list->FindObject("fHistPt");
 197       if(histPt){
 198         histPt->SetFillColor(9);
 199         histPt->Draw();
 200       }
 201
 202       cQA[iList][iCanvas]->cd(11);
 203       TH1F* histEta = (TH1F*)list->FindObject("fHistEta");
 204       if(histEta){
 205         histEta->SetFillColor(9);
 206         histEta->Draw();
 207       }
 208
 209       cQA[iList][iCanvas]->cd(12);
 210       TH1F* histPhi = (TH1F*)list->FindObject("fHistPhi");
 211       if(histPhi){
 212         histPhi->SetFillColor(9);
 213         histPhi->Draw();
 214       }
 215
 216       // centrality estimator QA
 217       cQAV0M->cd(iCanvas+1);
 218       cQAV0M->cd(iCanvas+1)->SetLogz();
 219       TH1F* histV0M = (TH1F*)list->FindObject("fHistV0M");
 220       if(histV0M){
 221         histV0M->Draw("colz");
 222       }
 223
 224       cQARef->cd(iCanvas+1);
 225       cQARef->cd(iCanvas+1)->SetLogz();
 226       TH1F* histRef = (TH1F*)list->FindObject("fHistRefTracks");
 227       if(histRef){
 228         histRef->Draw("colz");
 229       }
 230     }
 231     // ----------------------------------------------------
 232
 233     // ----------------------------------------------------
 234     // calculate and plot BF
 235     if(listName.Contains("BF_")&&listName.Contains(centEst.Data())){
 236
 237       for(iCanvas = 0; iCanvas < nrOfCentralities; iCanvas++){
 238
 239         cBF[iList][iCanvas] = new TCanvas(Form("cBF_%d_%d",iList,iCanvas),Form("cBF_%d_%d",iList,iCanvas));
 240         cBF[iList][iCanvas]->Divide(3,3);
 241
 242       }
 243
 244
 245       for(Int_t a = 0; a < 7; a++){
 246
 247         cout<<"ANALYSE "<<gBFAnalysisType[a]<<endl;
 248
 249         // create the BF object
 250         bf[iList][a]  = new AliBalance();
 251
 252         fHistP[iList][a] = (TH2D*)list->FindObject(Form("fHistP%s%s",gBFAnalysisType[a].Data(),centEst.Data()));
 253         fHistN[iList][a] = (TH2D*)list->FindObject(Form("fHistN%s%s",gBFAnalysisType[a].Data(),centEst.Data()));
 254         fHistPP[iList][a] = (TH2D*)list->FindObject(Form("fHistPP%s%s",gBFAnalysisType[a].Data(),centEst.Data()));
 255         fHistPN[iList][a] = (TH2D*)list->FindObject(Form("fHistPN%s%s",gBFAnalysisType[a].Data(),centEst.Data()));
 256         fHistNP[iList][a] = (TH2D*)list->FindObject(Form("fHistNP%s%s",gBFAnalysisType[a].Data(),centEst.Data()));
 257         fHistNN[iList][a] = (TH2D*)list->FindObject(Form("fHistNN%s%s",gBFAnalysisType[a].Data(),centEst.Data()));
 258
 259         // rebin histograms (be careful with divider!)
 260         if(a==6){
 261           fHistP[iList][a]->RebinY(fRebinPhi);
 262           fHistN[iList][a]->RebinY(fRebinPhi);
 263           fHistPP[iList][a]->RebinY(fRebinPhi);
 264           fHistPN[iList][a]->RebinY(fRebinPhi);
 265           fHistNP[iList][a]->RebinY(fRebinPhi);
 266           fHistNN[iList][a]->RebinY(fRebinPhi);
 267         }
 268         else{
 269           fHistP[iList][a]->RebinY(fRebin);
 270           fHistN[iList][a]->RebinY(fRebin);
 271           fHistPP[iList][a]->RebinY(fRebin);
 272           fHistPN[iList][a]->RebinY(fRebin);
 273           fHistNP[iList][a]->RebinY(fRebin);
 274           fHistNN[iList][a]->RebinY(fRebin);
 275         }
 276
 277         fHistP[iList][a]->SetName(Form("%s_%d_%d",fHistP[iList][a]->GetName(),iList,iCanvas));
 278         fHistN[iList][a]->SetName(Form("%s_%d_%d",fHistN[iList][a]->GetName(),iList,iCanvas));
 279         fHistPP[iList][a]->SetName(Form("%s_%d_%d",fHistPP[iList][a]->GetName(),iList,iCanvas));
 280         fHistPN[iList][a]->SetName(Form("%s_%d_%d",fHistPN[iList][a]->GetName(),iList,iCanvas));
 281         fHistNP[iList][a]->SetName(Form("%s_%d_%d",fHistNP[iList][a]->GetName(),iList,iCanvas));
 282         fHistNN[iList][a]->SetName(Form("%s_%d_%d",fHistNN[iList][a]->GetName(),iList,iCanvas));
 283
 284         // set histograms in AliBalance object
 285         bf[iList][a]->SetHistNp(a, fHistP[iList][a]);
 286         bf[iList][a]->SetHistNn(a, fHistN[iList][a]);
 287         bf[iList][a]->SetHistNpp(a, fHistPP[iList][a]);
 288         bf[iList][a]->SetHistNpn(a, fHistPN[iList][a]);
 289         bf[iList][a]->SetHistNnp(a, fHistNP[iList][a]);
 290         bf[iList][a]->SetHistNnn(a, fHistNN[iList][a]);
 291
 292         for(iCanvas = 0; iCanvas < nrOfCentralities; iCanvas++){
 293
 294           gbf[iList][iCanvas][a] = bf[iList][a]->GetBalanceFunctionHistogram(a,centralityArray[iCanvas],centralityArray[iCanvas+1],etaWindow);
 295           gbf[iList][iCanvas][a]->SetName(Form("BF_%s_Cent_%.0f_%.0f_%d",gBFAnalysisType[a].Data(),centralityArray[iCanvas],centralityArray[iCanvas+1],iList));
 296
 297           cBF[iList][iCanvas]->cd(a+1);
 298           gbf[iList][iCanvas][a]->SetMarkerStyle(20);
 299
 300           if(!bHistos){
 301             gbf[iList][iCanvas][a]->DrawCopy("AP");
 302             if(a==1){
 303               GetWeightedMean(gbf[iList][iCanvas][a],fStartBinBFWidth,WM[iList][iCanvas],WME[iList][iCanvas]); // for eta recalculate width
 304               GetIntegral(gbf[iList][iCanvas][a],integ[iList][iCanvas],integE[iList][iCanvas]);
 305             }
 306             else if(a==6){
 307               GetWeightedMean(gbf[iList][iCanvas][a],fStartBinBFWidthPhi,WMP[iList][iCanvas],WMPE[iList][iCanvas]); // for phi calculate width
 308               GetIntegral(gbf[iList][iCanvas][a],integP[iList][iCanvas],integPE[iList][iCanvas]);
 309             }
 310           }
 311           else{
 312             fHistPN[iList][a]->SetLineColor(2);
 313             fHistPN[iList][a]->ProjectionY(Form("pn%d",a))->DrawCopy();
 314             fHistPP[iList][a]->SetLineColor(1);
 315             fHistPP[iList][a]->ProjectionY(Form("pp%d",a))->DrawCopy("same");
 316             fHistNP[iList][a]->SetLineColor(4);
 317             fHistNP[iList][a]->ProjectionY(Form("np%d",a))->DrawCopy("same");
 318             fHistNN[iList][a]->SetLineColor(8);
 319             fHistNN[iList][a]->ProjectionY(Form("nn%d",a))->DrawCopy("same");
 320           }
 321         }
 322       }
 323     }
 324     // ----------------------------------------------------
 325
 326     // ----------------------------------------------------
 327     // calculate and plot BF (shuffled)
 328     if(listName.Contains("BFShuffled")&&listName.Contains(centEst.Data())&&listName.Contains(extraString.Data())){
 329
 330       for(iCanvas = 0; iCanvas < nrOfCentralities; iCanvas++){
 331
 332         cBFS[iList][iCanvas] = new TCanvas(Form("Shuffled_%d_%d",iList,iCanvas),Form("Shuffled_%d_%d",iList,iCanvas));
 333         cBFS[iList][iCanvas]->Divide(3,3);
 334
 335       }
 336
 337       for(Int_t a = 0; a < 7; a++){
 338
 339         // create the BF object
 340         bfs[iList][a]  = new AliBalance();
 341
 342         fHistPS[iList][a] = (TH2D*)list->FindObject(Form("fHistP%s_shuffle%s",gBFAnalysisType[a].Data(),centEst.Data()));
 343         fHistNS[iList][a] = (TH2D*)list->FindObject(Form("fHistN%s_shuffle%s",gBFAnalysisType[a].Data(),centEst.Data()));
 344         fHistPPS[iList][a] = (TH2D*)list->FindObject(Form("fHistPP%s_shuffle%s",gBFAnalysisType[a].Data(),centEst.Data()));
 345         fHistPNS[iList][a] = (TH2D*)list->FindObject(Form("fHistPN%s_shuffle%s",gBFAnalysisType[a].Data(),centEst.Data()));
 346         fHistNPS[iList][a] = (TH2D*)list->FindObject(Form("fHistNP%s_shuffle%s",gBFAnalysisType[a].Data(),centEst.Data()));
 347         fHistNNS[iList][a] = (TH2D*)list->FindObject(Form("fHistNN%s_shuffle%s",gBFAnalysisType[a].Data(),centEst.Data()));
 348
 349         // rebin histograms (be careful with divider!)
 350         if(a==6){
 351           fHistPS[iList][a]->RebinY(fRebinPhi);
 352           fHistNS[iList][a]->RebinY(fRebinPhi);
 353           fHistPPS[iList][a]->RebinY(fRebinPhi);
 354           fHistPNS[iList][a]->RebinY(fRebinPhi);
 355           fHistNPS[iList][a]->RebinY(fRebinPhi);
 356           fHistNNS[iList][a]->RebinY(fRebinPhi);
 357         }
 358         else{
 359           fHistPS[iList][a]->RebinY(fRebin);
 360           fHistNS[iList][a]->RebinY(fRebin);
 361           fHistPPS[iList][a]->RebinY(fRebin);
 362           fHistPNS[iList][a]->RebinY(fRebin);
 363           fHistNPS[iList][a]->RebinY(fRebin);
 364           fHistNNS[iList][a]->RebinY(fRebin);
 365         }
 366
 367         fHistPS[iList][a]->SetName(Form("%s_%d_%d",fHistPS[iList][a]->GetName(),iList,iCanvas));
 368         fHistNS[iList][a]->SetName(Form("%s_%d_%d",fHistNS[iList][a]->GetName(),iList,iCanvas));
 369         fHistPPS[iList][a]->SetName(Form("%s_%d_%d",fHistPPS[iList][a]->GetName(),iList,iCanvas));
 370         fHistPNS[iList][a]->SetName(Form("%s_%d_%d",fHistPNS[iList][a]->GetName(),iList,iCanvas));
 371         fHistNPS[iList][a]->SetName(Form("%s_%d_%d",fHistNPS[iList][a]->GetName(),iList,iCanvas));
 372         fHistNNS[iList][a]->SetName(Form("%s_%d_%d",fHistNNS[iList][a]->GetName(),iList,iCanvas));
 373
 374         // set histograms in AliBalance object
 375         bfs[iList][a]->SetHistNp(a, fHistPS[iList][a]);
 376         bfs[iList][a]->SetHistNn(a, fHistNS[iList][a]);
 377         bfs[iList][a]->SetHistNpp(a, fHistPPS[iList][a]);
 378         bfs[iList][a]->SetHistNpn(a, fHistPNS[iList][a]);
 379         bfs[iList][a]->SetHistNnp(a, fHistNPS[iList][a]);
 380         bfs[iList][a]->SetHistNnn(a, fHistNNS[iList][a]);
 381
 382         for(iCanvas = 0; iCanvas < nrOfCentralities; iCanvas++){
 383
 384           gbfs[iList][iCanvas][a] = bfs[iList][a]->GetBalanceFunctionHistogram(a,centralityArray[iCanvas],centralityArray[iCanvas+1],etaWindow);
 385           gbfs[iList][iCanvas][a]->SetName(Form("BFS_%s_Cent_%.0f_%.0f_%d",gBFAnalysisType[a].Data(),centralityArray[iCanvas],centralityArray[iCanvas+1],iList));
 386
 387           cBFS[iList][iCanvas]->cd(a+1);
 388           gbfs[iList][iCanvas][a]->SetMarkerStyle(20);
 389           if(!bHistos){
 390             gbfs[iList][iCanvas][a]->DrawCopy("AP");
 391             if(a==1){
 392               GetWeightedMean(gbfs[iList][iCanvas][a],fStartBinBFWidth,WMS[iList][iCanvas],WMSE[iList][iCanvas]);
 393               GetIntegral(gbfs[iList][iCanvas][a],integS[iList][iCanvas],integSE[iList][iCanvas]);
 394             }
 395             else if(a==6){
 396               GetWeightedMean(gbfs[iList][iCanvas][a],fStartBinBFWidthPhi,WMPS[iList][iCanvas],WMPSE[iList][iCanvas]); // for phi calculate width
 397               GetIntegral(gbfs[iList][iCanvas][a],integPS[iList][iCanvas],integPSE[iList][iCanvas]);
 398             }
 399           }
 400           else{
 401             fHistPNS[iList][a]->SetLineColor(2);
 402             fHistPNS[iList][a]->ProjectionY(Form("pns%d",a))->DrawCopy();
 403             fHistPPS[iList][a]->SetLineColor(1);
 404             fHistPPS[iList][a]->ProjectionY(Form("pps%d",a))->DrawCopy("same");
 405             fHistNPS[iList][a]->SetLineColor(4);
 406             fHistNPS[iList][a]->ProjectionY(Form("nps%d",a))->DrawCopy("same");
 407             fHistNNS[iList][a]->SetLineColor(8);
 408             fHistNNS[iList][a]->ProjectionY(Form("nns%d",a))->DrawCopy("same");
 409           }
 410         }
 411       }
 412     }
 413     // ----------------------------------------------------
 414   }
 415
 416   // for BF calculation create also graphs with weighted mean for eta
 417   if(iCanvas == 0) return;
 418
 419   TGraphErrors *gWM[13];
 420   TGraphErrors *gWMS[13];
 421   TGraphErrors *gWMP[13];
 422   TGraphErrors *gWMPS[13];
 423   TGraphErrors *ginteg[13];
 424   TGraphErrors *gintegS[13];
 425   TGraphErrors *gintegP[13];
 426   TGraphErrors *gintegPS[13];
 427   for(Int_t i = 0; i < 13; i++){
 428     gWM[i] = NULL;
 429     gWMS[i] = NULL;
 430     ginteg[i] = NULL;
 431     gintegS[i] = NULL;
 432     gWMP[i] = NULL;
 433     gWMPS[i] = NULL;
 434     gintegP[i] = NULL;
 435     gintegPS[i] = NULL;
 436   }
 437
 438
 439   if(!bHistos){
 440     for(Int_t i = 0; i < iList+1; i++){
 441       gWM[i] = new TGraphErrors(iCanvas,cent,WM[i],centE,WME[i]);
 442       if(iList==0) gWM[i]->SetName("gCentrality");
 443       else gWM[i]->SetName(Form("gCentrality_%d",i));
 444       gWMS[i] = new TGraphErrors(iCanvas,cent,WMS[i],centE,WMSE[i]);
 445       if(iList==0) gWMS[i]->SetName("gCentralityS");
 446       else gWMS[i]->SetName(Form("gCentralityS_%d",i));
 447
 448       gWMP[i] = new TGraphErrors(iCanvas,cent,WMP[i],centE,WMPE[i]);
 449       if(iList==0) gWMP[i]->SetName("gCentralityPhi");
 450       else gWMP[i]->SetName(Form("gCentralityPhi_%d",i));
 451       gWMPS[i] = new TGraphErrors(iCanvas,cent,WMPS[i],centE,WMPSE[i]);
 452       if(iList==0) gWMPS[i]->SetName("gCentralityPhiS");
 453       else gWMPS[i]->SetName(Form("gCentralityPhiS_%d",i));
 454
 455       ginteg[i] = new TGraphErrors(iCanvas,cent,integ[i],centE,integE[i]);
 456       if(iList==0) ginteg[i]->SetName("gIntegral");
 457       else ginteg[i]->SetName(Form("gIntegral_%d",i));
 458       gintegS[i] = new TGraphErrors(iCanvas,cent,integS[i],centE,integSE[i]);
 459       if(iList==0) gintegS[i]->SetName("gIntegralS");
 460       else gintegS[i]->SetName(Form("gIntegralS_%d",i));
 461
 462       gintegP[i] = new TGraphErrors(iCanvas,cent,integP[i],centE,integPE[i]);
 463       if(iList==0) gintegP[i]->SetName("gIntegraPhil");
 464       else gintegP[i]->SetName(Form("gIntegralPhi_%d",i));
 465       gintegPS[i] = new TGraphErrors(iCanvas,cent,integPS[i],centE,integPSE[i]);
 466       if(iList==0) gintegPS[i]->SetName("gIntegralPhiS");
 467       else gintegPS[i]->SetName(Form("gIntegralPhiS_%d",i));
 468     }
 469   }
 470
 471   TFile *fOut = TFile::Open(Form("Histograms_WMstart%d_rebin%d_WMstartPhi%d_rebinPhi%d_%s", fStartBinBFWidth, fRebin,fStartBinBFWidthPhi, fRebinPhi,inFile.Data()),"RECREATE");
 472   fOut->cd();
 473   for(Int_t i = 0; i < iList+1; i++){
 474     cout<<"PROCESS LIST "<<i<<" NOW!"<<endl;
 475     for(Int_t a = 0; a < 7; a++){
 476     cout<<"PROCESS VARIABLE "<<a<<" NOW!"<<endl;
 477
 478       if(fHistPN[i][a]){
 479         (fHistPN[i][a]->ProjectionY(Form("hPN_%s_%d",gBFAnalysisType[a].Data(),i)))->Write();
 480         (fHistPP[i][a]->ProjectionY(Form("hPP_%s_%d",gBFAnalysisType[a].Data(),i)))->Write();
 481         (fHistNP[i][a]->ProjectionY(Form("hNP_%s_%d",gBFAnalysisType[a].Data(),i)))->Write();
 482         (fHistNN[i][a]->ProjectionY(Form("hNN_%s_%d",gBFAnalysisType[a].Data(),i)))->Write();
 483         (fHistP[i][a]->ProjectionY(Form("hP_%s_%d",gBFAnalysisType[a].Data(),i)))->Write();
 484         (fHistN[i][a]->ProjectionY(Form("hN_%s_%d",gBFAnalysisType[a].Data(),i)))->Write();
 485       }
 486
 487       if(fHistPNS[i][a]){
 488         (fHistPNS[i][a]->ProjectionY(Form("hPNS_%s_%d",gBFAnalysisType[a].Data(),i)))->Write();
 489         (fHistPPS[i][a]->ProjectionY(Form("hPPS_%s_%d",gBFAnalysisType[a].Data(),i)))->Write();
 490         (fHistNPS[i][a]->ProjectionY(Form("hNPS_%s_%d",gBFAnalysisType[a].Data(),i)))->Write();
 491         (fHistNNS[i][a]->ProjectionY(Form("hNNS_%s_%d",gBFAnalysisType[a].Data(),i)))->Write();
 492         (fHistPS[i][a]->ProjectionY(Form("hPS_%s_%d",gBFAnalysisType[a].Data(),i)))->Write();
 493         (fHistNS[i][a]->ProjectionY(Form("hNS_%s_%d",gBFAnalysisType[a].Data(),i)))->Write();
 494       }
 495
 496       //printout in text format for delta eta
 497       for(Int_t j = 0; j < iCanvas; j++){
 498         cout<<"//=========================Centrality "<<centralityArray[j]<<"-"<<centralityArray[j+1]<<"%==================//"<<endl;
 499
 500         if(gbf[i][j][a]){
 501           if(a==1){
 502             cout<<"Double_t gALICEDataBalanceFunctionInDeltaEtaCentrality"<<centralityArray[j]<<"to"<<centralityArray[j+1]<<"[nBinsInDeltaEta] = {";
 503             for(Int_t k = 0; k < gbf[i][j][a]->GetNbinsX()-1;k++) cout<<gbf[i][j][a]->GetBinContent(k+1)<<", ";
 504             cout<<gbf[i][j][a]->GetBinContent(k+1)<<"};"<<endl;
 505             cout<<"Double_t gALICEDataBalanceFunctionInDeltaEtaCentrality"<<centralityArray[j]<<"to"<<centralityArray[j+1]<<"Error[nBinsInDeltaEta] = {";
 506             for(Int_t k = 0; k < gbf[i][j][a]->GetNbinsX()-1;k++) cout<<gbf[i][j][a]->GetBinError(k+1)<<", ";
 507             cout<<gbf[i][j][a]->GetBinError(k+1)<<"};"<<endl;
 508           }
 509           else if(a==6){
 510             cout<<"Double_t gALICEDataBalanceFunctionInDeltaPhiCentrality"<<centralityArray[j]<<"to"<<centralityArray[j+1]<<"[nBinsInDeltaPhi] = {";
 511             for(Int_t k = 0; k < gbf[i][j][a]->GetNbinsX()-1;k++) cout<<gbf[i][j][a]->GetBinContent(k+1)<<", ";
 512             cout<<gbf[i][j][a]->GetBinContent(k+1)<<"};"<<endl;
 513             cout<<"Double_t gALICEDataBalanceFunctionInDeltaPhiCentrality"<<centralityArray[j]<<"to"<<centralityArray[j+1]<<"Error[nBinsInDeltaPhi] = {";
 514             for(Int_t k = 0; k < gbf[i][j][a]->GetNbinsX()-1;k++) cout<<gbf[i][j][a]->GetBinError(k+1)<<", ";
 515             cout<<gbf[i][j][a]->GetBinError(k+1)<<"};"<<endl;
 516           }
 517           gbf[i][j][a]->Write();
 518           gbf[i][j][a]->Delete();
 519         }
 520         if(gbfs[i][j][a]){
 521           if(a==1){
 522             cout<<"Double_t gALICEShuffledDataBalanceFunctionInDeltaEtaCentrality"<<centralityArray[j]<<"to"<<centralityArray[j+1]<<"[nBinsInDeltaEta] = {";
 523             for(Int_t k = 0; k < gbfs[i][j][a]->GetNbinsX()-1;k++) cout<<gbfs[i][j][a]->GetBinContent(k+1)<<", ";
 524             cout<<gbfs[i][j][a]->GetBinContent(k+1)<<"};"<<endl;
 525             cout<<"Double_t gALICEShuffledDataBalanceFunctionInDeltaEtaCentrality"<<centralityArray[j]<<"to"<<centralityArray[j+1]<<"Error[nBinsInDeltaEta] = {";
 526             for(Int_t k = 0; k < gbfs[i][j][a]->GetNbinsX()-1;k++) cout<<gbfs[i][j][a]->GetBinError(k+1)<<", ";
 527             cout<<gbfs[i][j][a]->GetBinError(k+1)<<"};"<<endl;
 528           }
 529           else if(a==6){
 530             cout<<"Double_t gALICEShuffledDataBalanceFunctionInDeltaPhiCentrality"<<centralityArray[j]<<"to"<<centralityArray[j+1]<<"[nBinsInDeltaPhi] = {";
 531             for(Int_t k = 0; k < gbfs[i][j][a]->GetNbinsX()-1;k++) cout<<gbfs[i][j][a]->GetBinContent(k+1)<<", ";
 532             cout<<gbfs[i][j][a]->GetBinContent(k+1)<<"};"<<endl;
 533             cout<<"Double_t gALICEShuffledDataBalanceFunctionInDeltaPhiCentrality"<<centralityArray[j]<<"to"<<centralityArray[j+1]<<"Error[nBinsInDeltaPhi] = {";
 534             for(Int_t k = 0; k < gbfs[i][j][a]->GetNbinsX()-1;k++) cout<<gbfs[i][j][a]->GetBinError(k+1)<<", ";
 535             cout<<gbfs[i][j][a]->GetBinError(k+1)<<"};"<<endl;
 536           }
 537           gbfs[i][j][a]->Write();
 538           gbfs[i][j][a]->Delete();
 539         }
 540         cout<<"//=========================Centrality "<<centralityArray[j]<<"-"<<centralityArray[j+1]<<"%==================//"<<endl;
 541         cout<<endl;
 542       }
 543     }
 544
 545     Double_t x,y;
 546
 547     cout<<"//================================ALICE================================//"<<endl;
 548     if(gWM[i]){
 549       cout<<"Double_t gWeightedMeanInEtaAlice[nCentralityBins] = {";
 550       for(Int_t k = 0; k < gWM[i]->GetN()-1;k++){
 551         gWM[i]->GetPoint(k,x,y);
 552         cout<<y<<", ";
 553       }
 554       gWM[i]->GetPoint(k,x,y);
 555       cout<<y<<"};"<<endl;
 556
 557       cout<<"Double_t gWeightedMeanInEtaAliceError[nCentralityBins] = {";
 558       for(Int_t k = 0; k < gWM[i]->GetN()-1;k++){
 559         cout<<gWM[i]->GetErrorY(k)<<", ";
 560       }
 561       cout<<gWM[i]->GetErrorY(k)<<"};"<<endl;
 562       gWM[i]->Write();
 563     }
 564     if(gWMS[i]){
 565       cout<<"Double_t gShuffledWeightedMeanInEtaAlice[nCentralityBins] = {";
 566       for(Int_t k = 0; k < gWMS[i]->GetN()-1;k++){
 567         gWMS[i]->GetPoint(k,x,y);
 568         cout<<y<<", ";
 569       }
 570       gWMS[i]->GetPoint(k,x,y);
 571       cout<<y<<"};"<<endl;
 572
 573       cout<<"Double_t gShuffledWeightedMeanInEtaAliceError[nCentralityBins] = {";
 574       for(Int_t k = 0; k < gWMS[i]->GetN()-1;k++){
 575         cout<<gWMS[i]->GetErrorY(k)<<", ";
 576       }
 577       cout<<gWMS[i]->GetErrorY(k)<<"};"<<endl;
 578       cout<<endl;
 579       gWMS[i]->Write();
 580     }
 581
 582     if(gWMP[i]){
 583       cout<<"Double_t gWeightedMeanInPhiAlice[nCentralityBins] = {";
 584       for(Int_t k = 0; k < gWMP[i]->GetN()-1;k++){
 585         gWMP[i]->GetPoint(k,x,y);
 586         cout<<y<<", ";
 587       }
 588       gWMP[i]->GetPoint(k,x,y);
 589       cout<<y<<"};"<<endl;
 590
 591       cout<<"Double_t gWeightedMeanInPhiAliceError[nCentralityBins] = {";
 592       for(Int_t k = 0; k < gWMP[i]->GetN()-1;k++){
 593         cout<<gWMP[i]->GetErrorY(k)<<", ";
 594       }
 595       cout<<gWMP[i]->GetErrorY(k)<<"};"<<endl;
 596       gWMP[i]->Write();
 597     }
 598     if(gWMPS[i]){
 599       cout<<"Double_t gShuffledWeightedMeanInPhiAlice[nCentralityBins] = {";
 600       for(Int_t k = 0; k < gWMPS[i]->GetN()-1;k++){
 601         gWMPS[i]->GetPoint(k,x,y);
 602         cout<<y<<", ";
 603       }
 604       gWMPS[i]->GetPoint(k,x,y);
 605       cout<<y<<"};"<<endl;
 606
 607       cout<<"Double_t gShuffledWeightedMeanInPhiAliceError[nCentralityBins] = {";
 608       for(Int_t k = 0; k < gWMPS[i]->GetN()-1;k++){
 609         cout<<gWMPS[i]->GetErrorY(k)<<", ";
 610       }
 611       cout<<gWMPS[i]->GetErrorY(k)<<"};"<<endl;
 612       cout<<endl;
 613       gWMPS[i]->Write();
 614     }
 615     if(ginteg[i]){
 616       cout<<"Double_t gIntegralEtaAlice[nCentralityBins] = {";
 617       for(Int_t k = 0; k < ginteg[i]->GetN()-1;k++){
 618         ginteg[i]->GetPoint(k,x,y);
 619         cout<<y<<", ";
 620       }
 621       ginteg[i]->GetPoint(k,x,y);
 622       cout<<y<<"};"<<endl;
 623       cout<<endl;
 624
 625       cout<<"Double_t gIntegralErrorEtaAlice[nCentralityBins] = {";
 626       for(Int_t k = 0; k < ginteg[i]->GetN()-1;k++){
 627         cout<<ginteg[i]->GetErrorY(k)<<", ";
 628       }
 629       cout<<ginteg[i]->GetErrorY(k)<<"};"<<endl;
 630       cout<<endl;
 631
 632       ginteg[i]->Write();
 633
 634       cout<<"Double_t gIntegralPhiAlice[nCentralityBins] = {";
 635       for(Int_t k = 0; k < gintegP[i]->GetN()-1;k++){
 636         gintegP[i]->GetPoint(k,x,y);
 637         cout<<y<<", ";
 638       }
 639       gintegP[i]->GetPoint(k,x,y);
 640       cout<<y<<"};"<<endl;
 641       cout<<endl;
 642
 643       cout<<"Double_t gIntegralErrorPhiAlice[nCentralityBins] = {";
 644       for(Int_t k = 0; k < gintegP[i]->GetN()-1;k++){
 645         cout<<gintegP[i]->GetErrorY(k)<<", ";
 646       }
 647       cout<<gintegP[i]->GetErrorY(k)<<"};"<<endl;
 648       cout<<endl;
 649
 650       gintegP[i]->Write();
 651
 652       cout<<"Double_t gIntegralEtaShuffledAlice[nCentralityBins] = {";
 653       for(Int_t k = 0; k < gintegS[i]->GetN()-1;k++){
 654         gintegS[i]->GetPoint(k,x,y);
 655         cout<<y<<", ";
 656       }
 657       gintegS[i]->GetPoint(k,x,y);
 658       cout<<y<<"};"<<endl;
 659       cout<<endl;
 660
 661       cout<<"Double_t gIntegralErrorEtaShuffledAlice[nCentralityBins] = {";
 662       for(Int_t k = 0; k < gintegS[i]->GetN()-1;k++){
 663         cout<<gintegS[i]->GetErrorY(k)<<", ";
 664       }
 665       cout<<gintegS[i]->GetErrorY(k)<<"};"<<endl;
 666       cout<<endl;
 667
 668       gintegS[i]->Write();
 669
 670       cout<<"Double_t gIntegralPhiShuffledAlice[nCentralityBins] = {";
 671       for(Int_t k = 0; k < gintegPS[i]->GetN()-1;k++){
 672         cout<<y<<", ";
 673       }
 674       gintegPS[i]->GetPoint(k,x,y);
 675       cout<<y<<"};"<<endl;
 676       cout<<endl;
 677
 678       cout<<"Double_t gIntegralErrorPhiShuffledAlice[nCentralityBins] = {";
 679       for(Int_t k = 0; k < gintegPS[i]->GetN()-1;k++){
 680         cout<<gintegPS[i]->GetErrorY(k)<<", ";
 681       }
 682       cout<<gintegPS[i]->GetErrorY(k)<<"};"<<endl;
 683       cout<<endl;
 684
 685       gintegPS[i]->Write();
 686
 687     }
 688   }
 689   fOut->Close();
 690   f->Close();
 691   gROOT->Reset();
 692 }
 693
 694 //____________________________________________________________________//
 695 void GetWeightedMean(TH1D *gHistBalance, Int_t fStartBin = 1, Double_t &WM, Double_t &WME) {
 696
 697   //Prints the calculated width of the BF and its error
 698   Double_t gSumXi = 0.0, gSumBi = 0.0, gSumBiXi = 0.0;
 699   Double_t gSumBiXi2 = 0.0, gSumBi2Xi2 = 0.0;
 700   Double_t gSumDeltaBi2 = 0.0, gSumXi2DeltaBi2 = 0.0;
 701   Double_t deltaBalP2 = 0.0, integral = 0.0;
 702   Double_t deltaErrorNew = 0.0;
 703
 704   //Retrieve this variables from Histogram
 705   Int_t fNumberOfBins = gHistBalance->GetNbinsX();
 706   Double_t fP2Step    = gHistBalance->GetBinWidth(1); // assume equal binning!
 707
 708   cout<<"=================================================="<<endl;
 709   cout<<"RECALCULATION OF BF WIDTH (StartBin = "<<fStartBin<<")"<<endl;
 710   cout<<"HISTOGRAM has "<<fNumberOfBins<<" bins with bin size of "<<fP2Step<<endl;
 711   for(Int_t i = fStartBin; i <= fNumberOfBins; i++) {
 712     cout<<"B: "<<gHistBalance->GetBinContent(i)<<"\t Error: "<<gHistBalance->GetBinError(i)<<"\t bin: "<<gHistBalance->GetBinCenter(i)<<endl;
 713   }
 714   cout<<"=================================================="<<endl;
 715   for(Int_t i = fStartBin; i <= fNumberOfBins; i++) {
 716     gSumXi += gHistBalance->GetBinCenter(i);
 717     gSumBi += gHistBalance->GetBinContent(i);
 718     gSumBiXi += gHistBalance->GetBinContent(i)*gHistBalance->GetBinCenter(i);
 719     gSumBiXi2 += gHistBalance->GetBinContent(i)*TMath::Power(gHistBalance->GetBinCenter(i),2);
 720     gSumBi2Xi2 += TMath::Power(gHistBalance->GetBinContent(i),2)*TMath::Power(gHistBalance->GetBinCenter(i),2);
 721     gSumDeltaBi2 +=  TMath::Power(gHistBalance->GetBinError(i),2);
 722     gSumXi2DeltaBi2 += TMath::Power(gHistBalance->GetBinCenter(i),2) * TMath::Power(gHistBalance->GetBinError(i),2);
 723
 724     deltaBalP2 += fP2Step*TMath::Power(gHistBalance->GetBinError(i),2);
 725     integral += fP2Step*gHistBalance->GetBinContent(i);
 726   }
 727   for(Int_t i = fStartBin; i < fNumberOfBins; i++)
 728     deltaErrorNew += gHistBalance->GetBinError(i)*(gHistBalance->GetBinCenter(i)*gSumBi - gSumBiXi)/TMath::Power(gSumBi,2);
 729
 730   Double_t integralError = TMath::Sqrt(deltaBalP2);
 731
 732   Double_t delta = gSumBiXi / gSumBi;
 733   Double_t deltaError = (gSumBiXi / gSumBi) * TMath::Sqrt(TMath::Power((TMath::Sqrt(gSumXi2DeltaBi2)/gSumBiXi),2) + TMath::Power((gSumDeltaBi2/gSumBi),2) );
 734
 735   cout<<"Width: "<<delta<<"\t Error: "<<deltaError<<endl;
 736   cout<<"New error: "<<deltaErrorNew<<endl;
 737   cout<<"Integral: "<<integral<<"\t Error: "<<integralError<<endl;
 738   cout<<"=================================================="<<endl;
 739
 740   WM  = delta;
 741   WME = deltaError;
 742 }
 743
 744 //____________________________________________________________________//
 745 void GetIntegral(TH1D *gHistBalance, Double_t &integ, Double_t &integE) {
 746
 747   //always start with 1
 748   Int_t fStartBin = 1;
 749
 750   //Prints the calculated width of the BF and its error
 751   Double_t gSumXi = 0.0, gSumBi = 0.0, gSumBiXi = 0.0;
 752   Double_t gSumBiXi2 = 0.0, gSumBi2Xi2 = 0.0;
 753   Double_t gSumDeltaBi2 = 0.0, gSumXi2DeltaBi2 = 0.0;
 754   Double_t deltaBalP2 = 0.0, integral = 0.0;
 755   Double_t deltaErrorNew = 0.0;
 756
 757   //Retrieve this variables from Histogram
 758   Int_t fNumberOfBins = gHistBalance->GetNbinsX();
 759   Double_t fP2Step    = gHistBalance->GetBinWidth(1); // assume equal binning!
 760
 761   cout<<"=================================================="<<endl;
 762   cout<<"RECALCULATION OF BF WIDTH (StartBin = "<<fStartBin<<")"<<endl;
 763   cout<<"HISTOGRAM has "<<fNumberOfBins<<" bins with bin size of "<<fP2Step<<endl;
 764   for(Int_t i = fStartBin; i <= fNumberOfBins; i++) {
 765     cout<<"B: "<<gHistBalance->GetBinContent(i)<<"\t Error: "<<gHistBalance->GetBinError(i)<<"\t bin: "<<gHistBalance->GetBinCenter(i)<<endl;
 766   }
 767   cout<<"=================================================="<<endl;
 768   for(Int_t i = fStartBin; i <= fNumberOfBins; i++) {
 769     gSumXi += gHistBalance->GetBinCenter(i);
 770     gSumBi += gHistBalance->GetBinContent(i);
 771     gSumBiXi += gHistBalance->GetBinContent(i)*gHistBalance->GetBinCenter(i);
 772     gSumBiXi2 += gHistBalance->GetBinContent(i)*TMath::Power(gHistBalance->GetBinCenter(i),2);
 773     gSumBi2Xi2 += TMath::Power(gHistBalance->GetBinContent(i),2)*TMath::Power(gHistBalance->GetBinCenter(i),2);
 774     gSumDeltaBi2 +=  TMath::Power(gHistBalance->GetBinError(i),2);
 775     gSumXi2DeltaBi2 += TMath::Power(gHistBalance->GetBinCenter(i),2) * TMath::Power(gHistBalance->GetBinError(i),2);
 776
 777     deltaBalP2 += fP2Step*TMath::Power(gHistBalance->GetBinError(i),2);
 778     integral += fP2Step*gHistBalance->GetBinContent(i);
 779   }
 780   for(Int_t i = fStartBin; i < fNumberOfBins; i++)
 781     deltaErrorNew += gHistBalance->GetBinError(i)*(gHistBalance->GetBinCenter(i)*gSumBi - gSumBiXi)/TMath::Power(gSumBi,2);
 782
 783   Double_t integralError = TMath::Sqrt(deltaBalP2);
 784
 785   Double_t delta = gSumBiXi / gSumBi;
 786   Double_t deltaError = (gSumBiXi / gSumBi) * TMath::Sqrt(TMath::Power((TMath::Sqrt(gSumXi2DeltaBi2)/gSumBiXi),2) + TMath::Power((gSumDeltaBi2/gSumBi),2) );
 787
 788   cout<<"Width: "<<delta<<"\t Error: "<<deltaError<<endl;
 789   cout<<"New error: "<<deltaErrorNew<<endl;
 790   cout<<"Integral: "<<integral<<"\t Error: "<<integralError<<endl;
 791   cout<<"=================================================="<<endl;
 792
 793   integ  = integral;
 794   integE = integralError;
 795 }
 796
 797 //___________________________________________________________//
 798 // NOT USED any more
 799 //___________________________________________________________//
 800 void mergeOutput(const char* outputDir) {
 801   //Function to merge the output of the sub-jobs
 802   //Create a BF object
 803   AliBalance *bf = new AliBalance();
 804
 805   //connect to AliEn's API services
 806   TGrid::Connect("alien://");
 807
 808   //Getting the output dir from the env. variable
 809   //(JDL field: JDLVariables={"OutputDir"};)
 810   TGridResult* result = gGrid->Query(outputDir,"*/root_archive.zip","","-l 1000");
 811
 812   Int_t nEntries = result->GetEntries();
 813
 814   TString alienUrl;
 815   TDirectoryFile *dirSubJob;
 816
 817   TString gCutName[4] = {"Total","Offline trigger",
 818                          "Vertex","Analyzed"};
 819   TH1F *fHistEventStats = new TH1F("fHistEventStats",
 820                                    "Event statistics;;N_{events}",
 821                                    4,0.5,4.5);
 822   for(Int_t i = 1; i <= 4; i++)
 823     fHistEventStats->GetXaxis()->SetBinLabel(i,gCutName[i-1].Data());
 824
 825   AliESDtrackCuts *bfTrackCuts = new AliESDtrackCuts("bfTrackCuts");
 826   for(Int_t i = 0; i < nEntries; i++) {
 827     alienUrl = result->GetKey(i,"turl");
 828     alienUrl += "#AnalysisResults.root";
 829     Printf("Opening file: %s",alienUrl.Data());
 830     TFile *file = TFile::Open(alienUrl.Data());
 831     dirSubJob = dynamic_cast<TDirectoryFile *>(file->Get("PWGCFEbyE.outputBalanceFunctionAnalysis.root"));
 832
 833     //merge BF
 834     AliBalance *bfSubJob = dynamic_cast<AliBalance *>(dirSubJob->Get("AliBalance"));
 835     //bfSubJob->PrintResults();
 836     bf->Merge(bfSubJob);
 837     //delete bfSubJob;
 838
 839     //merge event stats
 840     TList *listSubJob = dynamic_cast<TList *>(dirSubJob->Get("listQA"));
 841     fHistEventStats->Add(dynamic_cast<TH1F *>(listSubJob->At(0)));
 842
 843     bfTrackCuts = dynamic_cast<AliESDtrackCuts *>(listSubJob->At(1));
 844     delete listSubJob;
 845   }
 846
 847   //Create the output file
 848   TString outputFile = "AnalysisResults.Merged.root";
 849   TFile *foutput = TFile::Open(outputFile.Data(),"recreate");
 850   TDirectoryFile *dirOutput = new TDirectoryFile();
 851   dirOutput->SetName("PWGCFEbyE.outputBalanceFunctionAnalysis.root");
 852   //dirOutput->cd();
 853   dirOutput->Add(bf);
 854   TList *list = new TList();
 855   list->SetName("listQA");
 856   list->Add(fHistEventStats);
 857   list->Add(bfTrackCuts);
 858   dirOutput->Add(list);
 859   dirOutput->Write();
 860   bf->Write();
 861   list->Write();
 862   foutput->Close();
 863
 864     //cout<<alienUrl.Data()<<endl;
 865 }