PWG4/JetTasks/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 "TH2F.h"
  33 #include "TH1F.h"
  34 #include "TH3F.h"
  35 #include "TMath.h"
  36
  37 ClassImp(AliUEHistograms)
  38
  39 const Int_t AliUEHistograms::fgkUEHists = 3;
  40
  41 AliUEHistograms::AliUEHistograms(const char* name, const char* histograms) :
  42   TNamed(name, name),
  43   fNumberDensitypT(0),
  44   fSumpT(0),
  45   fNumberDensityPhi(0),
  46   fCorrelationpT(0),
  47   fCorrelationEta(0),
  48   fCorrelationPhi(0),
  49   fCorrelationR(0),
  50   fCorrelationLeading2Phi(0),
  51   fCorrelationMultiplicity(0),
  52   fEventCount(0),
  53   fEventCountDifferential(0),
  54   fVertexContributors(0),
  55   fCentralityDistribution(0),
  56   fITSClusterMap(0),
  57   fSelectCharge(0),
  58   fRunNumber(0)
  59 {
  60   // Constructor
  61   //
  62   // the string histograms defines which histograms are created:
  63   //    1 = NumberDensitypT
  64   //    2 = SumpT
  65   //    3 = NumberDensityPhi
  66   //    4 = NumberDensityPhiCentrality (other multiplicity for Pb)
  67
  68   fTwoTrackDistancePt[0] = 0;
  69   fTwoTrackDistancePt[1] = 0;
  70
  71   TString histogramsStr(histograms);
  72
  73   if (histogramsStr.Contains("1"))
  74     fNumberDensitypT = new AliUEHist("NumberDensitypT");
  75   if (histogramsStr.Contains("2"))
  76     fSumpT = new AliUEHist("SumpT");
  77
  78   if (histogramsStr.Contains("3"))
  79     fNumberDensityPhi = new AliUEHist("NumberDensityPhi");
  80   else if (histogramsStr.Contains("4"))
  81     fNumberDensityPhi = new AliUEHist("NumberDensityPhiCentrality");
  82
  83   // do not add this hists to the directory
  84   Bool_t oldStatus = TH1::AddDirectoryStatus();
  85   TH1::AddDirectory(kFALSE);
  86
  87   if (!histogramsStr.Contains("4"))
  88   {
  89     fCorrelationpT  = new TH2F("fCorrelationpT", ";p_{T,lead} (MC);p_{T,lead} (RECO)", 200, 0, 50, 200, 0, 50);
  90     fCorrelationEta = new TH2F("fCorrelationEta", ";#eta_{lead} (MC);#eta_{T,lead} (RECO)", 200, -1, 1, 200, -1, 1);
  91     fCorrelationPhi = new TH2F("fCorrelationPhi", ";#varphi_{lead} (MC);#varphi_{T,lead} (RECO)", 200, 0, TMath::TwoPi(), 200, 0, TMath::TwoPi());
  92   }
  93   else
  94   {
  95     fCorrelationpT  = new TH2F("fCorrelationpT", ";Centrality;p_{T} (RECO)", 100, 0, 100.001, 200, 0, 50);
  96     fCorrelationEta = new TH2F("fCorrelationEta", ";Centrality;#eta (RECO)", 100, 0, 100.001, 200, -1, 1);
  97     fCorrelationPhi = new TH2F("fCorrelationPhi", ";Centrality;#varphi (RECO)", 100, 0, 100.001, 200, 0, TMath::TwoPi());
  98   }
  99
 100   fCorrelationR =   new TH2F("fCorrelationR", ";R;p_{T,lead} (MC)", 200, 0, 2, 200, 0, 50);
 101   fCorrelationLeading2Phi = new TH2F("fCorrelationLeading2Phi", ";#Delta #varphi;p_{T,lead} (MC)", 200, -TMath::Pi(), TMath::Pi(), 200, 0, 50);
 102   fCorrelationMultiplicity = new TH2F("fCorrelationMultiplicity", ";MC tracks;Reco tracks", 100, -0.5, 99.5, 100, -0.5, 99.5);
 103
 104   if (!histogramsStr.Contains("4"))
 105   {
 106     fEventCount = new TH2F("fEventCount", ";step;event type;count", AliUEHist::fgkCFSteps+2, -2.5, -0.5 + AliUEHist::fgkCFSteps, 3, -0.5, 2.5);
 107     fEventCount->GetYaxis()->SetBinLabel(1, "ND");
 108     fEventCount->GetYaxis()->SetBinLabel(2, "SD");
 109     fEventCount->GetYaxis()->SetBinLabel(3, "DD");
 110   }
 111   else
 112   {
 113     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());
 114   }
 115
 116   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);
 117   fEventCountDifferential->GetZaxis()->SetBinLabel(1, "ND");
 118   fEventCountDifferential->GetZaxis()->SetBinLabel(2, "SD");
 119   fEventCountDifferential->GetZaxis()->SetBinLabel(3, "DD");
 120
 121   fVertexContributors = new TH1F("fVertexContributors", ";contributors;count", 100, -0.5, 99.5);
 122
 123   if (fNumberDensityPhi)
 124     fCentralityDistribution = new TH1F("fCentralityDistribution", ";;count", fNumberDensityPhi->GetEventHist()->GetNBins(1), fNumberDensityPhi->GetEventHist()->GetAxis(1, 0)->GetXbins()->GetArray());
 125
 126   fITSClusterMap = new TH3F("fITSClusterMap", "; its cluster map; centrality; pT", 256, -0.5, 255.5, 20, 0, 100.001, 100, 0, 20);
 127
 128   TH1::AddDirectory(oldStatus);
 129 }
 130
 131 //_____________________________________________________________________________
 132 AliUEHistograms::AliUEHistograms(const AliUEHistograms &c) :
 133   TNamed(fName, fTitle),
 134   fNumberDensitypT(0),
 135   fSumpT(0),
 136   fNumberDensityPhi(0),
 137   fCorrelationpT(0),
 138   fCorrelationEta(0),
 139   fCorrelationPhi(0),
 140   fCorrelationR(0),
 141   fCorrelationLeading2Phi(0),
 142   fCorrelationMultiplicity(0),
 143   fEventCount(0),
 144   fEventCountDifferential(0),
 145   fVertexContributors(0),
 146   fCentralityDistribution(0),
 147   fITSClusterMap(0),
 148   fSelectCharge(0),
 149   fRunNumber(0)
 150 {
 151   //
 152   // AliUEHistograms copy constructor
 153   //
 154
 155   fTwoTrackDistancePt[0] = 0;
 156   fTwoTrackDistancePt[1] = 0;
 157
 158   ((AliUEHistograms &) c).Copy(*this);
 159 }
 160
 161 //____________________________________________________________________
 162 AliUEHistograms::~AliUEHistograms()
 163 {
 164   // Destructor
 165
 166   if (fNumberDensitypT)
 167   {
 168     delete fNumberDensitypT;
 169     fNumberDensitypT = 0;
 170   }
 171
 172   if (fSumpT)
 173   {
 174     delete fSumpT;
 175     fSumpT = 0;
 176   }
 177
 178   if (fNumberDensityPhi)
 179   {
 180     delete fNumberDensityPhi;
 181     fNumberDensityPhi = 0;
 182   }
 183
 184   if (fCorrelationpT)
 185   {
 186     delete fCorrelationpT;
 187     fCorrelationpT = 0;
 188   }
 189
 190   if (fCorrelationEta)
 191   {
 192     delete fCorrelationEta;
 193     fCorrelationEta = 0;
 194   }
 195
 196   if (fCorrelationPhi)
 197   {
 198     delete fCorrelationPhi;
 199     fCorrelationPhi = 0;
 200   }
 201
 202   if (fCorrelationR)
 203   {
 204     delete fCorrelationR;
 205     fCorrelationR = 0;
 206   }
 207
 208   if (fCorrelationLeading2Phi)
 209   {
 210     delete fCorrelationLeading2Phi;
 211     fCorrelationLeading2Phi = 0;
 212   }
 213
 214   if (fCorrelationMultiplicity)
 215   {
 216     delete fCorrelationMultiplicity;
 217     fCorrelationMultiplicity = 0;
 218   }
 219
 220   if (fEventCount)
 221   {
 222     delete fEventCount;
 223     fEventCount = 0;
 224   }
 225
 226   if (fEventCountDifferential)
 227   {
 228     delete fEventCountDifferential;
 229     fEventCountDifferential = 0;
 230   }
 231
 232   if (fVertexContributors)
 233   {
 234     delete fVertexContributors;
 235     fVertexContributors = 0;
 236   }
 237
 238   if (fCentralityDistribution)
 239   {
 240     delete fCentralityDistribution;
 241     fCentralityDistribution = 0;
 242   }
 243
 244   if (fITSClusterMap)
 245   {
 246     delete fITSClusterMap;
 247     fITSClusterMap = 0;
 248   }
 249
 250   for (Int_t i=0; i<2; i++)
 251     if (fTwoTrackDistancePt[i])
 252     {
 253       delete fTwoTrackDistancePt[i];
 254       fTwoTrackDistancePt[i] = 0;
 255     }
 256 }
 257
 258 AliUEHist* AliUEHistograms::GetUEHist(Int_t id)
 259 {
 260   // returns AliUEHist object, useful for loops
 261
 262   switch (id)
 263   {
 264     case 0: return fNumberDensitypT; break;
 265     case 1: return fSumpT; break;
 266     case 2: return fNumberDensityPhi; break;
 267   }
 268
 269   return 0;
 270 }
 271
 272 //____________________________________________________________________
 273 Int_t AliUEHistograms::CountParticles(TList* list, Float_t ptMin)
 274 {
 275   // counts the number of particles in the list with a pT above ptMin
 276   // TODO eta cut needed here?
 277
 278   Int_t count = 0;
 279   for (Int_t j=0; j<list->GetEntries(); j++)
 280     if (((AliVParticle*) list->At(j))->Pt() > ptMin)
 281       count++;
 282
 283   return count;
 284 }
 285
 286 //____________________________________________________________________
 287 void AliUEHistograms::Fill(Int_t eventType, Float_t zVtx, AliUEHist::CFStep step, AliVParticle* leading, TList* toward, TList* away, TList* min, TList* max)
 288 {
 289   // fills the UE event histograms
 290   //
 291   // 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
 292
 293   // if leading is not set, just fill event statistics
 294   if (leading)
 295   {
 296     Int_t multiplicity = 0;
 297
 298     // TODO configurable?
 299     Float_t ptMin = 0.15;
 300     if (leading->Pt() > ptMin)
 301       multiplicity++;
 302
 303     multiplicity += CountParticles(toward, ptMin);
 304     multiplicity += CountParticles(away, ptMin);
 305     multiplicity += CountParticles(min, ptMin);
 306     multiplicity += CountParticles(max, ptMin);
 307
 308     FillRegion(AliUEHist::kToward, zVtx, step, leading, toward, multiplicity);
 309     FillRegion(AliUEHist::kAway,   zVtx, step, leading, away, multiplicity);
 310     FillRegion(AliUEHist::kMin,    zVtx, step, leading, min, multiplicity);
 311     FillRegion(AliUEHist::kMax,    zVtx, step, leading, max, multiplicity);
 312
 313     Double_t vars[3];
 314     vars[0] = leading->Pt();
 315     vars[1] = multiplicity;
 316     vars[2] = zVtx;
 317     for (Int_t i=0; i<fgkUEHists; i++)
 318       if (GetUEHist(i))
 319         GetUEHist(i)->GetEventHist()->Fill(vars, step);
 320
 321     fEventCountDifferential->Fill(leading->Pt(), step, eventType);
 322   }
 323
 324   FillEvent(eventType, step);
 325 }
 326
 327 //____________________________________________________________________
 328 void AliUEHistograms::FillRegion(AliUEHist::Region region, Float_t zVtx, AliUEHist::CFStep step, AliVParticle* leading, TList* list, Int_t multiplicity)
 329 {
 330   // loops over AliVParticles in list and fills the given region at the given step
 331   //
 332   // See also Fill(...)
 333
 334   for (Int_t i=0; i<list->GetEntries(); i++)
 335   {
 336     AliVParticle* particle = (AliVParticle*) list->At(i);
 337
 338     Double_t vars[6];
 339     vars[0] = particle->Eta();
 340     vars[1] = particle->Pt();
 341     vars[2] = leading->Pt();
 342     vars[3] = multiplicity;
 343     vars[4] = leading->Phi() - particle->Phi();
 344     if (vars[4] > 1.5 * TMath::Pi())
 345       vars[4] -= TMath::TwoPi();
 346     if (vars[4] < -0.5 * TMath::Pi())
 347       vars[4] += TMath::TwoPi();
 348     vars[5] = zVtx;
 349
 350     if (fNumberDensitypT)
 351       fNumberDensitypT->GetTrackHist(region)->Fill(vars, step);
 352
 353     if (fSumpT)
 354       fSumpT->GetTrackHist(region)->Fill(vars, step, particle->Pt());
 355
 356     // fill all in toward region (is anyway as function of delta phi!)
 357     if (fNumberDensityPhi)
 358       fNumberDensityPhi->GetTrackHist(AliUEHist::kToward)->Fill(vars, step);
 359   }
 360 }
 361
 362 //____________________________________________________________________
 363 void AliUEHistograms::Fill(AliVParticle* leadingMC, AliVParticle* leadingReco)
 364 {
 365   // fills the correlation histograms
 366
 367   if (leadingMC)
 368   {
 369     fCorrelationpT->Fill(leadingMC->Pt(), leadingReco->Pt());
 370     if (leadingMC->Pt() > 0.5)
 371     {
 372       fCorrelationEta->Fill(leadingMC->Eta(), leadingReco->Eta());
 373       fCorrelationPhi->Fill(leadingMC->Phi(), leadingReco->Phi());
 374     }
 375
 376     Float_t phiDiff = leadingMC->Phi() - leadingReco->Phi();
 377     if (phiDiff > TMath::Pi())
 378       phiDiff -= TMath::TwoPi();
 379     if (phiDiff < -TMath::Pi())
 380       phiDiff += TMath::TwoPi();
 381
 382     Float_t etaDiff = leadingMC->Eta() - leadingReco->Eta();
 383
 384     fCorrelationR->Fill(TMath::Sqrt(phiDiff * phiDiff + etaDiff * etaDiff), leadingMC->Pt());
 385     fCorrelationLeading2Phi->Fill(phiDiff, leadingMC->Pt());
 386   }
 387   else
 388   {
 389     fCorrelationpT->Fill(1.0, leadingReco->Pt());
 390     if (leadingReco->Pt() > 0.5)
 391     {
 392       fCorrelationEta->Fill(0.0, leadingReco->Eta());
 393       fCorrelationPhi->Fill(0.0, leadingReco->Phi());
 394     }
 395   }
 396 }
 397
 398 //____________________________________________________________________
 399 void AliUEHistograms::FillCorrelations(Double_t centrality, Float_t zVtx, AliUEHist::CFStep step, TObjArray* particles, TObjArray* mixed, Float_t weight, Bool_t firstTime)
 400 {
 401   // fills the fNumberDensityPhi histogram
 402   //
 403   // this function need a list of AliVParticles which contain the particles/tracks to be filled
 404   //
 405   // if mixed is non-0, mixed events are filled, the trigger particle is from particles, the associated from mixed
 406   // if weight < 0, then the pt of the associated particle is filled as weight
 407
 408   Bool_t fillpT = kFALSE;
 409   if (weight < 0)
 410     fillpT = kTRUE;
 411
 412   // Eta() is extremely time consuming, therefore cache it for the inner loop here:
 413   TObjArray* input = (mixed) ? mixed : particles;
 414   TArrayF eta(input->GetEntriesFast());
 415   for (Int_t i=0; i<input->GetEntriesFast(); i++)
 416     eta[i] = ((AliVParticle*) input->At(i))->Eta();
 417
 418   // if particles is not set, just fill event statistics
 419   if (particles)
 420   {
 421     Int_t jMax = particles->GetEntriesFast();
 422     if (mixed)
 423       jMax = mixed->GetEntriesFast();
 424
 425     for (Int_t i=0; i<particles->GetEntriesFast(); i++)
 426     {
 427       AliVParticle* triggerParticle = (AliVParticle*) particles->At(i);
 428
 429       // some optimization
 430       Float_t triggerEta = triggerParticle->Eta();
 431
 432       if (!mixed)
 433       {
 434         // QA
 435         fCorrelationpT->Fill(centrality, triggerParticle->Pt());
 436         fCorrelationEta->Fill(centrality, triggerEta);
 437         fCorrelationPhi->Fill(centrality, triggerParticle->Phi());
 438 /*        if (dynamic_cast<AliAODTrack*>(triggerParticle))
 439           fITSClusterMap->Fill(((AliAODTrack*) triggerParticle)->GetITSClusterMap(), centrality, triggerParticle->Pt());*/
 440       }
 441
 442       for (Int_t j=0; j<jMax; j++)
 443       {
 444         if (!mixed && i == j)
 445           continue;
 446
 447         AliVParticle* particle = 0;
 448         if (!mixed)
 449           particle = (AliVParticle*) particles->At(j);
 450         else
 451           particle = (AliVParticle*) mixed->At(j);
 452
 453         // check if both particles point to the same element (does not occur for mixed events, but if subsets are mixed within the same event for cross-checks)
 454         if (mixed && triggerParticle == particle)
 455           continue;
 456
 457         if (particle->Pt() > triggerParticle->Pt())
 458           continue;
 459
 460         if (fSelectCharge > 0)
 461         {
 462           // skip like sign
 463           if (fSelectCharge == 1 && particle->Charge() * triggerParticle->Charge() > 0)
 464             continue;
 465
 466           // skip unlike sign
 467           if (fSelectCharge == 2 && particle->Charge() * triggerParticle->Charge() < 0)
 468             continue;
 469         }
 470
 471         Double_t vars[6];
 472         vars[0] = triggerEta - eta[j];
 473         vars[1] = particle->Pt();
 474         vars[2] = triggerParticle->Pt();
 475         vars[3] = centrality;
 476         vars[4] = triggerParticle->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 (fillpT)
 484           weight = particle->Pt();
 485
 486         // fill all in toward region and do not use the other regions
 487         fNumberDensityPhi->GetTrackHist(AliUEHist::kToward)->Fill(vars, step, weight);
 488       }
 489
 490       if (firstTime)
 491       {
 492         // once per trigger particle
 493         Double_t vars[3];
 494         vars[0] = triggerParticle->Pt();
 495         vars[1] = centrality;
 496         vars[2] = zVtx;
 497         fNumberDensityPhi->GetEventHist()->Fill(vars, step);
 498       }
 499     }
 500   }
 501
 502   fCentralityDistribution->Fill(centrality);
 503   FillEvent(centrality, step);
 504 }
 505
 506 //____________________________________________________________________
 507 void AliUEHistograms::FillTrackingEfficiency(TObjArray* mc, TObjArray* recoPrim, TObjArray* recoAll, Int_t particleType, Double_t centrality)
 508 {
 509   // fills the tracking efficiency objects
 510   //
 511   // mc: all primary MC particles
 512   // recoPrim: reconstructed primaries (again MC particles)
 513   // recoAll: reconstructed (again MC particles)
 514   // particleType is: 0 for pion, 1 for kaon, 2 for proton, 3 for others
 515
 516   for (Int_t step=0; step<3; step++)
 517   {
 518     TObjArray* list = mc;
 519     if (step == 1)
 520       list = recoPrim;
 521     else if (step == 2)
 522       list = recoAll;
 523
 524     for (Int_t i=0; i<list->GetEntriesFast(); i++)
 525     {
 526       AliVParticle* particle = (AliVParticle*) list->At(i);
 527       Double_t vars[4];
 528       vars[0] = particle->Eta();
 529       vars[1] = particle->Pt();
 530       vars[2] = particleType;
 531       vars[3] = centrality;
 532
 533       for (Int_t j=0; j<fgkUEHists; j++)
 534         if (GetUEHist(j))
 535           GetUEHist(j)->GetTrackHistEfficiency()->Fill(vars, step);
 536     }
 537   }
 538 }
 539
 540 //____________________________________________________________________
 541 void AliUEHistograms::FillEvent(Int_t eventType, Int_t step)
 542 {
 543   // fills the number of events at the given step and the given enty type
 544   //
 545   // WARNING: This function is called from Fill, so only call it for steps where Fill is not called
 546
 547   fEventCount->Fill(step, eventType);
 548 }
 549
 550 //____________________________________________________________________
 551 void AliUEHistograms::FillEvent(Double_t centrality, Int_t step)
 552 {
 553   // fills the number of events at the given step and the given centrality
 554   //
 555   // WARNING: This function is called from Fill, so only call it for steps where Fill is not called
 556
 557   fEventCount->Fill(step, centrality);
 558 }
 559
 560 //____________________________________________________________________
 561 void AliUEHistograms::SetEtaRange(Float_t etaMin, Float_t etaMax)
 562 {
 563   // sets eta min and max for all contained AliUEHist classes
 564
 565   for (Int_t i=0; i<fgkUEHists; i++)
 566     if (GetUEHist(i))
 567       GetUEHist(i)->SetEtaRange(etaMin, etaMax);
 568 }
 569
 570 //____________________________________________________________________
 571 void AliUEHistograms::SetPtRange(Float_t ptMin, Float_t ptMax)
 572 {
 573   // sets pT min and max for all contained AliUEHist classes
 574
 575   for (Int_t i=0; i<fgkUEHists; i++)
 576     if (GetUEHist(i))
 577       GetUEHist(i)->SetPtRange(ptMin, ptMax);
 578 }
 579
 580 //____________________________________________________________________
 581 void AliUEHistograms::SetZVtxRange(Float_t min, Float_t max)
 582 {
 583   // sets pT min and max for all contained AliUEHist classes
 584
 585   for (Int_t i=0; i<fgkUEHists; i++)
 586     if (GetUEHist(i))
 587       GetUEHist(i)->SetZVtxRange(min, max);
 588 }
 589
 590 //____________________________________________________________________
 591 void AliUEHistograms::SetContaminationEnhancement(TH1F* hist)
 592 {
 593   // sets the contamination enhancement histogram in all contained AliUEHist classes
 594
 595   for (Int_t i=0; i<fgkUEHists; i++)
 596     if (GetUEHist(i))
 597       GetUEHist(i)->SetContaminationEnhancement(hist);
 598 }
 599
 600 //____________________________________________________________________
 601 void AliUEHistograms::SetCombineMinMax(Bool_t flag)
 602 {
 603   // sets pT min and max for all contained AliUEHist classes
 604
 605   for (Int_t i=0; i<fgkUEHists; i++)
 606     if (GetUEHist(i))
 607       GetUEHist(i)->SetCombineMinMax(flag);
 608 }
 609
 610 //____________________________________________________________________
 611 void AliUEHistograms::Correct(AliUEHistograms* corrections)
 612 {
 613   // corrects the contained histograms by calling AliUEHist::Correct
 614
 615   for (Int_t i=0; i<fgkUEHists; i++)
 616     if (GetUEHist(i))
 617       GetUEHist(i)->Correct(corrections->GetUEHist(i));
 618 }
 619
 620 //____________________________________________________________________
 621 AliUEHistograms &AliUEHistograms::operator=(const AliUEHistograms &c)
 622 {
 623   // assigment operator
 624
 625   if (this != &c)
 626     ((AliUEHistograms &) c).Copy(*this);
 627
 628   return *this;
 629 }
 630
 631 //____________________________________________________________________
 632 void AliUEHistograms::Copy(TObject& c) const
 633 {
 634   // copy function
 635
 636   AliUEHistograms& target = (AliUEHistograms &) c;
 637
 638   if (fNumberDensitypT)
 639     target.fNumberDensitypT = dynamic_cast<AliUEHist*> (fNumberDensitypT->Clone());
 640
 641   if (fSumpT)
 642     target.fSumpT = dynamic_cast<AliUEHist*> (fSumpT->Clone());
 643
 644   if (fNumberDensityPhi)
 645     target.fNumberDensityPhi = dynamic_cast<AliUEHist*> (fNumberDensityPhi->Clone());
 646
 647   if (fCorrelationpT)
 648     target.fCorrelationpT = dynamic_cast<TH2F*> (fCorrelationpT->Clone());
 649
 650   if (fCorrelationEta)
 651     target.fCorrelationEta = dynamic_cast<TH2F*> (fCorrelationEta->Clone());
 652
 653   if (fCorrelationPhi)
 654     target.fCorrelationPhi = dynamic_cast<TH2F*> (fCorrelationPhi->Clone());
 655
 656   if (fCorrelationR)
 657     target.fCorrelationR = dynamic_cast<TH2F*> (fCorrelationR->Clone());
 658
 659   if (fCorrelationLeading2Phi)
 660     target.fCorrelationLeading2Phi = dynamic_cast<TH2F*> (fCorrelationLeading2Phi->Clone());
 661
 662   if (fCorrelationMultiplicity)
 663     target.fCorrelationMultiplicity = dynamic_cast<TH2F*> (fCorrelationMultiplicity->Clone());
 664
 665   if (fEventCount)
 666     target.fEventCount = dynamic_cast<TH2F*> (fEventCount->Clone());
 667
 668   if (fEventCountDifferential)
 669     target.fEventCountDifferential = dynamic_cast<TH3F*> (fEventCountDifferential->Clone());
 670
 671   if (fVertexContributors)
 672     target.fVertexContributors = dynamic_cast<TH1F*> (fVertexContributors->Clone());
 673
 674   if (fCentralityDistribution)
 675     target.fCentralityDistribution = dynamic_cast<TH1F*> (fCentralityDistribution->Clone());
 676
 677   if (fITSClusterMap)
 678     target.fITSClusterMap = dynamic_cast<TH3F*> (fITSClusterMap->Clone());
 679
 680   for (Int_t i=0; i<2; i++)
 681     if (fTwoTrackDistancePt[i])
 682       target.fTwoTrackDistancePt[i] = dynamic_cast<TH3F*> (fTwoTrackDistancePt[i]->Clone());
 683
 684   target.fSelectCharge = fSelectCharge;
 685   target.fRunNumber = fRunNumber;
 686 }
 687
 688 //____________________________________________________________________
 689 Long64_t AliUEHistograms::Merge(TCollection* list)
 690 {
 691   // Merge a list of AliUEHistograms objects with this (needed for
 692   // PROOF).
 693   // Returns the number of merged objects (including this).
 694
 695   if (!list)
 696     return 0;
 697
 698   if (list->IsEmpty())
 699     return 1;
 700
 701   TIterator* iter = list->MakeIterator();
 702   TObject* obj;
 703
 704   // collections of objects
 705   const Int_t kMaxLists = 16;
 706   TList* lists[kMaxLists];
 707
 708   for (Int_t i=0; i<kMaxLists; i++)
 709     lists[i] = new TList;
 710
 711   Int_t count = 0;
 712   while ((obj = iter->Next())) {
 713
 714     AliUEHistograms* entry = dynamic_cast<AliUEHistograms*> (obj);
 715     if (entry == 0)
 716       continue;
 717
 718     if (entry->fNumberDensitypT)
 719       lists[0]->Add(entry->fNumberDensitypT);
 720     if (entry->fSumpT)
 721       lists[1]->Add(entry->fSumpT);
 722     if (entry->fNumberDensityPhi)
 723       lists[2]->Add(entry->fNumberDensityPhi);
 724     lists[3]->Add(entry->fCorrelationpT);
 725     lists[4]->Add(entry->fCorrelationEta);
 726     lists[5]->Add(entry->fCorrelationPhi);
 727     lists[6]->Add(entry->fCorrelationR);
 728     lists[7]->Add(entry->fCorrelationLeading2Phi);
 729     lists[8]->Add(entry->fCorrelationMultiplicity);
 730     lists[9]->Add(entry->fEventCount);
 731     lists[10]->Add(entry->fEventCountDifferential);
 732     lists[11]->Add(entry->fVertexContributors);
 733     lists[12]->Add(entry->fCentralityDistribution);
 734     lists[13]->Add(entry->fITSClusterMap);
 735     if (fTwoTrackDistancePt[0])
 736       lists[14]->Add(entry->fTwoTrackDistancePt[0]);
 737     if (fTwoTrackDistancePt[1])
 738       lists[15]->Add(entry->fTwoTrackDistancePt[1]);
 739
 740     count++;
 741   }
 742
 743   if (fNumberDensitypT)
 744     fNumberDensitypT->Merge(lists[0]);
 745   if (fSumpT)
 746     fSumpT->Merge(lists[1]);
 747   if (fNumberDensityPhi)
 748     fNumberDensityPhi->Merge(lists[2]);
 749   fCorrelationpT->Merge(lists[3]);
 750   fCorrelationEta->Merge(lists[4]);
 751   fCorrelationPhi->Merge(lists[5]);
 752   fCorrelationR->Merge(lists[6]);
 753   fCorrelationLeading2Phi->Merge(lists[7]);
 754   fCorrelationMultiplicity->Merge(lists[8]);
 755   fEventCount->Merge(lists[9]);
 756   fEventCountDifferential->Merge(lists[10]);
 757   fVertexContributors->Merge(lists[11]);
 758   fCentralityDistribution->Merge(lists[12]);
 759   fITSClusterMap->Merge(lists[13]);
 760   if (fTwoTrackDistancePt[0])
 761     fTwoTrackDistancePt[0]->Merge(lists[14]);
 762   if (fTwoTrackDistancePt[1])
 763     fTwoTrackDistancePt[1]->Merge(lists[15]);
 764
 765   for (Int_t i=0; i<kMaxLists; i++)
 766     delete lists[i];
 767
 768   return count+1;
 769 }
 770
 771 void AliUEHistograms::CopyReconstructedData(AliUEHistograms* from)
 772 {
 773   // copies those histograms extracted from ESD to this object
 774
 775   for (Int_t i=0; i<fgkUEHists; i++)
 776     if (GetUEHist(i))
 777       GetUEHist(i)->CopyReconstructedData(from->GetUEHist(i));
 778 }
 779
 780 void AliUEHistograms::ExtendTrackingEfficiency(Bool_t verbose)
 781 {
 782   // delegates to AliUEHists
 783
 784   for (Int_t i=0; i<fgkUEHists; i++)
 785     if (GetUEHist(i))
 786       GetUEHist(i)->ExtendTrackingEfficiency(verbose);
 787 }
 788
 789 void AliUEHistograms::Scale(Double_t factor)
 790 {
 791   // scales all contained histograms by the given factor
 792
 793   for (Int_t i=0; i<fgkUEHists; i++)
 794     if (GetUEHist(i))
 795       GetUEHist(i)->Scale(factor);
 796
 797   TList list;
 798   list.Add(fCorrelationpT);
 799   list.Add(fCorrelationEta);
 800   list.Add(fCorrelationPhi);
 801   list.Add(fCorrelationR);
 802   list.Add(fCorrelationLeading2Phi);
 803   list.Add(fCorrelationMultiplicity);
 804   list.Add(fEventCount);
 805   list.Add(fEventCountDifferential);
 806   list.Add(fVertexContributors);
 807   list.Add(fCentralityDistribution);
 808   list.Add(fITSClusterMap);
 809   list.Add(fTwoTrackDistancePt[0]);
 810   list.Add(fTwoTrackDistancePt[1]);
 811
 812   for (Int_t i=0; i<list.GetEntries(); i++)
 813     ((TH1*) list.At(i))->Scale(factor);
 814 }
 815
 816 void AliUEHistograms::Reset()
 817 {
 818   // delegates to AliUEHists
 819
 820   for (Int_t i=0; i<fgkUEHists; i++)
 821     if (GetUEHist(i))
 822       GetUEHist(i)->Reset();
 823 }
 824
 825 TObjArray* AliUEHistograms::ApplyTwoTrackCut(TObjArray* tracks)
 826 {
 827   // takes the input list <tracks> and applies two-track efficiency cuts
 828   // returns the tracks which pass the cuts (if a pair fails the cut, both are removed)
 829   // while the cut is applied, control histograms are filled: fTwoTrackDistancePt[i] (i = 0 before, i = 1 after)
 830   // the cut has been developed by the HBT group and removes tracks which are spatially close inside the TPC volume
 831   // see https://indico.cern.ch/materialDisplay.py?contribId=36&sessionId=6&materialId=slides&confId=142700
 832
 833   if (!fTwoTrackDistancePt[0])
 834   {
 835     fTwoTrackDistancePt[0] = new TH3F("fTwoTrackDistancePt[0]", ";#Delta#eta;#Delta#varphi^{*}_{min};#Delta p_{T}", 100, -0.05, 0.05, 400, -0.2, 0.2, 20, 0, 10);
 836     fTwoTrackDistancePt[1] = (TH3F*) fTwoTrackDistancePt[0]->Clone("fTwoTrackDistancePt[1]");
 837   }
 838
 839   TObjArray* accepted = new TObjArray(*tracks);
 840
 841   // Eta() is extremely time consuming, therefore cache it for the inner loop here:
 842   TArrayF eta(tracks->GetEntriesFast());
 843   for (Int_t i=0; i<tracks->GetEntriesFast(); i++)
 844     eta[i] = ((AliVParticle*) tracks->At(i))->Eta();
 845
 846   for (Int_t i=0; i<tracks->GetEntriesFast(); i++)
 847   {
 848     AliVParticle* particle1 = (AliVParticle*) tracks->At(i);
 849     Float_t phi1 = particle1->Phi();
 850     Float_t pt1 = particle1->Pt();
 851     Float_t charge1 = particle1->Charge();
 852
 853     // analyze region for IAA paper
 854     if (pt1 < 8 || pt1 > 15)
 855       continue;
 856
 857     for (Int_t j=0; j<tracks->GetEntriesFast(); j++)
 858     {
 859       if (i == j)
 860         continue;
 861
 862       AliVParticle* particle2 = (AliVParticle*) tracks->At(j);
 863       Float_t phi2 = particle2->Phi();
 864       Float_t pt2 = particle2->Pt();
 865       Float_t charge2 = particle2->Charge();
 866
 867       if (pt2 > pt1)
 868         continue;
 869
 870 //       Double_t dpt = TMath::Abs(pt1 - pt2);
 871       Float_t deta = eta[i] - eta[j];
 872       Float_t detaabs = TMath::Abs(deta);
 873
 874       // optimization
 875 /*      if (detaabs > 0.05)
 876         continue;*/
 877
 878       Bool_t cutPassed = kTRUE;
 879
 880       Float_t dphistarmin = 1e5;
 881       Float_t dphistarminabs = 1e5;
 882
 883       for (Double_t rad=0.8; rad<2.51; rad+=0.01)
 884       {
 885         Float_t dphistar = phi1 - phi2 - TMath::ASin(charge1 * 0.075 * rad / pt1) + TMath::ASin(charge2 * 0.075 * rad / pt2);
 886         Float_t dphistarabs = TMath::Abs(dphistar);
 887
 888         if (dphistarabs < dphistarminabs)
 889         {
 890           dphistarmin = dphistar;
 891           dphistarminabs = dphistarabs;
 892         }
 893       }
 894
 895       // hardcoded cut values for the moment
 896       if (detaabs < 0.011 && dphistarminabs < 0.01)
 897       {
 898         cutPassed = kFALSE;
 899         //Printf("%d %d failed: %.3f %.3f; %.3f %.3f %.4f; %.2f %.2f (%p %p)", i, j, eta[i], eta[j], phi1, phi2, dphistarminabs, pt1, pt2, particle1, particle2);
 900       }
 901
 902       fTwoTrackDistancePt[0]->Fill(deta, dphistarmin, pt2);
 903       if (cutPassed)
 904         fTwoTrackDistancePt[1]->Fill(deta, dphistarmin, pt2);
 905       else
 906       {
 907         // remove tracks from list
 908         accepted->Remove(particle1);
 909         accepted->Remove(particle2);
 910       }
 911     }
 912   }
 913
 914   accepted->Compress();
 915
 916   //Printf("AliUEHistograms::ApplyTwoTrackCut: Accepted %d out of %d tracks", accepted->GetEntriesFast(), tracks->GetEntriesFast());
 917
 918   return accepted;
 919 }