PWG0/dNdEta/AliMultiplicityESDSelector.cxx

   1 /* $Id$ */
   2
   3 #include "AliMultiplicityESDSelector.h"
   4
   5 #include <TVector3.h>
   6 #include <TFile.h>
   7 #include <TH2F.h>
   8 #include <TH3F.h>
   9 #include <TTree.h>
  10
  11 #include <AliLog.h>
  12 #include <AliESD.h>
  13 #include <AliMultiplicity.h>
  14
  15 #include "esdTrackCuts/AliESDtrackCuts.h"
  16 #include "AliPWG0Helper.h"
  17 #include "dNdEta/AliMultiplicityCorrection.h"
  18
  19 //#define TPCMEASUREMENT
  20 #define ITSMEASUREMENT
  21
  22 ClassImp(AliMultiplicityESDSelector)
  23
  24 AliMultiplicityESDSelector::AliMultiplicityESDSelector() :
  25   AliSelector(),
  26   fMultiplicity(0),
  27   fEsdTrackCuts(0)
  28 {
  29   //
  30   // Constructor. Initialization of pointers
  31   //
  32 }
  33
  34 AliMultiplicityESDSelector::~AliMultiplicityESDSelector()
  35 {
  36   //
  37   // Destructor
  38   //
  39
  40   // histograms are in the output list and deleted when the output
  41   // list is deleted by the TSelector dtor
  42 }
  43
  44 void AliMultiplicityESDSelector::Begin(TTree* tree)
  45 {
  46   // Begin function
  47
  48   AliSelector::Begin(tree);
  49
  50   ReadUserObjects(tree);
  51 }
  52
  53 void AliMultiplicityESDSelector::ReadUserObjects(TTree* tree)
  54 {
  55   // read the user objects, called from slavebegin and begin
  56
  57   if (!fEsdTrackCuts && fInput)
  58     fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (fInput->FindObject("AliESDtrackCuts"));
  59
  60   if (!fEsdTrackCuts && tree)
  61     fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (tree->GetUserInfo()->FindObject("AliESDtrackCuts"));
  62
  63   if (!fEsdTrackCuts)
  64      AliDebug(AliLog::kError, "ERROR: Could not read EsdTrackCuts from input list.");
  65 }
  66
  67 void AliMultiplicityESDSelector::SlaveBegin(TTree* tree)
  68 {
  69   // The SlaveBegin() function is called after the Begin() function.
  70   // When running with PROOF SlaveBegin() is called on each slave server.
  71   // The tree argument is deprecated (on PROOF 0 is passed).
  72
  73   AliSelector::SlaveBegin(tree);
  74
  75   ReadUserObjects(tree);
  76
  77   fMultiplicity = new AliMultiplicityCorrection("Multiplicity", "Multiplicity");
  78 }
  79
  80 void AliMultiplicityESDSelector::Init(TTree* tree)
  81 {
  82   // read the user objects
  83
  84   AliSelector::Init(tree);
  85
  86   // enable only the needed branches
  87   if (tree)
  88   {
  89     tree->SetBranchStatus("*", 0);
  90     tree->SetBranchStatus("fTriggerMask", 1);
  91     tree->SetBranchStatus("fSPDVertex*", 1);
  92
  93     #ifdef ITSMEASUREMENT
  94       tree->SetBranchStatus("fSPDMult*", 1);
  95     #endif
  96
  97     #ifdef TPCMEASUREMENT
  98       AliESDtrackCuts::EnableNeededBranches(tree);
  99     #endif
 100   }
 101 }
 102
 103 Bool_t AliMultiplicityESDSelector::Process(Long64_t entry)
 104 {
 105   // The Process() function is called for each entry in the tree (or possibly
 106   // keyed object in the case of PROOF) to be processed. The entry argument
 107   // specifies which entry in the currently loaded tree is to be processed.
 108   // It can be passed to either TTree::GetEntry() or TBranch::GetEntry()
 109   // to read either all or the required parts of the data. When processing
 110   // keyed objects with PROOF, the object is already loaded and is available
 111   // via the fObject pointer.
 112   //
 113   // This function should contain the "body" of the analysis. It can contain
 114   // simple or elaborate selection criteria, run algorithms on the data
 115   // of the event and typically fill histograms.
 116
 117   // WARNING when a selector is used with a TChain, you must use
 118   //  the pointer to the current TTree to call GetEntry(entry).
 119   //  The entry is always the local entry number in the current tree.
 120   //  Assuming that fTree is the pointer to the TChain being processed,
 121   //  use fTree->GetTree()->GetEntry(entry).
 122
 123   if (AliSelector::Process(entry) == kFALSE)
 124     return kFALSE;
 125
 126   // Check prerequisites
 127   if (!fESD)
 128   {
 129     AliDebug(AliLog::kError, "ESD branch not available");
 130     return kFALSE;
 131   }
 132
 133   Bool_t eventTriggered = AliPWG0Helper::IsEventTriggered(fESD);
 134   Bool_t eventVertex = AliPWG0Helper::IsVertexReconstructed(fESD);
 135
 136   if (!eventTriggered || !eventVertex)
 137     return kTRUE;
 138
 139   // get the ESD vertex
 140   const AliESDVertex* vtxESD = fESD->GetVertex();
 141   Double_t vtx[3];
 142   vtxESD->GetXYZ(vtx);
 143
 144   Int_t nESDTracks05 = 0;
 145   Int_t nESDTracks10 = 0;
 146   Int_t nESDTracks15 = 0;
 147   Int_t nESDTracks20 = 0;
 148
 149 #ifdef ITSMEASUREMENT
 150   // get tracklets
 151   const AliMultiplicity* mult = fESD->GetMultiplicity();
 152   if (!mult)
 153   {
 154     AliDebug(AliLog::kError, "AliMultiplicity not available");
 155     return kFALSE;
 156   }
 157
 158   // get multiplicity from ITS tracklets
 159   for (Int_t i=0; i<mult->GetNumberOfTracklets(); ++i)
 160   {
 161     //printf("%d %f %f %f\n", i, mult->GetTheta(i), mult->GetPhi(i), mult->GetDeltaPhi(i));
 162
 163     // this removes non-tracklets. Very bad solution. SPD guys are working on better solution...
 164     if (mult->GetDeltaPhi(i) < -1000)
 165       continue;
 166
 167     Float_t theta = mult->GetTheta(i);
 168     Float_t eta   = -TMath::Log(TMath::Tan(theta/2.));
 169
 170     if (TMath::Abs(eta) < 0.5)
 171       nESDTracks05++;
 172
 173     if (TMath::Abs(eta) < 1.0)
 174       nESDTracks10++;
 175
 176     if (TMath::Abs(eta) < 1.5)
 177       nESDTracks15++;
 178
 179     if (TMath::Abs(eta) < 2.0)
 180       nESDTracks20++;
 181   }
 182 #endif
 183
 184 #ifdef TPCMEASUREMENT
 185   if (!fEsdTrackCuts)
 186   {
 187     AliDebug(AliLog::kError, "fESDTrackCuts not available");
 188     return kFALSE;
 189   }
 190
 191   // get multiplicity from ESD tracks
 192   TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD);
 193   Int_t nGoodTracks = list->GetEntries();
 194   // loop over esd tracks
 195   for (Int_t i=0; i<nGoodTracks; i++)
 196   {
 197     AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (list->At(i));
 198     if (!esdTrack)
 199     {
 200       AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", i));
 201       continue;
 202     }
 203
 204     Double_t p[3];
 205     esdTrack->GetConstrainedPxPyPz(p); // ### TODO should be okay because we have a vertex, however GetInnerPxPyPy / GetOuterPxPyPy also exist
 206     TVector3 vector(p);
 207
 208     Float_t theta = vector.Theta();
 209     Float_t eta   = -TMath::Log(TMath::Tan(theta/2.));
 210     Float_t pt = vector.Pt();
 211
 212     //if (pt < kPtCut)
 213     //  continue;
 214
 215     if (TMath::Abs(eta) < 0.5)
 216       nESDTracks05++;
 217
 218     if (TMath::Abs(eta) < 1.0)
 219       nESDTracks10++;
 220
 221     if (TMath::Abs(eta) < 1.5)
 222       nESDTracks15++;
 223
 224     if (TMath::Abs(eta) < 2.0)
 225       nESDTracks20++;
 226   }
 227 #endif
 228
 229   fMultiplicity->FillMeasured(vtx[2], nESDTracks05, nESDTracks10, nESDTracks15, nESDTracks20);
 230
 231   return kTRUE;
 232 }
 233
 234 void AliMultiplicityESDSelector::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   fOutput->Add(fMultiplicity);
 250 }
 251
 252 void AliMultiplicityESDSelector::Terminate()
 253 {
 254   // The Terminate() function is the last function to be called during
 255   // a query. It always runs on the client, it can be used to present
 256   // the results graphically or save the results to file.
 257
 258   AliSelector::Terminate();
 259
 260   fMultiplicity = dynamic_cast<AliMultiplicityCorrection*> (fOutput->FindObject("Multiplicity"));
 261
 262   if (!fMultiplicity)
 263   {
 264     AliDebug(AliLog::kError, Form("ERROR: Histograms not available %p", (void*) fMultiplicity));
 265     return;
 266   }
 267
 268   TFile* file = TFile::Open("multiplicityESD.root", "RECREATE");
 269
 270   fMultiplicity->SaveHistograms();
 271
 272   file->Close();
 273
 274   fMultiplicity->DrawHistograms();
 275 }