PWG0/dNdEta/AliMultiplicityMCSelector.cxx

   1 /* $Id$ */
   2
   3 #include "AliMultiplicityMCSelector.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 #include <TH2F.h>
  13 #include <TParticle.h>
  14
  15 #include <AliLog.h>
  16 #include <AliESD.h>
  17 #include <AliStack.h>
  18
  19 #include "esdTrackCuts/AliESDtrackCuts.h"
  20 #include "AliPWG0Helper.h"
  21
  22 ClassImp(AliMultiplicityMCSelector)
  23
  24 AliMultiplicityMCSelector::AliMultiplicityMCSelector() :
  25   AliSelectorRL(),
  26   fMultiplicityESD(0),
  27   fMultiplicityMC(0),
  28   fCorrelation(0),
  29   fEsdTrackCuts(0)
  30 {
  31   //
  32   // Constructor. Initialization of pointers
  33   //
  34 }
  35
  36 AliMultiplicityMCSelector::~AliMultiplicityMCSelector()
  37 {
  38   //
  39   // Destructor
  40   //
  41
  42   // histograms are in the output list and deleted when the output
  43   // list is deleted by the TSelector dtor
  44 }
  45
  46 void AliMultiplicityMCSelector::Begin(TTree* tree)
  47 {
  48   // Begin function
  49
  50   AliSelectorRL::Begin(tree);
  51
  52   ReadUserObjects(tree);
  53 }
  54
  55 void AliMultiplicityMCSelector::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 void AliMultiplicityMCSelector::SlaveBegin(TTree* tree)
  70 {
  71   // The SlaveBegin() function is called after the Begin() function.
  72   // When running with PROOF SlaveBegin() is called on each slave server.
  73   // The tree argument is deprecated (on PROOF 0 is passed).
  74
  75   AliSelector::SlaveBegin(tree);
  76
  77   ReadUserObjects(tree);
  78
  79   fMultiplicityESD = new TH1F("fMultiplicityESD", "fMultiplicityESD;Ntracks;Count", 201, -0.5, 200.5);
  80   fMultiplicityMC = new TH1F("fMultiplicityMC", "fMultiplicityMC;Npart;Count", 201, -0.5, 200.5);
  81
  82   fCorrelation = new TH2F("fCorrelation", "fCorrelation;Npart;Ntracks", 201, -0.5, 200.5, 201, -0.5, 200.5);
  83 }
  84
  85 Bool_t AliMultiplicityMCSelector::Process(Long64_t entry)
  86 {
  87   // The Process() function is called for each entry in the tree (or possibly
  88   // keyed object in the case of PROOF) to be processed. The entry argument
  89   // specifies which entry in the currently loaded tree is to be processed.
  90   // It can be passed to either TTree::GetEntry() or TBranch::GetEntry()
  91   // to read either all or the required parts of the data. When processing
  92   // keyed objects with PROOF, the object is already loaded and is available
  93   // via the fObject pointer.
  94   //
  95   // This function should contain the "body" of the analysis. It can contain
  96   // simple or elaborate selection criteria, run algorithms on the data
  97   // of the event and typically fill histograms.
  98
  99   // WARNING when a selector is used with a TChain, you must use
 100   //  the pointer to the current TTree to call GetEntry(entry).
 101   //  The entry is always the local entry number in the current tree.
 102   //  Assuming that fTree is the pointer to the TChain being processed,
 103   //  use fTree->GetTree()->GetEntry(entry).
 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   AliStack* stack = GetStack();
 122   if (!stack)
 123   {
 124     AliDebug(AliLog::kError, "Stack not available");
 125     return kFALSE;
 126   }
 127
 128   if (AliPWG0Helper::IsEventTriggered(fESD) == kFALSE)
 129     return kTRUE;
 130
 131   if (AliPWG0Helper::IsVertexReconstructed(fESD) == kFALSE)
 132     return kTRUE;
 133
 134   // TODO put pt cut
 135
 136   // get number of "good" tracks from ESD
 137   Int_t nESDTracks = fEsdTrackCuts->CountAcceptedTracks(fESD);
 138   fMultiplicityESD->Fill(nESDTracks);
 139
 140   // get number of tracks from MC
 141   // loop over mc particles
 142   Int_t nPrim  = stack->GetNprimary();
 143   Int_t nMCTracks = 0;
 144   for (Int_t iMc = 0; iMc < nPrim; ++iMc)
 145   {
 146     AliDebug(AliLog::kDebug+1, Form("MC Loop: Processing particle %d.", iMc));
 147
 148     TParticle* particle = stack->Particle(iMc);
 149
 150     if (!particle)
 151     {
 152       AliDebug(AliLog::kError, Form("UNEXPECTED: particle with label %d not found in stack (mc loop).", iMc));
 153       continue;
 154     }
 155
 156     if (AliPWG0Helper::IsPrimaryCharged(particle, nPrim) == kFALSE)
 157       continue;
 158
 159     Float_t eta = particle->Eta();
 160
 161     if (TMath::Abs(eta) > 0.9)
 162       continue;
 163
 164     nMCTracks++;
 165   }// end of mc particle
 166
 167   fMultiplicityMC->Fill(nMCTracks);
 168
 169   fCorrelation->Fill(nMCTracks, nESDTracks);
 170
 171   return kTRUE;
 172 }
 173
 174 void AliMultiplicityMCSelector::SlaveTerminate()
 175 {
 176   // The SlaveTerminate() function is called after all entries or objects
 177   // have been processed. When running with PROOF SlaveTerminate() is called
 178   // on each slave server.
 179
 180   AliSelector::SlaveTerminate();
 181
 182   // Add the histograms to the output on each slave server
 183   if (!fOutput)
 184   {
 185     AliDebug(AliLog::kError, Form("ERROR: Output list not initialized."));
 186     return;
 187   }
 188
 189   fOutput->Add(fMultiplicityESD);
 190   fOutput->Add(fMultiplicityMC);
 191   fOutput->Add(fCorrelation);
 192 }
 193
 194 void AliMultiplicityMCSelector::Terminate()
 195 {
 196   // The Terminate() function is the last function to be called during
 197   // a query. It always runs on the client, it can be used to present
 198   // the results graphically or save the results to file.
 199
 200   AliSelector::Terminate();
 201
 202   fOutput->Print();
 203
 204   fMultiplicityESD = dynamic_cast<TH1F*> (fOutput->FindObject("fMultiplicityESD"));
 205   fMultiplicityMC = dynamic_cast<TH1F*> (fOutput->FindObject("fMultiplicityMC"));
 206   fCorrelation = dynamic_cast<TH2F*> (fOutput->FindObject("fCorrelation"));
 207
 208   if (!fMultiplicityESD || !fMultiplicityMC || !fCorrelation)
 209   {
 210     AliDebug(AliLog::kError, Form("ERROR: Histograms not available %p %p %p", (void*) fMultiplicityESD, (void*) fMultiplicityMC, (void*) fCorrelation));
 211     return;
 212   }
 213
 214   TFile* file = TFile::Open("multiplicityMC.root", "RECREATE");
 215
 216   fMultiplicityMC->Write();
 217   fMultiplicityESD->Write();
 218   fCorrelation->Write();
 219
 220   file->Close();
 221 }