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