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 "dNdEta/AlidNdEtaCorrection.h"
  19 #include "AliPWG0Helper.h"
  20
  21 ClassImp(AlidNdEtaAnalysisESDSelector)
  22
  23 AlidNdEtaAnalysisESDSelector::AlidNdEtaAnalysisESDSelector() :
  24   AliSelector(),
  25   fdNdEtaAnalysisMBVtx(0),
  26   fdNdEtaAnalysisMB(0),
  27   fdNdEtaAnalysis(0),
  28   fEsdTrackCuts(0),
  29   fdNdEtaCorrection(0)
  30   {
  31   //
  32   // Constructor. Initialization of pointers
  33   //
  34
  35   AliLog::SetClassDebugLevel("AlidNdEtaAnalysisESDSelector", AliLog::kDebug);
  36 }
  37
  38 AlidNdEtaAnalysisESDSelector::~AlidNdEtaAnalysisESDSelector()
  39 {
  40   //
  41   // Destructor
  42   //
  43
  44   // histograms are in the output list and deleted when the output
  45   // list is deleted by the TSelector dtor
  46 }
  47
  48 void AlidNdEtaAnalysisESDSelector::Begin(TTree* tree)
  49 {
  50   // Begin function
  51
  52   ReadUserObjects(tree);
  53 }
  54
  55 void AlidNdEtaAnalysisESDSelector::ReadUserObjects(TTree* tree)
  56 {
  57   // read the user objects, called from slavebegin and begin
  58
  59   if (!fEsdTrackCuts && fInput)
  60     fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (fInput->FindObject("AliESDtrackCuts"));
  61
  62   if (!fEsdTrackCuts && tree)
  63     fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (tree->GetUserInfo()->FindObject("AliESDtrackCuts"));
  64
  65   if (!fEsdTrackCuts)
  66      AliDebug(AliLog::kError, "ERROR: Could not read EsdTrackCuts from input list.");
  67
  68
  69   if (!fdNdEtaCorrection && fInput)
  70     fdNdEtaCorrection = dynamic_cast<AlidNdEtaCorrection*> (fInput->FindObject("dndeta_correction"));
  71
  72   if (!fdNdEtaCorrection && tree)
  73     fdNdEtaCorrection = dynamic_cast<AlidNdEtaCorrection*> (tree->GetUserInfo()->FindObject("dndeta_correction"));
  74
  75   if (!fdNdEtaCorrection)
  76      AliDebug(AliLog::kError, "ERROR: Could not read dndeta_correction from input list.");
  77 }
  78
  79 void AlidNdEtaAnalysisESDSelector::SlaveBegin(TTree* tree)
  80 {
  81   // The SlaveBegin() function is called after the Begin() function.
  82   // When running with PROOF SlaveBegin() is called on each slave server.
  83   // The tree argument is deprecated (on PROOF 0 is passed).
  84
  85   AliSelector::SlaveBegin(tree);
  86
  87   ReadUserObjects(tree);
  88
  89   fdNdEtaAnalysisMBVtx = new dNdEtaAnalysis("dndeta_mbvtx", "dndeta_mbvtx");
  90   fdNdEtaAnalysisMB = new dNdEtaAnalysis("dndeta_mb", "dndeta_mb");
  91   fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta");
  92 }
  93
  94 void AlidNdEtaAnalysisESDSelector::Init(TTree* tree)
  95 {
  96   // read the user objects
  97
  98   AliSelector::Init(tree);
  99
 100   // Enable only the needed branches
 101   if (tree)
 102   {
 103     tree->SetBranchStatus("*", 0);
 104     tree->SetBranchStatus("fTriggerMask", 1);
 105     tree->SetBranchStatus("fSPDVertex*", 1);
 106
 107     AliESDtrackCuts::EnableNeededBranches(tree);
 108   }
 109 }
 110
 111 Bool_t AlidNdEtaAnalysisESDSelector::Process(Long64_t entry)
 112 {
 113   // The Process() function is called for each entry in the tree (or possibly
 114   // keyed object in the case of PROOF) to be processed. The entry argument
 115   // specifies which entry in the currently loaded tree is to be processed.
 116   // It can be passed to either TTree::GetEntry() or TBranch::GetEntry()
 117   // to read either all or the required parts of the data. When processing
 118   // keyed objects with PROOF, the object is already loaded and is available
 119   // via the fObject pointer.
 120   //
 121   // This function should contain the "body" of the analysis. It can contain
 122   // simple or elaborate selection criteria, run algorithms on the data
 123   // of the event and typically fill histograms.
 124
 125   // WARNING when a selector is used with a TChain, you must use
 126   //  the pointer to the current TTree to call GetEntry(entry).
 127   //  The entry is always the local entry number in the current tree.
 128   //  Assuming that fTree is the pointer to the TChain being processed,
 129   //  use fTree->GetTree()->GetEntry(entry).
 130
 131   if (AliSelector::Process(entry) == kFALSE)
 132     return kFALSE;
 133
 134   // Check prerequisites
 135   if (!fESD)
 136   {
 137     AliDebug(AliLog::kError, "ESD branch not available");
 138     return kFALSE;
 139   }
 140
 141   if (!fEsdTrackCuts)
 142   {
 143     AliDebug(AliLog::kError, "fESDTrackCuts not available");
 144     return kFALSE;
 145   }
 146
 147   if (!fdNdEtaCorrection)
 148   {
 149     AliDebug(AliLog::kError, "fdNdEtaCorrection not available");
 150     return kFALSE;
 151   }
 152
 153   if (AliPWG0Helper::IsEventTriggered(fESD) == kFALSE)
 154   {
 155     AliDebug(AliLog::kDebug+1, Form("Skipping event %d because it was not triggered", (Int_t) entry));
 156     return kTRUE;
 157   }
 158
 159   if (AliPWG0Helper::IsVertexReconstructed(fESD) == kFALSE)
 160   {
 161     AliDebug(AliLog::kDebug+1, Form("Skipping event %d because its vertex was not reconstructed", (Int_t) entry));
 162     return kTRUE;
 163   }
 164
 165   // ########################################################
 166   // get the EDS vertex
 167   const AliESDVertex* vtxESD = fESD->GetVertex();
 168   Double_t vtx[3];
 169   vtxESD->GetXYZ(vtx);
 170
 171   // get number of "good" tracks
 172   TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD);
 173   Int_t nGoodTracks = list->GetEntries();
 174
 175   Float_t vertexRecoCorr = fdNdEtaCorrection->GetVertexRecoCorrection(vtx[2], nGoodTracks);
 176   if (vertexRecoCorr <= 0)
 177   {
 178     AliDebug(AliLog::kError, Form("INFO: Skipping event because vertexRecoCorr is <= 0 (%f)", vertexRecoCorr));
 179     delete list;
 180     return kTRUE;
 181   }
 182
 183   Float_t triggerCorr = fdNdEtaCorrection->GetTriggerBiasCorrection(vtx[2], nGoodTracks);
 184   if (triggerCorr <= 0)
 185   {
 186     AliDebug(AliLog::kError, Form("INFO: Skipping event because triggerCorr is <= 0 (%f)", triggerCorr));
 187     delete list;
 188     return kTRUE;
 189   }
 190
 191   // loop over esd tracks
 192   for (Int_t t=0; t<nGoodTracks; t++)
 193   {
 194     AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (list->At(t));
 195     if (!esdTrack)
 196     {
 197       AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", t));
 198       continue;
 199     }
 200
 201     Double_t p[3];
 202     esdTrack->GetConstrainedPxPyPz(p); // ### TODO should be okay because we have a vertex, however GetInnerPxPyPy / GetOuterPxPyPy also exist
 203     TVector3 vector(p);
 204
 205     Float_t theta = vector.Theta();
 206     Float_t eta   = -TMath::Log(TMath::Tan(theta/2.));
 207     Float_t pt = vector.Pt();
 208
 209     Float_t track2particleCorr = fdNdEtaCorrection->GetTrack2ParticleCorrection(vtx[2], eta, pt);
 210
 211     Float_t weight = track2particleCorr * vertexRecoCorr * triggerCorr;
 212     if (weight <= 0)
 213     {
 214       AliDebug(AliLog::kError, Form("INFO: Skipping track because weight is <= 0 (track %d, weight %f) (vtx %f, eta %f, pt %f)", t, weight, vtx[2], eta, pt));
 215       continue;
 216     }
 217
 218     fdNdEtaAnalysisMBVtx->FillTrack(vtx[2], eta, pt, track2particleCorr);
 219     fdNdEtaAnalysisMB->FillTrack(vtx[2], eta, pt, track2particleCorr * vertexRecoCorr);
 220     fdNdEtaAnalysis->FillTrack(vtx[2], eta, pt, weight);
 221   } // end of track loop
 222
 223   delete list;
 224   list = 0;
 225
 226   // for event count per vertex
 227   fdNdEtaAnalysisMBVtx->FillEvent(vtx[2], 1);
 228   fdNdEtaAnalysisMB->FillEvent(vtx[2], vertexRecoCorr);
 229   fdNdEtaAnalysis->FillEvent(vtx[2], vertexRecoCorr * triggerCorr);
 230
 231   return kTRUE;
 232 }
 233
 234 void AlidNdEtaAnalysisESDSelector::SlaveTerminate()
 235 {
 236   // The SlaveTerminate() function is called after all entries or objects
 237   // have been processed. When running with PROOF SlaveTerminate() is called
 238   // on each slave server.
 239
 240   AliSelector::SlaveTerminate();
 241
 242   // Add the histograms to the output on each slave server
 243   if (!fOutput)
 244   {
 245     AliDebug(AliLog::kError, Form("ERROR: Output list not initialized."));
 246     return;
 247   }
 248
 249   // Add the objects to the output list and set them to 0, so that the destructor does not delete them.
 250
 251   fOutput->Add(fdNdEtaAnalysis);
 252   fdNdEtaAnalysis = 0;
 253
 254   fOutput->Add(fdNdEtaAnalysisMB);
 255   fdNdEtaAnalysisMB = 0;
 256
 257   fOutput->Add(fdNdEtaAnalysisMBVtx);
 258   fdNdEtaAnalysisMBVtx = 0;
 259 }
 260
 261 void AlidNdEtaAnalysisESDSelector::Terminate()
 262 {
 263   // The Terminate() function is the last function to be called during
 264   // a query. It always runs on the client, it can be used to present
 265   // the results graphically or save the results to file.
 266
 267   AliSelector::Terminate();
 268
 269   fdNdEtaAnalysis = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta"));
 270   fdNdEtaAnalysisMB = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta_mb"));
 271   fdNdEtaAnalysisMBVtx = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta_mbvtx"));
 272
 273   if (!fdNdEtaAnalysis || !fdNdEtaAnalysisMB || !fdNdEtaAnalysisMBVtx)
 274   {
 275     AliDebug(AliLog::kError, Form("ERROR: Histograms not available %p %p %p", (void*) fdNdEtaAnalysis, (void*) fdNdEtaAnalysisMB, (void*) fdNdEtaAnalysisMBVtx));
 276     return;
 277   }
 278
 279   if (fdNdEtaAnalysis)
 280     fdNdEtaAnalysis->Finish(fdNdEtaCorrection, 0.3);
 281
 282   if (fdNdEtaAnalysisMB)
 283     fdNdEtaAnalysisMB->Finish(fdNdEtaCorrection, 0.3);
 284
 285   if (fdNdEtaAnalysisMBVtx)
 286     fdNdEtaAnalysisMBVtx->Finish(fdNdEtaCorrection, 0.3);
 287
 288   TFile* fout = new TFile("analysis_esd.root","RECREATE");
 289
 290   if (fdNdEtaAnalysis)
 291     fdNdEtaAnalysis->SaveHistograms();
 292
 293   if (fdNdEtaAnalysisMB)
 294     fdNdEtaAnalysisMB->SaveHistograms();
 295
 296   if (fdNdEtaAnalysisMBVtx)
 297     fdNdEtaAnalysisMBVtx->SaveHistograms();
 298
 299   if (fEsdTrackCuts)
 300     fEsdTrackCuts->SaveHistograms("esd_tracks_cuts");
 301
 302   if (fdNdEtaCorrection)
 303     fdNdEtaCorrection->SaveHistograms();
 304
 305   fout->Write();
 306   fout->Close();
 307 }