PWGLF/SPECTRA/PiKaPr/TestAOD/AliSpectraAODEventCuts.cxx

   1
   2 /**************************************************************************
   3  * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *
   4  *                                                                        *
   5  * Author: The ALICE Off-line Project.                                    *
   6  * Contributors are mentioned in the code where appropriate.              *
   7  *                                                                        *
   8  * Permission to use, copy, modify and distribute this software and its   *
   9  * documentation strictly for non-commercial purposes is hereby granted   *
  10  * without fee, provided that the above copyright notice appears in all   *
  11  * copies and that both the copyright notice and this permission notice   *
  12  * appear in the supporting documentation. The authors make no claims     *
  13  * about the suitability of this software for any purpose. It is          *
  14  * provided "as is" without express or implied warranty.                  *
  15  **************************************************************************/
  16
  17 //-----------------------------------------------------------------
  18 //         AliSpectraAODEventCuts class
  19 //-----------------------------------------------------------------
  20
  21 #include "TChain.h"
  22 #include "TTree.h"
  23 #include "TLegend.h"
  24 #include "TH1F.h"
  25 #include "TH2F.h"
  26 #include "TCanvas.h"
  27 #include "AliAnalysisTask.h"
  28 #include "AliAnalysisManager.h"
  29 #include "AliAODTrack.h"
  30 #include "AliAODMCParticle.h"
  31 #include "AliAODEvent.h"
  32 #include "AliAODInputHandler.h"
  33 #include "AliAnalysisTaskESDfilter.h"
  34 #include "AliAnalysisDataContainer.h"
  35 #include "AliSpectraAODEventCuts.h"
  36 #include "AliSpectraAODTrackCuts.h"
  37 #include "AliSpectraAODHistoManager.h"
  38 #include <iostream>
  39
  40 using namespace std;
  41
  42 ClassImp(AliSpectraAODEventCuts)
  43
  44 AliSpectraAODEventCuts::AliSpectraAODEventCuts(const char *name) : TNamed(name, "AOD Event Cuts"), fAOD(0), fIsSelected(0), fCentralityCutMin(0), fCentralityCutMax(0), fHistoCuts(0),fHistoVtxBefSel(0),fHistoVtxAftSel(0),fHistoEtaBefSel(0),fHistoEtaAftSel(0),fHistoNChAftSel(0)
  45 {
  46   // Constructor
  47   fHistoCuts = new TH1I("fEventCuts", "Event Cuts", kNVtxCuts, -0.5, kNVtxCuts - 0.5);
  48   fHistoVtxBefSel = new TH1F("fHistoVtxBefSel", "Vtx distr before event selection",500,-15,15);
  49   fHistoVtxAftSel = new TH1F("fHistoVtxAftSel", "Vtx distr after event selection",500,-15,15);
  50   fHistoEtaBefSel = new TH1F("fHistoEtaBefSel", "Eta distr before event selection",500,-2,2);
  51   fHistoEtaAftSel = new TH1F("fHistoEtaAftSel", "Eta distr after event selection",500,-2,2);
  52   fHistoNChAftSel = new TH1F("fHistoNChAftSel", "NCh distr after event selection",3000,-0.5,2999.5);
  53   fCentralityCutMin = 0.0;      // default value of centrality cut minimum, 0 ~ no cut
  54   fCentralityCutMax = 10000.0;  // default value of centrality cut maximum,  ~ no cut
  55
  56 }
  57
  58 //______________________________________________________
  59 Bool_t AliSpectraAODEventCuts::IsSelected(AliAODEvent * aod,AliSpectraAODTrackCuts     *trackcuts)
  60 {
  61   // Returns true if Event Cuts are selected and applied
  62   fAOD = aod;
  63   fTrackCuts = trackcuts;
  64   fHistoCuts->Fill(kProcessedEvents);
  65   //loop on tracks, before event selection, filling QA histos
  66   AliAODVertex * vertex = fAOD->GetPrimaryVertex();//FIXME vertex is recreated
  67   if(vertex)fHistoVtxBefSel->Fill(vertex->GetZ());
  68   fIsSelected = (CheckVtxRange() && CheckCentralityCut());
  69   if(fIsSelected){
  70     fHistoCuts->Fill(kAcceptedEvents);
  71     if(vertex)fHistoVtxAftSel->Fill(vertex->GetZ());
  72   }
  73   Int_t Nch=0;
  74   for (Int_t iTracks = 0; iTracks < fAOD->GetNumberOfTracks(); iTracks++) {
  75     AliAODTrack* track = fAOD->GetTrack(iTracks);
  76     if (!fTrackCuts->IsSelected(track)) continue;
  77     fHistoEtaBefSel->Fill(track->Eta());
  78     if(fIsSelected){
  79       fHistoEtaAftSel->Fill(track->Eta());
  80       Nch++;
  81     }
  82   }
  83   if(fIsSelected)fHistoNChAftSel->Fill(Nch);
  84   return fIsSelected;
  85 }
  86
  87 //______________________________________________________
  88 Bool_t AliSpectraAODEventCuts::CheckVtxRange()
  89 {
  90   // reject events outside of range
  91   AliAODVertex * vertex = fAOD->GetPrimaryVertex();
  92   if (!vertex)
  93     {
  94       fHistoCuts->Fill(kVtxNoEvent);
  95       return kFALSE;
  96     }
  97   if (TMath::Abs(vertex->GetZ()) < 10)
  98     {
  99       return kTRUE;
 100     }
 101   fHistoCuts->Fill(kVtxRange);
 102   return kFALSE;
 103 }
 104
 105 //______________________________________________________
 106 Bool_t AliSpectraAODEventCuts::CheckCentralityCut()
 107 {
 108   // Check centrality cut
 109   if ( (fAOD->GetCentrality()->GetCentralityPercentile("V0M") <= fCentralityCutMax)  &&  (fAOD->GetCentrality()->GetCentralityPercentile("V0M") >= fCentralityCutMin) )  return kTRUE;
 110   fHistoCuts->Fill(kVtxCentral);
 111   return kFALSE;
 112 }
 113
 114 //______________________________________________________
 115 void AliSpectraAODEventCuts::PrintCuts()
 116 {
 117   // print info about event cuts
 118   cout << "Event Stats" << endl;
 119   cout << " > Number of accepted events: " << fHistoCuts->GetBinContent(kAcceptedEvents + 1) << endl;
 120   cout << " > Number of processed events: " << fHistoCuts->GetBinContent(kProcessedEvents + 1) << endl;
 121   cout << " > Vertex out of range: " << fHistoCuts->GetBinContent(kVtxRange + 1) << endl;
 122   cout << " > Events cut by centrality: " << fHistoCuts->GetBinContent(kVtxCentral + 1) << endl;
 123   cout << " > Events without vertex: " << fHistoCuts->GetBinContent(kVtxNoEvent + 1) << endl;
 124 }
 125 //______________________________________________________
 126
 127 Long64_t AliSpectraAODEventCuts::Merge(TCollection* list)
 128 {
 129   // Merge a list of AliSpectraAODEventCuts objects with this.
 130   // Returns the number of merged objects (including this).
 131
 132   //  AliInfo("Merging");
 133
 134   if (!list)
 135     return 0;
 136
 137   if (list->IsEmpty())
 138     return 1;
 139
 140   TIterator* iter = list->MakeIterator();
 141   TObject* obj;
 142
 143   // collections of all histograms
 144   TList collections;//FIXME we should use only 1 collection
 145   TList collections_histoVtxBefSel;
 146   TList collections_histoVtxAftSel;
 147   TList collections_histoEtaBefSel;
 148   TList collections_histoEtaAftSel;
 149   TList collections_histoNChAftSel;
 150
 151   Int_t count = 0;
 152
 153   while ((obj = iter->Next())) {
 154     AliSpectraAODEventCuts* entry = dynamic_cast<AliSpectraAODEventCuts*> (obj);
 155     if (entry == 0)
 156       continue;
 157
 158     TH1I * histo = entry->GetHistoCuts();
 159     collections.Add(histo);
 160     TH1F * histo_histoVtxBefSel = entry->GetHistoVtxBefSel();
 161     collections_histoVtxBefSel.Add(histo_histoVtxBefSel);
 162     TH1F * histo_histoVtxAftSel = entry->GetHistoVtxAftSel();
 163     collections_histoVtxAftSel.Add(histo_histoVtxAftSel);
 164     TH1F * histo_histoEtaBefSel = entry->GetHistoEtaBefSel();
 165     collections_histoEtaBefSel.Add(histo_histoEtaBefSel);
 166     TH1F * histo_histoEtaAftSel = entry->GetHistoEtaAftSel();
 167     collections_histoEtaAftSel.Add(histo_histoEtaAftSel);
 168     TH1F * histo_histoNChAftSel = entry->GetHistoNChAftSel();
 169     collections_histoNChAftSel.Add(histo_histoNChAftSel);
 170     count++;
 171   }
 172
 173   fHistoCuts->Merge(&collections);
 174   fHistoVtxBefSel->Merge(&collections_histoVtxBefSel);
 175   fHistoVtxAftSel->Merge(&collections_histoVtxAftSel);
 176   fHistoEtaBefSel->Merge(&collections_histoEtaBefSel);
 177   fHistoEtaAftSel->Merge(&collections_histoEtaAftSel);
 178   fHistoNChAftSel->Merge(&collections_histoNChAftSel);
 179
 180   delete iter;
 181
 182   return count+1;
 183 }
 184
 185 /// FIXME: Q vector
 186 // //Selection on QVector, before ANY other selection on the event
 187 // //Spectra MUST be normalized wrt events AFTER the selection on Qvector
 188 // // Can we include this in fEventCuts
 189 // Double_t Qx2EtaPos = 0, Qy2EtaPos = 0;
 190 // Double_t Qx2EtaNeg = 0, Qy2EtaNeg = 0;
 191 // Int_t multPos = 0;
 192 // Int_t multNeg = 0;
 193 // for(Int_t iT = 0; iT < fAOD->GetNumberOfTracks(); iT++) {
 194 //   AliAODTrack* aodTrack = fAOD->GetTrack(iT);
 195 //   if (!fTrackCuts->IsSelected(aodTrack)) continue;
 196 //   if (aodTrack->Eta() >= 0){
 197 //     multPos++;
 198 //     Qx2EtaPos += TMath::Cos(2*aodTrack->Phi());
 199 //     Qy2EtaPos += TMath::Sin(2*aodTrack->Phi());
 200 //   } else {
 201 //     multNeg++;
 202 //     Qx2EtaNeg += TMath::Cos(2*aodTrack->Phi());
 203 //     Qy2EtaNeg += TMath::Sin(2*aodTrack->Phi());
 204 //   }
 205 // }
 206 // Double_t qPos=-999;
 207 // if(multPos!=0)qPos= TMath::Sqrt((Qx2EtaPos*Qx2EtaPos + Qy2EtaPos*Qy2EtaPos)/multPos);
 208 // Double_t qNeg=-999;
 209 // if(multNeg!=0)qNeg= TMath::Sqrt((Qx2EtaNeg*Qx2EtaNeg + Qy2EtaNeg*Qy2EtaNeg)/multNeg);
 210
 211 // if((qPos>fTrackCuts->GetQvecMin() && qPos<fTrackCuts->GetQvecMax()) || (qNeg>fTrackCuts->GetQvecMin() && qNeg<fTrackCuts->GetQvecMax())){
 212
 213 //fill q distributions vs centrality, after all event selection
 214 // fHistMan->GetqVecHistogram(kHistqVecPos)->Fill(qPos,fAOD->GetCentrality()->GetCentralityPercentile("V0M"));  // qVector distribution
 215 // fHistMan->GetqVecHistogram(kHistqVecNeg)->Fill(qNeg,fAOD->GetCentrality()->GetCentralityPercentile("V0M"));  // qVector distribution