PWG0/dNdEta/AlidNdEtaAnalysisESDSelector.cxx

   1 /* $Id$ */
   2
   3 #include "AlidNdEtaAnalysisESDSelector.h"
   4
   5 #include <TStyle.h>
   6 #include <TSystem.h>
   7 #include <TCanvas.h>
   8 #include <TVector3.h>
   9 #include <TChain.h>
  10 #include <TFile.h>
  11 #include <TH1F.h>
  12
  13 #include <AliLog.h>
  14 #include <AliESDVertex.h>
  15 #include <AliESD.h>
  16
  17 #include "esdTrackCuts/AliESDtrackCuts.h"
  18 #include "dNdEta/dNdEtaAnalysis.h"
  19 #include "AliPWG0Helper.h"
  20
  21 #include "AliGenEventHeader.h"
  22 #include "AliHeader.h"
  23 #include "AliStack.h"
  24 #include "TParticle.h"
  25
  26 ClassImp(AlidNdEtaAnalysisESDSelector)
  27
  28 AlidNdEtaAnalysisESDSelector::AlidNdEtaAnalysisESDSelector() :
  29   AliSelectorRL(),
  30   fdNdEtaAnalysis(0),
  31   fMult(0),
  32   fEsdTrackCuts(0)
  33 {
  34   //
  35   // Constructor. Initialization of pointers
  36   //
  37
  38   AliLog::SetClassDebugLevel("AlidNdEtaAnalysisESDSelector", AliLog::kDebug);
  39 }
  40
  41 AlidNdEtaAnalysisESDSelector::~AlidNdEtaAnalysisESDSelector()
  42 {
  43   //
  44   // Destructor
  45   //
  46
  47   // histograms are in the output list and deleted when the output
  48   // list is deleted by the TSelector dtor
  49 }
  50
  51 void AlidNdEtaAnalysisESDSelector::Begin(TTree* tree)
  52 {
  53   // Begin function
  54
  55   ReadUserObjects(tree);
  56 }
  57
  58 void AlidNdEtaAnalysisESDSelector::ReadUserObjects(TTree* tree)
  59 {
  60   // read the user objects, called from slavebegin and begin
  61
  62   if (!fEsdTrackCuts && fInput)
  63     fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (fInput->FindObject("AliESDtrackCuts"));
  64
  65   if (!fEsdTrackCuts && tree)
  66     fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (tree->GetUserInfo()->FindObject("AliESDtrackCuts"));
  67
  68   if (!fEsdTrackCuts)
  69      AliDebug(AliLog::kError, "ERROR: Could not read EsdTrackCuts from input list.");
  70 }
  71
  72 void AlidNdEtaAnalysisESDSelector::SlaveBegin(TTree* tree)
  73 {
  74   // The SlaveBegin() function is called after the Begin() function.
  75   // When running with PROOF SlaveBegin() is called on each slave server.
  76   // The tree argument is deprecated (on PROOF 0 is passed).
  77
  78   AliSelectorRL::SlaveBegin(tree);
  79
  80   ReadUserObjects(tree);
  81
  82   fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta");
  83   fMult = new TH1F("fMult", "fMult;Ntracks;Count", 201, -0.5, 200.5);
  84 }
  85
  86 void AlidNdEtaAnalysisESDSelector::Init(TTree* tree)
  87 {
  88   // read the user objects
  89
  90   AliSelectorRL::Init(tree);
  91
  92   // Enable only the needed branches
  93   if (tree)
  94   {
  95     tree->SetBranchStatus("*", 0);
  96     tree->SetBranchStatus("fTriggerMask", 1);
  97     tree->SetBranchStatus("fSPDVertex*", 1);
  98     tree->SetBranchStatus("fTracks.fLabel", 1);
  99
 100     AliESDtrackCuts::EnableNeededBranches(tree);
 101   }
 102 }
 103
 104 Bool_t AlidNdEtaAnalysisESDSelector::Process(Long64_t entry)
 105 {
 106   // loop over all events
 107
 108   if (AliSelectorRL::Process(entry) == kFALSE)
 109     return kFALSE;
 110
 111   // Check prerequisites
 112   if (!fESD)
 113   {
 114     AliDebug(AliLog::kError, "ESD branch not available");
 115     return kFALSE;
 116   }
 117
 118   if (!fEsdTrackCuts)
 119   {
 120     AliDebug(AliLog::kError, "fESDTrackCuts not available");
 121     return kFALSE;
 122   }
 123
 124   if (AliPWG0Helper::IsEventTriggered(fESD) == kFALSE)
 125   {
 126     AliDebug(AliLog::kDebug+1, Form("Skipping event %d because it was not triggered", (Int_t) entry));
 127     return kTRUE;
 128   }
 129
 130   if (AliPWG0Helper::IsVertexReconstructed(fESD) == kFALSE)
 131   {
 132     AliDebug(AliLog::kDebug+1, Form("Skipping event %d because its vertex was not reconstructed", (Int_t) entry));
 133     return kTRUE;
 134   }
 135
 136   AliHeader* header = GetHeader();
 137   if (!header)
 138   {
 139     AliDebug(AliLog::kError, "Header not available");
 140     return kFALSE;
 141   }
 142
 143   AliStack* stack = GetStack();
 144   if (!stack)
 145   {
 146     AliDebug(AliLog::kError, "Stack not available");
 147     return kFALSE;
 148   }
 149
 150   // get the MC vertex
 151   AliGenEventHeader* genHeader = header->GenEventHeader();
 152
 153   TArrayF vtxMC(3);
 154   genHeader->PrimaryVertex(vtxMC);
 155
 156   // ########################################################
 157   // get the ESD vertex
 158   const AliESDVertex* vtxESD = fESD->GetVertex();
 159   Double_t vtx[3];
 160   vtxESD->GetXYZ(vtx);
 161
 162   //vtx[2] = vtxMC[2];
 163   //vtx[2] -= 2.951034e-03 + 6.859620e-04 * vtxMC[2];
 164
 165   // get number of "good" tracks
 166   TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD);
 167   Int_t nGoodTracks = list->GetEntries();
 168
 169   // FAKE test!
 170   //Int_t nContributors = vtxESD->GetNContributors();
 171
 172   // loop over esd tracks
 173   for (Int_t t=0; t<nGoodTracks; t++)
 174   {
 175     AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (list->At(t));
 176     if (!esdTrack)
 177     {
 178       AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", t));
 179       continue;
 180     }
 181
 182     Int_t label = TMath::Abs(esdTrack->GetLabel());
 183
 184     TParticle* particle = stack->Particle(label);
 185     if (!particle)
 186     {
 187       AliDebug(AliLog::kError, Form("ERROR: Could not retrieve particle %d.", esdTrack->GetLabel()));
 188       continue;
 189     }
 190
 191     Double_t p[3];
 192     esdTrack->GetConstrainedPxPyPz(p); // ### TODO should be okay because we have a vertex, however GetInnerPxPyPy / GetOuterPxPyPy also exist
 193     TVector3 vector(p);
 194
 195     Float_t theta = vector.Theta();
 196     Float_t eta   = -TMath::Log(TMath::Tan(theta/2.));
 197     Float_t pt = vector.Pt();
 198
 199     //eta = particle->Eta();
 200     //pt = particle->Pt();
 201
 202     fdNdEtaAnalysis->FillTrack(vtx[2], eta, pt);
 203   } // end of track loop
 204
 205   delete list;
 206   list = 0;
 207
 208   // for event count per vertex
 209   fdNdEtaAnalysis->FillEvent(vtx[2], nGoodTracks);
 210   fMult->Fill(nGoodTracks);
 211
 212   return kTRUE;
 213 }
 214
 215 void AlidNdEtaAnalysisESDSelector::SlaveTerminate()
 216 {
 217   // The SlaveTerminate() function is called after all entries or objects
 218   // have been processed. When running with PROOF SlaveTerminate() is called
 219   // on each slave server.
 220
 221   AliSelectorRL::SlaveTerminate();
 222
 223   // Add the histograms to the output on each slave server
 224   if (!fOutput)
 225   {
 226     AliDebug(AliLog::kError, Form("ERROR: Output list not initialized."));
 227     return;
 228   }
 229
 230   // Add the objects to the output list and set them to 0, so that the destructor does not delete them.
 231
 232   fOutput->Add(fdNdEtaAnalysis);
 233   fOutput->Add(fMult);
 234
 235   fdNdEtaAnalysis = 0;
 236 }
 237
 238 void AlidNdEtaAnalysisESDSelector::Terminate()
 239 {
 240   // The Terminate() function is the last function to be called during
 241   // a query. It always runs on the client, it can be used to present
 242   // the results graphically or save the results to file.
 243
 244   AliSelectorRL::Terminate();
 245
 246   fdNdEtaAnalysis = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta"));
 247   fMult = dynamic_cast<TH1F*> (fOutput->FindObject("fMult"));
 248
 249   if (!fdNdEtaAnalysis)
 250   {
 251     AliDebug(AliLog::kError, "ERROR: Histograms not available");
 252     return;
 253   }
 254
 255   TFile* fout = new TFile("analysis_esd_raw.root", "RECREATE");
 256
 257   if (fdNdEtaAnalysis)
 258     fdNdEtaAnalysis->SaveHistograms();
 259
 260   if (fEsdTrackCuts)
 261     fEsdTrackCuts->SaveHistograms("esd_tracks_cuts");
 262
 263   if (fMult)
 264     fMult->Write();
 265
 266   fout->Write();
 267   fout->Close();
 268 }