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