PWGCF/Correlations/Base/AliUEHistograms.h

   1 #ifndef AliUEHistograms_H
   2 #define AliUEHistograms_H
   3
   4 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   5  * See cxx source for full Copyright notice                               */
   6
   7 /* $Id: AliUEHistograms.h 20164 2007-08-14 15:31:50Z morsch $ */
   8
   9 // encapsulates several AliUEHist objects for a full UE analysis plus additional control histograms
  10
  11 #include "TNamed.h"
  12 #include "AliUEHist.h"
  13 #include "TMath.h"
  14 #include "THn.h" // in cxx file causes .../THn.h:257: error: conflicting declaration ‘typedef class THnT<float> THnF’
  15
  16 class AliVParticle;
  17
  18 class TList;
  19 class TSeqCollection;
  20 class TObjArray;
  21 class TH1F;
  22 class TH2F;
  23 class TH3F;
  24
  25 class AliUEHistograms : public TNamed
  26 {
  27  public:
  28   AliUEHistograms(const char* name = "AliUEHistograms", const char* histograms = "", const char* binning = 0);
  29   virtual ~AliUEHistograms();
  30
  31   void Fill(Int_t eventType, Float_t zVtx, AliUEHist::CFStep step, AliVParticle* leading, TList* toward, TList* away, TList* min, TList* max);
  32   void FillCorrelations(Double_t centrality, Float_t zVtx, AliUEHist::CFStep step, TObjArray* particles, TObjArray* mixed = 0, Float_t weight = 1, Bool_t firstTime = kTRUE, Bool_t twoTrackEfficiencyCut = kFALSE, Float_t bSign = 0, Float_t twoTrackEfficiencyCutValue = 0.02, Bool_t applyEfficiency = kFALSE);
  33   void Fill(AliVParticle* leadingMC, AliVParticle* leadingReco);
  34   void FillEvent(Int_t eventType, Int_t step);
  35   void FillEvent(Double_t centrality, Int_t step);
  36   void FillTrackingEfficiency(TObjArray* mc, TObjArray* recoPrim, TObjArray* recoAll, TObjArray* recoPrimPID, TObjArray* recoAllPID, TObjArray* fake, Int_t particleType, Double_t centrality = 0, Double_t zVtx = 0);
  37   void FillFakePt(TObjArray* fake, Double_t centrality);
  38
  39   void CopyReconstructedData(AliUEHistograms* from);
  40   void DeepCopy(AliUEHistograms* from);
  41
  42   AliUEHist* GetUEHist(Int_t id);
  43
  44   AliUEHist* GetNumberDensitypT() { return fNumberDensitypT; }
  45   AliUEHist* GetSumpT() { return fSumpT; }
  46   AliUEHist* GetNumberDensityPhi() { return fNumberDensityPhi; }
  47
  48   void SetNumberDensitypT(AliUEHist* obj) { fNumberDensitypT = obj; }
  49   void SetSumpT(AliUEHist* obj) { fSumpT = obj; }
  50   void SetNumberDensityPhi(AliUEHist* obj) { fNumberDensityPhi = obj; }
  51
  52   void SetRunNumber(Long64_t runNumber) { fRunNumber = runNumber; }
  53
  54   void SetEfficiencyCorrectionTriggers(THnF* hist)   { fEfficiencyCorrectionTriggers = hist;   }
  55   void SetEfficiencyCorrectionAssociated(THnF* hist) { fEfficiencyCorrectionAssociated = hist; }
  56
  57   TH2F* GetCorrelationpT()  { return fCorrelationpT; }
  58   TH2F* GetCorrelationEta() { return fCorrelationEta; }
  59   TH2F* GetCorrelationPhi() { return fCorrelationPhi; }
  60   TH2F* GetCorrelationR()   { return fCorrelationR; }
  61   TH2F* GetCorrelationLeading2Phi() { return fCorrelationLeading2Phi; }
  62   TH2F* GetCorrelationMultiplicity() { return fCorrelationMultiplicity; }
  63   TH3F* GetYield() { return fYields; }
  64   TH2F* GetInvYield() { return fInvYield2; }
  65
  66   TH2F* GetEventCount()     { return fEventCount; }
  67   TH3F* GetEventCountDifferential() { return fEventCountDifferential; }
  68   TH1F* GetVertexContributors() { return fVertexContributors; }
  69   TH1F* GetCentralityDistribution() { return fCentralityDistribution; }
  70   TH2F* GetCentralityCorrelation() { return fCentralityCorrelation; }
  71   Long64_t GetRunNumber() { return fRunNumber; }
  72   Int_t GetMergeCount() { return fMergeCount; }
  73   TH3F* GetTwoTrackDistance(Int_t i) { return fTwoTrackDistancePt[i]; }
  74   Bool_t GetWeightPerEvent() { return fWeightPerEvent; }
  75
  76   void Correct(AliUEHistograms* corrections);
  77
  78   void SetEtaRange(Float_t etaMin, Float_t etaMax);
  79   void SetPtRange(Float_t ptMin, Float_t ptMax);
  80   void SetPartSpecies(Int_t species);
  81   void SetZVtxRange(Float_t min, Float_t max);
  82   void SetContaminationEnhancement(TH1F* hist);
  83   void SetCombineMinMax(Bool_t flag);
  84   void SetTrackEtaCut(Float_t value);
  85   void SetWeightPerEvent(Bool_t flag);
  86   void SetSelectCharge(Int_t selectCharge) { fSelectCharge = selectCharge; }
  87   void SetSelectTriggerCharge(Int_t selectCharge) { fTriggerSelectCharge = selectCharge; }
  88   void SetSelectAssociatedCharge(Int_t selectCharge) { fAssociatedSelectCharge = selectCharge; }
  89   void SetTriggerRestrictEta(Float_t eta) { fTriggerRestrictEta = eta; }
  90   void SetEtaOrdering(Bool_t flag) { fEtaOrdering = flag; }
  91   void SetPairCuts(Bool_t conversions, Bool_t resonances) { fCutConversions = conversions; fCutResonances = resonances; }
  92   void SetRejectResonanceDaughters(Int_t value) { fRejectResonanceDaughters = value; }
  93   void SetOnlyOneEtaSide(Int_t flag)    { fOnlyOneEtaSide = flag; }
  94   void SetPtOrder(Bool_t flag) { fPtOrder = flag; }
  95   void SetTwoTrackCutMinRadius(Float_t min) { fTwoTrackCutMinRadius = min; }
  96
  97   void ExtendTrackingEfficiency(Bool_t verbose = kFALSE);
  98   void Reset();
  99
 100   AliUEHistograms(const AliUEHistograms &c);
 101   AliUEHistograms& operator=(const AliUEHistograms& c);
 102   virtual void Copy(TObject& c) const;
 103
 104   virtual Long64_t Merge(TCollection* list);
 105   void Scale(Double_t factor);
 106
 107 protected:
 108   void FillRegion(AliUEHist::Region region, Float_t zVtx, AliUEHist::CFStep step, AliVParticle* leading, TList* list, Int_t multiplicity);
 109   Int_t CountParticles(TList* list, Float_t ptMin);
 110   void DeleteContainers();
 111   inline Float_t GetInvMassSquared(Float_t pt1, Float_t eta1, Float_t phi1, Float_t pt2, Float_t eta2, Float_t phi2, Float_t m0_1, Float_t m0_2);
 112   inline Float_t GetInvMassSquaredCheap(Float_t pt1, Float_t eta1, Float_t phi1, Float_t pt2, Float_t eta2, Float_t phi2, Float_t m0_1, Float_t m0_2);
 113   inline Float_t GetDPhiStar(Float_t phi1, Float_t pt1, Float_t charge1, Float_t phi2, Float_t pt2, Float_t charge2, Float_t radius, Float_t bSign);
 114
 115   static const Int_t fgkUEHists; // number of histograms
 116
 117   AliUEHist* fNumberDensitypT;   // d^2N/dphideta vs pT,lead
 118   AliUEHist* fSumpT;             // d^2 sum(pT)/dphideta vs pT,lead
 119   AliUEHist* fNumberDensityPhi;  // d^2N/dphideta vs delta phi,lead (in pT,lead bins)
 120
 121   TH2F* fCorrelationpT;         // pT,lead: true vs reco
 122   TH2F* fCorrelationEta;        // #eta,lead; true vs reco
 123   TH2F* fCorrelationPhi;        // #phi,lead; true vs reco
 124   TH2F* fCorrelationR;          // R = sqrt(delta eta^2 + delta phi^2) (true vs reco) vs pT,lead,MC
 125   TH2F* fCorrelationLeading2Phi;// delta phi (true vs reco) vs pT,lead,MC
 126   TH2F* fCorrelationMultiplicity; // number of mc particls vs reco particles (for pT > 0.5 GeV/c)
 127   TH3F* fYields;                // centrality vs pT vs eta
 128   TH2F* fInvYield2;             // invariant yield as cross check of tracking
 129
 130   TH2F* fEventCount;            // event count as function of step, (for pp: event type (plus additional step -1 for all events without vertex range even in MC)) (for PbPb: centrality)
 131   TH3F* fEventCountDifferential;// event count as function of leading pT, step, event type
 132
 133   TH1F* fVertexContributors;    // number of contributors to the vertex
 134   TH1F* fCentralityDistribution; // distribution of the variable used for centrality selection
 135   TH2F* fCentralityCorrelation;  // centrality vs multiplicity
 136
 137   TH3F* fITSClusterMap;          // its cluster map vs centrality vs pT
 138
 139   TH3F* fTwoTrackDistancePt[2];    // control histograms for two-track efficiency study: dphi*_min vs deta (0 = before cut, 1 = after cut)
 140   TH2F* fControlConvResoncances; // control histograms for cuts on conversions and resonances
 141
 142   THnF* fEfficiencyCorrectionTriggers;   // if non-0 this efficiency correction is applied on the fly to the filling for trigger particles. The factor is multiplicative, i.e. should contain 1/efficiency
 143   THnF* fEfficiencyCorrectionAssociated;   // if non-0 this efficiency correction is applied on the fly to the filling for associated particles. The factor is multiplicative, i.e. should contain 1/efficiency
 144
 145   Int_t fSelectCharge;           // (un)like sign selection when building correlations: 0: no selection; 1: unlike sign; 2: like sign
 146   Int_t fTriggerSelectCharge;    // select charge of trigger particle
 147   Int_t fAssociatedSelectCharge; // select charge of associated particle
 148   Float_t fTriggerRestrictEta;   // restrict eta range for trigger particle (default: -1 [off])
 149   Bool_t fEtaOrdering;           // activate eta ordering to prevent shape distortions. see FillCorrelation for the details
 150   Bool_t fCutConversions;        // cut on conversions (inv mass)
 151   Bool_t fCutResonances;         // cut on resonances (inv mass)
 152   Int_t fRejectResonanceDaughters; // reject all daughters of all resonance candidates (1: test method (cut at m_inv=0.9); 2: k0; 3: lambda)
 153   Int_t fOnlyOneEtaSide;       // decides that only trigger particle from one eta side are considered (0 = all; -1 = negative, 1 = positive)
 154   Bool_t fWeightPerEvent;       // weight with the number of trigger particles per event
 155   Bool_t fPtOrder;              // apply pT,a < pT,t condition
 156   Float_t fTwoTrackCutMinRadius; // min radius for TTR cut
 157
 158   Long64_t fRunNumber;           // run number that has been processed
 159
 160   Int_t fMergeCount;            // counts how many objects have been merged together
 161
 162   ClassDef(AliUEHistograms, 29)  // underlying event histogram container
 163 };
 164
 165 Float_t AliUEHistograms::GetDPhiStar(Float_t phi1, Float_t pt1, Float_t charge1, Float_t phi2, Float_t pt2, Float_t charge2, Float_t radius, Float_t bSign)
 166 {
 167   //
 168   // calculates dphistar
 169   //
 170
 171   Float_t dphistar = phi1 - phi2 - charge1 * bSign * TMath::ASin(0.075 * radius / pt1) + charge2 * bSign * TMath::ASin(0.075 * radius / pt2);
 172
 173   static const Double_t kPi = TMath::Pi();
 174
 175   // circularity
 176 //   if (dphistar > 2 * kPi)
 177 //     dphistar -= 2 * kPi;
 178 //   if (dphistar < -2 * kPi)
 179 //     dphistar += 2 * kPi;
 180
 181   if (dphistar > kPi)
 182     dphistar = kPi * 2 - dphistar;
 183   if (dphistar < -kPi)
 184     dphistar = -kPi * 2 - dphistar;
 185   if (dphistar > kPi) // might look funny but is needed
 186     dphistar = kPi * 2 - dphistar;
 187
 188   return dphistar;
 189 }
 190
 191 Float_t AliUEHistograms::GetInvMassSquared(Float_t pt1, Float_t eta1, Float_t phi1, Float_t pt2, Float_t eta2, Float_t phi2, Float_t m0_1, Float_t m0_2)
 192 {
 193   // calculate inv mass squared
 194   // same can be achieved, but with more computing time with
 195   /*TLorentzVector photon, p1, p2;
 196   p1.SetPtEtaPhiM(triggerParticle->Pt(), triggerEta, triggerParticle->Phi(), 0.510e-3);
 197   p2.SetPtEtaPhiM(particle->Pt(), eta[j], particle->Phi(), 0.510e-3);
 198   photon = p1+p2;
 199   photon.M()*/
 200
 201   Float_t tantheta1 = 1e10;
 202
 203   if (eta1 < -1e-10 || eta1 > 1e-10)
 204   {
 205     Float_t expTmp = TMath::Exp(-eta1);
 206     tantheta1 = 2.0 * expTmp / ( 1.0 - expTmp*expTmp);
 207   }
 208
 209   Float_t tantheta2 = 1e10;
 210   if (eta2 < -1e-10 || eta2 > 1e-10)
 211   {
 212     Float_t expTmp = TMath::Exp(-eta2);
 213     tantheta2 = 2.0 * expTmp / ( 1.0 - expTmp*expTmp);
 214   }
 215
 216   Float_t e1squ = m0_1 * m0_1 + pt1 * pt1 * (1.0 + 1.0 / tantheta1 / tantheta1);
 217   Float_t e2squ = m0_2 * m0_2 + pt2 * pt2 * (1.0 + 1.0 / tantheta2 / tantheta2);
 218
 219   Float_t mass2 = m0_1 * m0_1 + m0_2 * m0_2 + 2 * ( TMath::Sqrt(e1squ * e2squ) - ( pt1 * pt2 * ( TMath::Cos(phi1 - phi2) + 1.0 / tantheta1 / tantheta2 ) ) );
 220
 221 //   Printf(Form("%f %f %f %f %f %f %f %f %f", pt1, eta1, phi1, pt2, eta2, phi2, m0_1, m0_2, mass2));
 222
 223   return mass2;
 224 }
 225
 226 Float_t AliUEHistograms::GetInvMassSquaredCheap(Float_t pt1, Float_t eta1, Float_t phi1, Float_t pt2, Float_t eta2, Float_t phi2, Float_t m0_1, Float_t m0_2)
 227 {
 228   // calculate inv mass squared approximately
 229
 230   Float_t tantheta1 = 1e10;
 231
 232   if (eta1 < -1e-10 || eta1 > 1e-10)
 233   {
 234     Float_t expTmp = 1.0-eta1+eta1*eta1/2-eta1*eta1*eta1/6+eta1*eta1*eta1*eta1/24;
 235     tantheta1 = 2.0 * expTmp / ( 1.0 - expTmp*expTmp);
 236   }
 237
 238   Float_t tantheta2 = 1e10;
 239   if (eta2 < -1e-10 || eta2 > 1e-10)
 240   {
 241     Float_t expTmp = 1.0-eta2+eta2*eta2/2-eta2*eta2*eta2/6+eta2*eta2*eta2*eta2/24;
 242     tantheta2 = 2.0 * expTmp / ( 1.0 - expTmp*expTmp);
 243   }
 244
 245   Float_t e1squ = m0_1 * m0_1 + pt1 * pt1 * (1.0 + 1.0 / tantheta1 / tantheta1);
 246   Float_t e2squ = m0_2 * m0_2 + pt2 * pt2 * (1.0 + 1.0 / tantheta2 / tantheta2);
 247
 248   // fold onto 0...pi
 249   Float_t deltaPhi = TMath::Abs(phi1 - phi2);
 250   while (deltaPhi > TMath::TwoPi())
 251     deltaPhi -= TMath::TwoPi();
 252   if (deltaPhi > TMath::Pi())
 253     deltaPhi = TMath::TwoPi() - deltaPhi;
 254
 255   Float_t cosDeltaPhi = 0;
 256   if (deltaPhi < TMath::Pi()/3)
 257     cosDeltaPhi = 1.0 - deltaPhi*deltaPhi/2 + deltaPhi*deltaPhi*deltaPhi*deltaPhi/24;
 258   else if (deltaPhi < 2*TMath::Pi()/3)
 259     cosDeltaPhi = -(deltaPhi - TMath::Pi()/2) + 1.0/6 * TMath::Power((deltaPhi - TMath::Pi()/2), 3);
 260   else
 261     cosDeltaPhi = -1.0 + 1.0/2.0*(deltaPhi - TMath::Pi())*(deltaPhi - TMath::Pi()) - 1.0/24.0 * TMath::Power(deltaPhi - TMath::Pi(), 4);
 262
 263   Float_t mass2 = m0_1 * m0_1 + m0_2 * m0_2 + 2 * ( TMath::Sqrt(e1squ * e2squ) - ( pt1 * pt2 * ( cosDeltaPhi + 1.0 / tantheta1 / tantheta2 ) ) );
 264
 265 //   Printf(Form("%f %f %f %f %f %f %f %f %f", pt1, eta1, phi1, pt2, eta2, phi2, m0_1, m0_2, mass2));
 266
 267   return mass2;
 268 }
 269
 270 #endif