PWG3/hfe/AliESDv0KineCuts.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, 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 /* $Id$ */
  17
  18 /*
  19  * author: M.Kalisky@gsi.de
  20  * 08/Dec/2010
  21  *
  22  * Description: This class allows with purely kinematical cuts
  23  * to select clean samples of electrons, pions and protons from the
  24  * V0 online finder ESD V0 candidates for PID and dectector resonse
  25  * studies.
  26  */
  27
  28 #include <TVector3.h>
  29 #include <TDatabasePDG.h>
  30
  31 #include "AliESDv0.h"
  32 #include "AliESDtrack.h"
  33 #include "AliESDEvent.h"
  34 #include "AliVEvent.h"
  35 #include "AliLog.h"
  36 #include "AliKFParticle.h"
  37 #include "AliVTrack.h"
  38 #include "AliKFVertex.h"
  39
  40 #include "AliESDv0KineCuts.h"
  41
  42 ClassImp(AliESDv0KineCuts)
  43
  44 //____________________________________________________________________
  45 AliESDv0KineCuts::AliESDv0KineCuts() :
  46   fEvent(0x0)
  47   , fPrimaryVertex(0x0)
  48   , fType(0)
  49   , fMode(0)
  50   , fTPCNcls(1)
  51   , fTPCrefit(kTRUE)
  52   , fTPCchi2perCls(4.0)
  53   , fTPCclsRatio(0.6)
  54   , fNoKinks(kTRUE)
  55   , fGcutChi2NDF(10)
  56   , fGcutInvMass(0.05)
  57   , fK0cutChi2NDF(10)
  58   , fLcutChi2NDF(10)
  59 {
  60   //
  61   // Default constructor
  62   //
  63
  64   // default single track cuts
  65   fTPCNcls = 1;                // minimal number of the TPC clusters
  66   fTPCrefit = kTRUE;           // TPC refit
  67   fTPCchi2perCls = 4.0;        // chi2 per TPC cluster
  68   fTPCclsRatio = 0.6;          // minimal foun/findable TPC cluster ratio
  69   fNoKinks = kTRUE;            // kinks - no [kTRUE] or do not care [kFalse]
  70
  71
  72   // default gamma cuts values
  73   fGcutChi2NDF = 10;           // Chi2NF cut value for the AliKFparticle gamma
  74   fGcutCosPoint[0] = 0;        // cos of the pointing angle [min, max]
  75   fGcutCosPoint[1] = 0.02;     // cos of the pointing angle [min, max]
  76   fGcutDCA[0] = 0.;            // DCA between the daughter tracks [min, max]
  77   fGcutDCA[1] = 0.25;          // DCA between the daughter tracks [min, max]
  78   fGcutVertexR[0] = 3.;        // radius of the conversion point [min, max]
  79   fGcutVertexR[1] = 90.;       // radius of the conversion point [min, max]
  80   fGcutPsiPair[0] = 0.;        // value of the psi pair cut [min, max]
  81   fGcutPsiPair[1] = 0.05;      // value of the psi pair cut [min, max]
  82   fGcutInvMass = 0.05;         // upper value on the gamma invariant mass
  83   // default K0 cuts
  84   fK0cutChi2NDF = 10;          // Chi2NF cut value for the AliKFparticle K0
  85   fK0cutCosPoint[0] = 0.;      // cos of the pointing angle [min, max]
  86   fK0cutCosPoint[1] = 0.02;    // cos of the pointing angle [min, max]
  87   fK0cutDCA[0] = 0.;           // DCA between the daughter tracks [min, max]
  88   fK0cutDCA[1] = 0.2;          // DCA between the daughter tracks [min, max]
  89   fK0cutVertexR[0] = 2.0;      // radius of the decay point [min, max]
  90   fK0cutVertexR[1] = 30.0;     // radius of the decay point [min, max]
  91   fK0cutInvMass[0] = 0.486;    // invariant mass window
  92   fK0cutInvMass[1] = 0.508;    // invariant mass window
  93   // Lambda & anti-Lambda cut values
  94   fLcutChi2NDF = 10;           // Chi2NF cut value for the AliKFparticle K0
  95   fLcutCosPoint[0] = 0.;       // cos of the pointing angle [min, max]
  96   fLcutCosPoint[1] = 0.02;     // cos of the pointing angle [min, max]
  97   fLcutDCA[0] = 0.;            // DCA between the daughter tracks [min, max]
  98   fLcutDCA[1] = 0.2;           // DCA between the daughter tracks [min, max]
  99   fLcutVertexR[0] = 2.0;       // radius of the decay point [min, max]
 100   fLcutVertexR[1] = 40.0;      // radius of the decay point [min, max]
 101   fLcutInvMass[0] = 1.11;      // invariant mass window
 102   fLcutInvMass[1] = 1.12;      // invariant mass window
 103
 104 }
 105 //____________________________________________________________________
 106 AliESDv0KineCuts::~AliESDv0KineCuts(){
 107   //
 108   // Destructor
 109   //
 110
 111
 112 }
 113 //____________________________________________________________________
 114 AliESDv0KineCuts::AliESDv0KineCuts(const AliESDv0KineCuts &ref):
 115   TObject(ref)
 116   , fEvent(0x0)
 117   , fPrimaryVertex(0x0)
 118   , fType(0)
 119   , fMode(0)
 120   , fTPCNcls(1)
 121   , fTPCrefit(kTRUE)
 122   , fTPCchi2perCls(4.0)
 123   , fTPCclsRatio(0.6)
 124   , fNoKinks(kTRUE)
 125   , fGcutChi2NDF(10)
 126   , fGcutInvMass(0.05)
 127   , fK0cutChi2NDF(10)
 128   , fLcutChi2NDF(10)
 129 {
 130   //
 131   // Copy operator
 132   //
 133
 134   ref.Copy(*this);
 135 }
 136 //____________________________________________________________________
 137 AliESDv0KineCuts &AliESDv0KineCuts::operator=(const AliESDv0KineCuts &ref){
 138   //
 139   // assignment operator
 140   //
 141   if(this != &ref)
 142     ref.Copy(*this);
 143   return *this;
 144 }
 145 //____________________________________________________________________
 146 void AliESDv0KineCuts::Copy(TObject &ref) const {
 147   //
 148   // Performs the copying of the object
 149   //
 150
 151   TObject::Copy(ref);
 152
 153   AliESDv0KineCuts &target = dynamic_cast<AliESDv0KineCuts &>(ref);
 154
 155   // default single track cuts
 156   target.fTPCNcls = fTPCNcls;
 157   target.fTPCrefit = fTPCrefit;
 158   target.fTPCchi2perCls = fTPCchi2perCls;
 159   target.fTPCclsRatio = fTPCclsRatio;
 160   target.fNoKinks = fNoKinks;
 161
 162
 163   // default gamma cuts values
 164   target.fGcutChi2NDF = fGcutChi2NDF;
 165   memcpy(target.fGcutCosPoint, fGcutCosPoint, sizeof(Float_t) * 2);
 166   memcpy(target.fGcutDCA, fGcutDCA, sizeof(Float_t) * 2);
 167   memcpy(target.fGcutVertexR, fGcutVertexR, sizeof(Float_t) * 2);
 168   memcpy(target.fGcutPsiPair, fGcutPsiPair, sizeof(Float_t) * 2);
 169   target.fGcutInvMass = fGcutInvMass;
 170   // default K0 cuts
 171   target.fK0cutChi2NDF = fK0cutChi2NDF;
 172   memcpy(target.fK0cutCosPoint, fK0cutCosPoint, sizeof(Float_t) * 2);
 173   memcpy(target.fK0cutDCA, fK0cutDCA, sizeof(Float_t) * 2);
 174   memcpy(target.fK0cutVertexR, fK0cutVertexR, sizeof(Float_t) * 2);
 175   memcpy(target.fK0cutInvMass, fK0cutInvMass, sizeof(Float_t) * 2);
 176   // Lambda & anti-Lambda cut values
 177   target.fLcutChi2NDF = fLcutChi2NDF;
 178   memcpy(target.fLcutCosPoint, fLcutCosPoint, sizeof(Float_t) * 2);
 179   memcpy(target.fLcutDCA, fLcutDCA, sizeof(Float_t) * 2);
 180   memcpy(target.fLcutVertexR, fLcutVertexR, sizeof(Float_t) * 2);
 181   memcpy(target.fLcutInvMass, fLcutInvMass, sizeof(Float_t) * 2);
 182
 183 }
 184 //____________________________________________________________________
 185 Bool_t AliESDv0KineCuts::ProcessV0(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN){
 186   //
 187   // main user function
 188   //
 189
 190   if(!v0) return kFALSE;
 191   if(!fEvent){
 192     AliErrorClass("No valid Event pointer available, provide it first");
 193     return kFALSE;
 194   }
 195
 196   if(!V0CutsCommon(v0)) return kFALSE;
 197
 198   const Int_t id = PreselectV0(v0);
 199
 200   if(!SingleTrackCuts(v0)) return kFALSE;
 201
 202   switch(id){
 203   case kUndef:
 204     return kFALSE;
 205   case kGamma:
 206     return CaseGamma(v0, pdgV0, pdgP, pdgN);
 207   case kK0:
 208     return CaseK0(v0, pdgV0, pdgP, pdgN);
 209   case kLambda:
 210     return CaseLambda(v0, pdgV0, pdgP, pdgN, 0);
 211   case kALambda:
 212     return CaseLambda(v0, pdgV0, pdgP, pdgN, 1);
 213   default:
 214     return kFALSE;
 215   }
 216
 217   return kFALSE;
 218 }
 219 //____________________________________________________________________
 220 Bool_t AliESDv0KineCuts::ProcessV0(AliESDv0* const v0, Int_t &pdgP, Int_t &pdgN){
 221   //
 222   // main user function, simplified if the V0 identity is not necessary
 223   //
 224
 225   if(!v0) return kFALSE;
 226   if(!fEvent){
 227     AliErrorClass("No valid Event pointer available, provide it first");
 228     return kFALSE;
 229   }
 230
 231   Int_t idV0 = -1;
 232   return ProcessV0(v0, idV0, pdgP, pdgN);
 233
 234 }
 235 //____________________________________________________________________
 236 Int_t  AliESDv0KineCuts::PreselectV0(AliESDv0* const v0){
 237   //
 238   // Make a preselection (exclusive) of the V0 cadidates based on
 239   // Armenteros plot
 240   // the armenteros cut values are currently fixed and user is not able to set them via
 241   // set funcions. The reason is that these cuts are optimized and furneter changes should
 242   // not be necessary. To prove otherwise please study in detail before changing the values
 243   //
 244
 245   Float_t ap[2] = {-1., -1.};
 246   Armenteros(v0, ap);
 247   // for clarity
 248   const Float_t alpha = ap[0];
 249   const Float_t qt = ap[1];
 250
 251   // selection cuts
 252   // - the reagions for different candidates must not overlap
 253
 254   // Gamma cuts
 255   const Double_t cutAlphaG = 0.35;
 256   const Double_t cutQTG = 0.05;
 257   const Double_t cutAlphaG2[2] = {0.6, 0.8};
 258   const Double_t cutQTG2 = 0.04;
 259
 260   // K0 cuts
 261   const Float_t cutQTK0[2] = {0.1075, 0.215};
 262   const Float_t cutAPK0[2] = {0.199, 0.8};   // parameters for curved QT cut
 263
 264   // Lambda & A-Lambda cuts
 265   const Float_t cutQTL = 0.03;
 266   const Float_t cutAlphaL[2] = {0.35, 0.7};
 267   const Float_t cutAlphaAL[2] = {-0.7,  -0.35};
 268   const Float_t cutAPL[3] = {0.107, -0.69, 0.5};  // parameters fir curved QT cut
 269
 270
 271   if(kPurity == fMode){
 272   // Check for Gamma candidates
 273     if(qt < cutQTG){
 274       if( (TMath::Abs(alpha) < cutAlphaG) ) return kGamma;
 275     }
 276     // additional region - should help high pT gammas
 277     if(qt < cutQTG2){
 278       if( (TMath::Abs(alpha) > cutAlphaG2[0]) &&  (TMath::Abs(alpha) < cutAlphaG2[1]) ) return kGamma;
 279     }
 280   }
 281   if(kEffGamma == fMode){
 282     if(qt < cutQTG) return kGamma;
 283   }
 284
 285
 286   // Check for K0 candidates
 287   Float_t q = cutAPK0[0] * TMath::Sqrt(TMath::Abs(1 - alpha*alpha/(cutAPK0[1]*cutAPK0[1])));
 288   if( (qt > cutQTK0[0]) && (qt < cutQTK0[1]) && (qt > q) ){
 289     return kK0;
 290   }
 291
 292   // Check for Lambda candidates
 293   q = cutAPL[0] * TMath::Sqrt(TMath::Abs(1 - ( (alpha + cutAPL[1]) * (alpha + cutAPL[1]) ) / (cutAPL[2]*cutAPL[2]) ));
 294   if( (alpha > cutAlphaL[0]) && (alpha < cutAlphaL[1]) && (qt > cutQTL) && (qt < q)  ){
 295     return kLambda;
 296   }
 297
 298   // Check for A-Lambda candidates
 299   q = cutAPL[0] * TMath::Sqrt(TMath::Abs(1 - ( (alpha - cutAPL[1]) * (alpha - cutAPL[1]) ) / (cutAPL[2]*cutAPL[2]) ));
 300   if( (alpha > cutAlphaAL[0]) && (alpha < cutAlphaAL[1]) && (qt > cutQTL) && (qt < q)  ){
 301     return kALambda;
 302   }
 303
 304   return kUndef;
 305 }
 306 //____________________________________________________________________
 307 Bool_t AliESDv0KineCuts::SingleTrackCuts(AliESDv0 * const v0){
 308   //
 309   // apply single track cuts
 310   // correct sign not relevat here
 311   //
 312
 313   if(!v0) return kFALSE;
 314
 315   Int_t pIndex = 0, nIndex = 0;
 316   pIndex = v0->GetPindex();
 317   nIndex = v0->GetNindex();
 318   AliESDtrack* d[2];
 319   d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
 320   d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
 321
 322   for(Int_t i=0; i<2; ++i){
 323     if(!d[i]) return kFALSE;
 324
 325     // status word
 326     ULong_t status = d[i]->GetStatus();
 327
 328     // No. of TPC clusters leave to the users
 329     if(d[i]->GetTPCNcls() < 1) return kFALSE;
 330
 331     // TPC refit
 332     if(!(status & AliESDtrack::kTPCrefit)) return kFALSE;
 333
 334     // Chi2 per TPC cluster
 335     Int_t nTPCclusters = d[i]->GetTPCNcls();
 336     Float_t chi2perTPCcluster = d[i]->GetTPCchi2()/Float_t(nTPCclusters);
 337     if(chi2perTPCcluster > 4) return kFALSE;
 338
 339     // TPC cluster ratio
 340     Float_t cRatioTPC = d[i]->GetTPCNclsF() > 0. ? static_cast<Float_t>(d[i]->GetTPCNcls())/static_cast<Float_t> (d[i]->GetTPCNclsF()) : 1.;
 341     if(cRatioTPC < 0.6) return kFALSE;
 342
 343     // kinks
 344     if(d[i]->GetKinkIndex(0) != 0) return kFALSE;
 345
 346   }
 347
 348   return kTRUE;
 349 }
 350 //____________________________________________________________________
 351 Bool_t  AliESDv0KineCuts::CaseGamma(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN){
 352   //
 353   // process the gamma conversion candidate
 354   //
 355
 356   if(!v0) return kFALSE;
 357
 358   AliVTrack* daughter[2];
 359   Int_t pIndex = 0, nIndex = 0;
 360
 361   Bool_t sign = CheckSigns(v0);
 362   if(sign){
 363     pIndex = v0->GetPindex();
 364     nIndex = v0->GetNindex();
 365   }
 366   else{
 367     pIndex = v0->GetNindex();
 368     nIndex = v0->GetPindex();
 369   }
 370   daughter[0] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(pIndex));
 371   daughter[1] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(nIndex));
 372   if(!daughter[0] || !daughter[1]) return kFALSE;
 373
 374   AliKFParticle *kfMother = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kElectron), TMath::Abs(kElectron));
 375   if(!kfMother) return kFALSE;
 376
 377   AliESDtrack* d[2];
 378   d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
 379   d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
 380
 381   Float_t iMass = v0->GetEffMass(0, 0);
 382
 383   // cos pointing angle
 384   Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
 385   cosPoint = TMath::ACos(cosPoint);
 386
 387   // DCA between daughters
 388   Double_t dca = v0->GetDcaV0Daughters();
 389
 390   // Production vertex
 391   Double_t x, y, z;
 392   v0->GetXYZ(x,y,z);
 393   Double_t r = TMath::Sqrt(x*x + y*y);
 394
 395   Double_t xy[2];
 396   Double_t r2 = -1.;
 397   if ( GetConvPosXY(d[0], d[1], xy) ){
 398     r2 = TMath::Sqrt(xy[0]*xy[0] + xy[1]*xy[1]);
 399   }
 400
 401   // psi pair
 402   Double_t psiPair = PsiPair(v0);
 403
 404   // V0 chi2/ndf
 405   Double_t chi2ndf = kfMother->GetChi2()/kfMother->GetNDF();
 406
 407   if(kfMother) delete kfMother;
 408
 409   // apply the cuts
 410
 411   if(iMass > fGcutInvMass) return kFALSE;
 412
 413   if(chi2ndf > fGcutChi2NDF) return kFALSE;
 414
 415   if(cosPoint < fGcutCosPoint[0] || cosPoint > fGcutCosPoint[1]) return kFALSE;
 416
 417   if(dca < fGcutDCA[0] || dca > fGcutDCA[1]) return kFALSE;
 418
 419   if(r < fGcutVertexR[0] || r > fGcutVertexR[1]) return kFALSE;
 420
 421   if(psiPair < fGcutPsiPair[0] || psiPair > fGcutPsiPair[1]) return kFALSE;
 422
 423   // all cuts passed
 424
 425   pdgV0 = 22;
 426   if(sign){
 427     pdgP = -11;
 428     pdgN = 11;
 429   }
 430   else{
 431     pdgP = 11;
 432     pdgN = -11;
 433   }
 434
 435   return kTRUE;
 436 }
 437 //____________________________________________________________________
 438 Bool_t  AliESDv0KineCuts::CaseK0(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN){
 439   //
 440   // process the K0 candidate
 441   //
 442
 443   if(!v0) return kFALSE;
 444
 445   AliVTrack* daughter[2];
 446   Int_t pIndex = 0, nIndex = 0;
 447   Bool_t sign = CheckSigns(v0);
 448   if(sign){
 449     pIndex = v0->GetPindex();
 450     nIndex = v0->GetNindex();
 451   }
 452   else{
 453     pIndex = v0->GetNindex();
 454     nIndex = v0->GetPindex();
 455   }
 456
 457   daughter[0] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(pIndex));
 458   daughter[1] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(nIndex));
 459   if(!daughter[0] || !daughter[1]) return kFALSE;
 460
 461   AliKFParticle *kfMother = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kPiPlus), TMath::Abs(kPiPlus));
 462   if(!kfMother) return kFALSE;
 463
 464   AliESDtrack* d[2];
 465   d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
 466   d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
 467
 468   Float_t iMass = v0->GetEffMass(2, 2);
 469
 470   // cos pointing angle
 471   Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
 472   cosPoint = TMath::ACos(cosPoint);
 473
 474   // DCA between daughters
 475   Double_t dca = v0->GetDcaV0Daughters();
 476
 477   // Production vertex
 478   Double_t x, y, z;
 479   v0->GetXYZ(x,y,z);
 480
 481   Double_t r = TMath::Sqrt(x*x + y*y);
 482
 483   // V0 chi2/ndf
 484   Double_t chi2ndf = kfMother->GetChi2()/kfMother->GetNDF();
 485
 486   if(kfMother) delete kfMother;
 487
 488   //
 489   // apply the cuts
 490   //
 491   if(iMass < fK0cutInvMass[0] || iMass > fK0cutInvMass[1]) return kFALSE;
 492
 493   if(chi2ndf > fK0cutChi2NDF) return kFALSE;
 494
 495   if(cosPoint < fK0cutCosPoint[0] || cosPoint > fK0cutCosPoint[1]) return kFALSE;
 496
 497   if(dca < fK0cutDCA[0] || dca > fK0cutDCA[1]) return kFALSE;
 498
 499   if(r < fK0cutVertexR[0] || r > fK0cutVertexR[1]) return kFALSE;
 500
 501   // all cuts passed
 502   pdgV0 = 310;
 503   if(sign){
 504     pdgP = 211;
 505     pdgN = -211;
 506   }
 507   else{
 508     pdgP = -211;
 509     pdgN = 211;
 510   }
 511
 512   return kTRUE;
 513 }
 514 //____________________________________________________________________
 515 Bool_t  AliESDv0KineCuts::CaseLambda(AliESDv0* const v0, Int_t &pdgV0, Int_t &pdgP, Int_t &pdgN, Int_t id){
 516   //
 517   // process teh Lambda and Anti-Lambda candidate
 518   //
 519
 520   if(!v0) return kFALSE;
 521
 522     const Double_t cL0mass=TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass();  // PDG lambda mass
 523
 524   AliVTrack* daughter[2];
 525   Int_t pIndex = 0, nIndex = 0;
 526   Float_t mMass[2] = {-1., -1.};
 527   Bool_t sign = CheckSigns(v0);
 528   if(sign){
 529     pIndex = v0->GetPindex();
 530     nIndex = v0->GetNindex();
 531     mMass[0] = v0->GetEffMass(4, 2);
 532     mMass[1] = v0->GetEffMass(2, 4);
 533   }
 534   else{
 535     pIndex = v0->GetNindex();
 536     nIndex = v0->GetPindex();
 537     mMass[0] = v0->GetEffMass(2, 4);
 538     mMass[1] = v0->GetEffMass(4, 2);
 539   }
 540
 541   daughter[0] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(pIndex));
 542   daughter[1] = dynamic_cast<AliVTrack *>(fEvent->GetTrack(nIndex));
 543   if(!daughter[0] || !daughter[1]) return kFALSE;
 544
 545   AliKFParticle *kfMother[2] = {0x0, 0x0};
 546   // Lambda
 547   kfMother[0] = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kProton), TMath::Abs(kPiPlus));
 548   if(!kfMother[0]) return kFALSE;
 549
 550   // Anti-Lambda
 551   kfMother[1] = CreateMotherParticle(daughter[0], daughter[1], TMath::Abs(kPiPlus), TMath::Abs(kProton));
 552   if(!kfMother[1]) return kFALSE;
 553
 554   Float_t dMass[2] = {TMath::Abs(mMass[0] - cL0mass), TMath::Abs(mMass[1] - cL0mass)};
 555
 556   AliESDtrack* d[2];
 557   d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
 558   d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
 559   if(!d[0] || !d[1])    return kFALSE;
 560
 561   Float_t p[2] = {d[0]->GetP(), d[1]->GetP()};
 562
 563   // check the 3 lambda - antilambda variables
 564   Int_t check[2] = {-1, -1};   // 0 : lambda, 1 : antilambda
 565   // 1) momentum of the daughter particles - proton is expected to have higher momentum than pion
 566   check[0] = (p[0] > p[1]) ? 0 : 1;
 567   // 2) mass of the mother particle
 568   check[1] = (dMass[0] < dMass[1]) ? 0 : 1;
 569
 570   // require positive correlation of (1) and (2)
 571   if(check[0] != check[1]){
 572     if(kfMother[0]) delete kfMother[0];
 573     if(kfMother[1]) delete kfMother[1];
 574     return kFALSE;
 575   }
 576
 577   // now that the check[0] == check[1]
 578   const Int_t type = check[0];
 579
 580   // require that the input armenteros preselection agree:
 581   if(type != id) return kFALSE;
 582
 583   Float_t iMass =0.;
 584   if(sign){
 585     iMass = (type == 0) ? v0->GetEffMass(4, 2) : v0->GetEffMass(2, 4);
 586   }
 587   else{
 588     iMass = (type == 0) ? v0->GetEffMass(2, 4) : v0->GetEffMass(4, 2);
 589   }
 590
 591   // cos pointing angle
 592   Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
 593   cosPoint = TMath::ACos(cosPoint);
 594
 595   // DCA between daughters
 596   Double_t dca = v0->GetDcaV0Daughters();
 597
 598   // Production vertex
 599   Double_t x, y, z;
 600   v0->GetXYZ(x,y,z);
 601   Double_t r = TMath::Sqrt(x*x + y*y);
 602
 603   // proton - pion indices
 604   Int_t ix[2] = {0, 1};
 605   if(1 == type){
 606     ix[0] = 1;
 607     ix[1] = 0;
 608   }
 609
 610   // V0 chi2/ndf
 611   Double_t chi2ndf = kfMother[type]->GetChi2()/kfMother[type]->GetNDF();
 612
 613   if(kfMother[0]) delete kfMother[0];
 614   if(kfMother[1]) delete kfMother[1];
 615
 616   //
 617   // apply the cuts
 618   //
 619
 620   if(iMass < fLcutInvMass[0] || iMass > fLcutInvMass[1]) return kFALSE;
 621
 622   if(chi2ndf > fLcutChi2NDF) return kFALSE;
 623
 624   if(cosPoint < fLcutCosPoint[0] || cosPoint > fLcutCosPoint[1]) return kFALSE;
 625
 626   if(dca < fLcutDCA[0] || dca > fLcutDCA[1]) return kFALSE;
 627
 628   if(r < fLcutVertexR[0] || r > fLcutVertexR[1]) return kFALSE;
 629
 630   // all cuts passed
 631
 632   if(0 == type){
 633     pdgV0 = 3122;
 634     if(sign){
 635       pdgP = 2212;
 636       pdgN = -211;
 637     }
 638     else{
 639       pdgP = -211;
 640       pdgN = 2212;
 641     }
 642   }
 643   else{
 644     pdgV0 = -3122;
 645     if(sign){
 646       pdgP = 211;
 647       pdgN = -2212;
 648     }
 649     else{
 650       pdgP = -2212;
 651       pdgN = 211;
 652     }
 653   }
 654
 655   return kTRUE;
 656 }
 657 //____________________________________________________________________
 658 Bool_t AliESDv0KineCuts::V0CutsCommon(const AliESDv0 * const v0){
 659   //
 660   // V0 cuts common to all V0s
 661   //
 662
 663   AliESDtrack* dN, *dP;
 664
 665   dP = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(v0->GetPindex()));
 666   dN = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(v0->GetNindex()));
 667
 668   if(!dN || !dP) return kFALSE;
 669
 670   Int_t qP = dP->Charge();
 671   Int_t qN = dN->Charge();
 672
 673   if((qP*qN) != -1) return kFALSE;
 674
 675   return kTRUE;
 676 }
 677 //____________________________________________________________________
 678 void   AliESDv0KineCuts::Armenteros(AliESDv0* const v0, Float_t val[2]){
 679   //
 680   // computes the Armenteros variables for given V0
 681   // fills the histogram
 682   // returns the values via "val"
 683   //
 684
 685   Double_t mn[3] = {0,0,0};
 686   Double_t mp[3] = {0,0,0};
 687   Double_t mm[3] = {0,0,0};
 688
 689   if(CheckSigns(v0)){
 690     v0->GetNPxPyPz(mn[0],mn[1],mn[2]); //reconstructed cartesian momentum components of negative daughter
 691     v0->GetPPxPyPz(mp[0],mp[1],mp[2]); //reconstructed cartesian momentum components of positive daughter
 692   }
 693   else{
 694     v0->GetPPxPyPz(mn[0],mn[1],mn[2]); //reconstructed cartesian momentum components of negative daughter
 695     v0->GetNPxPyPz(mp[0],mp[1],mp[2]); //reconstructed cartesian momentum components of positive daughter
 696   }
 697   v0->GetPxPyPz(mm[0],mm[1],mm[2]); //reconstructed cartesian momentum components of mother
 698
 699   TVector3 vecN(mn[0],mn[1],mn[2]);
 700   TVector3 vecP(mp[0],mp[1],mp[2]);
 701   TVector3 vecM(mm[0],mm[1],mm[2]);
 702
 703   Double_t thetaP = acos((vecP * vecM)/(vecP.Mag() * vecM.Mag()));
 704   Double_t thetaN = acos((vecN * vecM)/(vecN.Mag() * vecM.Mag()));
 705
 706   Double_t alfa = ((vecP.Mag())*cos(thetaP)-(vecN.Mag())*cos(thetaN))/
 707     ((vecP.Mag())*cos(thetaP)+(vecN.Mag())*cos(thetaN)) ;
 708   Double_t qt = vecP.Mag()*sin(thetaP);
 709
 710   val[0] = alfa;
 711   val[1] = qt;
 712 }
 713 //____________________________________________________________________
 714 Bool_t AliESDv0KineCuts::CheckSigns(AliESDv0* const v0){
 715   //
 716   // check wheter the sign was correctly applied to
 717   // V0 daughter tracks
 718   //
 719
 720   Bool_t correct = kFALSE;
 721
 722   Int_t pIndex = 0, nIndex = 0;
 723   pIndex = v0->GetPindex();
 724   nIndex = v0->GetNindex();
 725
 726   AliESDtrack* d[2];
 727   d[0] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(pIndex));
 728   d[1] = dynamic_cast<AliESDtrack*>(fEvent->GetTrack(nIndex));
 729
 730   Int_t sign[2];
 731   sign[0] = (int)d[0]->GetSign();
 732   sign[1] = (int)d[1]->GetSign();
 733
 734   if(-1 == sign[0] && 1 == sign[1]){
 735     correct = kFALSE;
 736   }
 737   else{
 738     correct = kTRUE;
 739   }
 740
 741   return correct;
 742 }
 743 //________________________________________________________________
 744 Double_t AliESDv0KineCuts::PsiPair(AliESDv0* const v0) {
 745   //
 746   // Angle between daughter momentum plane and plane
 747   //
 748
 749   if(!fEvent) return -1.;
 750
 751   Float_t magField = fEvent->GetMagneticField();
 752
 753   Int_t pIndex = -1;
 754   Int_t nIndex = -1;
 755   if(CheckSigns(v0)){
 756     pIndex = v0->GetPindex();
 757     nIndex = v0->GetNindex();
 758   }
 759   else{
 760     pIndex = v0->GetNindex();
 761     nIndex = v0->GetPindex();
 762   }
 763
 764
 765   AliESDtrack* daughter[2];
 766
 767   daughter[0] = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(pIndex));
 768   daughter[1] = dynamic_cast<AliESDtrack *>(fEvent->GetTrack(nIndex));
 769
 770   Double_t x, y, z;
 771   v0->GetXYZ(x,y,z);//Reconstructed coordinates of V0; to be replaced by Markus Rammler's method in case of conversions!
 772
 773   Double_t mn[3] = {0,0,0};
 774   Double_t mp[3] = {0,0,0};
 775
 776
 777   v0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;
 778   v0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter;
 779
 780
 781   Double_t deltat = 1.;
 782   deltat = TMath::ATan(mp[2]/(TMath::Sqrt(mp[0]*mp[0] + mp[1]*mp[1])+1.e-13)) -  TMath::ATan(mn[2]/(TMath::Sqrt(mn[0]*mn[0] + mn[1]*mn[1])+1.e-13));//difference of angles of the two daughter tracks with z-axis
 783
 784   Double_t radiussum = TMath::Sqrt(x*x + y*y) + 50;//radius to which tracks shall be propagated
 785
 786   Double_t momPosProp[3];
 787   Double_t momNegProp[3];
 788
 789   AliExternalTrackParam pt(*daughter[0]), nt(*daughter[1]);
 790
 791   Double_t psiPair = 4.;
 792
 793   if(nt.PropagateTo(radiussum,magField) == 0)//propagate tracks to the outside
 794     psiPair =  -5.;
 795   if(pt.PropagateTo(radiussum,magField) == 0)
 796     psiPair = -5.;
 797   pt.GetPxPyPz(momPosProp);//Get momentum vectors of tracks after propagation
 798   nt.GetPxPyPz(momNegProp);
 799
 800   Double_t pEle =
 801     TMath::Sqrt(momNegProp[0]*momNegProp[0]+momNegProp[1]*momNegProp[1]+momNegProp[2]*momNegProp[2]);//absolute momentum value of negative daughter
 802   Double_t pPos =
 803     TMath::Sqrt(momPosProp[0]*momPosProp[0]+momPosProp[1]*momPosProp[1]+momPosProp[2]*momPosProp[2]);//absolute momentum value of positive daughter
 804
 805   Double_t scalarproduct =
 806     momPosProp[0]*momNegProp[0]+momPosProp[1]*momNegProp[1]+momPosProp[2]*momNegProp[2];//scalar product of propagated positive and negative daughters' momenta
 807
 808   Double_t chipair = TMath::ACos(scalarproduct/(pEle*pPos));//Angle between propagated daughter tracks
 809
 810   psiPair =  TMath::Abs(TMath::ASin(deltat/chipair));
 811
 812   return psiPair;
 813 }
 814 //___________________________________________________________________
 815 Bool_t AliESDv0KineCuts::GetConvPosXY(AliESDtrack * const ptrack, AliESDtrack * const ntrack, Double_t convpos[2]){
 816   //
 817   // recalculate the gamma conversion XY postition
 818   //
 819
 820   const Double_t b = fEvent->GetMagneticField();
 821
 822   Double_t helixcenterpos[2];
 823   GetHelixCenter(ptrack,b,ptrack->Charge(),helixcenterpos);
 824
 825   Double_t helixcenterneg[2];
 826   GetHelixCenter(ntrack,b,ntrack->Charge(),helixcenterneg);
 827
 828   Double_t  poshelix[6];
 829   ptrack->GetHelixParameters(poshelix,b);
 830   Double_t posradius = TMath::Abs(1./poshelix[4]);
 831
 832   Double_t  neghelix[6];
 833   ntrack->GetHelixParameters(neghelix,b);
 834   Double_t negradius = TMath::Abs(1./neghelix[4]);
 835
 836   Double_t xpos = helixcenterpos[0];
 837   Double_t ypos = helixcenterpos[1];
 838   Double_t xneg = helixcenterneg[0];
 839   Double_t yneg = helixcenterneg[1];
 840
 841   convpos[0] = (xpos*negradius + xneg*posradius)/(negradius+posradius);
 842   convpos[1] = (ypos*negradius+  yneg*posradius)/(negradius+posradius);
 843
 844   return 1;
 845 }
 846 //___________________________________________________________________
 847 Bool_t AliESDv0KineCuts::GetHelixCenter(AliESDtrack * const track, Double_t b,Int_t charge, Double_t center[2]){
 848   //
 849   // computes the center of the track helix
 850   //
 851
 852   Double_t pi = TMath::Pi();
 853
 854   Double_t  helix[6];
 855   track->GetHelixParameters(helix,b);
 856
 857   Double_t xpos =  helix[5];
 858   Double_t ypos =  helix[0];
 859   Double_t radius = TMath::Abs(1./helix[4]);
 860   Double_t phi = helix[2];
 861
 862   if(phi < 0){
 863     phi = phi + 2*pi;
 864   }
 865
 866   phi -= pi/2.;
 867   Double_t xpoint =  radius * TMath::Cos(phi);
 868   Double_t ypoint =  radius * TMath::Sin(phi);
 869
 870   if(b<0){
 871     if(charge > 0){
 872       xpoint = - xpoint;
 873       ypoint = - ypoint;
 874     }
 875
 876     if(charge < 0){
 877       xpoint =  xpoint;
 878       ypoint =  ypoint;
 879     }
 880   }
 881   if(b>0){
 882     if(charge > 0){
 883       xpoint =  xpoint;
 884       ypoint =  ypoint;
 885     }
 886
 887     if(charge < 0){
 888       xpoint = - xpoint;
 889       ypoint = - ypoint;
 890     }
 891   }
 892   center[0] =  xpos + xpoint;
 893   center[1] =  ypos + ypoint;
 894
 895   return 1;
 896 }
 897 //___________________________________________________________________
 898 AliKFParticle *AliESDv0KineCuts::CreateMotherParticle(const AliVTrack* const pdaughter, const AliVTrack* const ndaughter, Int_t pspec, Int_t nspec){
 899   //
 900   // Creates a mother particle
 901   //
 902   AliKFParticle pkfdaughter(*pdaughter, pspec);
 903   AliKFParticle nkfdaughter(*ndaughter, nspec);
 904
 905
 906   // Create the mother particle
 907   AliKFParticle *m = new AliKFParticle(pkfdaughter, nkfdaughter);
 908   // DEBUG - testing
 909   if(TMath::Abs(kElectron) == pspec && TMath::Abs(kElectron) == nspec) m->SetMassConstraint(0, 0.001);
 910   else if(TMath::Abs(kPiPlus) == pspec && TMath::Abs(kPiPlus) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kK0Short)->Mass(), 0.);
 911   else if(TMath::Abs(kProton) == pspec && TMath::Abs(kPiPlus) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(), 0.);
 912   else if(TMath::Abs(kPiPlus) == pspec && TMath::Abs(kProton) == nspec) m->SetMassConstraint(TDatabasePDG::Instance()->GetParticle(kLambda0)->Mass(), 0.);
 913   else{
 914     AliErrorClass("Wrong daughter ID - mass constraint can not be set");
 915   }
 916
 917   AliKFVertex improvedVertex = *fPrimaryVertex;
 918   improvedVertex += *m;
 919   m->SetProductionVertex(improvedVertex);
 920
 921   // update 15/06/2010
 922   // mother particle will not be added to primary vertex but only to its copy
 923   // as this confilcts with calling
 924   // m->SetPrimaryVertex() function and
 925   // subsequently removing the mother particle afterwards
 926   // Source: Sergey Gorbunov
 927
 928   return m;
 929 }
 930 //____________________________________________________________________
 931 void  AliESDv0KineCuts::SetEvent(AliESDEvent* const event){
 932   //
 933   // direct setter of ESD event
 934   //
 935   fEvent = event;
 936   if(!fEvent){
 937     AliErrorClass("Invalid input event pointer");
 938     return;
 939   }
 940
 941 }
 942 //____________________________________________________________________
 943 void  AliESDv0KineCuts::SetEvent(AliVEvent* const event){
 944   //
 945   // direct setter of ESD event
 946   //
 947   if(event)
 948     fEvent = static_cast<AliESDEvent*>(event);
 949   if(!fEvent){
 950     AliErrorClass("Invalid input event pointer");
 951     return;
 952   }
 953
 954 }
 955 //________________________________________________________________
 956 void AliESDv0KineCuts::SetPrimaryVertex(AliKFVertex* const v){
 957   //
 958   // set the primary vertex of the event
 959   //
 960   fPrimaryVertex = v;
 961   if(!fPrimaryVertex){
 962     AliErrorClass("Failed to initialize the primary vertex");
 963     return;
 964   }
 965 }
 966 //___________________________________________________________________
 967 void AliESDv0KineCuts::SetMode(Int_t mode, Int_t type){
 968   //
 969   // this function allows the user to select (prior running the 'ProcessV0' function)
 970   // to select different approaches to V0 selection - the 'mode'
 971   // - and -
 972   // different systems (pp, PbPb) - 'type'
 973   //
 974   // To see the cut values for different modes please refer to the
 975   // function SetCuts()
 976   //
 977   // Important notice: based on the parameters particular sets of cuts will
 978   // be activated for teh V0 selection. If some additional changes to single
 979   // cuts are needed please us the SetXXXcut function (see the header file)
 980   //
 981
 982   switch(mode){
 983   case kPurity:
 984     fMode = kPurity;  // used to obtain highest purity possible - the efficiency may be low
 985   case kEffGamma:
 986     fMode = kEffGamma; // used to obtain highes efficiency possible - the purity may be worse
 987   default:
 988     AliError("V0 selection mode not recognozed, setting 'kPurity'");
 989     fMode = kPurity;
 990   }
 991
 992   switch(type){
 993   case kPP:
 994     fType = kPP;  // cuts optimized for low multiplicity
 995   case kPbPb:
 996     fType = kPbPb;  // cuts optimized for high multiplicity
 997   }
 998
 999   // setup the cut values for selected mode & type
1000   SetCuts();
1001
1002 }
1003 //___________________________________________________________________
1004 void AliESDv0KineCuts::SetMode(Int_t mode, const char* type){
1005   //
1006   // overloaded function - please see above
1007   //
1008
1009   Int_t t = -1;
1010
1011   if(!strcmp("pp", type)) t = kPP;
1012   else if(!(strcmp("PbPb", type))) t = kPbPb;
1013   else{
1014     AliError("data type not recognized, setting 'pp'");
1015     t = kPP;
1016   }
1017
1018   SetMode(mode, t);
1019
1020 }
1021 //___________________________________________________________________
1022 void AliESDv0KineCuts::SetCuts(){
1023   //
1024   // this funciton sets the default cut values based on the selected
1025   // fMode and fType.
1026   // please note that only the cuts that have different values than the default
1027   // cuts are updated here
1028   //
1029
1030   // last update: 14/02/2011
1031   // as a very preliminary  - the only change to default cuts is to apply
1032   // less restricting gamma conversion selection in PreselectV0() function
1033
1034
1035
1036 }