PWGPP/TRD/info/AliTRDv0Info.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: AliTRDv0Info.cxx 27496 2008-07-22 08:35:45Z cblume $ */
  17
  18 ////////////////////////////////////////////////////////////////////////////
  19 //                                                                        //
  20 //  Reconstruction QA                                                     //
  21 //                                                                        //
  22 //  Gathers all information necessary for reference data selection about  //
  23 //  the track and (in case) its corresponding V0.                         //
  24 //  Carries out the selection of electrons (from gamma conversions),      //
  25 //  pions (from K0s decays) and protons (from Lambda and Anti-Lambda      //
  26 //  decays) by cuts specific for the respective decay and particle        //
  27 //  species.                                                              //
  28 //  (M.Heide, 2009/10/06)                                                 //
  29 //                                                                        //
  30 //  Authors:                                                              //
  31 //   Alex Bercuci <A.Bercuci@gsi.de>                                      //
  32 //   Alex Wilk    <wilka@uni-muenster.de>                                 //
  33 //   Markus Heide <mheide@uni-muenster.de>                                //
  34 //                                                                        //
  35 ////////////////////////////////////////////////////////////////////////////
  36
  37 #include "TMath.h"
  38 #include "TDatabasePDG.h"
  39
  40 #include "AliESDtrack.h"
  41 #include "AliESDv0.h"
  42 #include "AliLog.h"
  43 #include "TVector3.h"
  44 #include "AliKFParticle.h"
  45 #include "AliKFVertex.h"
  46
  47 #include "AliTRDv0Info.h"
  48 #include "AliTRDtrackInfo.h"
  49 #include "AliTRDtrackInfo.h"
  50
  51 ClassImp(AliTRDv0Info)
  52
  53 //_________________________________________________
  54 AliTRDv0Info::AliTRDv0Info()
  55   : TObject()
  56   ,fQuality(0)
  57   ,fDCA(10)
  58   ,fPointingAngle(10)
  59   ,fOpenAngle(10)
  60   ,fPsiPair(99)
  61   ,fMagField(0)
  62   ,fRadius(0)
  63   ,fV0Momentum(0)
  64   ,fNindex(0)
  65   ,fPindex(0)
  66   ,fInputEvent(NULL)
  67   ,fPrimaryVertex(NULL)
  68   ,fTrackP(NULL)
  69   ,fTrackN(NULL)
  70 {
  71   //
  72   // Default constructor
  73   //
  74
  75   memset(fDetPID, 0, 2*kNDaughters*kNDetectors*AliPID::kSPECIES*sizeof(Float_t));
  76   memset(fComPID, 0, 2*kNDaughters*AliPID::kSPECIES*sizeof(Float_t));
  77   memset(fInvMass, 0, kNDecays*sizeof(Double_t));
  78   memset(fArmenteros, 0, kNDecays*sizeof(Bool_t));
  79   memset(fTPCdEdx, 0, kNDaughters*sizeof(Float_t));
  80   memset(fChi2ndf, 0, kNDecays*sizeof(Double_t));
  81   memset(fDownOpenAngle, 0, kNDecays*sizeof(Float_t));
  82   memset(fDownPsiPair, 0, kNDecays*sizeof(Float_t));
  83   /////////////////////////////////////////////////////////////////////////////
  84   //Set Cut values: First specify decay in brackets, then the actual cut value!
  85   /////////////////////////////////////////////////////////////////////////////
  86
  87   //Upper limit for distance of closest approach of two daughter tracks :
  88   fUpDCA[kGamma] = 1000.;
  89   fUpDCA[kK0s] = 0.08;
  90   fUpDCA[kLambda] = 0.2;
  91   fUpDCA[kAntiLambda] = 0.2;
  92
  93   //Upper limit for pointing angle (= angle between between vector from primary to secondary vertex and reconstructed momentum of V0 mother particle) :
  94   fUpPointingAngle[kGamma] = 0.03;
  95   fUpPointingAngle[kK0s] = 0.03;
  96   fUpPointingAngle[kLambda] = 0.04;
  97   fUpPointingAngle[kAntiLambda] = 0.04;
  98
  99   //Upper limit for invariant mass of V0 mother :
 100   fUpInvMass[kGamma][0] = 0.05;// second pair of brackets is for momentum bin: 0: below mother momentm of 2.5 GeV
 101   fUpInvMass[kGamma][1] = 0.07;//1: above 2.5 GeV
 102   fUpInvMass[kK0s][0] = fUpInvMass[kK0s][1] = 0.50265;
 103   fUpInvMass[kLambda][0] = fUpInvMass[kLambda][1] = 1.1207;
 104   fUpInvMass[kAntiLambda][0] = fUpInvMass[kAntiLambda][1] = 1.1207;
 105
 106   //Lower limit for invariant mass of V0 mother :
 107   fDownInvMass[kGamma] = -1.;
 108   fDownInvMass[kK0s] = 0.49265;
 109   fDownInvMass[kLambda] = 1.107;
 110   fDownInvMass[kAntiLambda] = 1.107;
 111
 112   //Upper limit for KF Chi2/NDF value;
 113   fUpChi2ndf[kGamma] = 10000.;//7.;
 114   fUpChi2ndf[kK0s] = 10000.;//5.;
 115   fUpChi2ndf[kLambda] = 10000.;//5.;
 116   fUpChi2ndf[kAntiLambda] = 10000.;//5.;
 117
 118   //Lower limit for distance from secondary vertex to primary vertex in x-y plane :
 119   fDownRadius[kGamma] = 6.;
 120   fDownRadius[kK0s] = 0.;
 121   fDownRadius[kLambda] = 0.;
 122   fDownRadius[kAntiLambda] = 0.;
 123
 124   //Upper limit for distance from secondary vertex to primary vertex in x-y plane :
 125   fUpRadius[kGamma] = 1000.;
 126   fUpRadius[kK0s] = 20.;
 127   fUpRadius[kLambda] = 1000.;
 128   fUpRadius[kAntiLambda] = 1000.;
 129
 130   //Upper limit for opening angle between two daughter tracks (characteristically near zero for conversions) :
 131   fUpOpenAngle[kGamma] = 0.1;
 132   fUpOpenAngle[kK0s] = 3.15;
 133   fUpOpenAngle[kLambda] = 3.15;
 134   fUpOpenAngle[kAntiLambda] = 3.15;
 135
 136   //Upper limit for angle between daughter momentum plane and plane perpendicular to magnetic field (characteristically around zero for conversions) :
 137   fUpPsiPair[kGamma] = 0.05;
 138   fUpPsiPair[kK0s] = 1.6;
 139   fUpPsiPair[kLambda] = 1.6;
 140   fUpPsiPair[kAntiLambda] = 1.6;
 141
 142   //Lower limit for likelihood value of TPC PID :
 143   fDownTPCPIDneg[AliPID::kElectron] = 0.;
 144   fDownTPCPIDpos[AliPID::kElectron] = 0.;
 145
 146   fDownTPCPIDneg[AliPID::kMuon] = 0.;
 147   fDownTPCPIDpos[AliPID::kMuon] = 0.;
 148
 149   fDownTPCPIDneg[AliPID::kPion] = 0.;
 150   fDownTPCPIDpos[AliPID::kPion] = 0.;
 151
 152   fDownTPCPIDneg[AliPID::kKaon] = 0.;
 153   fDownTPCPIDpos[AliPID::kKaon] = 0.;
 154
 155   fDownTPCPIDneg[AliPID::kProton] = 0.;
 156   fDownTPCPIDpos[AliPID::kProton] = 0.;
 157
 158  //Lower limit for likelihood value of combined PID :
 159   fDownComPIDneg[AliPID::kElectron] = 0.;
 160   fDownComPIDpos[AliPID::kElectron] = 0.;
 161
 162   fDownComPIDneg[AliPID::kMuon] = 0.;
 163   fDownComPIDpos[AliPID::kMuon] = 0.;
 164
 165   fDownComPIDneg[AliPID::kPion] = 0.;
 166   fDownComPIDpos[AliPID::kPion] = 0.;
 167
 168   fDownComPIDneg[AliPID::kKaon] = 0.;
 169   fDownComPIDpos[AliPID::kKaon] = 0.;
 170
 171   fDownComPIDneg[AliPID::kProton] = 0.;
 172   fDownComPIDpos[AliPID::kProton] = 0.;
 173
 174  //Lower limit for likelihood value of combined PID for daughter track which doesn't enter reference data (here: pion daughters from Lambda decays:
 175   fDownComPIDnegPart[AliPID::kElectron] = 0.;
 176   fDownComPIDposPart[AliPID::kElectron] = 0.;
 177
 178   fDownComPIDnegPart[AliPID::kMuon] = 0.;
 179   fDownComPIDposPart[AliPID::kMuon] = 0.;
 180
 181   fDownComPIDnegPart[AliPID::kPion] = 0.;
 182   fDownComPIDposPart[AliPID::kPion] = 0.;
 183
 184   fDownComPIDnegPart[AliPID::kKaon] = 0.;
 185   fDownComPIDposPart[AliPID::kKaon] = 0.;
 186
 187   fDownComPIDnegPart[AliPID::kProton] = 0.;
 188   fDownComPIDposPart[AliPID::kProton] = 0.;
 189
 190   //Parameters for data with well-calibrated PID (after usage of tender):
 191   /* //Lower limit for likelihood value of TPC PID :
 192   fDownTPCPIDneg[AliPID::kElectron] = 0.21;
 193   fDownTPCPIDpos[AliPID::kElectron] = 0.21;
 194
 195   fDownTPCPIDneg[AliPID::kMuon] = 0.21;
 196   fDownTPCPIDpos[AliPID::kMuon] = 0.21;
 197
 198   fDownTPCPIDneg[AliPID::kPion] = 0.21;
 199   fDownTPCPIDpos[AliPID::kPion] = 0.21;
 200
 201   fDownTPCPIDneg[AliPID::kKaon] = 0.21;
 202   fDownTPCPIDpos[AliPID::kKaon] = 0.21;
 203
 204   fDownTPCPIDneg[AliPID::kProton] = 0.21;
 205   fDownTPCPIDpos[AliPID::kProton] = 0.21;
 206
 207   //Lower limit for likelihood value of combined PID :
 208   fDownComPIDneg[AliPID::kElectron] = 0.21;
 209   fDownComPIDpos[AliPID::kElectron] = 0.21;
 210
 211   fDownComPIDneg[AliPID::kMuon] = 0.21;
 212   fDownComPIDpos[AliPID::kMuon] = 0.21;
 213
 214   fDownComPIDneg[AliPID::kPion] = 0.9;
 215   fDownComPIDpos[AliPID::kPion] = 0.9;
 216
 217   fDownComPIDneg[AliPID::kKaon] = 0.21;
 218   fDownComPIDpos[AliPID::kKaon] = 0.21;
 219
 220   fDownComPIDneg[AliPID::kProton] = 0.9;
 221   fDownComPIDpos[AliPID::kProton] = 0.9;
 222
 223  //Lower limit for likelihood value of combined PID for daughter track which doesn't enter reference data (here: pion daughters from Lambda decays:
 224   fDownComPIDnegPart[AliPID::kElectron] = 0.05;
 225   fDownComPIDposPart[AliPID::kElectron] = 0.05;
 226
 227   fDownComPIDnegPart[AliPID::kMuon] = 0.05;
 228   fDownComPIDposPart[AliPID::kMuon] = 0.05;
 229
 230   fDownComPIDnegPart[AliPID::kPion] = 0.05;
 231   fDownComPIDposPart[AliPID::kPion] = 0.05;
 232
 233   fDownComPIDnegPart[AliPID::kKaon] = 0.05;
 234   fDownComPIDposPart[AliPID::kKaon] = 0.05;
 235
 236   fDownComPIDnegPart[AliPID::kProton] = 0.05;
 237   fDownComPIDposPart[AliPID::kProton] = 0.05;*/
 238 }
 239
 240 //_________________________________________________
 241 AliTRDv0Info::AliTRDv0Info(const AliTRDv0Info &ref)
 242   : TObject((TObject&)ref)
 243   ,fQuality(ref.fQuality)
 244   ,fDCA(ref.fDCA)
 245   ,fPointingAngle(ref.fPointingAngle)
 246   ,fOpenAngle(ref.fOpenAngle)
 247   ,fPsiPair(ref.fPsiPair)
 248   ,fMagField(ref.fMagField)
 249   ,fRadius(ref.fRadius)
 250   ,fV0Momentum(ref.fV0Momentum)
 251   ,fNindex(ref.fNindex)
 252   ,fPindex(ref.fPindex)
 253   ,fInputEvent(ref.fInputEvent)
 254   ,fPrimaryVertex(ref.fPrimaryVertex)
 255   ,fTrackP(ref.fTrackP)
 256   ,fTrackN(ref.fTrackN)
 257 {
 258   //
 259   // Copy constructor
 260   //
 261
 262   memcpy(fDetPID, ref.fDetPID, 2*kNDaughters*kNDetectors*AliPID::kSPECIES*sizeof(Float_t));
 263   memcpy(fComPID, ref.fComPID, 2*kNDaughters*AliPID::kSPECIES*sizeof(Float_t));
 264   memcpy(fInvMass, ref.fInvMass, kNDecays*sizeof(Double_t));
 265   memcpy(fArmenteros, ref.fArmenteros, kNDecays*sizeof(Bool_t));
 266   memcpy(fChi2ndf, ref.fChi2ndf, kNDecays*sizeof(Double_t));
 267   memcpy(fTPCdEdx, ref.fTPCdEdx, kNDaughters*sizeof(Float_t));
 268
 269   //Upper limit for distance of closest approach of two daughter tracks :
 270   memcpy(fUpDCA, ref.fUpDCA, kNDecays*sizeof(Float_t));
 271   memcpy(fUpPointingAngle, ref.fUpPointingAngle, kNDecays*sizeof(Float_t));
 272   memcpy(fUpOpenAngle, ref.fUpOpenAngle, kNDecays*sizeof(Float_t));
 273   memcpy(fDownOpenAngle, ref.fDownOpenAngle, kNDecays*sizeof(Float_t));
 274   memcpy(fUpPsiPair, ref.fUpPsiPair, kNDecays*sizeof(Float_t));
 275   memcpy(fDownPsiPair, ref.fDownPsiPair, kNDecays*sizeof(Float_t));
 276   memcpy(fUpInvMass, ref.fUpInvMass, kNDecays*kNMomBins*sizeof(Double_t));
 277   memcpy(fDownInvMass, ref.fDownInvMass, kNDecays*sizeof(Double_t));
 278   memcpy(fUpChi2ndf, ref.fUpChi2ndf, kNDecays*sizeof(Double_t));
 279   memcpy(fUpRadius, ref.fUpRadius, kNDecays*sizeof(Float_t));
 280   memcpy(fDownRadius, ref.fDownRadius, kNDecays*sizeof(Float_t));
 281   memcpy(fDownTPCPIDneg, ref.fDownTPCPIDneg, AliPID::kSPECIES*sizeof(Float_t));
 282   memcpy(fDownTPCPIDpos, ref.fDownTPCPIDpos, AliPID::kSPECIES*sizeof(Float_t));
 283   memcpy(fDownComPIDneg, ref.fDownComPIDneg, AliPID::kSPECIES*sizeof(Float_t));
 284   memcpy(fDownComPIDpos, ref.fDownComPIDpos, AliPID::kSPECIES*sizeof(Float_t));
 285   memcpy(fDownComPIDnegPart, ref.fDownComPIDnegPart, AliPID::kSPECIES*sizeof(Float_t));
 286   memcpy(fDownComPIDposPart, ref.fDownComPIDposPart, AliPID::kSPECIES*sizeof(Float_t));
 287 }
 288
 289 //_________________________________________________
 290 void AliTRDv0Info::SetV0Info(const AliESDv0 *esdv0)
 291 {
 292   //Gets values of ESDv0 and daughter track properties
 293   //See header file for description of variables
 294
 295   fQuality = Quality(esdv0);//Attributes an Int_t to the V0 due to quality cuts (= 1 if V0 is accepted, other integers depending on cut which excludes the vertex)
 296
 297   fRadius = Radius(esdv0);//distance from secondary vertex to primary vertex in x-y plane
 298
 299   fDCA = esdv0->GetDcaV0Daughters();//distance of closest approach of two daughter tracks
 300
 301   fPointingAngle = TMath::ACos(esdv0->GetV0CosineOfPointingAngle());// pointing angle (= angle between between vector from primary to secondary vertex and reconstructed momentum of V0 mother particle)
 302
 303   fOpenAngle = OpenAngle(esdv0);//Opening angle between two daughter tracks
 304
 305   fPsiPair = PsiPair(esdv0);//Angle between daughter momentum plane and plane perpendicular to magnetic field
 306
 307   fV0Momentum = V0Momentum(esdv0);//Reconstructed momentum of the mother particle
 308
 309   //4 decay types : conversions, K0s, Lambda, Anti-Lambda
 310     //five particle types: electrons, muons, pions, kaons, protons (muons and kaons not involved)
 311   for(Int_t idecay(0), part1(-1), part2(-1); idecay < kNDecays; idecay++){
 312
 313     fArmenteros[idecay]=Armenteros(esdv0, idecay);//Attribute the Armenteros yes/no decision for every decay type
 314     switch(idecay){
 315     case kLambda: //protons and pions from Lambda
 316       part1 = AliPID::kProton;
 317       part2 = AliPID::kPion;
 318       break;
 319     case kAntiLambda: //antiprotons and pions from Anti-Lambda
 320       part1 = AliPID::kPion;
 321       part2 = AliPID::kProton;
 322       break;
 323     case kK0s: //pions from K0s
 324       part1 = part2 = AliPID::kPion;
 325       break;
 326     case kGamma://electrons from conversions
 327       part1 = part2 = AliPID::kElectron;
 328       break;
 329     }
 330     fInvMass[idecay] = InvMass(part1, part2, esdv0);//Calculate invariant mass for all of our four supposed decays
 331
 332     // Comment out until bug fix is provided
 333     // A.Bercuci 14. July 2010
 334     //fChi2ndf[idecay] = KFChi2ndf(part1, part2,idecay);
 335
 336   }
 337   //Gets all likelihood values from TPC, TOF and ITS PID for the fDetPID[kNDaughters][kNDetectors][AliPID::kSPECIES] array
 338   GetDetectorPID();
 339   //Bayesian combination of likelihoods from TPC and TOF
 340   CombinePID();
 341   //TPC dE/dx values for both tracks
 342   GetTPCdEdx();
 343
 344 }
 345 //_________________________________________________
 346 Float_t  AliTRDv0Info::V0Momentum(const AliESDv0 *esdv0) const
 347 {
 348   //
 349   // Reconstructed momentum of V0 mother particle
 350   //
 351
 352   Double_t mn[3] = {0,0,0};
 353   Double_t mp[3] = {0,0,0};
 354
 355
 356   esdv0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;
 357   esdv0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter;
 358
 359
 360   return TMath::Sqrt((mn[0]+mp[0])*(mn[0]+mp[0]) + (mn[1]+mp[1])*(mn[1]+mp[1])+(mn[2]+mp[2])*(mn[2]+mp[2]));
 361 }
 362
 363 //_________________________________________________
 364 Double_t AliTRDv0Info::InvMass(Int_t part1, Int_t part2, const AliESDv0 *esdv0) const
 365 {
 366   //
 367   // Invariant mass of reconstructed V0 mother
 368   //
 369
 370   const Double_t kpmass[5] = {AliPID::ParticleMass(AliPID::kElectron),AliPID::ParticleMass(AliPID::kMuon),AliPID::ParticleMass(AliPID::kPion),AliPID::ParticleMass(AliPID::kKaon),AliPID::ParticleMass(AliPID::kProton)};
 371   //Masses of electrons, muons, pions, kaons and protons, as implemented in ROOT
 372
 373
 374   Double_t mn[3] = {0,0,0};
 375   Double_t mp[3] = {0,0,0};
 376
 377   esdv0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;
 378   esdv0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter;
 379
 380   Double_t mass1 = kpmass[part1];//sets supposed rest masses for both daughters: positive
 381   Double_t mass2 = kpmass[part2];//negative
 382
 383   //Calculate daughters' energies :
 384   Double_t e1    = TMath::Sqrt(mass1*mass1+
 385             mp[0]*mp[0]+
 386             mp[1]*mp[1]+
 387             mp[2]*mp[2]);
 388   Double_t e2    = TMath::Sqrt(mass2*mass2+
 389             mn[0]*mn[0]+
 390             mn[1]*mn[1]+
 391             mn[2]*mn[2]);
 392
 393   //Sum of daughter momenta :
 394   Double_t momsum =
 395     (mn[0]+mp[0])*(mn[0]+mp[0])+
 396     (mn[1]+mp[1])*(mn[1]+mp[1])+
 397     (mn[2]+mp[2])*(mn[2]+mp[2]);
 398
 399   //invariant mass :
 400   Double_t mInv = TMath::Sqrt((e1+e2)*(e1+e2)-momsum);
 401
 402   return mInv;
 403
 404 }
 405 //_________________________________________________
 406 Float_t AliTRDv0Info::OpenAngle(const AliESDv0 *esdv0)
 407 {
 408   //Opening angle between two daughter tracks
 409   Double_t mn[3] = {0,0,0};
 410   Double_t mp[3] = {0,0,0};
 411
 412
 413   esdv0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;
 414   esdv0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter;
 415
 416
 417   fOpenAngle = TMath::ACos((mp[0]*mn[0] + mp[1]*mn[1] + mp[2]*mn[2])/(TMath::Sqrt(mp[0]*mp[0] + mp[1]*mp[1] + mp[2]*mp[2])*TMath::Sqrt(mn[0]*mn[0] + mn[1]*mn[1] + mn[2]*mn[2])));
 418
 419   return fOpenAngle;
 420 }
 421
 422 //_________________________________________________
 423 Float_t AliTRDv0Info::PsiPair(const AliESDv0 *esdv0)
 424 {
 425   //Angle between daughter momentum plane and plane perpendicular to magnetic field
 426   Double_t x, y, z;
 427   esdv0->GetXYZ(x,y,z);//Reconstructed coordinates of V0; to be replaced by Markus Rammler's method in case of conversions!
 428
 429   Double_t mn[3] = {0,0,0};
 430   Double_t mp[3] = {0,0,0};
 431
 432
 433   esdv0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;
 434   esdv0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter;
 435
 436
 437   Double_t deltat = 1.;
 438   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
 439
 440   Double_t radiussum = TMath::Sqrt(x*x + y*y) + 50;//radius to which tracks shall be propagated
 441
 442   Double_t momPosProp[3];
 443   Double_t momNegProp[3];
 444
 445   AliExternalTrackParam nt(*fTrackN), pt(*fTrackP);
 446
 447   fPsiPair = 4.;
 448
 449   if(nt.PropagateTo(radiussum,fMagField) == 0)//propagate tracks to the outside
 450     fPsiPair =  -5.;
 451   if(pt.PropagateTo(radiussum,fMagField) == 0)
 452     fPsiPair = -5.;
 453   pt.GetPxPyPz(momPosProp);//Get momentum vectors of tracks after propagation
 454   nt.GetPxPyPz(momNegProp);
 455
 456   Double_t pEle =
 457     TMath::Sqrt(momNegProp[0]*momNegProp[0]+momNegProp[1]*momNegProp[1]+momNegProp[2]*momNegProp[2]);//absolute momentum value of negative daughter
 458   Double_t pPos =
 459     TMath::Sqrt(momPosProp[0]*momPosProp[0]+momPosProp[1]*momPosProp[1]+momPosProp[2]*momPosProp[2]);//absolute momentum value of positive daughter
 460
 461   Double_t scalarproduct =
 462     momPosProp[0]*momNegProp[0]+momPosProp[1]*momNegProp[1]+momPosProp[2]*momNegProp[2];//scalar product of propagated positive and negative daughters' momenta
 463
 464   Double_t chipair = TMath::ACos(scalarproduct/(pEle*pPos));//Angle between propagated daughter tracks
 465
 466   fPsiPair =  TMath::Abs(TMath::ASin(deltat/chipair));
 467
 468   return fPsiPair;
 469
 470 }
 471 //_________________________________________________
 472 Double_t AliTRDv0Info::KFChi2ndf(Int_t part1, Int_t part2,Int_t decay){
 473   //Calculates Kalman filter Chi2/NDF
 474   Int_t mothers[4]={22,310,3122,3122};
 475
 476   const Double_t partMass=TDatabasePDG::Instance()->GetParticle(mothers[decay])->Mass();
 477   const Double_t massWidth[4] = {0.001, 0., 0., 0.};
 478
 479   AliKFParticle *kfMother = CreateMotherParticle(fTrackP, fTrackN, part1, part2);
 480
 481   // Lambda
 482   if(!kfMother) {
 483   return kFALSE;
 484   }
 485
 486   // production vertex is set in the 'CreateMotherParticle' function
 487   kfMother->SetMassConstraint(partMass, massWidth[decay]);
 488
 489   Double_t chi2ndf = (kfMother->GetChi2()/kfMother->GetNDF());
 490
 491   if(kfMother)delete kfMother;
 492   return chi2ndf;
 493 }
 494 //________________________________________________________________
 495 AliKFParticle *AliTRDv0Info::CreateMotherParticle(const AliESDtrack *pdaughter, const AliESDtrack *ndaughter, Int_t pspec, Int_t nspec){
 496   //
 497   // Creates a mother particle on the HEAP !! User code is responsible for its deletion
 498   //
 499   AliKFParticle pkfdaughter(*pdaughter, pspec);
 500   AliKFParticle nkfdaughter(*ndaughter, nspec);
 501
 502
 503   // Create the mother particle
 504   AliKFParticle *m = new AliKFParticle(pkfdaughter, nkfdaughter);
 505
 506   AliKFVertex improvedVertex = *fPrimaryVertex;
 507   improvedVertex += *m;
 508   m->SetProductionVertex(improvedVertex);
 509
 510
 511   return m;
 512 }
 513 //_________________________________________________
 514 Int_t AliTRDv0Info::HasTrack(const AliTRDtrackInfo * const track) const
 515 {
 516   //Checks if track is a secondary vertex daughter (due to V0 finder)
 517
 518   if(!track) return 0;
 519   if(!fTrackP->GetID()) return 0;
 520   if(!fTrackN->GetID()) return 0;
 521
 522   Int_t trackID(track->GetTrackId());//index of the track
 523   return HasTrack(trackID);
 524 }
 525
 526 //_________________________________________________
 527 Int_t AliTRDv0Info::HasTrack(Int_t trackID) const
 528 {
 529   //comparing index of track with indices of pos./neg. V0 daughter :
 530   if(fNindex==trackID) return -1;
 531   else if(fPindex==trackID) return 1;
 532   else return 0;
 533 }
 534
 535 //_________________________________________________
 536 void AliTRDv0Info::GetDetectorPID()
 537 {
 538   //PID likelihoods from TPC, TOF, and ITS, for all particle species
 539
 540   fTrackN->GetTPCpid(fDetPID[kNeg][kTPC]);
 541   fTrackP->GetTPCpid(fDetPID[kPos][kTPC]);
 542   fTrackN->GetTOFpid(fDetPID[kNeg][kTOF]);
 543   fTrackP->GetTOFpid(fDetPID[kPos][kTOF]);
 544   fTrackN->GetITSpid(fDetPID[kNeg][kITS]);
 545   fTrackP->GetITSpid(fDetPID[kPos][kITS]);
 546
 547   Long_t statusN = fTrackN->GetStatus();
 548   Long_t statusP = fTrackP->GetStatus();
 549
 550   if(!(statusN & AliESDtrack::kTPCpid)){
 551          for(Int_t iPart = 0; iPart < AliPID::kSPECIES; iPart++){
 552         fDetPID[kNeg][kTPC][iPart] = 0.2;
 553       }
 554   }
 555   if(!(statusN & AliESDtrack::kTOFpid)){
 556     for(Int_t iPart = 0; iPart < AliPID::kSPECIES; iPart++){
 557       fDetPID[kNeg][kTOF][iPart] = 0.2;
 558     }
 559
 560   }
 561   if(!(statusN & AliESDtrack::kITSpid)){
 562     for(Int_t iPart = 0; iPart < AliPID::kSPECIES; iPart++){
 563       fDetPID[kNeg][kITS][iPart] = 0.2;
 564     }
 565   }
 566   if(!(statusP & AliESDtrack::kTPCpid)){
 567     for(Int_t iPart = 0; iPart < AliPID::kSPECIES; iPart++){
 568       fDetPID[kPos][kTPC][iPart] = 0.2;
 569     }
 570   }
 571   if(!(statusP & AliESDtrack::kTOFpid)){
 572     for(Int_t iPart = 0; iPart < AliPID::kSPECIES; iPart++){
 573       fDetPID[kPos][kTOF][iPart] = 0.2;
 574     }
 575
 576   }
 577   if(!(statusP & AliESDtrack::kITSpid)){
 578     for(Int_t iPart = 0; iPart < AliPID::kSPECIES; iPart++){
 579       fDetPID[kPos][kITS][iPart] = 0.2;
 580     }
 581   }
 582
 583 }
 584 //____________________________________________________________________________________
 585 void AliTRDv0Info::CombinePID()
 586 {
 587   //combined bayesian PID from TPC and TOF
 588   Double_t partrat[AliPID::kSPECIES] = {0.208, 0.010, 0.662, 0.019, 0.101};
 589
 590   for(Int_t iSign = 0; iSign < kNDaughters; iSign++)
 591     {
 592     for(Int_t iPart = 0; iPart < AliPID::kSPECIES; iPart++)
 593       {
 594       fComPID[iSign][iPart] = (partrat[iPart]*fDetPID[iSign][kTPC][iPart]*fDetPID[iSign][kTOF][iPart])/((partrat[0]*fDetPID[iSign][kTPC][0]*fDetPID[iSign][kTOF][0])+(partrat[1]*fDetPID[iSign][kTPC][1]*fDetPID[iSign][kTOF][1])+(partrat[2]*fDetPID[iSign][kTPC][2]*fDetPID[iSign][kTOF][2])+(partrat[3]*fDetPID[iSign][kTPC][3]*fDetPID[iSign][kTOF][3])+(partrat[4]*fDetPID[iSign][kTPC][4]*fDetPID[iSign][kTOF][4]));
 595
 596       }
 597     }
 598 }
 599 //_________________________________________________
 600 Bool_t AliTRDv0Info::GetTPCdEdx()
 601 {
 602   //gets the TPC dE/dx for both daughter tracks
 603   if(!fTrackP->GetID()) return 0;
 604   if(!fTrackN->GetID()) return 0;
 605
 606   fTPCdEdx[kNeg] = fTrackN->GetTPCsignal();
 607   fTPCdEdx[kPos] = fTrackP->GetTPCsignal();
 608   return 1;
 609
 610 }
 611 //_________________________________________________
 612 Bool_t AliTRDv0Info::TPCdEdxCuts(Int_t part, const AliTRDtrackInfo * const track)
 613 {
 614   //applies cuts on TPC dE/dx according to particle species; cutting lines are drawn shifted to the Bethe-Bloch paremeterization
 615   if(!fTrackP->GetID()) return 0;
 616   if(!fTrackN->GetID()) return 0;
 617
 618   //Bethe-Bloch lines
 619   Double_t alephParameters[5];
 620
 621   // data
 622   alephParameters[0] = 0.0283086;
 623   alephParameters[1] = 2.63394e+01;
 624   alephParameters[2] = 5.04114e-11;
 625   alephParameters[3] = 2.12543e+00;
 626   alephParameters[4] = 4.88663e+00;
 627
 628
 629   Double_t deposit = 0;
 630   Float_t x = 0;
 631   if(HasTrack(track) == 1){
 632     x = fTrackP->P();
 633     deposit = fTPCdEdx[kPos];
 634   }
 635   else if(HasTrack(track) == -1){
 636     x = fTrackN->P();
 637     deposit = fTPCdEdx[kNeg];
 638   }
 639   else{
 640     printf("No track found");
 641     return 0;
 642   }
 643   if(x < 0.2)return 0;
 644
 645   Float_t upLimits[5]={85,1000,50*AliExternalTrackParam::BetheBlochAleph(x/0.13957, alephParameters[0], alephParameters[1], alephParameters[2], alephParameters[3],  alephParameters[4])+6,1000,50*AliExternalTrackParam::BetheBlochAleph(x/0.93827, alephParameters[0], alephParameters[1], alephParameters[2], alephParameters[3],  alephParameters[4])+10};
 646   Float_t downLimits[5]={62,40,50*AliExternalTrackParam::BetheBlochAleph(x/0.13957, alephParameters[0], alephParameters[1], alephParameters[2], alephParameters[3],  alephParameters[4])-6,40,50*AliExternalTrackParam::BetheBlochAleph(x/0.93827, alephParameters[0], alephParameters[1], alephParameters[2], alephParameters[3],  alephParameters[4])-11};
 647
 648
 649   if(x < 0.7){
 650     downLimits[4]=90;
 651   }
 652   if(x < 1.25){
 653     upLimits[0] = 85;
 654   }
 655   else{
 656     downLimits[0] = 64;
 657   }
 658
 659
 660   if(deposit < downLimits[part])
 661     return 0;
 662   if(deposit > upLimits[part])
 663     return 0;
 664
 665
 666   return 1;
 667
 668 }
 669 //_________________________________________________
 670 Float_t AliTRDv0Info::Radius(const AliESDv0 *esdv0)
 671 {
 672   //distance from secondary vertex to primary vertex in x-y plane
 673   Double_t x, y, z;
 674   esdv0->GetXYZ(x,y,z); //Reconstructed coordinates of V0
 675   fRadius = TMath::Sqrt(x*x + y*y);
 676   return fRadius;
 677
 678 }
 679
 680 //_________________________________________________
 681 Int_t AliTRDv0Info::Quality(const AliESDv0 *const esdv0)
 682 {
 683   //
 684   // Checking track and V0 quality status in order to exclude vertices based on poor information
 685   //
 686
 687   Float_t nClsN;
 688   nClsN = fTrackN->GetTPCNcls();//number of found clusters in TPC for negative track
 689   Float_t nClsFN;
 690   nClsFN = fTrackN->GetTPCNclsF();//number of findable clusters in TPC for negative track
 691   Float_t nClsP;
 692   nClsP = fTrackP->GetTPCNcls();//number of found clusters in TPC for positive track
 693   Float_t nClsFP;
 694   nClsFP = fTrackP->GetTPCNclsF();//number of findable clusters in TPC for positive track
 695
 696   fQuality = 0;
 697
 698
 699   if (!(esdv0->GetOnFlyStatus()))//accept only vertices from online V0 finder
 700     return -1;
 701
 702   Float_t clsRatioN;
 703   Float_t clsRatioP;
 704
 705   if((nClsFN < 80) || (nClsFP < 80)) return -2;//reject all V0s where at least one track has less than 80 TPC clusters
 706
 707  // Chi2 per TPC cluster
 708   Int_t nTPCclustersP = fTrackP->GetTPCclusters(0);
 709   Int_t nTPCclustersN = fTrackN->GetTPCclusters(0);
 710   Float_t chi2perTPCclusterP = fTrackP->GetTPCchi2()/Float_t(nTPCclustersP);
 711   Float_t chi2perTPCclusterN = fTrackN->GetTPCchi2()/Float_t(nTPCclustersN);
 712
 713   if((chi2perTPCclusterN > 3.5)||(chi2perTPCclusterP > 3.5)) return -3;//reject all V0s where at least one track has a chi2 above 3.5
 714
 715   clsRatioN = nClsN/nClsFN; //ratios of found to findable clusters in TPC
 716   clsRatioP = nClsP/nClsFP;
 717
 718   if((clsRatioN < 0.6)||(clsRatioP < 0.6))//exclude tracks with low ratio of found to findable TPC clusters
 719     return -4;
 720
 721   if (!((fTrackP->GetStatus() &
 722   AliESDtrack::kTPCrefit)))//accept only vertices in which both tracks have TPC refit
 723     return -5;
 724   if (!((fTrackN->GetStatus() &
 725   AliESDtrack::kTPCrefit)))
 726     return -6;
 727   if (fTrackP->GetKinkIndex(0)>0  ||
 728       fTrackN->GetKinkIndex(0)>0 )//exclude tracks with kinks
 729     return -7;
 730
 731   if(!(V0SignCheck()))
 732        return -8;
 733   fQuality = 1;
 734   return fQuality;
 735 }
 736
 737 //________________________________________________________________
 738 Bool_t AliTRDv0Info::V0SignCheck(){
 739   //
 740   // Check if v0 daughters really carry opposite charges
 741   //
 742
 743   Int_t qP = fTrackP->Charge();
 744   Int_t qN = fTrackN->Charge();
 745
 746   if((qP*qN) != -1) return kFALSE;
 747
 748   return kTRUE;
 749 }
 750
 751 //___________________________________________________________________
 752 Bool_t AliTRDv0Info::Armenteros(const AliESDv0 *esdv0, Int_t decay){
 753   //
 754   // computes the Armenteros variables for given V0
 755   //
 756   Double_t mn[3] = {0,0,0};
 757   Double_t mp[3] = {0,0,0};
 758   Double_t mm[3] = {0,0,0};
 759
 760   if(V0SignCheck()){
 761     esdv0->GetNPxPyPz(mn[0],mn[1],mn[2]); //reconstructed cartesian momentum components of negative daughter
 762     esdv0->GetPPxPyPz(mp[0],mp[1],mp[2]); //reconstructed cartesian momentum components of positive daughter
 763   }
 764   else{
 765     esdv0->GetPPxPyPz(mn[0],mn[1],mn[2]); //reconstructed cartesian momentum components of negative daughter
 766     esdv0->GetNPxPyPz(mp[0],mp[1],mp[2]); //reconstructed cartesian momentum components of positive daughter
 767   }
 768   esdv0->GetPxPyPz(mm[0],mm[1],mm[2]); //reconstructed cartesian momentum components of mother
 769
 770   TVector3 vecN(mn[0],mn[1],mn[2]);
 771   TVector3 vecP(mp[0],mp[1],mp[2]);
 772   TVector3 vecM(mm[0],mm[1],mm[2]);
 773
 774   Double_t thetaP = acos((vecP * vecM)/(vecP.Mag() * vecM.Mag()));
 775   Double_t thetaN = acos((vecN * vecM)/(vecN.Mag() * vecM.Mag()));
 776
 777   Double_t alfa = ((vecP.Mag())*cos(thetaP)-(vecN.Mag())*cos(thetaN))/
 778     ((vecP.Mag())*cos(thetaP)+(vecN.Mag())*cos(thetaN)) ;
 779   Double_t qt = vecP.Mag()*sin(thetaP);
 780
 781   Float_t ap[2];
 782   ap[0] = alfa;
 783   ap[1] = qt;
 784
 785   Double_t lCutAP[2];//Lambda/Anti-Lambda cuts
 786   if(decay == 0){
 787     // armenteros cuts
 788     const Double_t cutAlpha[2] = {0.35, 0.45};   // [0.35, 0.45]
 789     const Double_t cutQT = 0.015;
 790     if(TMath::Abs(ap[0]) > cutAlpha[0] && TMath::Abs(ap[0]) < cutAlpha[1]) return kFALSE;
 791
 792     if(ap[1] > cutQT) return kFALSE;
 793   }
 794
 795   else if(decay == 1){
 796     const Double_t cutQT = 0.1075;
 797     const Double_t cutAP = 0.22 * TMath::Sqrt( TMath::Abs( (1-ap[0]*ap[0]/(0.92*0.92)) ) );
 798     if(ap[1] < cutQT) return kFALSE;
 799     if(ap[1] > cutAP) return kFALSE;
 800   }
 801   else if(decay == 2){
 802     const Double_t cutQT = 0.03;
 803     const Double_t cutAlpha = 0.7;  // VERY strong - should supress the overlap with K0
 804     lCutAP[0] = 1.0 - (ap[0]-0.7 * ap[0]-0.7)*1.1 - 0.87;
 805     if(TMath::Abs(ap[0]) > cutAlpha) return kFALSE;
 806     if(ap[1] < cutQT) return kFALSE;
 807     if(ap[1] > lCutAP[0]) return kFALSE;
 808
 809   }
 810   else if(decay == 3){
 811     const Double_t cutQT = 0.03;
 812     const Double_t cutAlpha = 0.7;  // VERY strong - should supress the overlap with K0
 813     lCutAP[1] = 1.0 - (ap[0]+0.7 * ap[0]+0.7)*1.1 - 0.87;
 814     if(TMath::Abs(ap[0]) > cutAlpha) return kFALSE;
 815     if(ap[1] < cutQT) return kFALSE;
 816     if(ap[1] > lCutAP[1]) return kFALSE;
 817   }
 818   return kTRUE;
 819 }
 820
 821 //_________________________________________________
 822 Int_t AliTRDv0Info::GetPID(Int_t ipart, AliTRDtrackInfo *track)
 823 {
 824   // Decides if track is accepted for one of the reference data samples
 825
 826   Int_t cutCode = -99;
 827   if(!(track)) {
 828     AliError("No track info");
 829     return -1;
 830   }
 831   if(!HasTrack(track)){
 832     AliDebug(2, "Track not attached to v0.");
 833     return -2;
 834   }
 835
 836   //translate ipart to decay (Anti-Lambda will be treated separately)
 837   Int_t iDecay = -1;
 838   switch(ipart){
 839   case AliPID::kElectron: iDecay = kGamma; break;
 840   case AliPID::kPion: iDecay = kK0s; break;
 841   case AliPID::kProton: iDecay = kLambda; break;
 842   default:
 843     AliDebug(1, Form("Hypothesis \"ipart=%d\" not handled", ipart));
 844     return -3;
 845   }
 846
 847   //... it fulfills our quality criteria
 848   if(!(fQuality == 1)) return -4;
 849   //... distance of closest approach between daughters is reasonably small
 850   if((fDCA > fUpDCA[iDecay])) return -5;
 851   //... pointing angle between momentum of mother particle and vector from prim. to sec. vertex is small
 852   if((fPointingAngle > fUpPointingAngle[iDecay])) return -6;
 853   //... x-y plane distance of decay point to prim. vertex is bigger than a certain minimum value (for conversions)
 854   if((fRadius < fDownRadius[iDecay])) return -7;
 855   //...or smaller than a maximum value (for K0s)
 856   if((fRadius > fUpRadius[iDecay])) return -8;
 857   //... opening angle is close enough to zero (for conversions)
 858   if((fOpenAngle > fUpOpenAngle[iDecay])) return -9;
 859   //... Psi-pair angle is close enough to zero(for conversions)
 860   if((TMath::Abs(fPsiPair) > fUpPsiPair[iDecay])) return -10;
 861
 862
 863
 864   //Mother momentum slots above/below 2.5 GeV
 865   Int_t iPSlot(fV0Momentum > 2.5);
 866   Int_t trackID(track->GetTrackId());
 867
 868   //specific cut criteria :
 869   if(ipart == AliPID::kProton) {
 870     if((fInvMass[kK0s] < fUpInvMass[kK0s][iPSlot]) && (fInvMass[kK0s] > fDownInvMass[kK0s])) return -11;//explicit exclusion of K0s decays
 871     if(fOpenAngle < (0.3 - 0.2*fV0Momentum)) return -9;
 872
 873     //for proton sample: separate treatment of Lamba and Anti-Lambda decays:
 874     //for Anti-Lambda:
 875     //Combined PID likelihoods high enough for pi+ and anti-proton ; invariant mass calculated postulating these two particle species...
 876     //if((fComPID[kNeg][AliPID::kProton] > fDownComPIDneg[AliPID::kProton]) && (fComPID[kPos][AliPID::kPion] > fDownComPIDposPart[AliPID::kPion])) {
 877     //if((fDetPID[kNeg][kTPC][AliPID::kProton] > fDownTPCPIDneg[AliPID::kProton]) && (fDetPID[kPos][kTPC][AliPID::kPion] > fDownTPCPIDpos[AliPID::kPion])){
 878     if((TPCdEdxCuts(ipart, track))){//momentary solution: direct cut on TPC dE/dx
 879       if(fNindex == trackID) {//we're only interested in the anti-proton
 880         if(fArmenteros[kAntiLambda]){//Armenteros condition has to be fulfilled
 881           if(fChi2ndf[kAntiLambda] < fUpChi2ndf[kAntiLambda]){//Kalman filter Chi2/NDF not allowed to be too large
 882             if((fInvMass[kAntiLambda] < fUpInvMass[kAntiLambda][iPSlot]) && (fInvMass[kAntiLambda] > fDownInvMass[kAntiLambda])){
 883               return 1;
 884             } else cutCode = -15;
 885           } else cutCode =-14;
 886         } else cutCode = -13;
 887       }
 888     } else cutCode = -12;
 889     //for Lambda:
 890     //TPC PID likelihoods high enough for pi- and proton ; invariant mass calculated accordingly
 891     //if((fComPID[kNeg][AliPID::kPion] > fDownComPIDnegPart[AliPID::kPion]) && (fComPID[kPos][AliPID::kProton] > fDownComPIDpos[AliPID::kProton])) {
 892     //if((fDetPID[kNeg][kTPC][AliPID::kPion] > fDownTPCPIDneg[AliPID::kPion]) && (fDetPID[kPos][kTPC][AliPID::kProton] > fDownTPCPIDpos[AliPID::kProton])){
 893     if((TPCdEdxCuts(ipart, track))){//momentary solution: direct TPC dE/dx cuts
 894       if(fPindex == trackID) {
 895         if(fArmenteros[kLambda]){
 896           if(fChi2ndf[kLambda] < fUpChi2ndf[kLambda]){
 897             if((fInvMass[kLambda] < fUpInvMass[kLambda][iPSlot]) && (fInvMass[kLambda] > fDownInvMass[kLambda])){
 898               return 1;
 899             } else cutCode = -15;
 900           } else cutCode = -14;
 901         } else cutCode = -13;
 902       }
 903     } else cutCode = -12;
 904     return cutCode;
 905   }
 906
 907   //for K0s decays: equal TPC PID likelihood criteria for both daughters ; invariant mass calculated postulating two pions
 908   if(ipart == AliPID::kPion) {
 909     if(fOpenAngle < (1.0/(fV0Momentum + 0.3) - 0.1)) return -9;
 910     //explicit exclusion of Lambda decays
 911     if((fInvMass[kLambda] < fUpInvMass[kLambda][iPSlot]) && (fInvMass[kLambda] > fDownInvMass[kLambda])) return -11;
 912     //explicit exclusion of Anti-Lambda decays
 913     if((fInvMass[kAntiLambda] < fUpInvMass[kAntiLambda][iPSlot]) && (fInvMass[kAntiLambda] > fDownInvMass[kAntiLambda])) return -11;
 914     //if((fDetPID[kNeg][kTPC][ipart] < fDownTPCPIDneg[ipart]) || (fDetPID[kPos][kTPC][ipart] < fDownTPCPIDpos[ipart])) return -12;
 915     if(!(TPCdEdxCuts(ipart, track))){//momentary solution: direct TPC dE/dx cuts
 916       return -12;
 917     }
 918   }
 919
 920
 921   //for photon conversions: equal combined PID likelihood criteria for both daughters ; invariant mass calculated postulating two electrons
 922   //No Lambda/K0s exclusion is provided, since these contributions hardly ever interfere with gamma invariant mass!
 923   //Float_t momentum(track->GetESDinfo()->GetOuterParam()->P());
 924   if(ipart == AliPID::kElectron) {
 925     //if(momentum > 1.75) {//since combined PID performs a little worse in simulations than TPC standalone for higher momenta, ONLY TPC PID is used here
 926     //if((fDetPID[kNeg][kTPC][ipart] < fDownTPCPIDneg[ipart]) || (fDetPID[kPos][kTPC][ipart] < fDownTPCPIDpos[ipart])) return -12;
 927     //} else {//for low momenta, combined PID from TOF and TPC is used to get rid of proton contamination
 928     //if((fComPID[kNeg][ipart] > fDownComPIDneg[ipart]) && (fComPID[kPos][ipart] > fDownComPIDpos[ipart])) return 1;
 929     //}
 930     if(!(TPCdEdxCuts(ipart, track))){//momentary solution for direct TPC dE/dx cut
 931       return -12;
 932     }
 933   }
 934
 935   //Armenteros-Polanski cut
 936   if(!(fArmenteros[iDecay])) return -13;
 937   //Kalman filter Chi2/NDF cut
 938   if(fChi2ndf[iDecay] > fUpChi2ndf[iDecay]) return -14;
 939   //Invariant mass cut for K0s and photons, assuming two pions/two electrons as daughters:
 940   if((fInvMass[iDecay] > fUpInvMass[iDecay][iPSlot]) || (fInvMass[iDecay] < fDownInvMass[iDecay])) return -15;
 941
 942   return 1;
 943 }
 944
 945
 946 //_________________________________________________
 947 void AliTRDv0Info::Print(Option_t *opt) const
 948 {
 949   //prints text for debugging etc.
 950   printf("V0 legs :: %4d[+] %4d[-]\n", fPindex, fNindex);
 951   printf("  Decay :: Gamma[%c] K0s[%c] Lambda[%c] AntiLambda[%c]\n",
 952     IsDecay(kGamma)?'y':'n', IsDecay(kK0s)?'y':'n', IsDecay(kLambda)?'y':'n', IsDecay(kAntiLambda)?'y':'n');
 953   printf("  Kine  :: DCA[%5.3f] Radius[%5.3f]\n", fDCA, fRadius);
 954   printf("  Angle :: Pointing[%5.3f] Open[%5.3f] Psi[%5.3f]\n", fPointingAngle, fOpenAngle, fPsiPair);
 955   if(strcmp(opt, "a")!=0) return;
 956   printf("  PID   ::\n"
 957          "             ITS   TPC   TOF   COM\n");
 958   for(Int_t idt=0; idt<kNDaughters; idt++){
 959     for(Int_t is(0); is<AliPID::kSPECIES; is++){
 960       printf("      %s%c%s", AliPID::ParticleShortName(is), idt?'-':'+', (is==1||is==2)?"  ":"   ");
 961       for(Int_t id(0); id<kNDetectors; id++){
 962         printf("%5.1f ", 1.e2*fDetPID[idt][id][is]);
 963       }
 964       printf("%5.1f\n", 1.e2*fComPID[idt][is]);
 965     }
 966   }
 967 }
 968
 969 //_________________________________________________
 970 void AliTRDv0Info::SetV0tracks(AliESDtrack *p, AliESDtrack *n)
 971 {
 972   //sets the two daughter trex and their indices
 973   fTrackP = p; fPindex = p->GetID();
 974   fTrackN = n; fNindex = n->GetID();
 975 }
 976 //_________________________________________________
 977
 978 AliESDtrack *AliTRDv0Info::GetV0Daughter(Int_t sign)
 979 {
 980   //Gets positive of negative daughter of decay
 981   if(sign>0)
 982     return fTrackP;
 983   else if(sign < 0)
 984     return fTrackN;
 985
 986   return 0;
 987 }