PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskRhoMassBase.cxx

   1 // $Id$
   2 //
   3 // Base class for rho mass calculation.
   4 // Calculates parameterized rho mass for given centrality independent of input.
   5 //
   6 // Author: M. Verweij. Similar to AliAnalysisTaskRhoBase
   7
   8 #include <TF1.h>
   9 #include <TH1F.h>
  10 #include <TH2F.h>
  11 #include <TClonesArray.h>
  12
  13 #include "AliLog.h"
  14 #include "AliRhoParameter.h"
  15 #include "AliEmcalJet.h"
  16
  17 #include "AliAnalysisTaskRhoMassBase.h"
  18
  19 ClassImp(AliAnalysisTaskRhoMassBase)
  20
  21 //________________________________________________________________________
  22 AliAnalysisTaskRhoMassBase::AliAnalysisTaskRhoMassBase() :
  23   AliAnalysisTaskEmcalJet("AliAnalysisTaskRhoMassBase", kFALSE),
  24   fOutRhoMassName(),
  25   fOutRhoMassScaledName(),
  26   fCompareRhoMassName(),
  27   fCompareRhoMassScaledName(),
  28   fRhoMassFunction(0),
  29   fScaleFunction(0),
  30   fAttachToEvent(kTRUE),
  31   fOutRhoMass(0),
  32   fOutRhoMassScaled(0),
  33   fCompareRhoMass(0),
  34   fCompareRhoMassScaled(0),
  35   fHistJetMassvsCent(0),
  36   fHistRhoMassvsCent(0),
  37   fHistRhoMassScaledvsCent(0),
  38   fHistDeltaRhoMassvsCent(0),
  39   fHistDeltaRhoMassScalevsCent(0),
  40   fHistRhoMassvsNtrack(0),
  41   fHistRhoMassScaledvsNtrack(0),
  42   fHistDeltaRhoMassvsNtrack(0),
  43   fHistDeltaRhoMassScalevsNtrack(0),
  44   fHistRhoMassvsNcluster(0),
  45   fHistRhoMassScaledvsNcluster(0),
  46   fHistGammaVsNtrack(0)
  47 {
  48   // Constructor.
  49 }
  50
  51 //________________________________________________________________________
  52 AliAnalysisTaskRhoMassBase::AliAnalysisTaskRhoMassBase(const char *name, Bool_t histo) :
  53   AliAnalysisTaskEmcalJet(name, histo),
  54   fOutRhoMassName(),
  55   fOutRhoMassScaledName(),
  56   fCompareRhoMassName(),
  57   fCompareRhoMassScaledName(),
  58   fRhoMassFunction(0),
  59   fScaleFunction(0),
  60   fAttachToEvent(kTRUE),
  61   fOutRhoMass(0),
  62   fOutRhoMassScaled(0),
  63   fCompareRhoMass(0),
  64   fCompareRhoMassScaled(0),
  65   fHistJetMassvsCent(0),
  66   fHistRhoMassvsCent(0),
  67   fHistRhoMassScaledvsCent(0),
  68   fHistDeltaRhoMassvsCent(0),
  69   fHistDeltaRhoMassScalevsCent(0),
  70   fHistRhoMassvsNtrack(0),
  71   fHistRhoMassScaledvsNtrack(0),
  72   fHistDeltaRhoMassvsNtrack(0),
  73   fHistDeltaRhoMassScalevsNtrack(0),
  74   fHistRhoMassvsNcluster(0),
  75   fHistRhoMassScaledvsNcluster(0),
  76   fHistGammaVsNtrack(0)
  77 {
  78   // Constructor.
  79
  80   SetMakeGeneralHistograms(histo);
  81 }
  82
  83 //________________________________________________________________________
  84 void AliAnalysisTaskRhoMassBase::UserCreateOutputObjects()
  85 {
  86   // User create output objects, called at the beginning of the analysis.
  87
  88   if (!fCreateHisto)
  89     return;
  90
  91   AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
  92
  93   const Int_t nBinsRhom = 200;
  94   const Double_t minRhom = 0.;
  95   const Double_t maxRhom = 20.;
  96
  97   const Int_t nBinsM  = 100;
  98   const Double_t minM = -20.;
  99   const Double_t maxM = 80.;
 100
 101   fHistRhoMassvsCent = new TH2F("fHistRhoMassvsCent", "fHistRhoMassvsCent", 101, -1,  100, nBinsRhom,minRhom,maxRhom);
 102   fHistRhoMassvsCent->GetXaxis()->SetTitle("Centrality (%)");
 103   fHistRhoMassvsCent->GetYaxis()->SetTitle("#rho_{m} (GeV/c * rad^{-1})");
 104   fOutput->Add(fHistRhoMassvsCent);
 105
 106   if (fParticleCollArray.GetEntriesFast()>0) {
 107     fHistRhoMassvsNtrack = new TH2F("fHistRhoMassvsNtrack", "fHistRhoMassvsNtrack", 150, 0, 6000, nBinsRhom,minRhom,maxRhom);
 108     fHistRhoMassvsNtrack->GetXaxis()->SetTitle("No. of tracks");
 109     fHistRhoMassvsNtrack->GetYaxis()->SetTitle("#rho_{m} (GeV/c * rad^{-1})");
 110     fOutput->Add(fHistRhoMassvsNtrack);
 111
 112     //fHistGammaVsNtrack
 113     fHistGammaVsNtrack = new TH2F("fHistGammaVsNtrack", "fHistGammaVsNtrack", 150, 0, 6000, 100,0.,10.);
 114     fHistGammaVsNtrack->GetXaxis()->SetTitle("No. of tracks");
 115     fHistGammaVsNtrack->GetYaxis()->SetTitle("#gamma = #LT E #GT / #LT M #GT");
 116     fOutput->Add(fHistGammaVsNtrack);
 117   }
 118
 119   if (fClusterCollArray.GetEntriesFast()>0) {
 120     fHistRhoMassvsNcluster = new TH2F("fHistRhoMassvsNcluster", "fHistRhoMassvsNcluster", 50, 0, 1500, nBinsRhom,minRhom,maxRhom);
 121     fHistRhoMassvsNcluster->GetXaxis()->SetTitle("No. of tracks");
 122     fHistRhoMassvsNcluster->GetYaxis()->SetTitle("#rho_{m} (GeV/c * rad^{-1})");
 123     fOutput->Add(fHistRhoMassvsNcluster);
 124   }
 125
 126   if (fJetCollArray.GetEntriesFast()>0) {
 127     fHistJetMassvsCent = new TH2F("fHistJetMassvsCent", "fHistJetMassvsCent", 101, -1,  100, nBinsM,minM,maxM);
 128     fHistJetMassvsCent->GetXaxis()->SetTitle("Centrality (%)");
 129     fHistJetMassvsCent->GetYaxis()->SetTitle("#it{M}_{jet} (GeV/c)");
 130     fOutput->Add(fHistJetMassvsCent);
 131   }
 132
 133   if (!fCompareRhoMassName.IsNull()) {
 134     fHistDeltaRhoMassvsCent = new TH2F("fHistDeltaRhoMassvsCent", "fHistDeltaRhoMassvsCent", 101, -1, 100, fNbins, -fMaxBinPt, fMaxBinPt);
 135     fHistDeltaRhoMassvsCent->GetXaxis()->SetTitle("Centrality (%)");
 136     fHistDeltaRhoMassvsCent->GetYaxis()->SetTitle("#Delta#rho (GeV/c * rad^{-1})");
 137     fOutput->Add(fHistDeltaRhoMassvsCent);
 138
 139     if (fParticleCollArray.GetEntriesFast()>0) {
 140       fHistDeltaRhoMassvsNtrack = new TH2F("fHistDeltaRhoMassvsNtrack", "fHistDeltaRhoMassvsNtrack", 150, 0, 6000, fNbins, -fMaxBinPt, fMaxBinPt);
 141       fHistDeltaRhoMassvsNtrack->GetXaxis()->SetTitle("No. of tracks");
 142       fHistDeltaRhoMassvsNtrack->GetYaxis()->SetTitle("#Delta#rho (GeV/c * rad^{-1})");
 143       fOutput->Add(fHistDeltaRhoMassvsNtrack);
 144     }
 145   }
 146
 147   if (fScaleFunction) {
 148     fHistRhoMassScaledvsCent = new TH2F("fHistRhoMassScaledvsCent", "fHistRhoMassScaledvsCent", 101, -1, 100, fNbins, fMinBinPt , fMaxBinPt*2);
 149     fHistRhoMassScaledvsCent->GetXaxis()->SetTitle("Centrality (%)");
 150     fHistRhoMassScaledvsCent->GetYaxis()->SetTitle("#rho_{m,scaled} (GeV/c * rad^{-1})");
 151     fOutput->Add(fHistRhoMassScaledvsCent);
 152
 153     if (fParticleCollArray.GetEntriesFast()>0) {
 154       fHistRhoMassScaledvsNtrack = new TH2F("fHistRhoMassScaledvsNtrack", "fHistRhoMassScaledvsNtrack", 150, 0, 6000, fNbins, fMinBinPt, fMaxBinPt*2);
 155       fHistRhoMassScaledvsNtrack->GetXaxis()->SetTitle("No. of tracks");
 156       fHistRhoMassScaledvsNtrack->GetYaxis()->SetTitle("#rho_{m,scaled} (GeV/c * rad^{-1})");
 157       fOutput->Add(fHistRhoMassScaledvsNtrack);
 158     }
 159
 160     if (fClusterCollArray.GetEntriesFast()>0) {
 161       fHistRhoMassScaledvsNcluster = new TH2F("fHistRhoMassScaledvsNcluster", "fHistRhoMassScaledvsNcluster", 50, 0, 1500, fNbins, fMinBinPt, fMaxBinPt*2);
 162       fHistRhoMassScaledvsNcluster->GetXaxis()->SetTitle("No. of clusters");
 163       fHistRhoMassScaledvsNcluster->GetYaxis()->SetTitle("#rho_{m,scaled} (GeV/c * rad^{-1})");
 164       fOutput->Add(fHistRhoMassScaledvsNcluster);
 165     }
 166
 167     if (!fCompareRhoMassScaledName.IsNull()) {
 168       fHistDeltaRhoMassScalevsCent = new TH2F("fHistDeltaRhoMassScalevsCent", "fHistDeltaRhoMassScalevsCent", 101, -1, 100, fNbins, -fMaxBinPt, fMaxBinPt);
 169       fHistDeltaRhoMassScalevsCent->GetXaxis()->SetTitle("Centrality (%)");
 170       fHistDeltaRhoMassScalevsCent->GetYaxis()->SetTitle("#Delta#rho_{m,scaled} (GeV/c * rad^{-1})");
 171       fOutput->Add(fHistDeltaRhoMassScalevsCent);
 172
 173       if (fParticleCollArray.GetEntriesFast()>0) {
 174         fHistDeltaRhoMassScalevsNtrack = new TH2F("fHistDeltaRhoMassScalevsNtrack", "fHistDeltaRhoMassScalevsNtrack", 150, 0, 6000, fNbins, -fMaxBinPt, fMaxBinPt);
 175         fHistDeltaRhoMassScalevsNtrack->GetXaxis()->SetTitle("No. of tracks");
 176         fHistDeltaRhoMassScalevsNtrack->GetYaxis()->SetTitle("#Delta#rho_{m,scaled} (GeV/c * rad^{-1})");
 177         fOutput->Add(fHistDeltaRhoMassScalevsNtrack);
 178       }
 179     }
 180   }
 181 }
 182
 183 //________________________________________________________________________
 184 Bool_t AliAnalysisTaskRhoMassBase::Run()
 185 {
 186   // Run the analysis.
 187
 188   Double_t rhom = GetRhoMassFactor(fCent);
 189   fOutRhoMass->SetVal(rhom);
 190
 191   if (fScaleFunction) {
 192     Double_t rhomScaled = rhom * GetScaleFactor(fCent);
 193     fOutRhoMassScaled->SetVal(rhomScaled);
 194   }
 195
 196   return kTRUE;
 197 }
 198
 199 //________________________________________________________________________
 200 Bool_t AliAnalysisTaskRhoMassBase::FillHistograms()
 201 {
 202   // Fill histograms.
 203
 204   Int_t Ntracks   = 0;
 205   Int_t Nclusters = 0;
 206
 207   if (fTracks)
 208     Ntracks = fTracks->GetEntries();
 209   if (fCaloClusters)
 210     Nclusters = fCaloClusters->GetEntries();
 211
 212   if (fJets) {
 213     Int_t    Njets         = fJets->GetEntries();
 214     Int_t    NjetAcc       = 0;
 215
 216     for (Int_t i = 0; i < Njets; ++i) {
 217       AliEmcalJet *jet = static_cast<AliEmcalJet*>(fJets->At(i));
 218       if (!jet) {
 219         AliError(Form("%s: Could not receive jet %d", GetName(), i));
 220         continue;
 221       }
 222
 223       if (!AcceptJet(jet))
 224         continue;
 225
 226       fHistJetMassvsCent->Fill(fCent, jet->M());
 227       NjetAcc++;
 228     }
 229   }
 230
 231   fHistRhoMassvsCent->Fill(fCent, fOutRhoMass->GetVal());
 232
 233   if (fTracks)
 234     fHistRhoMassvsNtrack->Fill(Ntracks, fOutRhoMass->GetVal());
 235   if (fCaloClusters)
 236     fHistRhoMassvsNcluster->Fill(Nclusters, fOutRhoMass->GetVal());
 237   if (fCompareRhoMass) {
 238     fHistDeltaRhoMassvsCent->Fill(fCent, fOutRhoMass->GetVal() - fCompareRhoMass->GetVal());
 239     if (fTracks)
 240       fHistDeltaRhoMassvsNtrack->Fill(Ntracks, fOutRhoMass->GetVal() - fCompareRhoMass->GetVal());
 241   }
 242
 243   if (fOutRhoMassScaled) {
 244     fHistRhoMassScaledvsCent->Fill(fCent, fOutRhoMassScaled->GetVal());
 245     if (fTracks)
 246       fHistRhoMassScaledvsNtrack->Fill(Ntracks, fOutRhoMassScaled->GetVal());
 247     if (fCaloClusters)
 248       fHistRhoMassScaledvsNcluster->Fill(Nclusters,  fOutRhoMassScaled->GetVal());
 249     if (fCompareRhoMassScaled) {
 250       fHistDeltaRhoMassScalevsCent->Fill(fCent, fOutRhoMassScaled->GetVal() - fCompareRhoMassScaled->GetVal());
 251       if (fTracks)
 252         fHistDeltaRhoMassScalevsNtrack->Fill(Ntracks, fOutRhoMassScaled->GetVal() - fCompareRhoMassScaled->GetVal());
 253     }
 254   }
 255
 256   return kTRUE;
 257 }
 258
 259
 260 //________________________________________________________________________
 261 void AliAnalysisTaskRhoMassBase::ExecOnce()
 262 {
 263   // Init the analysis.
 264
 265   if (!fOutRhoMass) {
 266     fOutRhoMass = new AliRhoParameter(fOutRhoMassName, 0);
 267
 268     if (fAttachToEvent) {
 269       if (!(InputEvent()->FindListObject(fOutRhoMassName))) {
 270         InputEvent()->AddObject(fOutRhoMass);
 271       } else {
 272         AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fOutRhoMassName.Data()));
 273         return;
 274       }
 275     }
 276   }
 277
 278   if (fScaleFunction && !fOutRhoMassScaled) {
 279     fOutRhoMassScaled = new AliRhoParameter(fOutRhoMassScaledName, 0);
 280
 281     if (fAttachToEvent) {
 282       if (!(InputEvent()->FindListObject(fOutRhoMassScaledName))) {
 283         InputEvent()->AddObject(fOutRhoMassScaled);
 284       } else {
 285         AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fOutRhoMassScaledName.Data()));
 286         return;
 287       }
 288     }
 289   }
 290
 291   if (!fCompareRhoMassName.IsNull() && !fCompareRhoMass) {
 292     fCompareRhoMass = dynamic_cast<AliRhoParameter*>(InputEvent()->FindListObject(fCompareRhoMassName));
 293     if (!fCompareRhoMass) {
 294       AliWarning(Form("%s: Could not retrieve rho %s!", GetName(), fCompareRhoMassName.Data()));
 295     }
 296   }
 297
 298   if (!fCompareRhoMassScaledName.IsNull() && !fCompareRhoMassScaled) {
 299     fCompareRhoMassScaled = dynamic_cast<AliRhoParameter*>(InputEvent()->FindListObject(fCompareRhoMassScaledName));
 300     if (!fCompareRhoMassScaled) {
 301       AliWarning(Form("%s: Could not retrieve rho %s!", GetName(), fCompareRhoMassScaledName.Data()));
 302     }
 303   }
 304
 305   AliAnalysisTaskEmcalJet::ExecOnce();
 306 }
 307
 308 //________________________________________________________________________
 309 Double_t AliAnalysisTaskRhoMassBase::GetRhoMassFactor(Double_t cent)
 310 {
 311   // Return rho per centrality.
 312
 313   Double_t rho = 0;
 314   if (fRhoMassFunction)
 315     rho = fRhoMassFunction->Eval(cent);
 316   return rho;
 317 }
 318
 319 //________________________________________________________________________
 320 Double_t AliAnalysisTaskRhoMassBase::GetScaleFactor(Double_t cent)
 321 {
 322   // Get scale factor.
 323
 324   Double_t scale = 1;
 325   if (fScaleFunction)
 326     scale = fScaleFunction->Eval(cent);
 327   return scale;
 328 }