PWGDQ/dielectron/AliDielectronV0Cuts.cxx

   1 /*************************************************************************
   2 * Copyright(c) 1998-2009, 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 //   Cut class providing cuts to V0 candidates                           //
  18 //   Selection or deselection of V0 candiates can be done.               //
  19 //                                                                       //
  20 // Authors:                                                              //
  21 //   Julian Book <Julian.Book@cern.ch>                                   //
  22 /*
  23
  24 Class to provide some V0 selection using exactley the same cuts in ESDs as in AODs.
  25 It implements the PID cut class AliDielectronPID and the standard AliDielectronVarCuts for
  26 the configuration of leg respective pair cuts. These pair cuts are applied on the KFparticle
  27 build by the legs.
  28
  29 Some QA cuts for the tracks are applied before the V0 pair is build. These cuts are:
  30 AliDielectronVarCuts dauQAcuts1;
  31 dauQAcuts1.AddCut(AliDielectronVarManager::kPt,            0.05, 100.0);
  32 dauQAcuts1.AddCut(AliDielectronVarManager::kEta,          -0.9,    0.9);
  33 dauQAcuts1.AddCut(AliDielectronVarManager::kNclsTPC,      50.0,  160.0);
  34 dauQAcuts1.AddCut(AliDielectronVarManager::kTPCchi2Cl,     0.0,    4.0);
  35 AliDielectronTrackCuts dauQAcuts2;
  36 dauQAcuts2.SetRequireTPCRefit(kTRUE);
  37
  38
  39
  40 Example configuration:
  41
  42   AliDielectronV0Cuts *gammaV0Cuts = new AliDielectronV0Cuts("IsGamma","IsGamma");
  43
  44   // which V0 finder you want to use
  45   gammaV0Cuts->SetV0finder(kOnTheFly);  // kAll(default), kOffline or kOnTheFly
  46
  47   // add some pdg codes (they are used then by the KF package and important for gamma conversions)
  48   gammaV0Cuts->SetPdgCodes(22,11,11); // mother, daughter1 and 2
  49
  50   // add default PID cuts (defined in AliDielectronPID)
  51   gammaV0Cuts->SetDefaultPID(16);
  52
  53   // add the pair cuts for V0 candidates
  54   gammaV0Cuts->AddCut(AliDielectronVarManager::kCosPointingAngle, TMath::Cos(0.02),   1.0, kFALSE);
  55   gammaV0Cuts->AddCut(AliDielectronVarManager::kChi2NDF,                       0.0,  10.0, kFALSE);
  56   gammaV0Cuts->AddCut(AliDielectronVarManager::kLegDist,                       0.0,   0.25, kFALSE);
  57   gammaV0Cuts->AddCut(AliDielectronVarManager::kR,                             3.0,  90.0, kFALSE);
  58   gammaV0Cuts->AddCut(AliDielectronVarManager::kPsiPair,                       0.0,   0.05, kFALSE);
  59   gammaV0Cuts->AddCut(AliDielectronVarManager::kM,                             0.0,   0.05, kFALSE);
  60   gammaV0Cuts->AddCut(AliDielectronVarManager::kArmPt,                         0.0,   0.05, kFALSE);
  61
  62   // selection or rejection of V0 tracks
  63   gammaV0Cuts->SetExcludeTracks(kTRUE);
  64
  65   // add the V0cuts directly to the track filter or to some cut group of it
  66
  67 */
  68 //                                                                       //
  69 ///////////////////////////////////////////////////////////////////////////
  70
  71
  72 #include "AliDielectronV0Cuts.h"
  73 #include "AliDielectronVarManager.h"
  74 #include "AliDielectronTrackCuts.h"
  75 #include "AliDielectronPID.h"
  76 #include "AliESDv0.h"
  77 #include "AliAODv0.h"
  78
  79 ClassImp(AliDielectronV0Cuts)
  80
  81
  82 AliDielectronV0Cuts::AliDielectronV0Cuts() :
  83   AliDielectronVarCuts(),
  84   fV0TrackArr(0),
  85   fExcludeTracks(kTRUE),
  86   fV0finder(kAll),
  87   fMotherPdg(0),
  88   fNegPdg(0),
  89   fPosPdg(0),
  90   fPID(-1),
  91   fOrbit(0),
  92   fPeriod(0),
  93   fBunchCross(0)
  94 {
  95   //
  96   // Default costructor
  97   //
  98 }
  99
 100 //________________________________________________________________________
 101 AliDielectronV0Cuts::AliDielectronV0Cuts(const char* name, const char* title) :
 102   AliDielectronVarCuts(name,title),
 103   fV0TrackArr(0),
 104   fExcludeTracks(kTRUE),
 105   fV0finder(kAll),
 106   fMotherPdg(0),
 107   fNegPdg(0),
 108   fPosPdg(0),
 109   fPID(-1),
 110   fOrbit(0),
 111   fPeriod(0),
 112   fBunchCross(0)
 113 {
 114   //
 115   // Named contructor
 116   //
 117 }
 118
 119 //________________________________________________________________________
 120 AliDielectronV0Cuts::~AliDielectronV0Cuts()
 121 {
 122   //
 123   // Destructor
 124   //
 125
 126 }
 127
 128 //________________________________________________________________________
 129 void AliDielectronV0Cuts::InitEvent(AliVTrack *trk)
 130 {
 131   //
 132   // Init the V0 candidates
 133   //
 134
 135   // take current event from the track
 136   // TODO: this should be simplyfied by AliVTrack::GetEvent() as soon as implemented
 137   const AliVEvent *ev=0;
 138   if(trk->IsA() == AliAODTrack::Class())
 139     ev=static_cast<const AliVEvent*>((static_cast<const AliAODTrack*>(trk))->GetAODEvent());
 140   else if(trk->IsA() == AliESDtrack::Class())
 141     ev=static_cast<const AliVEvent*>((static_cast<const AliESDtrack*>(trk))->GetESDEvent());
 142   else
 143     return;
 144
 145
 146   // IsNewEvent
 147   if(!ev) return;
 148   if(!IsNewEvent(ev)) return;
 149   //  printf(" Build V0 candidates according to the applied cuts \n");
 150
 151   // rest booleans
 152   fV0TrackArr.ResetAllBits();
 153
 154   // basic quality cut, /*at least one*/ both of the V0 daughters has to fullfill
 155   AliDielectronVarCuts dauQAcuts1;
 156   dauQAcuts1.AddCut(AliDielectronVarManager::kPt,            0.05, 100.0);
 157   dauQAcuts1.AddCut(AliDielectronVarManager::kEta,          -0.9,    0.9);
 158   dauQAcuts1.AddCut(AliDielectronVarManager::kNclsTPC,      50.0,  160.0);
 159   dauQAcuts1.AddCut(AliDielectronVarManager::kTPCchi2Cl,     0.0,    4.0);
 160   AliDielectronTrackCuts dauQAcuts2;
 161   //  dauQAcuts2.SetRequireITSRefit(kTRUE);
 162   dauQAcuts2.SetRequireTPCRefit(kTRUE);
 163   AliDielectronPID dauPIDcuts;
 164   if(fPID>=0) dauPIDcuts.SetDefaults(fPID);
 165
 166   Int_t nV0s      = 0;
 167   Int_t nV0stored = 0;
 168   AliDielectronPair candidate;
 169   candidate.SetPdgCode(fMotherPdg);
 170
 171   // ESD or AOD event
 172   if(ev->IsA() == AliESDEvent::Class()) {
 173     const AliESDEvent *esdev = static_cast<const AliESDEvent*>(ev);
 174
 175     //printf("there are %d V0s in the event \n",esdev->GetNumberOfV0s());
 176     // loop over V0s
 177     for (Int_t iv=0; iv<esdev->GetNumberOfV0s(); ++iv){
 178       AliESDv0 *v = esdev->GetV0(iv);
 179       if(!v) continue;
 180
 181       // check the v0 finder
 182       if( v->GetOnFlyStatus() && fV0finder==AliDielectronV0Cuts::kOffline  ) continue;
 183       if(!v->GetOnFlyStatus() && fV0finder==AliDielectronV0Cuts::kOnTheFly ) continue;
 184
 185       // should we make use of AliESDv0Cuts::GetPdgCode() to preselect candiadtes, e.g.:
 186       // if(fMotherPdg!=v->GetPdgCode()) continue;
 187
 188       AliESDtrack *trNeg=esdev->GetTrack(v->GetIndex(0));
 189       AliESDtrack *trPos=esdev->GetTrack(v->GetIndex(1));
 190       if(!trNeg || !trPos){
 191         printf("Error: Couldn't get V0 daughter: %p - %p\n",trNeg,trPos);
 192         continue;
 193       }
 194
 195       // protection against LS v0s
 196       if(trNeg->Charge() == trPos->Charge()) continue;
 197
 198       // PID default cuts
 199       if(fPID>=0) {
 200         if( !dauPIDcuts.IsSelected(trNeg) ) continue;
 201         if( !dauPIDcuts.IsSelected(trPos) ) continue;
 202       }
 203
 204       // basic track cuts
 205       if( !dauQAcuts2.IsSelected(trNeg) ) continue;
 206       if( !dauQAcuts2.IsSelected(trPos) ) continue;
 207       if( !dauQAcuts1.IsSelected(trNeg) ) continue;
 208       if( !dauQAcuts1.IsSelected(trPos) ) continue;
 209
 210       if(fMotherPdg==22) candidate.SetGammaTracks(trNeg, 11, trPos, 11);
 211       else candidate.SetTracks(trNeg, (trNeg->Charge()<0?fNegPdg:fPosPdg), trPos, (trPos->Charge()<0?fNegPdg:fPosPdg));
 212       // eventually take the external trackparam and build the KFparticles by hand (see AliESDv0::GetKFInfo)
 213       // the folowing is not needed, because the daughters where used in the v0 vertex fit (I guess)
 214       //      AliKFVertex v0vtx = *v;
 215       //      candidate.SetProductionVertex(v0vtx);
 216
 217       if(AliDielectronVarCuts::IsSelected(&candidate)) {
 218         nV0s++;
 219         //printf(" gamma found for vtx %p dau1id %d dau2id %d \n",v,trNeg->GetID(),trPos->GetID());
 220         fV0TrackArr.SetBitNumber(trNeg->GetID());
 221         fV0TrackArr.SetBitNumber(trPos->GetID());
 222       }
 223     }
 224
 225   }
 226   else if(ev->IsA() == AliAODEvent::Class()) {
 227     const AliAODEvent *aodEv = static_cast<const AliAODEvent*>(ev);
 228     if(!aodEv->GetV0s()) return; // protection for nano AODs
 229
 230     // loop over vertices
 231     for (Int_t ivertex=0; ivertex<aodEv->GetNumberOfV0s(); ++ivertex){
 232       AliAODv0 *v=aodEv->GetV0(ivertex);
 233       if(!v) continue;
 234
 235       // check the v0 finder
 236       if( v->GetOnFlyStatus() && fV0finder==AliDielectronV0Cuts::kOffline  ) continue;
 237       if(!v->GetOnFlyStatus() && fV0finder==AliDielectronV0Cuts::kOnTheFly ) continue;
 238
 239       AliAODTrack *trNeg=dynamic_cast<AliAODTrack*>(v->GetDaughter(0));
 240       AliAODTrack *trPos=dynamic_cast<AliAODTrack*>(v->GetDaughter(1));
 241       if(!trNeg || !trPos){
 242         printf("Error: Couldn't get V0 daughter: %p - %p\n",trNeg,trPos);
 243         continue;
 244       }
 245       nV0stored++;
 246
 247       // protection against LS v0s
 248       if(trNeg->Charge() == trPos->Charge()) continue;
 249
 250       // PID default cuts
 251       if(fPID>=0) {
 252         if( !dauPIDcuts.IsSelected(trNeg) ) continue;
 253         if( !dauPIDcuts.IsSelected(trPos) ) continue;
 254       }
 255
 256       // basic track cuts
 257       if( !dauQAcuts2.IsSelected(trNeg) ) continue;
 258       if( !dauQAcuts2.IsSelected(trPos) ) continue;
 259       if( !dauQAcuts1.IsSelected(trNeg) ) continue;
 260       if( !dauQAcuts1.IsSelected(trPos) ) continue;
 261
 262       AliKFVertex v0vtx = *(v->GetSecondaryVtx());
 263       if(fMotherPdg==22) candidate.SetGammaTracks(trNeg, 11, trPos, 11);
 264       else candidate.SetTracks(trNeg, (trNeg->Charge()<0?fNegPdg:fPosPdg), trPos, (trPos->Charge()<0?fNegPdg:fPosPdg));
 265       candidate.SetProductionVertex(v0vtx);
 266
 267       if(AliDielectronVarCuts::IsSelected(&candidate)) {
 268         nV0s++;
 269         //printf(" gamma found for vtx %p dau1id %d dau2id %d \n",v,trNeg->GetID(),trPos->GetID());
 270         fV0TrackArr.SetBitNumber(trNeg->GetID());
 271         fV0TrackArr.SetBitNumber(trPos->GetID());
 272       }
 273     }
 274     //printf("there are %d V0s in the event \n",nV0stored);
 275   }
 276   else
 277     return;
 278
 279   //  printf(" Number of V0s candiates found %d/%d \n",nV0s,nV0stored);
 280
 281 }
 282 //________________________________________________________________________
 283 Bool_t AliDielectronV0Cuts::IsSelected(TObject* track)
 284 {
 285   //
 286   // Make cut decision
 287   //
 288
 289   if(!track) return kFALSE;
 290
 291   AliVTrack *vtrack = static_cast<AliVTrack*>(track);
 292   InitEvent(vtrack);
 293   //printf(" track ID %d selected result %d %d \n",vtrack->GetID(),(fV0TrackArr.TestBitNumber(vtrack->GetID())),fExcludeTracks);
 294   return ( (fV0TrackArr.TestBitNumber(vtrack->GetID()))^fExcludeTracks );
 295 }
 296
 297 //________________________________________________________________________
 298 Bool_t AliDielectronV0Cuts::IsNewEvent(const AliVEvent *ev)
 299 {
 300   //
 301   // check weather we process a new event
 302   //
 303
 304   //  printf(" current ev %d %d %d \n",fBunchCross, fOrbit, fPeriod);
 305   //  printf(" new event %p %d %d %d \n",ev, ev->GetBunchCrossNumber(), ev->GetOrbitNumber(), ev->GetPeriodNumber());
 306
 307   if( fBunchCross == ev->GetBunchCrossNumber() ) {
 308     if( fOrbit == ev->GetOrbitNumber() )         {
 309       if( fPeriod == ev->GetPeriodNumber() )     {
 310         return kFALSE;
 311       }
 312     }
 313   }
 314
 315   fBunchCross = ev->GetBunchCrossNumber();
 316   fOrbit      = ev->GetOrbitNumber();
 317   fPeriod     = ev->GetPeriodNumber();
 318   return kTRUE;
 319 }