PWGCF/Correlations/DPhi/DiHadronPID/AliAODEventCutsDiHadronPID.cxx

   1 /*************************************************************************
   2 * Copyright(c) 1998-2008, 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 //  Event Cut class for the DiHadronPID analysis.
  18 // -----------------------------------------------------------------------
  19 //  Author: Misha Veldhoen (misha.veldhoen@cern.ch)
  20
  21 #include <iostream>
  22 #include "TList.h"
  23 #include "TMath.h"
  24 #include "TH1F.h"
  25 #include "TH2F.h"
  26 #include "TNamed.h"
  27 #include "TIterator.h"
  28 #include "AliAODEvent.h"
  29 #include "AliAODHeader.h"
  30 #include "AliAODVertex.h"
  31 #include "AliCentrality.h"
  32 #include "AliAnalysisManager.h"
  33 #include "AliInputEventHandler.h"
  34 #include "AliAODEventCutsDiHadronPID.h"
  35
  36 using namespace std;
  37
  38 ClassImp(AliAODEventCutsDiHadronPID);
  39
  40 // -----------------------------------------------------------------------
  41 AliAODEventCutsDiHadronPID::AliAODEventCutsDiHadronPID():
  42         TNamed(),
  43         fIsPbPb(kTRUE),
  44         fIsMC(kFALSE),
  45         fTrigger(AliVEvent::kMB),
  46         fMinCentrality(5.),
  47         fMaxCentrality(0.),
  48         fCentralityEstimator("V0M"),
  49         fMaxVertexZ(10.),
  50         fMinRefMult(0),
  51         fTestTrigger(kFALSE),
  52         fTestCentrality(kFALSE),
  53         fTestVertexZ(kFALSE),
  54         fTestMinRefMult(kFALSE),
  55         fSelectedEventQAHistos(0x0),
  56         fAllEventQAHistos(0x0),
  57         fHistTrigger(0x0),
  58         fHistRefMultiplicity(0x0),
  59         fHistCentrality(0x0),
  60         fHistCentralityQuality(0x0),
  61         fHistVertexZ(0x0),
  62         fDebug(0)
  63
  64 {
  65
  66         //
  67         // Default Constructor
  68         //
  69
  70         cout<<"AliAODEventCutsDiHadronPID Default Constructor Called."<<endl;
  71         if (fDebug > 1) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
  72
  73 }
  74
  75 // -----------------------------------------------------------------------
  76 AliAODEventCutsDiHadronPID::AliAODEventCutsDiHadronPID(const char* name):
  77         TNamed(name,"AOD Event Cuts"),
  78         fIsPbPb(kTRUE),
  79         fIsMC(kFALSE),
  80         fTrigger(AliVEvent::kMB),
  81         fMinCentrality(5.),
  82         fMaxCentrality(0.),
  83         fCentralityEstimator("V0M"),
  84         fMaxVertexZ(10.),
  85         fMinRefMult(0),
  86         fTestTrigger(kFALSE),
  87         fTestCentrality(kFALSE),
  88         fTestVertexZ(kFALSE),
  89         fTestMinRefMult(kFALSE),
  90         fSelectedEventQAHistos(0x0),
  91         fAllEventQAHistos(0x0),
  92         fHistTrigger(0x0),
  93         fHistRefMultiplicity(0x0),
  94         fHistCentrality(0x0),
  95         fHistCentralityQuality(0x0),
  96         fHistVertexZ(0x0),
  97         fDebug(0)
  98
  99 {
 100
 101         //
 102         // Named Constructor
 103         //
 104
 105         cout<<"AliAODEventCutsDiHadronPID Named Constructor Called."<<endl;
 106         if (fDebug > 1) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
 107
 108 }
 109
 110 // -----------------------------------------------------------------------
 111 AliAODEventCutsDiHadronPID::~AliAODEventCutsDiHadronPID() {
 112
 113         //
 114         // Destructor
 115         //
 116
 117         cout<<"AliAODEventCutsDiHadronPID Destructor Called."<<endl;
 118         if (fDebug > 1) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
 119
 120         if (fSelectedEventQAHistos) delete fSelectedEventQAHistos;
 121         fSelectedEventQAHistos = 0x0;
 122         if (fAllEventQAHistos) delete fAllEventQAHistos;
 123         fAllEventQAHistos = 0x0;
 124
 125 }
 126
 127 // -----------------------------------------------------------------------
 128 Long64_t AliAODEventCutsDiHadronPID::Merge(TCollection* list) {
 129
 130         //
 131         // Merger.
 132         //
 133
 134         cout<<"AliAODEventCutsDiHadronPID Merger Called."<<endl;
 135         if (fDebug > 1) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
 136
 137         if (!list) return 0;
 138         if (list->IsEmpty()) return 1;
 139
 140         if (!fSelectedEventQAHistos||!fAllEventQAHistos) {
 141                 cout<<"AliAODEventCutsDiHadronPID::Merge() - Warning, current object's histograms are missing... Generating."<<endl;
 142                 CreateHistos();
 143         }
 144
 145         TIterator* iter = list->MakeIterator();
 146         TObject* obj;
 147
 148         // List of collections
 149         TList collection_fSelectedEventQAHistos;
 150         TList collection_fAllEventQAHistos;
 151
 152         Int_t count = 0;
 153
 154         while ((obj = iter->Next())) {
 155         AliAODEventCutsDiHadronPID* entry = dynamic_cast<AliAODEventCutsDiHadronPID*> (obj);
 156         if (entry == 0) continue;
 157
 158         // Check if the object to be merged really has the same name! (FIXME!)
 159
 160         // Getting the lists from obj.
 161         TList* list_fSelectedEventQAHistos = entry->GetListOfSelectedEventQAHistos();
 162         TList* list_fAllEventQAHistos = entry->GetListOfAllEventQAHistos();
 163
 164         // Adding the retrieved lists to the collection.
 165         if (list_fSelectedEventQAHistos) collection_fSelectedEventQAHistos.Add(list_fSelectedEventQAHistos);
 166         if (list_fAllEventQAHistos) collection_fAllEventQAHistos.Add(list_fAllEventQAHistos);
 167
 168         count++;
 169     }
 170
 171     // Merging. Note that we require the original list to exist.
 172     //  * Assume that if the collection happens to be empty, then nothing will happen.
 173     //  * All other variables are taken from the original object.
 174     if (fSelectedEventQAHistos) fSelectedEventQAHistos->Merge(&collection_fSelectedEventQAHistos);
 175     if (fAllEventQAHistos) fAllEventQAHistos->Merge(&collection_fAllEventQAHistos);
 176
 177     delete iter;
 178
 179         return count+1;
 180
 181 }
 182
 183 // -----------------------------------------------------------------------
 184 void AliAODEventCutsDiHadronPID::CreateHistos() {
 185
 186         cout<<"AliAODEventCutsDiHadronPID - Creating histograms"<<endl;
 187         if (fDebug > 1) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
 188
 189         // Create list of Event related QA histograms (selected events).
 190         fSelectedEventQAHistos = new TList();
 191         fSelectedEventQAHistos->SetName("SelectedEventQAHistos");
 192         fSelectedEventQAHistos->SetOwner(kTRUE);
 193
 194         // The same, but for all events.
 195         fAllEventQAHistos = new TList();
 196         fAllEventQAHistos->SetName("AllEventQAHistos");
 197         fAllEventQAHistos->SetOwner(kTRUE);
 198
 199         // Creating arrays of pointers to the QA histos.
 200         fHistTrigger = new TH1F*[2];
 201         fHistRefMultiplicity = new TH1F*[2];
 202         fHistCentrality = new TH1F*[2];
 203         fHistCentralityQuality = new TH1F*[2];
 204         fHistVertexZ = new TH1F*[2];
 205
 206         const char* HistType[2] = {"Selected","All"};
 207
 208         for (Int_t iHistType = 0; iHistType < 2; iHistType++) {
 209
 210                 // Trigger Histogram.
 211                 fHistTrigger[iHistType] = new TH1F(Form("fHistTrigger%s",HistType[iHistType]),"Trigger;;Count",2,-0.5,1.5);
 212                 fHistTrigger[iHistType]->GetXaxis()->SetBinLabel(1,"kMB");
 213                 fHistTrigger[iHistType]->GetXaxis()->SetBinLabel(2,"Other");
 214                 if (iHistType == 0) fSelectedEventQAHistos->Add(fHistTrigger[iHistType]);
 215                 else fAllEventQAHistos->Add(fHistTrigger[iHistType]);
 216
 217                 // Ref Multiplicity Histogram.
 218                 fHistRefMultiplicity[iHistType] = new TH1F(Form("fHistRefMultiplicity%s",HistType[iHistType]),"Reference Multiplicity;N_{tracks};Count",100,0.,10000.);
 219                 if (iHistType == 0) fSelectedEventQAHistos->Add(fHistRefMultiplicity[iHistType]);
 220                 else fAllEventQAHistos->Add(fHistRefMultiplicity[iHistType]);
 221
 222                 // Centrality Histogram.
 223                 fHistCentrality[iHistType] = new TH1F(Form("fHistCentrality%s",HistType[iHistType]),"Centrality;Centrality;Count",20,0,100);
 224                 if (iHistType == 0) fSelectedEventQAHistos->Add(fHistCentrality[iHistType]);
 225                 else fAllEventQAHistos->Add(fHistCentrality[iHistType]);
 226
 227                 // Centrality Quality.
 228                 fHistCentralityQuality[iHistType] = new TH1F(Form("fHistCentralityQuality%s",HistType[iHistType]),"Centrality Quality;Quality;Count",2,-0.5,1.5);
 229                 fHistCentralityQuality[iHistType]->GetXaxis()->SetBinLabel(1,"0");
 230                 fHistCentralityQuality[iHistType]->GetXaxis()->SetBinLabel(2,"Other");
 231                 if (iHistType == 0) fSelectedEventQAHistos->Add(fHistCentralityQuality[iHistType]);
 232                 else fAllEventQAHistos->Add(fHistCentralityQuality[iHistType]);
 233
 234                 // VertexZ Histogram.
 235                 fHistVertexZ[iHistType] = new TH1F(Form("fHistVertexZ%s",HistType[iHistType]),"VertexZ;z (cm);Count",60,-15.,15.);
 236                 if (iHistType == 0) fSelectedEventQAHistos->Add(fHistVertexZ[iHistType]);
 237                 else fAllEventQAHistos->Add(fHistVertexZ[iHistType]);
 238
 239         }
 240
 241 }
 242
 243 // -----------------------------------------------------------------------
 244 Bool_t AliAODEventCutsDiHadronPID::IsSelected(AliAODEvent* event) {
 245
 246         if (fDebug > 1) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
 247
 248         if (!event) return kFALSE;
 249
 250         if (!fAllEventQAHistos||!fSelectedEventQAHistos) {cout<<"AliAODEventCutsDiHadronPID - Histograms were not created, you should have called CreateHistos()..."<<endl;}
 251
 252         // Input the event handler.
 253         AliInputEventHandler* InputHandler = (AliInputEventHandler*)((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
 254         if (!InputHandler) return kFALSE;
 255
 256         Bool_t select = kTRUE;
 257
 258         // Test Trigger.
 259         UInt_t trigger = InputHandler->IsEventSelected();
 260         Int_t triggerselect = 0; // 0 = selected.
 261     if (fTestTrigger) {
 262                 if (!(trigger & fTrigger)) {
 263                         select = kFALSE;
 264                         triggerselect = 1;      // 1 = not selected.
 265                 }
 266         }
 267
 268         AliCentrality* CurrentCentrality = 0x0;
 269         Int_t CurrentCentralityQuality = -999;
 270         Float_t percentile = -999.;
 271
 272         if (fIsPbPb) {
 273
 274                 // Get the centrality object.
 275             CurrentCentrality = event->GetCentrality();
 276             if (!CurrentCentrality) {select = kFALSE; return select;}
 277
 278             // Check the quality of the centrality estimation.
 279             // If 0 then quality is OK, c.f. TOF/PbPb276/macros/TOFmatchEff.C
 280             CurrentCentralityQuality = CurrentCentrality->GetQuality();
 281             //cout<<"Centrality: "<<CurrentCentrality->GetCentralityPercentile(fCentralityEstimator.Data())<<" Quality: "<<CurrentCentrality->GetQuality()<<endl;
 282             if (CurrentCentralityQuality) select = kFALSE;
 283
 284                 // Test Centrality.
 285             percentile = CurrentCentrality->GetCentralityPercentile(fCentralityEstimator.Data());
 286                 if (fTestCentrality) {
 287                 if ((percentile < fMaxCentrality)||(percentile > fMinCentrality)) select = kFALSE;
 288                 }
 289
 290         }
 291
 292         // Get the primary vertex.
 293         AliAODVertex* CurrentPrimaryVertex = event->GetPrimaryVertex();
 294     if (!CurrentPrimaryVertex) {select = kFALSE; return select;}
 295
 296         // Test Vertex Z.
 297     Double_t vtxz = CurrentPrimaryVertex->GetZ();
 298         if (fTestVertexZ) {
 299         if (TMath::Abs(vtxz) > fMaxVertexZ) select = kFALSE;
 300         }
 301
 302         // Get the event header.
 303         AliAODHeader* CurrentHeader = event->GetHeader();
 304
 305         // Test minimum reference multiplicity.
 306         Int_t CurrentRefMultiplicity = CurrentHeader->GetRefMultiplicity();
 307         if (fTestMinRefMult) {
 308                 if (CurrentRefMultiplicity < fMinRefMult) select = kFALSE;
 309         }
 310
 311         // Fill the histograms for selected events.
 312         if (select) {
 313                 fHistTrigger[0]->Fill(triggerselect);
 314                 fHistRefMultiplicity[0]->Fill(CurrentHeader->GetRefMultiplicity());
 315                 if (fIsPbPb) fHistCentrality[0]->Fill(percentile);
 316                 if (fIsPbPb) fHistCentralityQuality[0]->Fill(CurrentCentralityQuality);
 317                 fHistVertexZ[0]->Fill(vtxz);
 318         }
 319
 320         // Fill the histograms for all events.
 321         fHistTrigger[1]->Fill(triggerselect);
 322         fHistRefMultiplicity[1]->Fill(CurrentHeader->GetRefMultiplicity());
 323         if (fIsPbPb) fHistCentrality[1]->Fill(percentile);
 324         if (fIsPbPb) fHistCentralityQuality[1]->Fill(CurrentCentralityQuality);
 325         fHistVertexZ[1]->Fill(vtxz);
 326
 327         //cout<<"Event Selected: "<<select<<endl;
 328
 329         return select;
 330
 331 }
 332
 333 // -----------------------------------------------------------------------
 334 void AliAODEventCutsDiHadronPID::PrintCuts() {
 335
 336         // NOT IMPLEMENTED.
 337         if (fDebug > 1) {cout << Form("File: %s, Line: %i, Function: %s",__FILE__,__LINE__,__func__) << endl;}
 338
 339         return;
 340
 341 }