PWGCF/Correlations/Base/AliUEHistograms.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id: AliUEHistograms.cxx 20164 2007-08-14 15:31:50Z morsch $ */
  17
  18 //
  19 //
  20 // encapsulates several AliUEHist objects for a full UE analysis plus additional control histograms
  21 //
  22 //
  23 // Author: Jan Fiete Grosse-Oetringhaus, Sara Vallero
  24
  25 #include "AliUEHistograms.h"
  26
  27 #include "AliCFContainer.h"
  28 #include "AliVParticle.h"
  29 #include "AliAODTrack.h"
  30
  31 #include "TList.h"
  32 #include "TCanvas.h"
  33 #include "TH2F.h"
  34 #include "TH1F.h"
  35 #include "TH3F.h"
  36 #include "TMath.h"
  37 #include "TLorentzVector.h"
  38
  39 ClassImp(AliUEHistograms)
  40
  41 const Int_t AliUEHistograms::fgkUEHists = 3;
  42
  43 AliUEHistograms::AliUEHistograms(const char* name, const char* histograms, const char* binning) :
  44   TNamed(name, name),
  45   fNumberDensitypT(0),
  46   fSumpT(0),
  47   fNumberDensityPhi(0),
  48   fCorrelationpT(0),
  49   fCorrelationEta(0),
  50   fCorrelationPhi(0),
  51   fCorrelationR(0),
  52   fCorrelationLeading2Phi(0),
  53   fCorrelationMultiplicity(0),
  54   fYields(0),
  55   fEventCount(0),
  56   fEventCountDifferential(0),
  57   fVertexContributors(0),
  58   fCentralityDistribution(0),
  59   fCentralityCorrelation(0),
  60   fITSClusterMap(0),
  61   fControlConvResoncances(0),
  62   fEfficiencyCorrectionTriggers(0),
  63   fEfficiencyCorrectionAssociated(0),
  64   fSelectCharge(0),
  65   fTriggerSelectCharge(0),
  66   fAssociatedSelectCharge(0),
  67   fTriggerRestrictEta(-1),
  68   fEtaOrdering(kFALSE),
  69   fCutConversions(kFALSE),
  70   fCutResonances(kFALSE),
  71   fRejectResonanceDaughters(-1),
  72   fOnlyOneEtaSide(0),
  73   fWeightPerEvent(kFALSE),
  74   fPtOrder(kTRUE),
  75   fRunNumber(0),
  76   fMergeCount(1)
  77 {
  78   // Constructor
  79   //
  80   // the string histograms defines which histograms are created:
  81   //    1 = NumberDensitypT
  82   //    2 = SumpT
  83   //    3 = NumberDensityPhi
  84   //    4 = NumberDensityPhiCentrality (other multiplicity for Pb)
  85
  86   AliLog::SetClassDebugLevel("AliCFContainer", -1);
  87   AliLog::SetClassDebugLevel("AliCFGridSparse", -3);
  88
  89   fTwoTrackDistancePt[0] = 0;
  90   fTwoTrackDistancePt[1] = 0;
  91
  92   TString histogramsStr(histograms);
  93
  94   TString defaultBinningStr;
  95   defaultBinningStr = "eta: -1.0, -0.9, -0.8, -0.7, -0.6, -0.5, -0.4, -0.3, -0.2, -0.1, 0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0\n"
  96     "p_t_assoc: 0.5, 0.75, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 8.0\n"
  97     "p_t_leading: 0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0, 5.5, 6.0, 6.5, 7.0, 7.5, 8.0, 8.5, 9.0, 9.5, 10.0, 10.5, 11.0, 11.5, 12.0, 12.5, 13.0, 13.5, 14.0, 14.5, 15.0, 15.5, 16.0, 16.5, 17.0, 17.5, 18.0, 18.5, 19.0, 19.5, 20.0, 20.5, 21.0, 21.5, 22.0, 22.5, 23.0, 23.5, 24.0, 24.5, 25.0, 25.5, 26.0, 26.5, 27.0, 27.5, 28.0, 28.5, 29.0, 29.5, 30.0, 30.5, 31.0, 31.5, 32.0, 32.5, 33.0, 33.5, 34.0, 34.5, 35.0, 35.5, 36.0, 36.5, 37.0, 37.5, 38.0, 38.5, 39.0, 39.5, 40.0, 40.5, 41.0, 41.5, 42.0, 42.5, 43.0, 43.5, 44.0, 44.5, 45.0, 45.5, 46.0, 46.5, 47.0, 47.5, 48.0, 48.5, 49.0, 49.5, 50.0\n"
  98     "p_t_leading_course: 0.5, 1.0, 2.0, 3.0, 4.0, 6.0, 8.0\n"
  99     "p_t_eff: 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.25, 1.5, 1.75, 2.0, 2.25, 2.5, 2.75, 3.0, 3.25, 3.5, 3.75, 4.0, 4.5, 5.0, 6.0, 7.0, 8.0\n"
 100     "vertex_eff: -10, -8, -6, -4, -2, 0, 2, 4, 6, 8, 10\n"
 101   ;
 102
 103   if (histogramsStr.Contains("4") || histogramsStr.Contains("5") || histogramsStr.Contains("6")) // Dphi Corr
 104   {
 105     if (histogramsStr.Contains("C"))
 106       defaultBinningStr += "multiplicity: 0, 20, 40, 60, 80, 100.1\n"; // course
 107     else
 108       defaultBinningStr += "multiplicity: 0, 1, 2, 3, 4, 5, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100.1\n";
 109
 110     if (histogramsStr.Contains("5"))
 111       defaultBinningStr += "vertex: -7, -5, -3, -1, 1, 3, 5, 7\n";
 112     else
 113       defaultBinningStr += "vertex: -10, -8, -6, -4, -2, 0, 2, 4, 6, 8, 10\n";
 114
 115     if (histogramsStr.Contains("R"))
 116       defaultBinningStr += "delta_phi: -1.570796, -1.483530, -1.396263, -1.308997, -1.221730, -1.134464, -1.047198, -0.959931, -0.872665, -0.785398, -0.698132, -0.610865, -0.523599, -0.436332, -0.349066, -0.261799, -0.174533, -0.087266, 0.0, 0.087266, 0.174533, 0.261799, 0.349066, 0.436332, 0.523599, 0.610865, 0.698132, 0.785398, 0.872665, 0.959931, 1.047198, 1.134464, 1.221730, 1.308997, 1.396263, 1.483530, 1.570796, 1.658063, 1.745329, 1.832596, 1.919862, 2.007129, 2.094395, 2.181662, 2.268928, 2.356194, 2.443461, 2.530727, 2.617994, 2.705260, 2.792527, 2.879793, 2.967060, 3.054326, 3.141593, 3.228859, 3.316126, 3.403392, 3.490659, 3.577925, 3.665191, 3.752458, 3.839724, 3.926991, 4.014257, 4.101524, 4.188790, 4.276057, 4.363323, 4.450590, 4.537856, 4.625123, 4.712389\n" // this binning starts at -pi/2 and is modulo 3
 117         "delta_eta: -2.4, -2.3, -2.2, -2.1, -2.0, -1.9, -1.8, -1.7, -1.6, -1.5, -1.4, -1.3, -1.2, -1.1, -1.0, -0.9, -0.8, -0.7, -0.6, -0.5, -0.4, -0.3, -0.2, -0.1, 0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0,2.1, 2.2, 2.3, 2.4\n"
 118       ;
 119     else // for TTR studies
 120       defaultBinningStr += "delta_phi: -1.570796, -1.483530, -1.396263, -1.308997, -1.221730, -1.134464, -1.047198, -0.959931, -0.872665, -0.785398, -0.698132, -0.610865, -0.523599, -0.436332, -0.349066, -0.261799, -0.174533, -0.087266, -0.043633, -0.021817, 0.0, 0.021817, 0.043633, 0.087266, 0.174533, 0.261799, 0.349066, 0.436332, 0.523599, 0.610865, 0.698132, 0.785398, 0.872665, 0.959931, 1.047198, 1.134464, 1.221730, 1.308997, 1.396263, 1.483530, 1.570796, 1.658063, 1.745329, 1.832596, 1.919862, 2.007129, 2.094395, 2.181662, 2.268928, 2.356194, 2.443461, 2.530727, 2.617994, 2.705260, 2.792527, 2.879793, 2.967060, 3.054326, 3.141593, 3.228859, 3.316126, 3.403392, 3.490659, 3.577925, 3.665191, 3.752458, 3.839724, 3.926991, 4.014257, 4.101524, 4.188790, 4.276057, 4.363323, 4.450590, 4.537856, 4.625123, 4.712389\n"
 121         "delta_eta: -2.0, -1.9, -1.8, -1.7, -1.6, -1.5, -1.4, -1.3, -1.2, -1.1, -1.0, -0.9, -0.8, -0.7, -0.6, -0.5, -0.4, -0.3, -0.2, -0.1, -0.05, -0.025, 0, 0.025, 0.05, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0\n"
 122       ;
 123   }
 124   else // UE
 125     defaultBinningStr += "multiplicity: -0.5, 0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5, 9.5, 10.5, 11.5, 12.5, 13.5, 200.0\n"
 126       "delta_phi: -1.570796, -1.483530, -1.396263, -1.308997, -1.221730, -1.134464, -1.047198, -0.959931, -0.872665, -0.785398, -0.698132, -0.610865, -0.523599, -0.436332, -0.349066, -0.261799, -0.174533, -0.087266, 0.0, 0.087266, 0.174533, 0.261799, 0.349066, 0.436332, 0.523599, 0.610865, 0.698132, 0.785398, 0.872665, 0.959931, 1.047198, 1.134464, 1.221730, 1.308997, 1.396263, 1.483530, 1.570796, 1.658063, 1.745329, 1.832596, 1.919862, 2.007129, 2.094395, 2.181662, 2.268928, 2.356194, 2.443461, 2.530727, 2.617994, 2.705260, 2.792527, 2.879793, 2.967060, 3.054326, 3.141593, 3.228859, 3.316126, 3.403392, 3.490659, 3.577925, 3.665191, 3.752458, 3.839724, 3.926991, 4.014257, 4.101524, 4.188790, 4.276057, 4.363323, 4.450590, 4.537856, 4.625123, 4.712389\n"
 127       "vertex: -10, -8, -6, -4, -2, 0, 2, 4, 6, 8, 10\n"
 128     ;
 129
 130   // combine customBinning with defaultBinningStr -> use customBinning where available and otherwise defaultBinningStr
 131   TString customBinning(binning);
 132   TString binningStr;
 133
 134   TObjArray* lines = defaultBinningStr.Tokenize("\n");
 135   for (Int_t i=0; i<lines->GetEntriesFast(); i++)
 136   {
 137     TString line(lines->At(i)->GetName());
 138     TString tag = line(0, line.Index(":")+1);
 139     if (!customBinning.BeginsWith(tag) && !customBinning.Contains(TString("\n") + tag))
 140       binningStr += line + "\n";
 141     else
 142       Printf("Using custom binning for %s", tag.Data());
 143   }
 144   delete lines;
 145   binningStr += customBinning;
 146
 147   if (histogramsStr.Contains("1"))
 148     fNumberDensitypT = new AliUEHist("NumberDensitypT", binningStr);
 149   if (histogramsStr.Contains("2"))
 150     fSumpT = new AliUEHist("SumpT", binningStr);
 151
 152   if (histogramsStr.Contains("3"))
 153     fNumberDensityPhi = new AliUEHist("NumberDensityPhi", binningStr);
 154   else if (histogramsStr.Contains("4"))
 155     fNumberDensityPhi = new AliUEHist("NumberDensityPhiCentrality", binningStr);
 156   else if (histogramsStr.Contains("5") || histogramsStr.Contains("6"))
 157     fNumberDensityPhi = new AliUEHist("NumberDensityPhiCentralityVtx", binningStr);
 158
 159   // do not add this hists to the directory
 160   Bool_t oldStatus = TH1::AddDirectoryStatus();
 161   TH1::AddDirectory(kFALSE);
 162
 163   if (!histogramsStr.Contains("4") && !histogramsStr.Contains("5") && !histogramsStr.Contains("6"))
 164   {
 165     fCorrelationpT  = new TH2F("fCorrelationpT", ";p_{T,lead} (MC);p_{T,lead} (RECO)", 200, 0, 50, 200, 0, 50);
 166     fCorrelationEta = new TH2F("fCorrelationEta", ";#eta_{lead} (MC);#eta_{T,lead} (RECO)", 200, -1, 1, 200, -1, 1);
 167     fCorrelationPhi = new TH2F("fCorrelationPhi", ";#varphi_{lead} (MC);#varphi_{T,lead} (RECO)", 200, 0, TMath::TwoPi(), 200, 0, TMath::TwoPi());
 168   }
 169   else
 170   {
 171     fCorrelationpT  = new TH2F("fCorrelationpT", ";Centrality;p_{T} (RECO)", 100, 0, 100.001, 200, 0, 50);
 172     fCorrelationEta = new TH2F("fCorrelationEta", ";Centrality;#eta (RECO)", 100, 0, 100.001, 200, -5, 5);
 173     fCorrelationPhi = new TH2F("fCorrelationPhi", ";Centrality;#varphi (RECO)", 100, 0, 100.001, 200, 0, TMath::TwoPi());
 174   }
 175
 176   fCorrelationR =   new TH2F("fCorrelationR", ";R;p_{T,lead} (MC)", 200, 0, 2, 200, 0, 50);
 177   fCorrelationLeading2Phi = new TH2F("fCorrelationLeading2Phi", ";#Delta #varphi;p_{T,lead} (MC)", 200, -TMath::Pi(), TMath::Pi(), 200, 0, 50);
 178   fCorrelationMultiplicity = new TH2F("fCorrelationMultiplicity", ";MC tracks;Reco tracks", 100, -0.5, 99.5, 100, -0.5, 99.5);
 179   fYields = new TH3F("fYields", ";centrality;pT;eta", 100, 0, 100, 40, 0, 20, 100, -1, 1);
 180
 181   if (!histogramsStr.Contains("4") && !histogramsStr.Contains("5") && !histogramsStr.Contains("6"))
 182   {
 183     fEventCount = new TH2F("fEventCount", ";step;event type;count", AliUEHist::fgkCFSteps+2, -2.5, -0.5 + AliUEHist::fgkCFSteps, 3, -0.5, 2.5);
 184     fEventCount->GetYaxis()->SetBinLabel(1, "ND");
 185     fEventCount->GetYaxis()->SetBinLabel(2, "SD");
 186     fEventCount->GetYaxis()->SetBinLabel(3, "DD");
 187   }
 188   else
 189   {
 190     fEventCount = new TH2F("fEventCount", ";step;centrality;count", AliUEHist::fgkCFSteps+2, -2.5, -0.5 + AliUEHist::fgkCFSteps, fNumberDensityPhi->GetEventHist()->GetNBins(1), fNumberDensityPhi->GetEventHist()->GetAxis(1, 0)->GetXbins()->GetArray());
 191   }
 192
 193   fEventCountDifferential = new TH3F("fEventCountDifferential", ";p_{T,lead};step;event type", 100, 0, 50, AliUEHist::fgkCFSteps, -0.5, -0.5 + AliUEHist::fgkCFSteps, 3, -0.5, 2.5);
 194   fEventCountDifferential->GetZaxis()->SetBinLabel(1, "ND");
 195   fEventCountDifferential->GetZaxis()->SetBinLabel(2, "SD");
 196   fEventCountDifferential->GetZaxis()->SetBinLabel(3, "DD");
 197
 198   fVertexContributors = new TH1F("fVertexContributors", ";contributors;count", 100, -0.5, 99.5);
 199
 200   if (fNumberDensityPhi)
 201   {
 202     fCentralityDistribution = new TH1F("fCentralityDistribution", ";centrality;count", fNumberDensityPhi->GetEventHist()->GetNBins(1), fNumberDensityPhi->GetEventHist()->GetAxis(1, 0)->GetXbins()->GetArray());
 203     fCentralityCorrelation = new TH2F("fCentralityCorrelation", ";centrality;multiplicity", 101, 0, 101, 200, 0, 4000);
 204   }
 205
 206   fITSClusterMap = new TH3F("fITSClusterMap", "; its cluster map; centrality; pT", 256, -0.5, 255.5, 20, 0, 100.001, 100, 0, 20);
 207
 208   fControlConvResoncances = new TH2F("fControlConvResoncances", ";id;delta mass", 3, -0.5, 2.5, 100, -0.1, 0.1);
 209
 210   TH1::AddDirectory(oldStatus);
 211 }
 212
 213 //_____________________________________________________________________________
 214 AliUEHistograms::AliUEHistograms(const AliUEHistograms &c) :
 215   TNamed(fName, fTitle),
 216   fNumberDensitypT(0),
 217   fSumpT(0),
 218   fNumberDensityPhi(0),
 219   fCorrelationpT(0),
 220   fCorrelationEta(0),
 221   fCorrelationPhi(0),
 222   fCorrelationR(0),
 223   fCorrelationLeading2Phi(0),
 224   fCorrelationMultiplicity(0),
 225   fYields(0),
 226   fEventCount(0),
 227   fEventCountDifferential(0),
 228   fVertexContributors(0),
 229   fCentralityDistribution(0),
 230   fCentralityCorrelation(0),
 231   fITSClusterMap(0),
 232   fControlConvResoncances(0),
 233   fEfficiencyCorrectionTriggers(0),
 234   fEfficiencyCorrectionAssociated(0),
 235   fSelectCharge(0),
 236   fTriggerSelectCharge(0),
 237   fAssociatedSelectCharge(0),
 238   fTriggerRestrictEta(-1),
 239   fEtaOrdering(kFALSE),
 240   fCutConversions(kFALSE),
 241   fCutResonances(kFALSE),
 242   fRejectResonanceDaughters(-1),
 243   fOnlyOneEtaSide(0),
 244   fWeightPerEvent(kFALSE),
 245   fPtOrder(kTRUE),
 246   fRunNumber(0),
 247   fMergeCount(1)
 248 {
 249   //
 250   // AliUEHistograms copy constructor
 251   //
 252
 253   fTwoTrackDistancePt[0] = 0;
 254   fTwoTrackDistancePt[1] = 0;
 255
 256   ((AliUEHistograms &) c).Copy(*this);
 257 }
 258
 259 //____________________________________________________________________
 260 AliUEHistograms::~AliUEHistograms()
 261 {
 262   // Destructor
 263
 264   DeleteContainers();
 265 }
 266
 267 void AliUEHistograms::DeleteContainers()
 268 {
 269   if (fNumberDensitypT)
 270   {
 271     delete fNumberDensitypT;
 272     fNumberDensitypT = 0;
 273   }
 274
 275   if (fSumpT)
 276   {
 277     delete fSumpT;
 278     fSumpT = 0;
 279   }
 280
 281   if (fNumberDensityPhi)
 282   {
 283     delete fNumberDensityPhi;
 284     fNumberDensityPhi = 0;
 285   }
 286
 287   if (fCorrelationpT)
 288   {
 289     delete fCorrelationpT;
 290     fCorrelationpT = 0;
 291   }
 292
 293   if (fCorrelationEta)
 294   {
 295     delete fCorrelationEta;
 296     fCorrelationEta = 0;
 297   }
 298
 299   if (fCorrelationPhi)
 300   {
 301     delete fCorrelationPhi;
 302     fCorrelationPhi = 0;
 303   }
 304
 305   if (fCorrelationR)
 306   {
 307     delete fCorrelationR;
 308     fCorrelationR = 0;
 309   }
 310
 311   if (fCorrelationLeading2Phi)
 312   {
 313     delete fCorrelationLeading2Phi;
 314     fCorrelationLeading2Phi = 0;
 315   }
 316
 317   if (fCorrelationMultiplicity)
 318   {
 319     delete fCorrelationMultiplicity;
 320     fCorrelationMultiplicity = 0;
 321   }
 322
 323   if (fYields)
 324   {
 325     delete fYields;
 326     fYields = 0;
 327   }
 328
 329   if (fEventCount)
 330   {
 331     delete fEventCount;
 332     fEventCount = 0;
 333   }
 334
 335   if (fEventCountDifferential)
 336   {
 337     delete fEventCountDifferential;
 338     fEventCountDifferential = 0;
 339   }
 340
 341   if (fVertexContributors)
 342   {
 343     delete fVertexContributors;
 344     fVertexContributors = 0;
 345   }
 346
 347   if (fCentralityDistribution)
 348   {
 349     delete fCentralityDistribution;
 350     fCentralityDistribution = 0;
 351   }
 352
 353   if (fCentralityCorrelation)
 354   {
 355     delete fCentralityCorrelation;
 356     fCentralityCorrelation = 0;
 357   }
 358
 359   if (fITSClusterMap)
 360   {
 361     delete fITSClusterMap;
 362     fITSClusterMap = 0;
 363   }
 364
 365   for (Int_t i=0; i<2; i++)
 366     if (fTwoTrackDistancePt[i])
 367     {
 368       delete fTwoTrackDistancePt[i];
 369       fTwoTrackDistancePt[i] = 0;
 370     }
 371
 372   if (fControlConvResoncances)
 373   {
 374     delete fControlConvResoncances;
 375     fControlConvResoncances = 0;
 376   }
 377
 378   if (fEfficiencyCorrectionTriggers)
 379   {
 380     delete fEfficiencyCorrectionTriggers;
 381     fEfficiencyCorrectionTriggers = 0;
 382   }
 383
 384   if (fEfficiencyCorrectionAssociated)
 385   {
 386     delete fEfficiencyCorrectionAssociated;
 387     fEfficiencyCorrectionAssociated = 0;
 388   }
 389 }
 390
 391 AliUEHist* AliUEHistograms::GetUEHist(Int_t id)
 392 {
 393   // returns AliUEHist object, useful for loops
 394
 395   switch (id)
 396   {
 397     case 0: return fNumberDensitypT; break;
 398     case 1: return fSumpT; break;
 399     case 2: return fNumberDensityPhi; break;
 400   }
 401
 402   return 0;
 403 }
 404
 405 //____________________________________________________________________
 406 Int_t AliUEHistograms::CountParticles(TList* list, Float_t ptMin)
 407 {
 408   // counts the number of particles in the list with a pT above ptMin
 409   // TODO eta cut needed here?
 410
 411   Int_t count = 0;
 412   for (Int_t j=0; j<list->GetEntries(); j++)
 413     if (((AliVParticle*) list->At(j))->Pt() > ptMin)
 414       count++;
 415
 416   return count;
 417 }
 418
 419 //____________________________________________________________________
 420 void AliUEHistograms::Fill(Int_t eventType, Float_t zVtx, AliUEHist::CFStep step, AliVParticle* leading, TList* toward, TList* away, TList* min, TList* max)
 421 {
 422   // fills the UE event histograms
 423   //
 424   // this function needs the leading (track or jet or ...) and four lists of AliVParticles which contain the particles/tracks to be filled in the four regions
 425
 426   // if leading is not set, just fill event statistics
 427   if (leading)
 428   {
 429     Int_t multiplicity = 0;
 430
 431     // TODO configurable?
 432     Float_t ptMin = 0.15;
 433     if (leading->Pt() > ptMin)
 434       multiplicity++;
 435
 436     multiplicity += CountParticles(toward, ptMin);
 437     multiplicity += CountParticles(away, ptMin);
 438     multiplicity += CountParticles(min, ptMin);
 439     multiplicity += CountParticles(max, ptMin);
 440
 441     FillRegion(AliUEHist::kToward, zVtx, step, leading, toward, multiplicity);
 442     FillRegion(AliUEHist::kAway,   zVtx, step, leading, away, multiplicity);
 443     FillRegion(AliUEHist::kMin,    zVtx, step, leading, min, multiplicity);
 444     FillRegion(AliUEHist::kMax,    zVtx, step, leading, max, multiplicity);
 445
 446     Double_t vars[3];
 447     vars[0] = leading->Pt();
 448     vars[1] = multiplicity;
 449     vars[2] = zVtx;
 450     for (Int_t i=0; i<fgkUEHists; i++)
 451       if (GetUEHist(i))
 452         GetUEHist(i)->GetEventHist()->Fill(vars, step);
 453
 454     fEventCountDifferential->Fill(leading->Pt(), step, eventType);
 455   }
 456
 457   FillEvent(eventType, step);
 458 }
 459
 460 //____________________________________________________________________
 461 void AliUEHistograms::FillRegion(AliUEHist::Region region, Float_t zVtx, AliUEHist::CFStep step, AliVParticle* leading, TList* list, Int_t multiplicity)
 462 {
 463   // loops over AliVParticles in list and fills the given region at the given step
 464   //
 465   // See also Fill(...)
 466
 467   for (Int_t i=0; i<list->GetEntries(); i++)
 468   {
 469     AliVParticle* particle = (AliVParticle*) list->At(i);
 470
 471     Double_t vars[6];
 472     vars[0] = particle->Eta();
 473     vars[1] = particle->Pt();
 474     vars[2] = leading->Pt();
 475     vars[3] = multiplicity;
 476     vars[4] = leading->Phi() - particle->Phi();
 477     if (vars[4] > 1.5 * TMath::Pi())
 478       vars[4] -= TMath::TwoPi();
 479     if (vars[4] < -0.5 * TMath::Pi())
 480       vars[4] += TMath::TwoPi();
 481     vars[5] = zVtx;
 482
 483     if (fNumberDensitypT)
 484       fNumberDensitypT->GetTrackHist(region)->Fill(vars, step);
 485
 486     if (fSumpT)
 487       fSumpT->GetTrackHist(region)->Fill(vars, step, particle->Pt());
 488
 489     // fill all in toward region (is anyway as function of delta phi!)
 490     if (fNumberDensityPhi)
 491       fNumberDensityPhi->GetTrackHist(AliUEHist::kToward)->Fill(vars, step);
 492   }
 493 }
 494
 495 //____________________________________________________________________
 496 void AliUEHistograms::Fill(AliVParticle* leadingMC, AliVParticle* leadingReco)
 497 {
 498   // fills the correlation histograms
 499
 500   if (leadingMC)
 501   {
 502     fCorrelationpT->Fill(leadingMC->Pt(), leadingReco->Pt());
 503     if (leadingMC->Pt() > 0.5)
 504     {
 505       fCorrelationEta->Fill(leadingMC->Eta(), leadingReco->Eta());
 506       fCorrelationPhi->Fill(leadingMC->Phi(), leadingReco->Phi());
 507     }
 508
 509     Float_t phiDiff = leadingMC->Phi() - leadingReco->Phi();
 510     if (phiDiff > TMath::Pi())
 511       phiDiff -= TMath::TwoPi();
 512     if (phiDiff < -TMath::Pi())
 513       phiDiff += TMath::TwoPi();
 514
 515     Float_t etaDiff = leadingMC->Eta() - leadingReco->Eta();
 516
 517     fCorrelationR->Fill(TMath::Sqrt(phiDiff * phiDiff + etaDiff * etaDiff), leadingMC->Pt());
 518     fCorrelationLeading2Phi->Fill(phiDiff, leadingMC->Pt());
 519   }
 520   else
 521   {
 522     fCorrelationpT->Fill(1.0, leadingReco->Pt());
 523     if (leadingReco->Pt() > 0.5)
 524     {
 525       fCorrelationEta->Fill(0.0, leadingReco->Eta());
 526       fCorrelationPhi->Fill(0.0, leadingReco->Phi());
 527     }
 528   }
 529 }
 530
 531 //____________________________________________________________________
 532 void AliUEHistograms::FillCorrelations(Double_t centrality, Float_t zVtx, AliUEHist::CFStep step, TObjArray* particles, TObjArray* mixed, Float_t weight, Bool_t firstTime, Bool_t twoTrackEfficiencyCut, Float_t bSign, Float_t twoTrackEfficiencyCutValue, Bool_t applyEfficiency)
 533 {
 534   // fills the fNumberDensityPhi histogram
 535   //
 536   // this function need a list of AliVParticles which contain the particles/tracks to be filled
 537   //
 538   // if mixed is non-0, mixed events are filled, the trigger particle is from particles, the associated from mixed
 539   // if weight < 0, then the pt of the associated particle is filled as weight
 540
 541   Bool_t fillpT = kFALSE;
 542   if (weight < 0)
 543     fillpT = kTRUE;
 544
 545   if (twoTrackEfficiencyCut && !fTwoTrackDistancePt[0])
 546   {
 547     // do not add this hists to the directory
 548     Bool_t oldStatus = TH1::AddDirectoryStatus();
 549     TH1::AddDirectory(kFALSE);
 550
 551     fTwoTrackDistancePt[0] = new TH3F("fTwoTrackDistancePt[0]", ";#Delta#eta;#Delta#varphi^{*}_{min};#Delta p_{T}", 100, -0.15, 0.15, 100, -0.05, 0.05, 20, 0, 10);
 552     fTwoTrackDistancePt[1] = (TH3F*) fTwoTrackDistancePt[0]->Clone("fTwoTrackDistancePt[1]");
 553
 554     TH1::AddDirectory(oldStatus);
 555   }
 556
 557   // Eta() is extremely time consuming, therefore cache it for the inner loop here:
 558   TObjArray* input = (mixed) ? mixed : particles;
 559   TArrayF eta(input->GetEntriesFast());
 560   for (Int_t i=0; i<input->GetEntriesFast(); i++)
 561     eta[i] = ((AliVParticle*) input->UncheckedAt(i))->Eta();
 562
 563   // if particles is not set, just fill event statistics
 564   if (particles)
 565   {
 566     Int_t jMax = particles->GetEntriesFast();
 567     if (mixed)
 568       jMax = mixed->GetEntriesFast();
 569
 570     TH1* triggerWeighting = 0;
 571     if (fWeightPerEvent)
 572     {
 573       TAxis* axis = fNumberDensityPhi->GetTrackHist(AliUEHist::kToward)->GetGrid(0)->GetGrid()->GetAxis(2);
 574       triggerWeighting = new TH1F("triggerWeighting", "", axis->GetNbins(), axis->GetXbins()->GetArray());
 575
 576       for (Int_t i=0; i<particles->GetEntriesFast(); i++)
 577       {
 578         AliVParticle* triggerParticle = (AliVParticle*) particles->UncheckedAt(i);
 579
 580         // some optimization
 581         Float_t triggerEta = triggerParticle->Eta();
 582
 583         if (fTriggerRestrictEta > 0 && TMath::Abs(triggerEta) > fTriggerRestrictEta)
 584           continue;
 585
 586         if (fOnlyOneEtaSide != 0)
 587         {
 588           if (fOnlyOneEtaSide * triggerEta < 0)
 589             continue;
 590         }
 591
 592         if (fTriggerSelectCharge != 0)
 593           if (triggerParticle->Charge() * fTriggerSelectCharge < 0)
 594             continue;
 595
 596         triggerWeighting->Fill(triggerParticle->Pt());
 597       }
 598     }
 599
 600     // identify K, Lambda candidates and flag those particles
 601     // a TObject bit is used for this
 602     const UInt_t kResonanceDaughterFlag = 1 << 14;
 603     if (fRejectResonanceDaughters > 0)
 604     {
 605       Double_t resonanceMass = -1;
 606       Double_t massDaughter1 = -1;
 607       Double_t massDaughter2 = -1;
 608       const Double_t interval = 0.02;
 609
 610       switch (fRejectResonanceDaughters)
 611       {
 612         case 1: resonanceMass = 1.2; massDaughter1 = 0.1396; massDaughter2 = 0.9383; break; // method test
 613         case 2: resonanceMass = 0.4976; massDaughter1 = 0.1396; massDaughter2 = massDaughter1; break; // k0
 614         case 3: resonanceMass = 1.115; massDaughter1 = 0.1396; massDaughter2 = 0.9383; break; // lambda
 615         default: AliFatal(Form("Invalid setting %d", fRejectResonanceDaughters));
 616       }
 617
 618       for (Int_t i=0; i<particles->GetEntriesFast(); i++)
 619         particles->UncheckedAt(i)->ResetBit(kResonanceDaughterFlag);
 620       if (mixed)
 621         for (Int_t i=0; i<jMax; i++)
 622           mixed->UncheckedAt(i)->ResetBit(kResonanceDaughterFlag);
 623
 624       for (Int_t i=0; i<particles->GetEntriesFast(); i++)
 625       {
 626         AliVParticle* triggerParticle = (AliVParticle*) particles->UncheckedAt(i);
 627
 628         for (Int_t j=0; j<jMax; j++)
 629         {
 630           if (!mixed && i == j)
 631             continue;
 632
 633           AliVParticle* particle = 0;
 634           if (!mixed)
 635             particle = (AliVParticle*) particles->UncheckedAt(j);
 636           else
 637             particle = (AliVParticle*) mixed->UncheckedAt(j);
 638
 639           // check if both particles point to the same element (does not occur for mixed events, but if subsets are mixed within the same event)
 640           if (mixed && triggerParticle->IsEqual(particle))
 641             continue;
 642
 643           if (triggerParticle->Charge() * particle->Charge() > 0)
 644             continue;
 645
 646           Float_t mass = GetInvMassSquaredCheap(triggerParticle->Pt(), triggerParticle->Eta(), triggerParticle->Phi(), particle->Pt(), particle->Eta(), particle->Phi(), massDaughter1, massDaughter2);
 647
 648           if (TMath::Abs(mass - resonanceMass*resonanceMass) < interval*5)
 649           {
 650             mass = GetInvMassSquared(triggerParticle->Pt(), triggerParticle->Eta(), triggerParticle->Phi(), particle->Pt(), particle->Eta(), particle->Phi(), massDaughter1, massDaughter2);
 651
 652             if (mass > (resonanceMass-interval)*(resonanceMass-interval) && mass < (resonanceMass+interval)*(resonanceMass+interval))
 653             {
 654               triggerParticle->SetBit(kResonanceDaughterFlag);
 655               particle->SetBit(kResonanceDaughterFlag);
 656
 657 //            Printf("Flagged %d %d %f", i, j, TMath::Sqrt(mass));
 658             }
 659           }
 660         }
 661       }
 662     }
 663
 664     for (Int_t i=0; i<particles->GetEntriesFast(); i++)
 665     {
 666       AliVParticle* triggerParticle = (AliVParticle*) particles->UncheckedAt(i);
 667
 668       // some optimization
 669       Float_t triggerEta = triggerParticle->Eta();
 670
 671       if (fTriggerRestrictEta > 0 && TMath::Abs(triggerEta) > fTriggerRestrictEta)
 672         continue;
 673
 674       if (fOnlyOneEtaSide != 0)
 675       {
 676         if (fOnlyOneEtaSide * triggerEta < 0)
 677           continue;
 678       }
 679
 680       if (fTriggerSelectCharge != 0)
 681         if (triggerParticle->Charge() * fTriggerSelectCharge < 0)
 682           continue;
 683
 684       if (fRejectResonanceDaughters > 0)
 685         if (triggerParticle->TestBit(kResonanceDaughterFlag))
 686         {
 687 //        Printf("Skipped i=%d", i);
 688           continue;
 689         }
 690
 691       for (Int_t j=0; j<jMax; j++)
 692       {
 693         if (!mixed && i == j)
 694           continue;
 695
 696         AliVParticle* particle = 0;
 697         if (!mixed)
 698           particle = (AliVParticle*) particles->UncheckedAt(j);
 699         else
 700           particle = (AliVParticle*) mixed->UncheckedAt(j);
 701
 702         // check if both particles point to the same element (does not occur for mixed events, but if subsets are mixed within the same event)
 703         if (mixed && triggerParticle->IsEqual(particle))
 704           continue;
 705
 706         if (fPtOrder)
 707           if (particle->Pt() >= triggerParticle->Pt())
 708             continue;
 709
 710         if (fAssociatedSelectCharge != 0)
 711           if (particle->Charge() * fAssociatedSelectCharge < 0)
 712             continue;
 713
 714         if (fSelectCharge > 0)
 715         {
 716           // skip like sign
 717           if (fSelectCharge == 1 && particle->Charge() * triggerParticle->Charge() > 0)
 718             continue;
 719
 720           // skip unlike sign
 721           if (fSelectCharge == 2 && particle->Charge() * triggerParticle->Charge() < 0)
 722             continue;
 723         }
 724
 725         if (fEtaOrdering)
 726         {
 727           if (triggerEta < 0 && eta[j] < triggerEta)
 728             continue;
 729           if (triggerEta > 0 && eta[j] > triggerEta)
 730             continue;
 731         }
 732
 733         if (fRejectResonanceDaughters > 0)
 734           if (particle->TestBit(kResonanceDaughterFlag))
 735           {
 736 //          Printf("Skipped j=%d", j);
 737             continue;
 738           }
 739
 740         // conversions
 741         if (fCutConversions && particle->Charge() * triggerParticle->Charge() < 0)
 742         {
 743           Float_t mass = GetInvMassSquaredCheap(triggerParticle->Pt(), triggerEta, triggerParticle->Phi(), particle->Pt(), eta[j], particle->Phi(), 0.510e-3, 0.510e-3);
 744
 745           if (mass < 0.1)
 746           {
 747             mass = GetInvMassSquared(triggerParticle->Pt(), triggerEta, triggerParticle->Phi(), particle->Pt(), eta[j], particle->Phi(), 0.510e-3, 0.510e-3);
 748
 749             fControlConvResoncances->Fill(0.0, mass);
 750
 751             if (mass < 0.04*0.04)
 752               continue;
 753           }
 754         }
 755
 756         // K0s
 757         if (fCutResonances && particle->Charge() * triggerParticle->Charge() < 0)
 758         {
 759           Float_t mass = GetInvMassSquaredCheap(triggerParticle->Pt(), triggerEta, triggerParticle->Phi(), particle->Pt(), eta[j], particle->Phi(), 0.1396, 0.1396);
 760
 761           const Float_t kK0smass = 0.4976;
 762
 763           if (TMath::Abs(mass - kK0smass*kK0smass) < 0.1)
 764           {
 765             mass = GetInvMassSquared(triggerParticle->Pt(), triggerEta, triggerParticle->Phi(), particle->Pt(), eta[j], particle->Phi(), 0.1396, 0.1396);
 766
 767             fControlConvResoncances->Fill(1, mass - kK0smass*kK0smass);
 768
 769             if (mass > (kK0smass-0.02)*(kK0smass-0.02) && mass < (kK0smass+0.02)*(kK0smass+0.02))
 770               continue;
 771           }
 772         }
 773
 774         // Lambda
 775         if (fCutResonances && particle->Charge() * triggerParticle->Charge() < 0)
 776         {
 777           Float_t mass1 = GetInvMassSquaredCheap(triggerParticle->Pt(), triggerEta, triggerParticle->Phi(), particle->Pt(), eta[j], particle->Phi(), 0.1396, 0.9383);
 778           Float_t mass2 = GetInvMassSquaredCheap(triggerParticle->Pt(), triggerEta, triggerParticle->Phi(), particle->Pt(), eta[j], particle->Phi(), 0.9383, 0.1396);
 779
 780           const Float_t kLambdaMass = 1.115;
 781
 782           if (TMath::Abs(mass1 - kLambdaMass*kLambdaMass) < 0.1)
 783           {
 784             mass1 = GetInvMassSquared(triggerParticle->Pt(), triggerEta, triggerParticle->Phi(), particle->Pt(), eta[j], particle->Phi(), 0.1396, 0.9383);
 785
 786             fControlConvResoncances->Fill(2, mass1 - kLambdaMass*kLambdaMass);
 787
 788             if (mass1 > (kLambdaMass-0.02)*(kLambdaMass-0.02) && mass1 < (kLambdaMass+0.02)*(kLambdaMass+0.02))
 789               continue;
 790           }
 791           if (TMath::Abs(mass2 - kLambdaMass*kLambdaMass) < 0.1)
 792           {
 793             mass2 = GetInvMassSquared(triggerParticle->Pt(), triggerEta, triggerParticle->Phi(), particle->Pt(), eta[j], particle->Phi(), 0.9383, 0.1396);
 794
 795             fControlConvResoncances->Fill(2, mass2 - kLambdaMass*kLambdaMass);
 796
 797             if (mass2 > (kLambdaMass-0.02)*(kLambdaMass-0.02) && mass2 < (kLambdaMass+0.02)*(kLambdaMass+0.02))
 798               continue;
 799           }
 800         }
 801
 802         if (twoTrackEfficiencyCut)
 803         {
 804           // the variables & cuthave been developed by the HBT group
 805           // see e.g. https://indico.cern.ch/materialDisplay.py?contribId=36&sessionId=6&materialId=slides&confId=142700
 806
 807           Float_t phi1 = triggerParticle->Phi();
 808           Float_t pt1 = triggerParticle->Pt();
 809           Float_t charge1 = triggerParticle->Charge();
 810
 811           Float_t phi2 = particle->Phi();
 812           Float_t pt2 = particle->Pt();
 813           Float_t charge2 = particle->Charge();
 814
 815           Float_t deta = triggerEta - eta[j];
 816
 817           // optimization
 818           if (TMath::Abs(deta) < twoTrackEfficiencyCutValue * 2.5 * 3)
 819           {
 820             // check first boundaries to see if is worth to loop and find the minimum
 821             Float_t dphistar1 = GetDPhiStar(phi1, pt1, charge1, phi2, pt2, charge2, 0.8, bSign);
 822             Float_t dphistar2 = GetDPhiStar(phi1, pt1, charge1, phi2, pt2, charge2, 2.5, bSign);
 823
 824             const Float_t kLimit = twoTrackEfficiencyCutValue * 3;
 825
 826             Float_t dphistarminabs = 1e5;
 827             Float_t dphistarmin = 1e5;
 828             if (TMath::Abs(dphistar1) < kLimit || TMath::Abs(dphistar2) < kLimit || dphistar1 * dphistar2 < 0)
 829             {
 830               for (Double_t rad=0.8; rad<2.51; rad+=0.01)
 831               {
 832                 Float_t dphistar = GetDPhiStar(phi1, pt1, charge1, phi2, pt2, charge2, rad, bSign);
 833
 834                 Float_t dphistarabs = TMath::Abs(dphistar);
 835
 836                 if (dphistarabs < dphistarminabs)
 837                 {
 838                   dphistarmin = dphistar;
 839                   dphistarminabs = dphistarabs;
 840                 }
 841               }
 842
 843               fTwoTrackDistancePt[0]->Fill(deta, dphistarmin, TMath::Abs(pt1 - pt2));
 844
 845               if (dphistarminabs < twoTrackEfficiencyCutValue && TMath::Abs(deta) < twoTrackEfficiencyCutValue)
 846               {
 847 //              Printf("Removed track pair %d %d with %f %f %f %f %f %f %f %f %f", i, j, deta, dphistarminabs, phi1, pt1, charge1, phi2, pt2, charge2, bSign);
 848                 continue;
 849               }
 850
 851               fTwoTrackDistancePt[1]->Fill(deta, dphistarmin, TMath::Abs(pt1 - pt2));
 852             }
 853           }
 854         }
 855
 856         Double_t vars[6];
 857         vars[0] = triggerEta - eta[j];
 858         vars[1] = particle->Pt();
 859         vars[2] = triggerParticle->Pt();
 860         vars[3] = centrality;
 861         vars[4] = triggerParticle->Phi() - particle->Phi();
 862         if (vars[4] > 1.5 * TMath::Pi())
 863           vars[4] -= TMath::TwoPi();
 864         if (vars[4] < -0.5 * TMath::Pi())
 865           vars[4] += TMath::TwoPi();
 866         vars[5] = zVtx;
 867
 868         if (fillpT)
 869           weight = particle->Pt();
 870
 871         Double_t useWeight = weight;
 872         if (applyEfficiency)
 873         {
 874           if (fEfficiencyCorrectionAssociated)
 875           {
 876             Int_t effVars[4];
 877             // associated particle
 878             effVars[0] = fEfficiencyCorrectionAssociated->GetAxis(0)->FindBin(eta[j]);
 879             effVars[1] = fEfficiencyCorrectionAssociated->GetAxis(1)->FindBin(vars[1]); //pt
 880             effVars[2] = fEfficiencyCorrectionAssociated->GetAxis(2)->FindBin(vars[3]); //centrality
 881             effVars[3] = fEfficiencyCorrectionAssociated->GetAxis(3)->FindBin(vars[5]); //zVtx
 882
 883             //    Printf("%d %d %d %d %f", effVars[0], effVars[1], effVars[2], effVars[3], fEfficiencyCorrectionAssociated->GetBinContent(effVars));
 884
 885             useWeight *= fEfficiencyCorrectionAssociated->GetBinContent(effVars);
 886           }
 887           if (fEfficiencyCorrectionTriggers)
 888           {
 889             Int_t effVars[4];
 890
 891             effVars[0] = fEfficiencyCorrectionTriggers->GetAxis(0)->FindBin(triggerEta);
 892             effVars[1] = fEfficiencyCorrectionTriggers->GetAxis(1)->FindBin(vars[2]); //pt
 893             effVars[2] = fEfficiencyCorrectionTriggers->GetAxis(2)->FindBin(vars[3]); //centrality
 894             effVars[3] = fEfficiencyCorrectionTriggers->GetAxis(3)->FindBin(vars[5]); //zVtx
 895             useWeight *= fEfficiencyCorrectionTriggers->GetBinContent(effVars);
 896           }
 897         }
 898
 899         if (fWeightPerEvent)
 900         {
 901           Int_t weightBin = triggerWeighting->GetXaxis()->FindBin(vars[2]);
 902 //        Printf("Using weight %f", triggerWeighting->GetBinContent(weightBin));
 903           useWeight /= triggerWeighting->GetBinContent(weightBin);
 904         }
 905
 906         // fill all in toward region and do not use the other regions
 907         fNumberDensityPhi->GetTrackHist(AliUEHist::kToward)->Fill(vars, step, useWeight);
 908
 909 //      Printf("%.2f %.2f --> %.2f", triggerEta, eta[j], vars[0]);
 910       }
 911
 912       if (firstTime)
 913       {
 914         // once per trigger particle
 915         Double_t vars[3];
 916         vars[0] = triggerParticle->Pt();
 917         vars[1] = centrality;
 918         vars[2] = zVtx;
 919
 920         Double_t useWeight = 1;
 921         if (fEfficiencyCorrectionTriggers && applyEfficiency)
 922         {
 923           Int_t effVars[4];
 924
 925           // trigger particle
 926           effVars[0] = fEfficiencyCorrectionTriggers->GetAxis(0)->FindBin(triggerEta);
 927           effVars[1] = fEfficiencyCorrectionTriggers->GetAxis(1)->FindBin(vars[0]); //pt
 928           effVars[2] = fEfficiencyCorrectionTriggers->GetAxis(2)->FindBin(vars[1]); //centrality
 929           effVars[3] = fEfficiencyCorrectionTriggers->GetAxis(3)->FindBin(vars[2]); //zVtx
 930           useWeight *= fEfficiencyCorrectionTriggers->GetBinContent(effVars);
 931         }
 932         if (fWeightPerEvent)
 933         {
 934           // leads effectively to a filling of one entry per filled trigger particle pT bin
 935           Int_t weightBin = triggerWeighting->GetXaxis()->FindBin(vars[0]);
 936 //        Printf("Using weight %f", triggerWeighting->GetBinContent(weightBin));
 937           useWeight /= triggerWeighting->GetBinContent(weightBin);
 938         }
 939
 940         fNumberDensityPhi->GetEventHist()->Fill(vars, step, useWeight);
 941
 942         // QA
 943         fCorrelationpT->Fill(centrality, triggerParticle->Pt());
 944         fCorrelationEta->Fill(centrality, triggerEta);
 945         fCorrelationPhi->Fill(centrality, triggerParticle->Phi());
 946         fYields->Fill(centrality, triggerParticle->Pt(), triggerEta);
 947
 948 /*        if (dynamic_cast<AliAODTrack*>(triggerParticle))
 949           fITSClusterMap->Fill(((AliAODTrack*) triggerParticle)->GetITSClusterMap(), centrality, triggerParticle->Pt());*/
 950       }
 951     }
 952
 953     if (triggerWeighting)
 954     {
 955       delete triggerWeighting;
 956       triggerWeighting = 0;
 957     }
 958   }
 959
 960   fCentralityDistribution->Fill(centrality);
 961   fCentralityCorrelation->Fill(centrality, particles->GetEntriesFast());
 962   FillEvent(centrality, step);
 963 }
 964
 965 //____________________________________________________________________
 966 void AliUEHistograms::FillTrackingEfficiency(TObjArray* mc, TObjArray* recoPrim, TObjArray* recoAll, TObjArray* recoPrimPID, TObjArray* recoAllPID, TObjArray* fake, Int_t particleType, Double_t centrality, Double_t zVtx)
 967 {
 968   // fills the tracking efficiency objects
 969   //
 970   // mc: all primary MC particles
 971   // recoPrim: reconstructed primaries (again MC particles)
 972   // recoAll: reconstructed (again MC particles)
 973   // recoPrim: reconstructed primaries with checks on PID (again MC particles)
 974   // recoAll: reconstructed with checks on PID (again MC particles)
 975   // particleType is: 0 for pion, 1 for kaon, 2 for proton, 3 for others
 976
 977   for (Int_t step=0; step<6; step++)
 978   {
 979     TObjArray* list = mc;
 980     if (step == 1)
 981       list = recoPrim;
 982     else if (step == 2)
 983       list = recoAll;
 984     else if (step == 3)
 985       list = recoPrimPID;
 986     else if (step == 4)
 987       list = recoAllPID;
 988     else if (step == 5)
 989       list = fake;
 990
 991     if (!list)
 992       continue;
 993
 994     for (Int_t i=0; i<list->GetEntriesFast(); i++)
 995     {
 996       AliVParticle* particle = (AliVParticle*) list->UncheckedAt(i);
 997       Double_t vars[5];
 998       vars[0] = particle->Eta();
 999       vars[1] = particle->Pt();
1000       vars[2] = particleType;
1001       vars[3] = centrality;
1002       vars[4] = zVtx;
1003
1004       for (Int_t j=0; j<fgkUEHists; j++)
1005         if (GetUEHist(j))
1006           GetUEHist(j)->GetTrackHistEfficiency()->Fill(vars, step);
1007     }
1008   }
1009 }
1010
1011 //____________________________________________________________________
1012 void AliUEHistograms::FillFakePt(TObjArray* fake, Double_t centrality)
1013 {
1014   TObjArray* tracksReco = (TObjArray*) fake->At(0);
1015   TObjArray* tracksMC = (TObjArray*) fake->At(1);
1016
1017   if (tracksReco->GetEntriesFast() != tracksMC->GetEntriesFast())
1018     AliFatal(Form("Inconsistent arrays: %d vs %d", tracksReco->GetEntriesFast(), tracksMC->GetEntriesFast()));
1019
1020   for (Int_t i=0; i<tracksReco->GetEntriesFast(); i++)
1021   {
1022     AliVParticle* particle1 = (AliVParticle*) tracksReco->At(i);
1023     AliVParticle* particle2 = (AliVParticle*) tracksMC->At(i);
1024     Double_t vars[3];
1025     vars[0] = particle1->Pt();
1026     vars[1] = particle2->Pt();
1027     vars[2] = centrality;
1028     for (Int_t j=0; j<fgkUEHists; j++)
1029       if (GetUEHist(j))
1030         GetUEHist(j)->GetMCRecoPtCorrelation()->Fill(vars[0],vars[1],vars[2]);
1031   }
1032 }
1033
1034 //____________________________________________________________________
1035 void AliUEHistograms::FillEvent(Int_t eventType, Int_t step)
1036 {
1037   // fills the number of events at the given step and the given enty type
1038   //
1039   // WARNING: This function is called from Fill, so only call it for steps where Fill is not called
1040
1041   fEventCount->Fill(step, eventType);
1042 }
1043
1044 //____________________________________________________________________
1045 void AliUEHistograms::FillEvent(Double_t centrality, Int_t step)
1046 {
1047   // fills the number of events at the given step and the given centrality
1048   //
1049   // WARNING: This function is called from Fill, so only call it for steps where Fill is not called
1050
1051   fEventCount->Fill(step, centrality);
1052 }
1053
1054 //____________________________________________________________________
1055 void AliUEHistograms::SetEtaRange(Float_t etaMin, Float_t etaMax)
1056 {
1057   // sets eta min and max for all contained AliUEHist classes
1058
1059   for (Int_t i=0; i<fgkUEHists; i++)
1060     if (GetUEHist(i))
1061       GetUEHist(i)->SetEtaRange(etaMin, etaMax);
1062 }
1063
1064 //____________________________________________________________________
1065 void AliUEHistograms::SetPtRange(Float_t ptMin, Float_t ptMax)
1066 {
1067   // sets pT min and max for all contained AliUEHist classes
1068
1069   for (Int_t i=0; i<fgkUEHists; i++)
1070     if (GetUEHist(i))
1071       GetUEHist(i)->SetPtRange(ptMin, ptMax);
1072 }
1073
1074 //____________________________________________________________________
1075 void AliUEHistograms::SetPartSpecies(Int_t species)
1076 {
1077   // sets PartSpecie for all contained AliUEHist classes
1078
1079   for (Int_t i=0; i<fgkUEHists; i++)
1080     if (GetUEHist(i))
1081       GetUEHist(i)->SetPartSpecies(species);
1082 }
1083
1084 //____________________________________________________________________
1085 void AliUEHistograms::SetZVtxRange(Float_t min, Float_t max)
1086 {
1087   // sets pT min and max for all contained AliUEHist classes
1088
1089   for (Int_t i=0; i<fgkUEHists; i++)
1090     if (GetUEHist(i))
1091       GetUEHist(i)->SetZVtxRange(min, max);
1092 }
1093
1094 //____________________________________________________________________
1095 void AliUEHistograms::SetContaminationEnhancement(TH1F* hist)
1096 {
1097   // sets the contamination enhancement histogram in all contained AliUEHist classes
1098
1099   for (Int_t i=0; i<fgkUEHists; i++)
1100     if (GetUEHist(i))
1101       GetUEHist(i)->SetContaminationEnhancement(hist);
1102 }
1103
1104 //____________________________________________________________________
1105 void AliUEHistograms::SetCombineMinMax(Bool_t flag)
1106 {
1107   // sets pT min and max for all contained AliUEHist classes
1108
1109   for (Int_t i=0; i<fgkUEHists; i++)
1110     if (GetUEHist(i))
1111       GetUEHist(i)->SetCombineMinMax(flag);
1112 }
1113
1114 //____________________________________________________________________
1115 void AliUEHistograms::SetTrackEtaCut(Float_t value)
1116 {
1117   // sets track eta cut for all contained AliUEHist classes
1118
1119   for (Int_t i=0; i<fgkUEHists; i++)
1120     if (GetUEHist(i))
1121       GetUEHist(i)->SetTrackEtaCut(value);
1122 }
1123
1124 //____________________________________________________________________
1125 void AliUEHistograms::SetWeightPerEvent(Bool_t flag)
1126 {
1127   // sets fWeightPerEvent for all contained AliUEHist classes
1128
1129   fWeightPerEvent = flag;
1130
1131   for (Int_t i=0; i<fgkUEHists; i++)
1132     if (GetUEHist(i))
1133       GetUEHist(i)->SetWeightPerEvent(fWeightPerEvent);
1134 }
1135
1136 //____________________________________________________________________
1137 void AliUEHistograms::Correct(AliUEHistograms* corrections)
1138 {
1139   // corrects the contained histograms by calling AliUEHist::Correct
1140
1141   for (Int_t i=0; i<fgkUEHists; i++)
1142     if (GetUEHist(i))
1143       GetUEHist(i)->Correct(corrections->GetUEHist(i));
1144 }
1145
1146 //____________________________________________________________________
1147 AliUEHistograms &AliUEHistograms::operator=(const AliUEHistograms &c)
1148 {
1149   // assigment operator
1150
1151   DeleteContainers();
1152
1153   if (this != &c)
1154     ((AliUEHistograms &) c).Copy(*this);
1155
1156   return *this;
1157 }
1158
1159 //____________________________________________________________________
1160 void AliUEHistograms::Copy(TObject& c) const
1161 {
1162   // copy function
1163
1164   AliUEHistograms& target = (AliUEHistograms &) c;
1165
1166   if (fNumberDensitypT)
1167     target.fNumberDensitypT = dynamic_cast<AliUEHist*> (fNumberDensitypT->Clone());
1168
1169   if (fSumpT)
1170     target.fSumpT = dynamic_cast<AliUEHist*> (fSumpT->Clone());
1171
1172   if (fNumberDensityPhi)
1173     target.fNumberDensityPhi = dynamic_cast<AliUEHist*> (fNumberDensityPhi->Clone());
1174
1175   if (fCorrelationpT)
1176     target.fCorrelationpT = dynamic_cast<TH2F*> (fCorrelationpT->Clone());
1177
1178   if (fCorrelationEta)
1179     target.fCorrelationEta = dynamic_cast<TH2F*> (fCorrelationEta->Clone());
1180
1181   if (fCorrelationPhi)
1182     target.fCorrelationPhi = dynamic_cast<TH2F*> (fCorrelationPhi->Clone());
1183
1184   if (fCorrelationR)
1185     target.fCorrelationR = dynamic_cast<TH2F*> (fCorrelationR->Clone());
1186
1187   if (fCorrelationLeading2Phi)
1188     target.fCorrelationLeading2Phi = dynamic_cast<TH2F*> (fCorrelationLeading2Phi->Clone());
1189
1190   if (fCorrelationMultiplicity)
1191     target.fCorrelationMultiplicity = dynamic_cast<TH2F*> (fCorrelationMultiplicity->Clone());
1192
1193   if (fYields)
1194     target.fYields = dynamic_cast<TH3F*> (fYields->Clone());
1195
1196   if (fEventCount)
1197     target.fEventCount = dynamic_cast<TH2F*> (fEventCount->Clone());
1198
1199   if (fEventCountDifferential)
1200     target.fEventCountDifferential = dynamic_cast<TH3F*> (fEventCountDifferential->Clone());
1201
1202   if (fVertexContributors)
1203     target.fVertexContributors = dynamic_cast<TH1F*> (fVertexContributors->Clone());
1204
1205   if (fCentralityDistribution)
1206     target.fCentralityDistribution = dynamic_cast<TH1F*> (fCentralityDistribution->Clone());
1207
1208   if (fCentralityCorrelation)
1209     target.fCentralityCorrelation = dynamic_cast<TH2F*> (fCentralityCorrelation->Clone());
1210
1211   if (fITSClusterMap)
1212     target.fITSClusterMap = dynamic_cast<TH3F*> (fITSClusterMap->Clone());
1213
1214   if (fControlConvResoncances)
1215     target.fControlConvResoncances = dynamic_cast<TH2F*> (fControlConvResoncances->Clone());
1216
1217   for (Int_t i=0; i<2; i++)
1218     if (fTwoTrackDistancePt[i])
1219       target.fTwoTrackDistancePt[i] = dynamic_cast<TH3F*> (fTwoTrackDistancePt[i]->Clone());
1220
1221   if (fEfficiencyCorrectionTriggers)
1222     target.fEfficiencyCorrectionTriggers = dynamic_cast<THnF*> (fEfficiencyCorrectionTriggers->Clone());
1223
1224  if (fEfficiencyCorrectionAssociated)
1225     target.fEfficiencyCorrectionAssociated = dynamic_cast<THnF*> (fEfficiencyCorrectionAssociated->Clone());
1226
1227   target.fSelectCharge = fSelectCharge;
1228   target.fTriggerSelectCharge = fTriggerSelectCharge;
1229   target.fAssociatedSelectCharge = fAssociatedSelectCharge;
1230   target.fTriggerRestrictEta = fTriggerRestrictEta;
1231   target.fEtaOrdering = fEtaOrdering;
1232   target.fCutConversions = fCutConversions;
1233   target.fCutResonances = fCutResonances;
1234   target.fOnlyOneEtaSide = fOnlyOneEtaSide;
1235   target.fWeightPerEvent = fWeightPerEvent;
1236   target.fRunNumber = fRunNumber;
1237   target.fMergeCount = fMergeCount;
1238   target.fWeightPerEvent = fWeightPerEvent;
1239   target.fPtOrder = fPtOrder;
1240 }
1241
1242 //____________________________________________________________________
1243 Long64_t AliUEHistograms::Merge(TCollection* list)
1244 {
1245   // Merge a list of AliUEHistograms objects with this (needed for
1246   // PROOF).
1247   // Returns the number of merged objects (including this).
1248
1249   if (!list)
1250     return 0;
1251
1252   if (list->IsEmpty())
1253     return 1;
1254
1255   TIterator* iter = list->MakeIterator();
1256   TObject* obj;
1257
1258   // collections of objects
1259   const Int_t kMaxLists = 19;
1260   TList* lists[kMaxLists];
1261
1262   for (Int_t i=0; i<kMaxLists; i++)
1263     lists[i] = new TList;
1264
1265   Int_t count = 0;
1266   while ((obj = iter->Next())) {
1267
1268     AliUEHistograms* entry = dynamic_cast<AliUEHistograms*> (obj);
1269     if (entry == 0)
1270       continue;
1271
1272     if (entry->fNumberDensitypT)
1273       lists[0]->Add(entry->fNumberDensitypT);
1274     if (entry->fSumpT)
1275       lists[1]->Add(entry->fSumpT);
1276     if (entry->fNumberDensityPhi)
1277       lists[2]->Add(entry->fNumberDensityPhi);
1278     lists[3]->Add(entry->fCorrelationpT);
1279     lists[4]->Add(entry->fCorrelationEta);
1280     lists[5]->Add(entry->fCorrelationPhi);
1281     lists[6]->Add(entry->fCorrelationR);
1282     lists[7]->Add(entry->fCorrelationLeading2Phi);
1283     lists[8]->Add(entry->fCorrelationMultiplicity);
1284     lists[9]->Add(entry->fEventCount);
1285     lists[10]->Add(entry->fEventCountDifferential);
1286     lists[11]->Add(entry->fVertexContributors);
1287     lists[12]->Add(entry->fCentralityDistribution);
1288     lists[13]->Add(entry->fITSClusterMap);
1289     if (entry->fTwoTrackDistancePt[0])
1290       lists[14]->Add(entry->fTwoTrackDistancePt[0]);
1291     if (entry->fTwoTrackDistancePt[1])
1292       lists[15]->Add(entry->fTwoTrackDistancePt[1]);
1293     if (entry->fCentralityCorrelation)
1294       lists[16]->Add(entry->fCentralityCorrelation);
1295     if (entry->fYields)
1296       lists[17]->Add(entry->fYields);
1297     if (entry->fControlConvResoncances)
1298       lists[18]->Add(entry->fControlConvResoncances);
1299
1300     fMergeCount += entry->fMergeCount;
1301
1302     count++;
1303   }
1304
1305   if (fNumberDensitypT)
1306     fNumberDensitypT->Merge(lists[0]);
1307   if (fSumpT)
1308     fSumpT->Merge(lists[1]);
1309   if (fNumberDensityPhi)
1310     fNumberDensityPhi->Merge(lists[2]);
1311   fCorrelationpT->Merge(lists[3]);
1312   fCorrelationEta->Merge(lists[4]);
1313   fCorrelationPhi->Merge(lists[5]);
1314   fCorrelationR->Merge(lists[6]);
1315   fCorrelationLeading2Phi->Merge(lists[7]);
1316   fCorrelationMultiplicity->Merge(lists[8]);
1317   fEventCount->Merge(lists[9]);
1318   fEventCountDifferential->Merge(lists[10]);
1319   fVertexContributors->Merge(lists[11]);
1320   fCentralityDistribution->Merge(lists[12]);
1321   fITSClusterMap->Merge(lists[13]);
1322   if (fTwoTrackDistancePt[0] && lists[14]->GetEntries() > 0)
1323     fTwoTrackDistancePt[0]->Merge(lists[14]);
1324   if (fTwoTrackDistancePt[1] && lists[15]->GetEntries() > 0)
1325     fTwoTrackDistancePt[1]->Merge(lists[15]);
1326   if (fCentralityCorrelation)
1327     fCentralityCorrelation->Merge(lists[16]);
1328   if (fYields && lists[17]->GetEntries() > 0)
1329     fYields->Merge(lists[17]);
1330   if (fControlConvResoncances && lists[18]->GetEntries() > 0)
1331     fControlConvResoncances->Merge(lists[18]);
1332
1333   for (Int_t i=0; i<kMaxLists; i++)
1334     delete lists[i];
1335
1336   return count+1;
1337 }
1338
1339 void AliUEHistograms::CopyReconstructedData(AliUEHistograms* from)
1340 {
1341   // copies those histograms extracted from ESD to this object
1342
1343   for (Int_t i=0; i<fgkUEHists; i++)
1344     if (GetUEHist(i))
1345       GetUEHist(i)->CopyReconstructedData(from->GetUEHist(i));
1346 }
1347
1348 void AliUEHistograms::DeepCopy(AliUEHistograms* from)
1349 {
1350   // copies the entries of this object's members from the object <from> to this object
1351
1352   for (Int_t i=0; i<fgkUEHists; i++)
1353     if (GetUEHist(i) && from->GetUEHist(i))
1354       GetUEHist(i)->DeepCopy(from->GetUEHist(i));
1355 }
1356
1357 void AliUEHistograms::ExtendTrackingEfficiency(Bool_t verbose)
1358 {
1359   // delegates to AliUEHists
1360
1361   for (Int_t i=0; i<fgkUEHists; i++)
1362     if (GetUEHist(i))
1363       GetUEHist(i)->ExtendTrackingEfficiency(verbose);
1364 }
1365
1366 void AliUEHistograms::Scale(Double_t factor)
1367 {
1368   // scales all contained histograms by the given factor
1369
1370   for (Int_t i=0; i<fgkUEHists; i++)
1371     if (GetUEHist(i))
1372       GetUEHist(i)->Scale(factor);
1373
1374   TList list;
1375   list.Add(fCorrelationpT);
1376   list.Add(fCorrelationEta);
1377   list.Add(fCorrelationPhi);
1378   list.Add(fCorrelationR);
1379   list.Add(fCorrelationLeading2Phi);
1380   list.Add(fCorrelationMultiplicity);
1381   list.Add(fYields);
1382   list.Add(fEventCount);
1383   list.Add(fEventCountDifferential);
1384   list.Add(fVertexContributors);
1385   list.Add(fCentralityDistribution);
1386   list.Add(fCentralityCorrelation);
1387   list.Add(fITSClusterMap);
1388   list.Add(fTwoTrackDistancePt[0]);
1389   list.Add(fTwoTrackDistancePt[1]);
1390   list.Add(fControlConvResoncances);
1391
1392   for (Int_t i=0; i<list.GetEntries(); i++)
1393     ((TH1*) list.At(i))->Scale(factor);
1394 }
1395
1396 void AliUEHistograms::Reset()
1397 {
1398   // delegates to AliUEHists
1399
1400   for (Int_t i=0; i<fgkUEHists; i++)
1401     if (GetUEHist(i))
1402       GetUEHist(i)->Reset();
1403 }