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