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