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