PWG4/JetTasks/AliAnaESDSpectraQA.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 /* $Id: $ */
  17 #include "AliAnaESDSpectraQA.h"
  18
  19 #include "TH1.h"
  20 #include "TH2.h"
  21 #include "TH3.h"
  22 #include "TList.h"
  23 #include "TChain.h"
  24 #include "TDirectory.h"
  25
  26 #include "AliAnalysisManager.h"
  27 #include "AliLog.h"
  28 #include "AliESDInputHandler.h"
  29 #include "AliESDEvent.h"
  30 #include "AliESDVertex.h"
  31 #include "AliESDtrack.h"
  32 #include "AliESDtrackCuts.h"
  33
  34
  35 const Int_t AliAnaESDSpectraQA::fgkNPtBins=38;
  36 const Float_t AliAnaESDSpectraQA::fgkPtMin=2;
  37 const Float_t AliAnaESDSpectraQA::fgkPtMax=40;
  38 const Int_t AliAnaESDSpectraQA::fgkNPhiBins=18;
  39
  40 ClassImp(AliAnaESDSpectraQA)
  41
  42 AliAnaESDSpectraQA::AliAnaESDSpectraQA(): AliAnalysisTask("AliAnaESDSpectraQA", ""),
  43   fESD(0),
  44   fTrackCuts(0),
  45   fNEvent(0), // just to avoid warnings, inititialized in InitPointers too
  46   fPtAll(0),  //
  47   fPtSel(0),  //
  48   fHistList(0)
  49 {
  50   InitHistPointers();
  51 }
  52
  53 AliAnaESDSpectraQA::AliAnaESDSpectraQA(const char *name):
  54   AliAnalysisTask(name, ""),
  55   fESD(0),
  56   fNEvent(0), // just to avoid warnings, inititialized in InitPointers too
  57   fPtAll(0),  //
  58   fPtSel(0),  //
  59   fHistList(0) {
  60   // Input slot #0 works with a TChain ESD
  61   DefineInput(0, TChain::Class());
  62   // Output slot #0 writes into a TList
  63   DefineOutput(0, TList::Class());
  64   InitHistPointers();
  65   fTrackCuts = new AliESDtrackCuts;
  66   fTrackCuts->SetAcceptKingDaughters(kFALSE);
  67   fTrackCuts->SetRequireTPCRefit(kTRUE);
  68   fTrackCuts->SetEtaRange(-1,1);
  69   // Add chisq criterium to reject 'bad' tracks that might not make it as prim
  70 }
  71
  72 //________________________________________________________________________
  73 void AliAnaESDSpectraQA::ConnectInputData(Option_t *)
  74 {
  75   // Connect ESD here
  76   // Called once
  77   TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
  78
  79   if( handler && handler->InheritsFrom("AliESDInputHandler") ) {
  80      fESD  =  ((AliESDInputHandler*)handler)->GetEvent();
  81   }
  82   else {
  83     AliFatal("I can't get any ESD Event Handler");
  84   }
  85 }
  86
  87 void AliAnaESDSpectraQA::InitHistPointers() {
  88   fNEvent = 0;
  89   fPtAll = 0;
  90   fPtSel = 0;
  91   for (Int_t i=0; i< 4; i++) {
  92     fHists[i].PhiPtNPointTPC = 0;
  93     fHists[i].PhiPtNPointITS = 0;
  94     fHists[i].PhiPtChisqC = 0;
  95     fHists[i].PhiPtChisqTPC = 0;
  96     fHists[i].PhiPtDCAR = 0;
  97     fHists[i].PhiPtDCAZ = 0;
  98   }
  99 }
 100
 101 void AliAnaESDSpectraQA::CreateOutputObjects() {
 102   fHistList = new TList();
 103   //fDirectory = new TDirectory("trig_hists","Trigger histos");
 104   //fDirectory->cd();
 105   TString labels[4];
 106   labels[kNegA]="NegA";
 107   labels[kPosA]="PosA";
 108   labels[kNegC]="NegC";
 109   labels[kPosC]="PosC";
 110
 111   static const Float_t kMinPhi = 0;
 112   static const Float_t kMaxPhi = 2*TMath::Pi();
 113   fNEvent = new TH1F("NEvent","NEvent",1,-0.5,0.5);
 114   fHistList->Add(fNEvent);
 115   fPtAll = new TH1F("PtAll","PtAll",fgkNPtBins, fgkPtMin, fgkPtMax);
 116   fHistList->Add(fPtAll);
 117   fPtSel = new TH1F("PtSel","PtSel",fgkNPtBins, fgkPtMin, fgkPtMax);
 118   fHistList->Add(fPtSel);
 119   for (Int_t iSide = 0; iSide < 4; iSide++) {
 120     TString hname="PhiPtNpointTPC";
 121     hname += labels[iSide];
 122     fHists[iSide].PhiPtNPointTPC = new TH3F(hname,hname+";#phi;p_{T} (GeV);N_{point,TPC}",fgkNPhiBins,kMinPhi,kMaxPhi,fgkNPtBins,fgkPtMin,fgkPtMax,160,0.5,160.5);
 123     fHistList->Add(fHists[iSide].PhiPtNPointTPC);
 124
 125     hname="PhiPtNpointITS";
 126     hname += labels[iSide];
 127     fHists[iSide].PhiPtNPointITS = new TH3F(hname,hname+";#phi;p_{T} (GeV);N_{point,ITS}",fgkNPhiBins,kMinPhi,kMaxPhi,fgkNPtBins,fgkPtMin,fgkPtMax,9,-0.5,8.5);
 128     fHistList->Add(fHists[iSide].PhiPtNPointITS);
 129
 130     hname="PhiPtChisqC";
 131     hname += labels[iSide];
 132     fHists[iSide].PhiPtChisqC = new TH3F(hname,hname+";#phi;p_{T} (GeV);#Chi^{2}/NDF",fgkNPhiBins,kMinPhi,kMaxPhi,fgkNPtBins,fgkPtMin,fgkPtMax,160,0,80);
 133     fHistList->Add(fHists[iSide].PhiPtChisqC);
 134
 135     hname="PhiPtChisqTPC";
 136     hname += labels[iSide];
 137     fHists[iSide].PhiPtChisqTPC = new TH3F(hname,hname+";#phi;p_{T} (GeV);#Chi^{2}/NDF",fgkNPhiBins,kMinPhi,kMaxPhi,fgkNPtBins,fgkPtMin,fgkPtMax,50,0,5);
 138     fHistList->Add(fHists[iSide].PhiPtChisqTPC);
 139
 140     hname="PhiPtDCAR";
 141     hname += labels[iSide];
 142     fHists[iSide].PhiPtDCAR = new TH3F(hname,hname+";#phi;p_{T} (GeV);DCAR",fgkNPhiBins,kMinPhi,kMaxPhi,fgkNPtBins,fgkPtMin,fgkPtMax,200,-1,1);
 143     fHistList->Add(fHists[iSide].PhiPtDCAR);
 144
 145     hname="PhiPtDCAZ";
 146     hname += labels[iSide];
 147     fHists[iSide].PhiPtDCAZ = new TH3F(hname,hname+";#phi;p_{T} (GeV);DCAZ",fgkNPhiBins,kMinPhi,kMaxPhi,fgkNPtBins,fgkPtMin,fgkPtMax,200,-2,2);
 148     fHistList->Add(fHists[iSide].PhiPtDCAZ);
 149
 150     hname="PhiPtSigmaToVertex";
 151     hname += labels[iSide];
 152     fHists[iSide].PhiPtSigmaToVertex = new TH3F(hname,hname+";#phi;p_{T} (GeV);n#sigma to vtx",fgkNPhiBins,kMinPhi,kMaxPhi,fgkNPtBins,fgkPtMin,fgkPtMax,50,0,5);
 153     fHistList->Add(fHists[iSide].PhiPtSigmaToVertex);
 154   }
 155 }
 156
 157 void AliAnaESDSpectraQA::Exec(Option_t *) {
 158   // Main loop
 159   // Called for each event
 160
 161   if (!fESD) {
 162     Printf("ERROR: fESD not available");
 163     return;
 164   }
 165
 166   const AliESDVertex *vtx = fESD->GetPrimaryVertex();
 167
 168   // Need vertex cut
 169   if (vtx->GetNContributors() < 2)
 170     return;
 171
 172   printf("Vertex title %s, status %d, nCont %d\n",vtx->GetTitle(), vtx->GetStatus(), vtx->GetNContributors());
 173   // Need to keep track of evts without vertex
 174
 175   Int_t nTracks = fESD->GetNumberOfTracks();
 176   AliDebug(2,Form("nTracks %d", nTracks));
 177   printf("nTracks %d\n", nTracks);
 178   static Int_t Mult = 0;
 179   Mult = 0;   // Need extra init bc of static
 180   for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
 181     AliESDtrack *track = fESD->GetTrack(iTrack);
 182     hists *curTypeHists = 0;
 183
 184     if (fTrackCuts->AcceptTrack(track)) {
 185       Mult++;
 186
 187       Float_t dca2D, dcaZ;
 188       track->GetImpactParameters(dca2D,dcaZ);
 189       fPtAll->Fill(track->Pt());
 190
 191       if (track->Eta() > 0) { // A side (crude for now)
 192         if (track->Charge() > 0)
 193           curTypeHists = &fHists[kPosA];
 194         else
 195           curTypeHists = &fHists[kNegA];
 196       }
 197       else { // C side
 198         if (track->Charge() > 0)
 199           curTypeHists = &fHists[kPosC];
 200         else
 201           curTypeHists = &fHists[kNegC];
 202       }
 203
 204       Float_t phi = track->Phi();
 205       Float_t pt = track->Pt();
 206
 207       UChar_t itsMap = track->GetITSClusterMap();
 208       Int_t nPointITS = 0;
 209       for (Int_t i=0; i < 6; i++) {
 210         if (itsMap & (1 << i))
 211           nPointITS ++;
 212       }
 213
 214       Float_t sigToVertex = fTrackCuts->GetSigmaToVertex(track);
 215       Float_t chisqC = 1000;
 216       if (track->GetConstrainedParam())
 217         chisqC = track->GetConstrainedChi2();
 218
 219       curTypeHists->PhiPtNPointTPC->Fill(phi,pt,track->GetTPCNcls());
 220       curTypeHists->PhiPtNPointITS->Fill(phi,pt,nPointITS);
 221       curTypeHists->PhiPtChisqC->Fill(phi,pt,chisqC);
 222       curTypeHists->PhiPtChisqTPC->Fill(phi,pt,track->GetTPCchi2()/(track->GetTPCNclsF()-5));
 223       curTypeHists->PhiPtDCAR->Fill(phi,pt,dca2D);
 224       curTypeHists->PhiPtDCAZ->Fill(phi,pt,dcaZ);
 225       curTypeHists->PhiPtSigmaToVertex->Fill(phi,pt,sigToVertex);
 226     }
 227   }
 228
 229   // Post output data
 230   PostData(0, fHistList);
 231 }