PWGLF/SPECTRA/MultEvShape/AliAnalysisCentralCutEvtESD.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 //  ESD event level cuts for azimuthal isotropic
  18 //  expansion in highly central collisions analysis
  19 //  author: Cristian Andrei
  20 //           acristian@niham.nipne.ro
  21 //  --------------------------------------------------
  22
  23 #include "TMath.h"
  24 #include <TObjArray.h>
  25 #include <TList.h>
  26
  27 #include "AliESDEvent.h"
  28 #include "AliESDtrack.h"
  29 #include "AliESDtrackCuts.h"
  30  #include "AliMultiplicity.h"
  31
  32 #include "AliAnalysisCentralCutESD.h"
  33 #include "AliAnalysisCentralCutEvtESD.h"
  34
  35 class TObject;
  36
  37
  38 //____________________________________________________________________
  39 ClassImp(AliAnalysisCentralCutEvtESD)
  40
  41 //____________________________________________________________________
  42 AliAnalysisCentralCutEvtESD::AliAnalysisCentralCutEvtESD(const Char_t* name, const Char_t* title)
  43     :AliAnalysisCuts(name,title)
  44     ,fReqMult(kFALSE)
  45     ,fReqDir(kFALSE)
  46     ,fReqSPDMult(kFALSE)
  47     ,fReqSPDDir(kFALSE)
  48     ,fMultMin(0)
  49     ,fMultMax(0)
  50     ,fDirMin(0)
  51     ,fDirMax(0)
  52     ,fSPDMultMin(0)
  53     ,fSPDMultMax(0)
  54     ,fSPDDirMin(0)
  55     ,fSPDDirMax(0)
  56
  57 {
  58 //constructor
  59     for(Int_t i=0; i<10; i++){
  60                 fCutsList[i] = 0;
  61     }
  62
  63         InitCuts();
  64
  65 }
  66
  67 AliAnalysisCentralCutEvtESD::~AliAnalysisCentralCutEvtESD() {
  68 //Destructor
  69   for (Int_t i=0; i<10; i++)
  70     delete fCutsList[i];
  71
  72 }
  73
  74
  75
  76 Bool_t AliAnalysisCentralCutEvtESD::IsSelected(TObject *obj){
  77 // check whether the event passes the cuts
  78     AliESDEvent *esdEvent = dynamic_cast<AliESDEvent *>(obj);
  79     if(!esdEvent){
  80                 printf("AliAnalysisCentralCutEvtESD:IsSelected ->Can't get ESD event!\n");
  81                 return kFALSE;
  82     }
  83
  84         if(fReqMult){
  85                 Int_t mult = CalcMult(esdEvent);
  86                 if((mult<fMultMin)||(mult>fMultMax)){
  87                         return kFALSE;
  88                 }
  89     }
  90
  91     if(fReqDir){
  92                 Double_t dir = CalcDir(esdEvent);
  93                 if((dir<fDirMin)||(dir>fDirMax)){
  94                         return kFALSE;
  95                 }
  96     }
  97
  98     if(fReqSPDMult){
  99                 Double_t spdMult = CalcSPDMult(esdEvent);
 100                 if((spdMult<fSPDMultMin)||(spdMult>fSPDMultMax)){
 101                         return kFALSE;
 102                 }
 103     }
 104
 105     if(fReqSPDDir){
 106                 Double_t spdDir = CalcSPDDir(esdEvent);
 107                 if((spdDir<fSPDDirMin)||(spdDir>fSPDDirMax)){
 108                         return kFALSE;
 109                 }
 110     }
 111
 112     return kTRUE;
 113 }
 114
 115 //_________________________________________________________________________
 116 void AliAnalysisCentralCutEvtESD::InitCuts(){
 117 //Initialize internal cuts
 118
 119 ////////////////  FULL ALICE ///////////////////
 120 //------------General ESD Cuts------------------
 121     AliESDtrackCuts *esdCutsGen = new AliESDtrackCuts("AliESDtrackCuts", "CutEvtESDInternal");
 122     esdCutsGen->SetMinNClustersTPC(50);
 123     esdCutsGen->SetMaxChi2PerClusterTPC(2.2);
 124     esdCutsGen->SetMaxCovDiagonalElements(0.5,0.5,0.5,0.5,0.5);
 125     esdCutsGen->SetRequireTPCRefit(kTRUE);
 126     esdCutsGen->SetAcceptKinkDaughters(kFALSE);
 127     esdCutsGen->SetMaxNsigmaToVertex(2.0);
 128     esdCutsGen->SetRequireSigmaToVertex(kTRUE);
 129
 130     AliAnalysisCentralCutESD *esdCutsGen1 = new AliAnalysisCentralCutESD("AliAnalysisCentralCutESD","NIHAM");
 131     esdCutsGen1->SetReqIsCharged();
 132
 133 //-------------Specific ESD Cuts------------------
 134     AliESDtrackCuts *esdCutsMult = new AliESDtrackCuts("AliESDCutsMult", "midrapidity");
 135     esdCutsMult->SetEtaRange(-0.5,0.5);
 136     AliESDtrackCuts *esdCutsDir = new AliESDtrackCuts("AliESDCutsDir", "SPD coverage");
 137     esdCutsDir->SetEtaRange(0.0,1.9);
 138
 139 ///////////////////  SPD ONLY  ////////////////////
 140     AliESDtrackCuts *esdCutsSPD = new AliESDtrackCuts("AliESDtrackCuts", "CutEvtESDInternal");
 141     esdCutsSPD->SetAcceptKinkDaughters(kFALSE);
 142     esdCutsSPD->SetMaxNsigmaToVertex(2.0);
 143     esdCutsSPD->SetRequireSigmaToVertex(kTRUE);
 144
 145 //--------------set the cuts ----------------------
 146     TObjArray* esdListMult = new TObjArray();
 147     esdListMult->AddLast(esdCutsGen);
 148     esdListMult->AddLast(esdCutsGen1);
 149     esdListMult->AddLast(esdCutsMult);
 150
 151     TObjArray* esdListDir = new TObjArray();
 152     esdListDir->AddLast(esdCutsGen);
 153     esdListDir->AddLast(esdCutsGen1);
 154     esdListDir->AddLast(esdCutsDir);
 155
 156     TObjArray* esdListSPD = new TObjArray();
 157     esdListSPD->AddLast(esdCutsSPD);
 158
 159     fCutsList[0]=esdListDir;
 160     fCutsList[1]=esdListMult;
 161     fCutsList[2]=esdListSPD;
 162
 163 }
 164
 165
 166 //__________________________________________________________________________
 167 Bool_t AliAnalysisCentralCutEvtESD::CheckIntCuts(Int_t no, TObject *obj) const{
 168 // Check if the particle passes the internal cuts
 169
 170     if(no > 9){
 171                 printf("\n AliAnalysisCentralCutEvtESD::CheckIntCuts -> Cut number is not ok! \n");
 172                 return kFALSE;
 173     }
 174
 175     if(!fCutsList[no]){
 176                 printf("AliAnalysisCentralCutEvtESD::CheckIntCuts -> cuts list problem! \n");
 177                 return kFALSE;
 178     }
 179
 180
 181     TObjArrayIter iter(fCutsList[no]);
 182     AliAnalysisCuts *cut = 0;
 183
 184
 185     while((cut = (AliAnalysisCuts*)iter.Next())){
 186
 187                 if(!cut->IsSelected(obj)) return kFALSE;
 188     }
 189
 190     return kTRUE;
 191 }
 192
 193
 194 //__________________________________________________________________________
 195 Double_t AliAnalysisCentralCutEvtESD::CalcDir(AliESDEvent* const esdEv) {
 196
 197 //Compute the directivity - FULL ALICE
 198
 199     Double_t dir;
 200     Double_t px,py;
 201     Double_t sumaPt = 0;
 202     Double_t sumaPx = 0;
 203     Double_t sumaPy = 0;
 204
 205     Double_t pt;
 206
 207     if (!esdEv){
 208                 printf("NULL tree\n");
 209                 return -1;
 210     }
 211
 212     Int_t totTracks=esdEv->GetNumberOfTracks();
 213
 214     for(Int_t itrack = 0; itrack < totTracks; itrack++){//first loop->compute event directivity
 215
 216                 AliESDtrack* track = esdEv->GetTrack(itrack);
 217
 218                 if (!track) {
 219                                 Printf("ERROR: Could not receive track %d", itrack);
 220                         continue;
 221                 }
 222
 223                 if(!CheckIntCuts(0, track)) continue;
 224
 225                 px = track->Px();
 226                 py = track->Py();
 227                 pt = track->Pt();
 228
 229                 sumaPx = sumaPx + px;
 230                 sumaPy = sumaPy + py;
 231
 232                 sumaPt = sumaPt + pt;
 233
 234     }//end track loop
 235
 236
 237         if(sumaPt < 0.0000001){
 238                 return -1;
 239     }
 240
 241     dir = (sqrt(pow(sumaPx,2)+pow(sumaPy,2)))/sumaPt;
 242
 243     return dir;
 244 }
 245
 246 //__________________________________________________________________________
 247 Int_t AliAnalysisCentralCutEvtESD::CalcMult(AliESDEvent* const esdEv) {
 248
 249 //Compute multiplicity - FULL ALICE
 250
 251     Int_t charged = 0;
 252
 253     if (!esdEv){
 254                 printf("NULL tree\n");
 255                 return -1;
 256     }
 257
 258     Int_t totTracks=esdEv->GetNumberOfTracks();
 259
 260     for(Int_t iTrack = 0; iTrack < totTracks; iTrack++){//second track loop -> compute event multiplicity
 261
 262                 AliESDtrack* track = esdEv->GetTrack(iTrack);
 263
 264                 if (!track) {
 265                         Printf("ERROR: Could not receive track %d", iTrack);
 266                         continue;
 267                 }
 268
 269                 if(!CheckIntCuts(1, track)) continue;
 270
 271
 272                 charged++; //multiplicity
 273
 274     }//end second track loop
 275
 276
 277     return charged;
 278 }
 279
 280
 281 //__________________________________________________________________________
 282 Double_t AliAnalysisCentralCutEvtESD::CalcSPDDir(AliESDEvent* const esdEv) {
 283
 284 //Compute directivity - SPD ONLY
 285
 286     Double_t dirU;
 287     Double_t pxU,pyU;
 288     Double_t sumaPxU = 0;
 289     Double_t sumaPyU = 0;
 290
 291     Double_t goodtrack = 0;
 292     Double_t spdEta = 0.;
 293
 294     const AliMultiplicity *spdMult=esdEv->GetMultiplicity();
 295     if(!spdMult){
 296         printf("Unable to get multiplicity! \n");
 297         return -1;
 298     }
 299
 300     Int_t spdTracks = spdMult->GetNumberOfTracklets();  //SPD multiplicity
 301
 302     for(Int_t iTrack = 0; iTrack< spdTracks; iTrack++){  //SPD track loop -> Directivity
 303
 304                 AliESDtrack* track = esdEv->GetTrack(iTrack);
 305
 306                 if (!track) {
 307                                 Printf("ERROR: Could not receive track %d", iTrack);
 308                                 continue;
 309                 }
 310
 311                 if(!CheckIntCuts(2, track)) continue;
 312
 313                 spdEta = spdMult->GetEta(iTrack);
 314
 315                 if((spdEta<0.0)||(spdEta>1.9)) continue;
 316
 317                 Double_t phi = spdMult->GetPhi(iTrack);
 318
 319                 pxU = TMath::Cos(phi);
 320                 pyU = TMath::Sin(phi);
 321
 322                 sumaPxU = sumaPxU + pxU;
 323                 sumaPyU = sumaPyU + pyU;
 324
 325                 goodtrack++;
 326
 327     }//end SPD track loop
 328
 329     if(goodtrack < 1.) return -1;
 330
 331     dirU = (sqrt(pow(sumaPxU,2)+pow(sumaPyU,2)))/goodtrack;
 332
 333     return dirU;
 334 }
 335
 336 //__________________________________________________________________________
 337
 338 Int_t AliAnalysisCentralCutEvtESD::CalcSPDMult(AliESDEvent* const esdEv) {
 339
 340         //Compute multiplicity - SPD ONLY
 341
 342     const AliMultiplicity *spdMult=esdEv->GetMultiplicity();
 343     if(!spdMult){
 344         printf("Unable to get multiplicity! \n");
 345         return -1;
 346     }
 347
 348     Int_t spdTracks = spdMult->GetNumberOfTracklets();  //SPD multiplicity
 349
 350     Double_t spdEta;
 351     Int_t charged = 0;
 352
 353     for(Int_t iTrack = 0; iTrack< spdTracks; iTrack++){  //second SPD track loop -> Multiplicity
 354
 355                 AliESDtrack* track = esdEv->GetTrack(iTrack);
 356
 357                 if (!track) {
 358                                 Printf("ERROR: Could not receive track %d", iTrack);
 359                                 continue;
 360                 }
 361
 362                 if(!CheckIntCuts(2, track)) continue;
 363
 364                 spdEta = spdMult->GetEta(iTrack);
 365
 366                 if((spdEta<-0.5)||(spdEta>0.5)) continue;
 367
 368                 charged++;
 369
 370     }//end second SPD track loop
 371
 372     return charged;
 373 }