PWG4/JetTasks/AliPWG4HighPtSpectra.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 //-----------------------------------------------------------------------
  17 // Example of task (running locally, on AliEn and CAF),
  18 // which provides standard way of calculating acceptance and efficiency
  19 // between different steps of the procedure.
  20 // The ouptut of the task is a AliCFContainer from which the efficiencies
  21 // can be calculated
  22 //-----------------------------------------------------------------------
  23 // Author : Marta Verweij - UU
  24 //-----------------------------------------------------------------------
  25
  26
  27 #ifndef ALIPWG4HighPtSpectra_CXX
  28 #define ALIPWG4HighPtSpectra_CXX
  29
  30 #include "AliPWG4HighPtSpectra.h"
  31
  32 #include "AliStack.h"
  33 #include "TParticle.h"
  34 #include "TH1I.h"
  35 #include "AliMCEvent.h"
  36 #include "AliMCEventHandler.h"
  37 #include "AliAnalysisManager.h"
  38 #include "AliCFContainer.h"
  39 #include "TChain.h"
  40 #include "AliESDtrack.h"
  41 #include "AliESDtrackCuts.h"
  42 #include "AliExternalTrackParam.h"
  43 #include "AliESDInputHandler.h"
  44 #include "AliAnalysisHelperJetTasks.h"
  45
  46 using namespace std; //required for resolving the 'cout' symbol
  47
  48 ClassImp(AliPWG4HighPtSpectra)
  49
  50 //__________________________________________________________________________
  51 AliPWG4HighPtSpectra::AliPWG4HighPtSpectra() : AliAnalysisTask("AliPWG4HighPtSpectra", ""),
  52   fReadAODData(0),
  53   fCFManager(0x0),
  54   fESD(0),
  55   fTrackCuts(0),
  56   fTrigger(0),
  57   fHistList(0),
  58   fNEventAll(0),
  59   fNEventSel(0)
  60   //  fHistEventsProcessed(0x0)
  61 {
  62   //
  63   //Default ctor
  64   //
  65 }
  66 //___________________________________________________________________________
  67 AliPWG4HighPtSpectra::AliPWG4HighPtSpectra(const Char_t* name) :
  68   AliAnalysisTask(name,""),
  69   fReadAODData(0),
  70   fCFManager(0x0),
  71   fESD(0),
  72   fTrackCuts(),//new AliESDtrackCuts),
  73   fTrigger(0),
  74   fHistList(0),
  75   fNEventAll(0),
  76   fNEventSel(0)
  77   //  fHistEventsProcessed(0x0)
  78 {
  79   //
  80   // Constructor. Initialization of Inputs and Outputs
  81   //
  82   AliDebug(2,Form("AliPWG4HighPtQAMC","Calling Constructor"));
  83   // Input slot #0 works with a TChain ESD
  84   DefineInput(0, TChain::Class());
  85   DefineOutput(0,TList::Class());
  86   DefineOutput(1,AliCFContainer::Class());
  87 }
  88
  89 //___________________________________________________________________________
  90 AliPWG4HighPtSpectra& AliPWG4HighPtSpectra::operator=(const AliPWG4HighPtSpectra& c)
  91 {
  92   //
  93   // Assignment operator
  94   //
  95   if (this!=&c) {
  96     AliAnalysisTask::operator=(c) ;
  97     fReadAODData = c.fReadAODData ;
  98     fCFManager  = c.fCFManager;
  99     fHistList = c.fHistList;
 100     fNEventAll = c.fNEventAll;
 101     fNEventSel = c.fNEventSel;
 102     //    fHistEventsProcessed = c.fHistEventsProcessed;
 103   }
 104   return *this;
 105 }
 106
 107 //___________________________________________________________________________
 108 AliPWG4HighPtSpectra::AliPWG4HighPtSpectra(const AliPWG4HighPtSpectra& c) :
 109   AliAnalysisTask(c),
 110   fReadAODData(c.fReadAODData),
 111   fCFManager(c.fCFManager),
 112   fESD(c.fESD),
 113   fTrackCuts(c.fTrackCuts),
 114   fTrigger(c.fTrigger),
 115   fHistList(c.fHistList),
 116   fNEventAll(c.fNEventAll),
 117   fNEventSel(c.fNEventSel)
 118   //  fHistEventsProcessed(c.fHistEventsProcessed)
 119 {
 120   //
 121   // Copy Constructor
 122   //
 123 }
 124
 125 //___________________________________________________________________________
 126 AliPWG4HighPtSpectra::~AliPWG4HighPtSpectra() {
 127   //
 128   //destructor
 129   //
 130   Info("~AliPWG4HighPtSpectra","Calling Destructor");
 131   if (fCFManager)           delete fCFManager ;
 132   if (fNEventAll) delete fNEventAll ;
 133   if (fNEventSel) delete fNEventSel ;
 134 }
 135 //________________________________________________________________________
 136 void AliPWG4HighPtSpectra::ConnectInputData(Option_t *)
 137 {
 138   // Connect ESD here
 139   // Called once
 140   AliDebug(2,Form(">> AliPWG4HighPtSpectra::ConnectInputData \n"));
 141   //  printf(">> AliPWG4HighPtSpectra::ConnectInputData \n");
 142
 143   TTree* tree = dynamic_cast<TTree*> (GetInputData(0));
 144   if (!tree) {
 145     AliDebug(2,Form("ERROR: Could not read chain from input slot 0"));
 146   } else {
 147
 148     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
 149
 150     if (!esdH) {
 151       AliDebug(2,Form("ERROR: Could not get ESDInputHandler"));
 152     } else
 153       fESD = esdH->GetEvent();
 154   }
 155
 156 }
 157 //_________________________________________________
 158 void AliPWG4HighPtSpectra::Exec(Option_t *)//UserExec(Option_t *)
 159 {
 160   //
 161   // Main loop function
 162   //
 163   AliDebug(2,Form(">> AliPWG4HighPtSpectra::Exec \n"));
 164
 165   // All events without selection
 166   fNEventAll->Fill(0.);
 167
 168   if (!fESD) {
 169     AliDebug(2,Form("ERROR: fESD not available"));
 170     return;
 171   }
 172
 173   //Trigger selection
 174   AliAnalysisHelperJetTasks::Trigger trig;
 175   trig = (const enum AliAnalysisHelperJetTasks::Trigger)fTrigger;
 176   if (AliAnalysisHelperJetTasks::IsTriggerFired(fESD,trig)){
 177     AliDebug(2,Form(" Trigger Selection: event ACCEPTED ... "));
 178   }else{
 179     AliDebug(2,Form(" Trigger Selection: event REJECTED ... "));
 180     return;
 181   }
 182   //  if(!fESD->IsTriggerClassFired("CINT1B-ABCE-NOPF-ALL") || !fESD->IsTriggerClassFired("CSMBB-ABCE-NOPF-ALL")) return;
 183
 184   // Process MC truth, therefore we receive the AliAnalysisManager and ask it for the AliMCEventHandler
 185   // This handler can return the current MC event
 186   AliMCEventHandler* eventHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
 187   AliStack* stack = 0x0;
 188   AliMCEvent* mcEvent = 0x0;
 189
 190   if(eventHandler) {
 191     mcEvent = eventHandler->MCEvent();
 192     if (!mcEvent) {
 193       AliDebug(2,Form("ERROR: Could not retrieve MC event"));
 194       return;
 195     }
 196
 197     AliDebug(2,Form("MC particles: %d", mcEvent->GetNumberOfTracks()));
 198
 199     stack = mcEvent->Stack();                //Particles Stack
 200
 201     AliDebug(2,Form("MC particles stack: %d", stack->GetNtrack()));
 202   }
 203
 204   const AliESDVertex *vtx = fESD->GetPrimaryVertex();
 205   // Need vertex cut
 206   if (vtx->GetNContributors() < 2) return;
 207   double primVtx[3];
 208   vtx->GetXYZ(primVtx);
 209   //  printf("primVtx: %g  %g  %g \n",primVtx[0],primVtx[1],primVtx[2]);
 210   if(primVtx[0]>1. || primVtx[1]>1. || primVtx[2]>10.) return;
 211
 212   AliDebug(2,Form("Vertex title %s, status %d, nCont %d\n",vtx->GetTitle(), vtx->GetStatus(), vtx->GetNContributors()));
 213
 214   Int_t nTracks = fESD->GetNumberOfTracks();
 215   AliDebug(2,Form("nTracks %d", nTracks));
 216
 217   // Selected events for analysis
 218   fNEventSel->Fill(0.);
 219
 220   Double_t containerInputRec[1] ;
 221   Double_t containerInputTPConly[1] ;
 222   Double_t containerInputMC[1] ;
 223   //Now go to rec level
 224   for (Int_t iTrack = 0; iTrack<nTracks; iTrack++)
 225     {
 226       if(!fESD->GetTrack(iTrack) ) continue;
 227       AliESDtrack* track = fESD->GetTrack(iTrack);
 228       if(!(AliExternalTrackParam *)track->GetTPCInnerParam()) continue;
 229       AliExternalTrackParam *trackTPC = (AliExternalTrackParam *)track->GetTPCInnerParam();
 230       if(!track || !trackTPC) continue;
 231
 232
 233       //fill the container
 234       containerInputRec[0] = track->Pt();
 235       containerInputTPConly[0] = trackTPC->Pt();
 236
 237       if (fTrackCuts->AcceptTrack(track)) {
 238         fCFManager->GetParticleContainer()->Fill(containerInputRec,kStepReconstructed);
 239         fCFManager->GetParticleContainer()->Fill(containerInputTPConly,kStepReconstructedTPCOnly);
 240
 241         //Only fill the secondary particle container if MC information is available
 242         if(eventHandler) {
 243           Int_t label = TMath::Abs(track->GetLabel());
 244           TParticle *particle = stack->Particle(label) ;
 245           if(!particle) continue;
 246           containerInputMC[0] = particle->Pt();
 247           fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepReconstructedMC);
 248           if (!stack->IsPhysicalPrimary(label) ) {
 249             fCFManager->GetParticleContainer()->Fill(containerInputRec,kStepSecondaries);
 250           }
 251         }
 252       }
 253
 254     }
 255
 256   if(eventHandler) {
 257     for(int iPart = 1; iPart<(mcEvent->GetNumberOfTracks()); iPart++)//stack->GetNprimary();
 258       {
 259         AliMCParticle *mcPart  = (AliMCParticle*)mcEvent->GetTrack(iPart);
 260         if(!mcPart) continue;
 261
 262         //fill the container
 263         containerInputMC[0] = mcPart->Pt();
 264
 265         if (!fCFManager->CheckParticleCuts(3,mcPart)) continue ;
 266
 267         int counter;
 268
 269         Float_t trackLengthTPC = mcPart->GetTPCTrackLength(fESD->GetMagneticField(),0.1,counter,3.0);
 270
 271         if(trackLengthTPC>80.) fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepMCtrackable) ;
 272
 273       }
 274   }
 275
 276    PostData(0,fHistList);
 277    PostData(1,fCFManager->GetParticleContainer());
 278
 279 }
 280
 281
 282 //___________________________________________________________________________
 283 void AliPWG4HighPtSpectra::Terminate(Option_t*)
 284 {
 285   // The Terminate() function is the last function to be called during
 286   // a query. It always runs on the client, it can be used to present
 287   // the results graphically or save the results to file.
 288
 289
 290 }
 291
 292 //___________________________________________________________________________
 293 void AliPWG4HighPtSpectra::CreateOutputObjects() {
 294   //HERE ONE CAN CREATE OUTPUT OBJECTS, IN PARTICULAR IF THE OBJECT PARAMETERS DON'T NEED
 295   //TO BE SET BEFORE THE EXECUTION OF THE TASK
 296   //
 297   AliDebug(2,Form("CreateOutputObjects","CreateOutputObjects of task %s", GetName()));
 298
 299   //slot #1
 300   OpenFile(0);
 301   fHistList = new TList();
 302   fNEventAll = new TH1F("fNEventAll","NEventAll",1,-0.5,0.5);
 303   fHistList->Add(fNEventAll);
 304   fNEventSel = new TH1F("fNEventSel","NEvent Selected for analysis",1,-0.5,0.5);
 305   fHistList->Add(fNEventSel);
 306
 307 }
 308
 309 #endif