PWGCF/EBYE/AliAnalysisTaskChargeFluctuations.cxx

   1 #include "TChain.h"
   2 #include "TList.h"
   3 #include "TCanvas.h"
   4 #include "TLorentzVector.h"
   5 #include "TGraphErrors.h"
   6 #include "TH1F.h"
   7
   8 #include "AliAnalysisTaskSE.h"
   9 #include "AliAnalysisManager.h"
  10
  11 #include "AliESDVertex.h"
  12 #include "AliESDEvent.h"
  13 #include "AliESDInputHandler.h"
  14 #include "AliAODEvent.h"
  15 #include "AliAODTrack.h"
  16 #include "AliAODInputHandler.h"
  17 #include "AliMCEventHandler.h"
  18 #include "AliMCEvent.h"
  19 #include "AliStack.h"
  20 #include "AliESDtrackCuts.h"
  21
  22 #include "AliAnalysisTaskChargeFluctuations.h"
  23
  24 // Analysis task for the ChargeFluctuations code
  25 // Authors: Panos Cristakoglou@cern.ch
  26
  27 ClassImp(AliAnalysisTaskChargeFluctuations)
  28
  29 //________________________________________________________________________
  30 AliAnalysisTaskChargeFluctuations::AliAnalysisTaskChargeFluctuations(const char *name)
  31 : AliAnalysisTaskSE(name),
  32   fList(0),
  33   fHistEventStats(0),
  34   fHistVx(0),
  35   fHistVy(0),
  36   fHistVz(0),
  37   fESDtrackCuts(0),
  38   fUseOfflineTrigger(kFALSE),
  39   fVxMax(0.3),
  40   fVyMax(0.3),
  41   fVzMax(10.) {
  42   // Constructor
  43
  44   // Define input and output slots here
  45   // Input slot #0 works with a TChain
  46   DefineInput(0, TChain::Class());
  47   // Output slot #0 writes into a TH1 container
  48   DefineOutput(1, TList::Class());
  49 }
  50
  51 //________________________________________________________________________
  52 void AliAnalysisTaskChargeFluctuations::UserCreateOutputObjects() {
  53   // Create histograms
  54   // Called once
  55   fList = new TList();
  56   fList->SetName("outputList");
  57
  58   //Event stats.
  59   TString gCutName[4] = {"Total","Offline trigger",
  60                          "Vertex","Analyzed"};
  61   fHistEventStats = new TH1F("fHistEventStats",
  62                              "Event statistics;;N_{events}",
  63                              4,0.5,4.5);
  64   for(Int_t i = 1; i <= 4; i++)
  65     fHistEventStats->GetXaxis()->SetBinLabel(i,gCutName[i-1].Data());
  66   fList->Add(fHistEventStats);
  67
  68   //Vertex distributions
  69   fHistVx = new TH1F("fHistVx","Primary vertex distribution - x coordinate;V_{x} (cm);Entries",100,-0.5,0.5);
  70   fList->Add(fHistVx);
  71   fHistVy = new TH1F("fHistVy","Primary vertex distribution - y coordinate;V_{y} (cm);Entries",100,-0.5,0.5);
  72   fList->Add(fHistVy);
  73   fHistVz = new TH1F("fHistVz","Primary vertex distribution - z coordinate;V_{z} (cm);Entries",100,-20.,20.);
  74   fList->Add(fHistVz);
  75
  76   if(fESDtrackCuts) fList->Add(fESDtrackCuts);
  77
  78   // Post output data.
  79   PostData(1, fList);
  80 }
  81
  82 //________________________________________________________________________
  83 void AliAnalysisTaskChargeFluctuations::UserExec(Option_t *) {
  84   // Main loop
  85   // Called for each event
  86   TString gAnalysisLevel = "ESD";
  87
  88   //ESD analysis
  89   if(gAnalysisLevel == "ESD") {
  90     AliESDEvent* gESD = dynamic_cast<AliESDEvent*>(InputEvent()); // from TaskSE
  91     if (!gESD) {
  92       Printf("ERROR: gESD not available");
  93       return;
  94     }
  95
  96     fHistEventStats->Fill(1); //all events
  97     Bool_t isSelected = kTRUE;
  98     if(fUseOfflineTrigger)
  99       isSelected = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
 100     if(isSelected) {
 101       fHistEventStats->Fill(2); //triggered events
 102
 103       const AliESDVertex *vertex = gESD->GetPrimaryVertex();
 104       if(vertex) {
 105         if(vertex->GetNContributors() > 0) {
 106           if(vertex->GetZRes() != 0) {
 107             fHistEventStats->Fill(3); //events with a proper vertex
 108             if(TMath::Abs(vertex->GetXv()) < fVxMax) {
 109               if(TMath::Abs(vertex->GetYv()) < fVyMax) {
 110                 if(TMath::Abs(vertex->GetZv()) < fVzMax) {
 111                   fHistEventStats->Fill(4); //analayzed events
 112                   fHistVx->Fill(vertex->GetXv());
 113                   fHistVy->Fill(vertex->GetYv());
 114                   fHistVz->Fill(vertex->GetZv());
 115
 116                   //Printf("There are %d tracks in this event", gESD->GetNumberOfTracks());
 117                   for (Int_t iTracks = 0; iTracks < gESD->GetNumberOfTracks(); iTracks++) {
 118                     AliESDtrack* track = gESD->GetTrack(iTracks);
 119                     if (!track) {
 120                       Printf("ERROR: Could not receive track %d", iTracks);
 121                       continue;
 122                     }
 123
 124                     //ESD track cuts
 125                     if(fESDtrackCuts)
 126                       if(!fESDtrackCuts->AcceptTrack(track)) continue;
 127                   } //track loop
 128                 }//Vz cut
 129               }//Vy cut
 130             }//Vx cut
 131           }//proper vertex resolution
 132         }//proper number of contributors
 133       }//vertex object valid
 134     }//triggered event
 135   }//ESD analysis
 136   //AOD analysis
 137   else if(gAnalysisLevel == "AOD") {
 138     AliAODEvent* gAOD = dynamic_cast<AliAODEvent*>(InputEvent()); // from TaskSE
 139     if(!gAOD) {
 140       Printf("ERROR: gAOD not available");
 141       return;
 142     }
 143
 144     Printf("There are %d tracks in this event", gAOD->GetNumberOfTracks());
 145     for (Int_t iTracks = 0; iTracks < gAOD->GetNumberOfTracks(); iTracks++) {
 146       AliAODTrack* track = gAOD->GetTrack(iTracks);
 147       if (!track) {
 148         Printf("ERROR: Could not receive track %d", iTracks);
 149         continue;
 150       }
 151     } //track loop
 152   }//AOD analysis
 153   //MC analysis
 154   else if(gAnalysisLevel == "MC") {
 155
 156     AliMCEvent*  mcEvent = MCEvent();
 157     if (!mcEvent) {
 158       Printf("ERROR: mcEvent not available");
 159       return;
 160     }
 161
 162     Printf("There are %d tracks in this event", mcEvent->GetNumberOfPrimaries());
 163     for (Int_t iTracks = 0; iTracks < mcEvent->GetNumberOfPrimaries(); iTracks++) {
 164       AliMCParticle* track = dynamic_cast<AliMCParticle *>(mcEvent->GetTrack(iTracks));
 165       if (!track) {
 166         Printf("ERROR: Could not receive particle %d", iTracks);
 167         continue;
 168       }
 169     } //track loop
 170   }//MC analysis
 171 }
 172
 173 //________________________________________________________________________
 174 void AliAnalysisTaskChargeFluctuations::Terminate(Option_t *) {
 175   // Draw result to the screen
 176   // Called once at the end of the query
 177   fList = dynamic_cast<TList*> (GetOutputData(1));
 178   if (!fList) {
 179     Printf("ERROR: fList not available");
 180     return;
 181   }
 182 }