PWG2/FLOW/AliFlowTasks/AliFlowTrackCuts.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 // AliFlowTrackCuts:
  19 // ESD track cuts for flow framework
  20 //
  21 // origin: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
  22 //
  23 // This class gurantees consistency of cut methods, trackparameter
  24 // selection (global tracks, TPC only, etc..) and parameter mixing
  25 // in the flow framework. Transparently handles different input types:
  26 // ESD, MC, AOD.
  27 // This class works in 2 steps: first the requested track parameters are
  28 // constructed (to be set by SetParamType() ), then cuts are applied.
  29 // the constructed track can be requested AFTER checking the cuts by
  30 // calling GetTrack(), in this case the cut object stays in control,
  31 // caller does not have to delete the track.
  32 // Additionally caller can request an AliFlowTrack object to be constructed
  33 // according the parameter mixing scenario requested by SetParamMix().
  34 // AliFlowTrack is made using MakeFlowTrack() method, its an 'object factory'
  35 // so caller needs to take care of the freshly created object.
  36
  37 #include <limits.h>
  38 #include <float.h>
  39 #include <TMatrix.h>
  40 #include "TMCProcess.h"
  41 #include "TParticle.h"
  42 #include "TH2F.h"
  43 #include "AliStack.h"
  44 #include "TBrowser.h"
  45 #include "AliMCEvent.h"
  46 #include "AliESDEvent.h"
  47 #include "AliAODEvent.h"
  48 #include "AliVParticle.h"
  49 #include "AliMCParticle.h"
  50 #include "AliESDtrack.h"
  51 #include "AliMultiplicity.h"
  52 #include "AliAODTrack.h"
  53 #include "AliFlowTrackSimple.h"
  54 #include "AliFlowTrack.h"
  55 #include "AliFlowTrackCuts.h"
  56 #include "AliLog.h"
  57 #include "AliESDpid.h"
  58 #include "AliESDPmdTrack.h"
  59 #include "AliVVZERO.h"
  60
  61 ClassImp(AliFlowTrackCuts)
  62
  63 //-----------------------------------------------------------------------
  64 AliFlowTrackCuts::AliFlowTrackCuts():
  65   AliFlowTrackSimpleCuts(),
  66   fAliESDtrackCuts(NULL),
  67   fQA(NULL),
  68   fCutMC(kFALSE),
  69   fCutMChasTrackReferences(kFALSE),
  70   fCutMCprocessType(kFALSE),
  71   fMCprocessType(kPNoProcess),
  72   fCutMCPID(kFALSE),
  73   fMCPID(0),
  74   fIgnoreSignInMCPID(kFALSE),
  75   fCutMCisPrimary(kFALSE),
  76   fRequireTransportBitForPrimaries(kTRUE),
  77   fMCisPrimary(kFALSE),
  78   fRequireCharge(kFALSE),
  79   fFakesAreOK(kTRUE),
  80   fCutSPDtrackletDeltaPhi(kFALSE),
  81   fSPDtrackletDeltaPhiMax(FLT_MAX),
  82   fSPDtrackletDeltaPhiMin(-FLT_MAX),
  83   fIgnoreTPCzRange(kFALSE),
  84   fIgnoreTPCzRangeMax(FLT_MAX),
  85   fIgnoreTPCzRangeMin(-FLT_MAX),
  86   fCutChi2PerClusterTPC(kFALSE),
  87   fMaxChi2PerClusterTPC(FLT_MAX),
  88   fMinChi2PerClusterTPC(-FLT_MAX),
  89   fCutNClustersTPC(kFALSE),
  90   fNClustersTPCMax(INT_MAX),
  91   fNClustersTPCMin(INT_MIN),
  92   fCutNClustersITS(kFALSE),
  93   fNClustersITSMax(INT_MAX),
  94   fNClustersITSMin(INT_MIN),
  95   fUseAODFilterBit(kFALSE),
  96   fAODFilterBit(0),
  97   fCutDCAToVertexXY(kFALSE),
  98   fCutDCAToVertexZ(kFALSE),
  99   fCutMinimalTPCdedx(kFALSE),
 100   fMinimalTPCdedx(0.),
 101   fCutPmdDet(kFALSE),
 102   fPmdDet(0),
 103   fCutPmdAdc(kFALSE),
 104   fPmdAdc(0.),
 105   fCutPmdNcell(kFALSE),
 106   fPmdNcell(0.),
 107   fParamType(kGlobal),
 108   fParamMix(kPure),
 109   fTrack(NULL),
 110   fTrackPhi(0.),
 111   fTrackEta(0.),
 112   fTrackWeight(0.),
 113   fTrackLabel(INT_MIN),
 114   fMCevent(NULL),
 115   fMCparticle(NULL),
 116   fEvent(NULL),
 117   fTPCtrack(),
 118   fFlowTagType(AliFlowTrackSimple::kInvalid),
 119   fESDpid(),
 120   fPIDsource(kTOFpid),
 121   fTPCpidCuts(NULL),
 122   fTOFpidCuts(NULL),
 123   fParticleID(AliPID::kUnknown),
 124   fParticleProbability(.9),
 125   fAllowTOFmismatchFlag(kFALSE),
 126   fRequireStrictTOFTPCagreement(kFALSE),
 127   fCutRejectElectronsWithTPCpid(kFALSE)
 128 {
 129   //io constructor
 130   for ( Int_t i=0; i<5; i++ ) { fProbBayes[i]=0.0; }
 131   SetPriors(); //init arrays
 132 }
 133
 134 //-----------------------------------------------------------------------
 135 AliFlowTrackCuts::AliFlowTrackCuts(const char* name):
 136   AliFlowTrackSimpleCuts(),
 137   fAliESDtrackCuts(NULL),
 138   fQA(NULL),
 139   fCutMC(kFALSE),
 140   fCutMChasTrackReferences(kFALSE),
 141   fCutMCprocessType(kFALSE),
 142   fMCprocessType(kPNoProcess),
 143   fCutMCPID(kFALSE),
 144   fMCPID(0),
 145   fIgnoreSignInMCPID(kFALSE),
 146   fCutMCisPrimary(kFALSE),
 147   fRequireTransportBitForPrimaries(kTRUE),
 148   fMCisPrimary(kFALSE),
 149   fRequireCharge(kFALSE),
 150   fFakesAreOK(kTRUE),
 151   fCutSPDtrackletDeltaPhi(kFALSE),
 152   fSPDtrackletDeltaPhiMax(FLT_MAX),
 153   fSPDtrackletDeltaPhiMin(-FLT_MAX),
 154   fIgnoreTPCzRange(kFALSE),
 155   fIgnoreTPCzRangeMax(FLT_MAX),
 156   fIgnoreTPCzRangeMin(-FLT_MAX),
 157   fCutChi2PerClusterTPC(kFALSE),
 158   fMaxChi2PerClusterTPC(FLT_MAX),
 159   fMinChi2PerClusterTPC(-FLT_MAX),
 160   fCutNClustersTPC(kFALSE),
 161   fNClustersTPCMax(INT_MAX),
 162   fNClustersTPCMin(INT_MIN),
 163   fCutNClustersITS(kFALSE),
 164   fNClustersITSMax(INT_MAX),
 165   fNClustersITSMin(INT_MIN),
 166   fUseAODFilterBit(kFALSE),
 167   fAODFilterBit(0),
 168   fCutDCAToVertexXY(kFALSE),
 169   fCutDCAToVertexZ(kFALSE),
 170   fCutMinimalTPCdedx(kFALSE),
 171   fMinimalTPCdedx(0.),
 172   fCutPmdDet(kFALSE),
 173   fPmdDet(0),
 174   fCutPmdAdc(kFALSE),
 175   fPmdAdc(0.),
 176   fCutPmdNcell(kFALSE),
 177   fPmdNcell(0.),
 178   fParamType(kGlobal),
 179   fParamMix(kPure),
 180   fTrack(NULL),
 181   fTrackPhi(0.),
 182   fTrackEta(0.),
 183   fTrackWeight(0.),
 184   fTrackLabel(INT_MIN),
 185   fMCevent(NULL),
 186   fMCparticle(NULL),
 187   fEvent(NULL),
 188   fTPCtrack(),
 189   fFlowTagType(AliFlowTrackSimple::kInvalid),
 190   fESDpid(),
 191   fPIDsource(kTOFpid),
 192   fTPCpidCuts(NULL),
 193   fTOFpidCuts(NULL),
 194   fParticleID(AliPID::kUnknown),
 195   fParticleProbability(.9),
 196   fAllowTOFmismatchFlag(kFALSE),
 197   fRequireStrictTOFTPCagreement(kFALSE),
 198   fCutRejectElectronsWithTPCpid(kFALSE)
 199 {
 200   //constructor
 201   SetName(name);
 202   SetTitle("AliFlowTrackCuts");
 203   fESDpid.GetTPCResponse().SetBetheBlochParameters( 0.0283086,
 204                                                     2.63394e+01,
 205                                                     5.04114e-11,
 206                                                     2.12543e+00,
 207                                                     4.88663e+00 );
 208   for ( Int_t i=0; i<5; i++ ) { fProbBayes[i]=0.0; }
 209   SetPriors(); //init arrays
 210
 211 }
 212
 213 //-----------------------------------------------------------------------
 214 AliFlowTrackCuts::AliFlowTrackCuts(const AliFlowTrackCuts& that):
 215   AliFlowTrackSimpleCuts(that),
 216   fAliESDtrackCuts(NULL),
 217   fQA(NULL),
 218   fCutMC(that.fCutMC),
 219   fCutMChasTrackReferences(that.fCutMChasTrackReferences),
 220   fCutMCprocessType(that.fCutMCprocessType),
 221   fMCprocessType(that.fMCprocessType),
 222   fCutMCPID(that.fCutMCPID),
 223   fMCPID(that.fMCPID),
 224   fIgnoreSignInMCPID(that.fIgnoreSignInMCPID),
 225   fCutMCisPrimary(that.fCutMCisPrimary),
 226   fRequireTransportBitForPrimaries(that.fRequireTransportBitForPrimaries),
 227   fMCisPrimary(that.fMCisPrimary),
 228   fRequireCharge(that.fRequireCharge),
 229   fFakesAreOK(that.fFakesAreOK),
 230   fCutSPDtrackletDeltaPhi(that.fCutSPDtrackletDeltaPhi),
 231   fSPDtrackletDeltaPhiMax(that.fSPDtrackletDeltaPhiMax),
 232   fSPDtrackletDeltaPhiMin(that.fSPDtrackletDeltaPhiMin),
 233   fIgnoreTPCzRange(that.fIgnoreTPCzRange),
 234   fIgnoreTPCzRangeMax(that.fIgnoreTPCzRangeMax),
 235   fIgnoreTPCzRangeMin(that.fIgnoreTPCzRangeMin),
 236   fCutChi2PerClusterTPC(that.fCutChi2PerClusterTPC),
 237   fMaxChi2PerClusterTPC(that.fMaxChi2PerClusterTPC),
 238   fMinChi2PerClusterTPC(that.fMinChi2PerClusterTPC),
 239   fCutNClustersTPC(that.fCutNClustersTPC),
 240   fNClustersTPCMax(that.fNClustersTPCMax),
 241   fNClustersTPCMin(that.fNClustersTPCMin),
 242   fCutNClustersITS(that.fCutNClustersITS),
 243   fNClustersITSMax(that.fNClustersITSMax),
 244   fNClustersITSMin(that.fNClustersITSMin),
 245   fUseAODFilterBit(that.fUseAODFilterBit),
 246   fAODFilterBit(that.fAODFilterBit),
 247   fCutDCAToVertexXY(that.fCutDCAToVertexXY),
 248   fCutDCAToVertexZ(that.fCutDCAToVertexZ),
 249   fCutMinimalTPCdedx(that.fCutMinimalTPCdedx),
 250   fMinimalTPCdedx(that.fMinimalTPCdedx),
 251   fCutPmdDet(that.fCutPmdDet),
 252   fPmdDet(that.fPmdDet),
 253   fCutPmdAdc(that.fCutPmdAdc),
 254   fPmdAdc(that.fPmdAdc),
 255   fCutPmdNcell(that.fCutPmdNcell),
 256   fPmdNcell(that.fPmdNcell),
 257   fParamType(that.fParamType),
 258   fParamMix(that.fParamMix),
 259   fTrack(NULL),
 260   fTrackPhi(0.),
 261   fTrackEta(0.),
 262   fTrackWeight(0.),
 263   fTrackLabel(INT_MIN),
 264   fMCevent(NULL),
 265   fMCparticle(NULL),
 266   fEvent(NULL),
 267   fTPCtrack(),
 268   fFlowTagType(that.fFlowTagType),
 269   fESDpid(that.fESDpid),
 270   fPIDsource(that.fPIDsource),
 271   fTPCpidCuts(NULL),
 272   fTOFpidCuts(NULL),
 273   fParticleID(that.fParticleID),
 274   fParticleProbability(that.fParticleProbability),
 275   fAllowTOFmismatchFlag(that.fAllowTOFmismatchFlag),
 276   fRequireStrictTOFTPCagreement(that.fRequireStrictTOFTPCagreement),
 277   fCutRejectElectronsWithTPCpid(that.fCutRejectElectronsWithTPCpid)
 278 {
 279   //copy constructor
 280   if (that.fTPCpidCuts) fTPCpidCuts = new TMatrixF(*(that.fTPCpidCuts));
 281   if (that.fTOFpidCuts) fTOFpidCuts = new TMatrixF(*(that.fTOFpidCuts));
 282   if (that.fAliESDtrackCuts) fAliESDtrackCuts = new AliESDtrackCuts(*(that.fAliESDtrackCuts));
 283   memcpy(fProbBayes,that.fProbBayes,sizeof(fProbBayes));
 284   SetPriors(); //init arrays
 285   if (that.fQA) DefineHistograms();
 286 }
 287
 288 //-----------------------------------------------------------------------
 289 AliFlowTrackCuts& AliFlowTrackCuts::operator=(const AliFlowTrackCuts& that)
 290 {
 291   //assignment
 292   if (this==&that) return *this;
 293
 294   AliFlowTrackSimpleCuts::operator=(that);
 295   //the following may seem excessive but if AliESDtrackCuts properly does copy and clone
 296   //this approach is better memory-fragmentation-wise in some cases
 297   if (that.fAliESDtrackCuts && fAliESDtrackCuts) *fAliESDtrackCuts=*(that.fAliESDtrackCuts);
 298   if (that.fAliESDtrackCuts && !fAliESDtrackCuts) fAliESDtrackCuts=new AliESDtrackCuts(*(that.fAliESDtrackCuts));
 299   if (!that.fAliESDtrackCuts) delete fAliESDtrackCuts; fAliESDtrackCuts=NULL;
 300   //these guys we don't need to copy, just reinit
 301   if (that.fQA) {fQA->Delete(); delete fQA; fQA=NULL; DefineHistograms();}
 302   fCutMC=that.fCutMC;
 303   fCutMChasTrackReferences=that.fCutMChasTrackReferences;
 304   fCutMCprocessType=that.fCutMCprocessType;
 305   fMCprocessType=that.fMCprocessType;
 306   fCutMCPID=that.fCutMCPID;
 307   fMCPID=that.fMCPID;
 308   fIgnoreSignInMCPID=that.fIgnoreSignInMCPID,
 309   fCutMCisPrimary=that.fCutMCisPrimary;
 310   fRequireTransportBitForPrimaries=that.fRequireTransportBitForPrimaries;
 311   fMCisPrimary=that.fMCisPrimary;
 312   fRequireCharge=that.fRequireCharge;
 313   fFakesAreOK=that.fFakesAreOK;
 314   fCutSPDtrackletDeltaPhi=that.fCutSPDtrackletDeltaPhi;
 315   fSPDtrackletDeltaPhiMax=that.fSPDtrackletDeltaPhiMax;
 316   fSPDtrackletDeltaPhiMin=that.fSPDtrackletDeltaPhiMin;
 317   fIgnoreTPCzRange=that.fIgnoreTPCzRange;
 318   fIgnoreTPCzRangeMax=that.fIgnoreTPCzRangeMax;
 319   fIgnoreTPCzRangeMin=that.fIgnoreTPCzRangeMin;
 320   fCutChi2PerClusterTPC=that.fCutChi2PerClusterTPC;
 321   fMaxChi2PerClusterTPC=that.fMaxChi2PerClusterTPC;
 322   fMinChi2PerClusterTPC=that.fMinChi2PerClusterTPC;
 323   fCutNClustersTPC=that.fCutNClustersTPC;
 324   fNClustersTPCMax=that.fNClustersTPCMax;
 325   fNClustersTPCMin=that.fNClustersTPCMin;
 326   fCutNClustersITS=that.fCutNClustersITS;
 327   fNClustersITSMax=that.fNClustersITSMax;
 328   fNClustersITSMin=that.fNClustersITSMin;
 329   fUseAODFilterBit=that.fUseAODFilterBit;
 330   fAODFilterBit=that.fAODFilterBit;
 331   fCutDCAToVertexXY=that.fCutDCAToVertexXY;
 332   fCutDCAToVertexZ=that.fCutDCAToVertexZ;
 333   fCutMinimalTPCdedx=that.fCutMinimalTPCdedx;
 334   fMinimalTPCdedx=that.fMinimalTPCdedx;
 335   fCutPmdDet=that.fCutPmdDet;
 336   fPmdDet=that.fPmdDet;
 337   fCutPmdAdc=that.fCutPmdAdc;
 338   fPmdAdc=that.fPmdAdc;
 339   fCutPmdNcell=that.fCutPmdNcell;
 340   fPmdNcell=that.fPmdNcell;
 341   fFlowTagType=that.fFlowTagType;
 342
 343   fParamType=that.fParamType;
 344   fParamMix=that.fParamMix;
 345
 346   fTrack=NULL;
 347   fTrackEta=0.;
 348   fTrackPhi=0.;
 349   fTrackWeight=0.;
 350   fTrackLabel=INT_MIN;
 351   fMCevent=NULL;
 352   fMCparticle=NULL;
 353   fEvent=NULL;
 354
 355   fESDpid = that.fESDpid;
 356   fPIDsource = that.fPIDsource;
 357
 358   delete fTPCpidCuts;
 359   delete fTOFpidCuts;
 360   if (that.fTPCpidCuts) fTPCpidCuts = new TMatrixF(*(that.fTPCpidCuts));
 361   if (that.fTOFpidCuts) fTOFpidCuts = new TMatrixF(*(that.fTOFpidCuts));
 362
 363   fParticleID=that.fParticleID;
 364   fParticleProbability=that.fParticleProbability;
 365   fAllowTOFmismatchFlag=that.fAllowTOFmismatchFlag;
 366   fRequireStrictTOFTPCagreement=that.fRequireStrictTOFTPCagreement;
 367   fCutRejectElectronsWithTPCpid=that.fCutRejectElectronsWithTPCpid;
 368   memcpy(fProbBayes,that.fProbBayes,sizeof(fProbBayes));
 369
 370   return *this;
 371 }
 372
 373 //-----------------------------------------------------------------------
 374 AliFlowTrackCuts::~AliFlowTrackCuts()
 375 {
 376   //dtor
 377   delete fAliESDtrackCuts;
 378   delete fTPCpidCuts;
 379   delete fTOFpidCuts;
 380   if (fQA) { fQA->SetOwner(); fQA->Delete(); delete fQA; }
 381 }
 382
 383 //-----------------------------------------------------------------------
 384 void AliFlowTrackCuts::SetEvent(AliVEvent* event, AliMCEvent* mcEvent)
 385 {
 386   //set the event
 387   Clear();
 388   fEvent=event;
 389   fMCevent=mcEvent;
 390
 391   //do the magic for ESD
 392   AliESDEvent* myESD = dynamic_cast<AliESDEvent*>(event);
 393   if (fCutPID && myESD)
 394   {
 395     //TODO: maybe call it only for the TOF options?
 396     // Added by F. Noferini for TOF PID
 397     fESDpid.SetTOFResponse(myESD,AliESDpid::kTOF_T0);
 398     fESDpid.MakePID(myESD,kFALSE);
 399     // End F. Noferini added part
 400   }
 401
 402   //TODO: AOD
 403 }
 404
 405 //-----------------------------------------------------------------------
 406 void AliFlowTrackCuts::SetCutMC( Bool_t b )
 407 {
 408   //will we be cutting on MC information?
 409   fCutMC=b;
 410
 411   //if we cut on MC info then also the Bethe Bloch should be the one tuned for MC
 412   if (fCutMC)
 413   {
 414     fESDpid.GetTPCResponse().SetBetheBlochParameters( 2.15898e+00/50.,
 415                                                        1.75295e+01,
 416                                                        3.40030e-09,
 417                                                        1.96178e+00,
 418                                                        3.91720e+00);
 419   }
 420 }
 421
 422 //-----------------------------------------------------------------------
 423 Bool_t AliFlowTrackCuts::IsSelected(TObject* obj, Int_t id)
 424 {
 425   //check cuts
 426   AliVParticle* vparticle = dynamic_cast<AliVParticle*>(obj);
 427   if (vparticle) return PassesCuts(vparticle);
 428   AliFlowTrackSimple* flowtrack = dynamic_cast<AliFlowTrackSimple*>(obj);
 429   if (flowtrack) return PassesCuts(flowtrack);
 430   AliMultiplicity* tracklets = dynamic_cast<AliMultiplicity*>(obj);
 431   if (tracklets) return PassesCuts(tracklets,id);
 432   AliESDPmdTrack* pmdtrack = dynamic_cast<AliESDPmdTrack*>(obj);
 433   if (pmdtrack) return PassesPMDcuts(pmdtrack);
 434   AliVVZERO* vvzero = dynamic_cast<AliVVZERO*>(obj);
 435   if (vvzero) return PassesV0cuts(vvzero,id);
 436   return kFALSE;  //default when passed wrong type of object
 437 }
 438
 439 //-----------------------------------------------------------------------
 440 Bool_t AliFlowTrackCuts::IsSelectedMCtruth(TObject* obj, Int_t id)
 441 {
 442   //check cuts
 443   AliVParticle* vparticle = dynamic_cast<AliVParticle*>(obj);
 444   if (vparticle)
 445   {
 446     return PassesMCcuts(fMCevent,vparticle->GetLabel());
 447   }
 448   AliMultiplicity* tracklets = dynamic_cast<AliMultiplicity*>(obj);
 449   if (tracklets)
 450   {
 451     Int_t label0 = tracklets->GetLabel(id,0);
 452     Int_t label1 = tracklets->GetLabel(id,1);
 453     Int_t label = (label0==label1)?tracklets->GetLabel(id,1):-666;
 454     return PassesMCcuts(fMCevent,label);
 455   }
 456   return kFALSE;  //default when passed wrong type of object
 457 }
 458
 459 //-----------------------------------------------------------------------
 460 Bool_t AliFlowTrackCuts::PassesCuts(const AliFlowTrackSimple* track)
 461 {
 462   //check cuts on a flowtracksimple
 463
 464   //clean up from last iteration
 465   fTrack = NULL;
 466   return AliFlowTrackSimpleCuts::PassesCuts(track);
 467 }
 468
 469 //-----------------------------------------------------------------------
 470 Bool_t AliFlowTrackCuts::PassesCuts(const AliMultiplicity* tracklet, Int_t id)
 471 {
 472   //check cuts on a tracklets
 473
 474   //clean up from last iteration, and init label
 475   fTrack = NULL;
 476   fMCparticle=NULL;
 477   fTrackLabel=-1;
 478
 479   fTrackPhi = tracklet->GetPhi(id);
 480   fTrackEta = tracklet->GetEta(id);
 481   fTrackWeight = 1.0;
 482   if (fCutEta) {if (  fTrackEta < fEtaMin || fTrackEta >= fEtaMax ) return kFALSE;}
 483   if (fCutPhi) {if ( fTrackPhi < fPhiMin || fTrackPhi >= fPhiMax ) return kFALSE;}
 484
 485   //check MC info if available
 486   //if the 2 clusters have different label track cannot be good
 487   //and should therefore not pass the mc cuts
 488   Int_t label0 = tracklet->GetLabel(id,0);
 489   Int_t label1 = tracklet->GetLabel(id,1);
 490   //if possible get label and mcparticle
 491   fTrackLabel = (label0==label1)?tracklet->GetLabel(id,1):-1;
 492   if (!fFakesAreOK && fTrackLabel<0) return kFALSE;
 493   if (fTrackLabel>=0 && fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
 494   //check MC cuts
 495   if (fCutMC && !PassesMCcuts()) return kFALSE;
 496   return kTRUE;
 497 }
 498
 499 //-----------------------------------------------------------------------
 500 Bool_t AliFlowTrackCuts::PassesMCcuts(AliMCEvent* mcEvent, Int_t label)
 501 {
 502   //check the MC info
 503   if (!mcEvent) return kFALSE;
 504   if (label<0) return kFALSE;//otherwise AliCMevent prints a warning before returning NULL
 505   AliMCParticle* mcparticle = static_cast<AliMCParticle*>(mcEvent->GetTrack(label));
 506   if (!mcparticle) {AliError("no MC track"); return kFALSE;}
 507
 508   if (fCutMCisPrimary)
 509   {
 510     if (IsPhysicalPrimary(mcEvent,label,fRequireTransportBitForPrimaries) != fMCisPrimary) return kFALSE;
 511   }
 512   if (fCutMCPID)
 513   {
 514     Int_t pdgCode = mcparticle->PdgCode();
 515     if (fIgnoreSignInMCPID)
 516     {
 517       if (TMath::Abs(fMCPID) != TMath::Abs(pdgCode)) return kFALSE;
 518     }
 519     else
 520     {
 521       if (fMCPID != pdgCode) return kFALSE;
 522     }
 523   }
 524   if ( fCutMCprocessType )
 525   {
 526     TParticle* particle = mcparticle->Particle();
 527     Int_t processID = particle->GetUniqueID();
 528     if (processID != fMCprocessType ) return kFALSE;
 529   }
 530   if (fCutMChasTrackReferences)
 531   {
 532     if (mcparticle->GetNumberOfTrackReferences()<1) return kFALSE;
 533   }
 534   return kTRUE;
 535 }
 536
 537 //-----------------------------------------------------------------------
 538 Bool_t AliFlowTrackCuts::PassesMCcuts()
 539 {
 540   //check MC info
 541   if (!fMCevent) return kFALSE;
 542   if (fTrackLabel<0) return kFALSE;//otherwise AliCMevent prints a warning before returning NULL
 543   fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
 544   return PassesMCcuts(fMCevent,fTrackLabel);
 545 }
 546
 547 //-----------------------------------------------------------------------
 548 Bool_t AliFlowTrackCuts::PassesCuts(AliVParticle* vparticle)
 549 {
 550   //check cuts for an ESD vparticle
 551
 552   ////////////////////////////////////////////////////////////////
 553   //  start by preparing the track parameters to cut on //////////
 554   ////////////////////////////////////////////////////////////////
 555   //clean up from last iteration
 556   fTrack=NULL;
 557
 558   //get the label and the mc particle
 559   fTrackLabel = (fFakesAreOK)?TMath::Abs(vparticle->GetLabel()):vparticle->GetLabel();
 560   if (fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
 561   else fMCparticle=NULL;
 562
 563   Bool_t isMCparticle = kFALSE; //some things are different for MC particles, check!
 564   AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*>(vparticle);
 565   AliAODTrack* aodTrack = NULL;
 566   if (esdTrack)
 567   {
 568     //for an ESD track we do some magic sometimes like constructing TPC only parameters
 569     //or doing some hybrid, handle that here
 570     HandleESDtrack(esdTrack);
 571   }
 572   else
 573   {
 574     HandleVParticle(vparticle);
 575     //now check if produced particle is MC
 576     isMCparticle = (dynamic_cast<AliMCParticle*>(fTrack))!=NULL;
 577     aodTrack = dynamic_cast<AliAODTrack*>(vparticle); //keep the additional dynamic cast out of the way for ESDs
 578   }
 579   ////////////////////////////////////////////////////////////////
 580   ////////////////////////////////////////////////////////////////
 581
 582   if (!fTrack) return kFALSE;
 583   //because it may be different from global, not needed for aodTrack because we dont do anything funky there
 584   if (esdTrack) esdTrack = static_cast<AliESDtrack*>(fTrack);
 585
 586   Bool_t pass=kTRUE;
 587   //check the common cuts for the current particle fTrack (MC,AOD,ESD)
 588   Double_t pt = fTrack->Pt();
 589   Double_t p = fTrack->P();
 590   if (!fFakesAreOK) {if (fTrackLabel<0) pass=kFALSE;}
 591   if (fCutPt) {if (pt < fPtMin || pt >= fPtMax ) pass=kFALSE;}
 592   if (fCutEta) {if (fTrack->Eta() < fEtaMin || fTrack->Eta() >= fEtaMax ) pass=kFALSE;}
 593   if (fCutPhi) {if (fTrack->Phi() < fPhiMin || fTrack->Phi() >= fPhiMax ) pass=kFALSE;}
 594   if (fRequireCharge) {if (fTrack->Charge() == 0) pass=kFALSE;}
 595   if (fCutCharge && !isMCparticle) {if (fTrack->Charge() != fCharge) pass=kFALSE;}
 596   if (fCutCharge && isMCparticle)
 597   {
 598     //in case of an MC particle the charge is stored in units of 1/3|e|
 599     Int_t charge = TMath::Nint(fTrack->Charge()/3.0); //mc particles have charge in units of 1/3e
 600     if (charge!=fCharge) pass=kFALSE;
 601   }
 602   //if(fCutPID) {if (fTrack->PID() != fPID) pass=kFALSE;}
 603
 604   //when additionally MC info is required
 605   if (fCutMC && !PassesMCcuts()) pass=kFALSE;
 606
 607   //the case of ESD or AOD
 608   if (esdTrack) { if (!PassesESDcuts(esdTrack)) { pass=kFALSE; } }
 609   if (aodTrack) { if (!PassesAODcuts(aodTrack)) { pass=kFALSE; } }
 610
 611   if (fQA)
 612   {
 613     if (fMCparticle)
 614     {
 615       TParticle* tparticle=fMCparticle->Particle();
 616       Int_t processID = tparticle->GetUniqueID();
 617       Int_t firstMotherLabel = tparticle->GetFirstMother();
 618       //TLorentzVector v;
 619       //mcparticle->Particle()->ProductionVertex(v);
 620       //Double_t prodvtxX = v.X();
 621       //Double_t prodvtxY = v.Y();
 622
 623       Float_t pdg = 0;
 624       Int_t pdgcode = fMCparticle->PdgCode();
 625       Float_t pdgFirstMother = 0;
 626       Int_t pdgcodeFirstMother = fMCevent->GetTrack(firstMotherLabel)->PdgCode();
 627
 628       switch (TMath::Abs(pdgcode))
 629       {
 630         case 11:
 631           pdg = AliPID::kElectron + 0.5; break;
 632         case 13:
 633           pdg = AliPID::kMuon + 0.5; break;
 634         case 211:
 635           pdg = AliPID::kPion + 0.5; break;
 636         case 321:
 637           pdg = AliPID::kKaon + 0.5; break;
 638         case 2212:
 639           pdg = AliPID::kProton + 0.5; break;
 640         default:
 641           pdg = AliPID::kUnknown + 0.5; break;
 642       }
 643       pdg = TMath::Sign(pdg,static_cast<Float_t>(pdgcode));
 644
 645       switch (TMath::Abs(pdgcodeFirstMother))
 646       {
 647         case 3122: //lambda
 648           pdgFirstMother = 0.5; break;
 649         case 3222: case 3212: case 3112: //sigma+ sigma0 sigma-
 650           pdgFirstMother = 1.5; break;
 651         case 3322: case 3312: //xi0 xi+
 652           pdgFirstMother = 2.5;
 653         case 3332: //omega-
 654           pdgFirstMother = 3.5;
 655         default:
 656           pdgFirstMother = 1e10; break;
 657       }
 658       pdgFirstMother = TMath::Sign(pdgFirstMother,static_cast<Float_t>(pdgcodeFirstMother));
 659
 660       QAbefore(2)->Fill(p,pdg);
 661       QAbefore(3)->Fill(p,IsPhysicalPrimary()?0.5:-0.5);
 662       QAbefore(4)->Fill(p,static_cast<Float_t>(processID));
 663       QAbefore(7)->Fill(p,pdgFirstMother);
 664       if (pass) QAafter(2)->Fill(p,pdg);
 665       if (pass) QAafter(3)->Fill(p,IsPhysicalPrimary()?0.5:-0.5);
 666       if (pass) QAafter(4)->Fill(p,static_cast<Float_t>(processID));
 667       if (pass) QAafter(7)->Fill(p,pdgFirstMother);
 668     }
 669   }
 670
 671   //true by default, if we didn't set any cuts
 672   return pass;
 673 }
 674
 675 //_______________________________________________________________________
 676 Bool_t AliFlowTrackCuts::PassesAODcuts(const AliAODTrack* track)
 677 {
 678   //check cuts for AOD
 679   Bool_t pass = kTRUE;
 680
 681   if (fCutNClustersTPC)
 682   {
 683     Int_t ntpccls = track->GetTPCNcls();
 684     if (ntpccls < fNClustersTPCMin || ntpccls > fNClustersTPCMax) pass=kFALSE;
 685   }
 686
 687   if (fCutNClustersITS)
 688   {
 689     Int_t nitscls = track->GetITSNcls();
 690     if (nitscls < fNClustersITSMin || nitscls > fNClustersITSMax) pass=kFALSE;
 691   }
 692
 693    if (fCutChi2PerClusterTPC)
 694   {
 695     Double_t chi2tpc = track->Chi2perNDF();
 696     if (chi2tpc < fMinChi2PerClusterTPC || chi2tpc > fMaxChi2PerClusterTPC) pass=kFALSE;
 697   }
 698
 699   if (GetRequireTPCRefit() && !(track->GetStatus() & AliESDtrack::kTPCrefit) ) pass=kFALSE;
 700   if (GetRequireITSRefit() && !(track->GetStatus() & AliESDtrack::kITSrefit) ) pass=kFALSE;
 701
 702   if (fUseAODFilterBit && !track->TestFilterBit(fAODFilterBit)) pass=kFALSE;
 703
 704   if (fCutDCAToVertexXY && track->DCA()>GetMaxDCAToVertexXY()) pass=kFALSE;
 705
 706   if (fCutDCAToVertexZ && track->ZAtDCA()>GetMaxDCAToVertexZ()) pass=kFALSE;
 707
 708   if (track->GetTPCsignal() < fMinimalTPCdedx) pass=kFALSE;
 709
 710   return pass;
 711 }
 712
 713 //_______________________________________________________________________
 714 Bool_t AliFlowTrackCuts::PassesESDcuts(AliESDtrack* track)
 715 {
 716   //check cuts on ESD tracks
 717   Bool_t pass=kTRUE;
 718   Float_t dcaxy=0.0;
 719   Float_t dcaz=0.0;
 720   track->GetImpactParameters(dcaxy,dcaz);
 721   const AliExternalTrackParam* pout = track->GetOuterParam();
 722   const AliExternalTrackParam* pin = track->GetInnerParam();
 723   if (fIgnoreTPCzRange)
 724   {
 725     if (pin&&pout)
 726     {
 727       Double_t zin = pin->GetZ();
 728       Double_t zout = pout->GetZ();
 729       if (zin*zout<0) pass=kFALSE;   //reject if cross the membrane
 730       if (zin < fIgnoreTPCzRangeMin || zin > fIgnoreTPCzRangeMax) pass=kFALSE;
 731       if (zout < fIgnoreTPCzRangeMin || zout > fIgnoreTPCzRangeMax) pass=kFALSE;
 732     }
 733   }
 734
 735   Int_t ntpccls = ( fParamType==kTPCstandalone )?
 736                     track->GetTPCNclsIter1():track->GetTPCNcls();
 737   if (fCutChi2PerClusterTPC)
 738   {
 739     Float_t tpcchi2 = (fParamType==kTPCstandalone)?
 740                        track->GetTPCchi2Iter1():track->GetTPCchi2();
 741     tpcchi2 = (ntpccls>0)?tpcchi2/ntpccls:-FLT_MAX;
 742     if (tpcchi2<fMinChi2PerClusterTPC || tpcchi2 >=fMaxChi2PerClusterTPC)
 743       pass=kFALSE;
 744   }
 745
 746   if (fCutMinimalTPCdedx)
 747   {
 748     if (track->GetTPCsignal() < fMinimalTPCdedx) pass=kFALSE;
 749   }
 750
 751   if (fCutNClustersTPC)
 752   {
 753     if (ntpccls < fNClustersTPCMin || ntpccls > fNClustersTPCMax) pass=kFALSE;
 754   }
 755
 756   Int_t nitscls = track->GetNcls(0);
 757   if (fCutNClustersITS)
 758   {
 759     if (nitscls < fNClustersITSMin || nitscls > fNClustersITSMax) pass=kFALSE;
 760   }
 761
 762   //some stuff is still handled by AliESDtrackCuts class - delegate
 763   if (fAliESDtrackCuts)
 764   {
 765     if (!fAliESDtrackCuts->IsSelected(track)) pass=kFALSE;
 766   }
 767
 768   //PID part with pid QA
 769   Double_t beta = GetBeta(track);
 770   Double_t dedx = Getdedx(track);
 771   if (fQA)
 772   {
 773     if (pass) QAbefore(0)->Fill(track->GetP(),beta);
 774     if (pass) QAbefore(1)->Fill(pin->GetP(),dedx);
 775   }
 776   if (fCutPID && (fParticleID!=AliPID::kUnknown)) //if kUnknown don't cut on PID
 777   {
 778     switch (fPIDsource)
 779     {
 780       case kTPCpid:
 781         if (!PassesTPCpidCut(track)) pass=kFALSE;
 782         break;
 783       case kTPCdedx:
 784         if (!PassesTPCdedxCut(track)) pass=kFALSE;
 785         break;
 786       case kTOFpid:
 787         if (!PassesTOFpidCut(track)) pass=kFALSE;
 788         break;
 789       case kTOFbeta:
 790         if (!PassesTOFbetaCut(track)) pass=kFALSE;
 791         break;
 792       case kTOFbetaSimple:
 793         if (!PassesTOFbetaSimpleCut(track)) pass=kFALSE;
 794         break;
 795       case kTPCbayesian:
 796         if (!PassesTPCbayesianCut(track)) pass=kFALSE;
 797         break;
 798         // part added by F. Noferini
 799       case kTOFbayesian:
 800         if (!PassesTOFbayesianCut(track)) pass=kFALSE;
 801         break;
 802         // end part added by F. Noferini
 803
 804         //part added by Natasha
 805       case kTPCNuclei:
 806         if (!PassesNucleiSelection(track)) pass=kFALSE;
 807         break;
 808         //end part added by Natasha
 809
 810       default:
 811         printf("AliFlowTrackCuts::PassesCuts() this should never be called!\n");
 812         pass=kFALSE;
 813         break;
 814     }
 815   }
 816   if (fCutRejectElectronsWithTPCpid)
 817   {
 818     //reject electrons using standard TPC pid
 819     Double_t pidTPC[AliPID::kSPECIES];
 820     track->GetTPCpid(pidTPC);
 821     if (pidTPC[AliPID::kElectron]<fParticleProbability) pass=kFALSE;
 822   }
 823   if (fQA)
 824   {
 825     if (pass) QAafter(0)->Fill(track->GetP(),beta);
 826     if (pass) QAafter(1)->Fill(pin->GetP(),dedx);
 827   }
 828   //end pid part with qa
 829   if (fQA)
 830   {
 831     Double_t pt = track->Pt();
 832     QAbefore(5)->Fill(pt,dcaxy);
 833     QAbefore(6)->Fill(pt,dcaz);
 834     if (pass) QAafter(5)->Fill(pt,dcaxy);
 835     if (pass) QAafter(6)->Fill(pt,dcaz);
 836   }
 837
 838   return pass;
 839 }
 840
 841 //-----------------------------------------------------------------------
 842 void AliFlowTrackCuts::HandleVParticle(AliVParticle* track)
 843 {
 844   //handle the general case
 845   switch (fParamType)
 846   {
 847     default:
 848       fTrack = track;
 849       break;
 850   }
 851 }
 852
 853 //-----------------------------------------------------------------------
 854 void AliFlowTrackCuts::HandleESDtrack(AliESDtrack* track)
 855 {
 856   //handle esd track
 857   switch (fParamType)
 858   {
 859     case kGlobal:
 860       fTrack = track;
 861       break;
 862     case kTPCstandalone:
 863       if (!track->FillTPCOnlyTrack(fTPCtrack))
 864       {
 865         fTrack=NULL;
 866         fMCparticle=NULL;
 867         fTrackLabel=-1;
 868         return;
 869       }
 870       fTrack = &fTPCtrack;
 871       //recalculate the label and mc particle, they may differ as TPClabel != global label
 872       fTrackLabel = (fFakesAreOK)?TMath::Abs(fTrack->GetLabel()):fTrack->GetLabel();
 873       if (fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
 874       else fMCparticle=NULL;
 875       break;
 876     default:
 877       fTrack = track;
 878       break;
 879   }
 880 }
 881
 882 //-----------------------------------------------------------------------
 883 Int_t AliFlowTrackCuts::Count(AliVEvent* event)
 884 {
 885   //calculate the number of track in given event.
 886   //if argument is NULL(default) take the event attached
 887   //by SetEvent()
 888   Int_t multiplicity = 0;
 889   if (!event)
 890   {
 891     for (Int_t i=0; i<GetNumberOfInputObjects(); i++)
 892     {
 893       if (IsSelected(GetInputObject(i))) multiplicity++;
 894     }
 895   }
 896   else
 897   {
 898     for (Int_t i=0; i<event->GetNumberOfTracks(); i++)
 899     {
 900       if (IsSelected(event->GetTrack(i))) multiplicity++;
 901     }
 902   }
 903   return multiplicity;
 904 }
 905
 906 //-----------------------------------------------------------------------
 907 AliFlowTrackCuts* AliFlowTrackCuts::GetStandardVZEROOnlyTrackCuts()
 908 {
 909   AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard vzero flow cuts");
 910   cuts->SetParamType(kV0);
 911   cuts->SetEtaRange( -10, +10 );
 912   cuts->SetPhiMin( 0 );
 913   cuts->SetPhiMax( TMath::TwoPi() );
 914   return cuts;
 915 }
 916
 917 //-----------------------------------------------------------------------
 918 AliFlowTrackCuts* AliFlowTrackCuts::GetStandardGlobalTrackCuts2010()
 919 {
 920   //get standard cuts
 921   AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard Global tracks");
 922   cuts->SetParamType(kGlobal);
 923   cuts->SetPtRange(0.2,5.);
 924   cuts->SetEtaRange(-0.8,0.8);
 925   cuts->SetMinNClustersTPC(70);
 926   cuts->SetMinChi2PerClusterTPC(0.1);
 927   cuts->SetMaxChi2PerClusterTPC(4.0);
 928   cuts->SetMinNClustersITS(2);
 929   cuts->SetRequireITSRefit(kTRUE);
 930   cuts->SetRequireTPCRefit(kTRUE);
 931   cuts->SetMaxDCAToVertexXY(0.3);
 932   cuts->SetMaxDCAToVertexZ(0.3);
 933   cuts->SetAcceptKinkDaughters(kFALSE);
 934   cuts->SetMinimalTPCdedx(10.);
 935
 936   return cuts;
 937 }
 938
 939 //-----------------------------------------------------------------------
 940 AliFlowTrackCuts* AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts2010()
 941 {
 942   //get standard cuts
 943   AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard TPC standalone 2010");
 944   cuts->SetParamType(kTPCstandalone);
 945   cuts->SetPtRange(0.2,5.);
 946   cuts->SetEtaRange(-0.8,0.8);
 947   cuts->SetMinNClustersTPC(70);
 948   cuts->SetMinChi2PerClusterTPC(0.2);
 949   cuts->SetMaxChi2PerClusterTPC(4.0);
 950   cuts->SetMaxDCAToVertexXY(3.0);
 951   cuts->SetMaxDCAToVertexZ(3.0);
 952   cuts->SetDCAToVertex2D(kTRUE);
 953   cuts->SetAcceptKinkDaughters(kFALSE);
 954   cuts->SetMinimalTPCdedx(10.);
 955
 956   return cuts;
 957 }
 958
 959 //-----------------------------------------------------------------------
 960 AliFlowTrackCuts* AliFlowTrackCuts::GetStandardTPCStandaloneTrackCuts()
 961 {
 962   //get standard cuts
 963   AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard TPC standalone");
 964   cuts->SetParamType(kTPCstandalone);
 965   cuts->SetPtRange(0.2,5.);
 966   cuts->SetEtaRange(-0.8,0.8);
 967   cuts->SetMinNClustersTPC(70);
 968   cuts->SetMinChi2PerClusterTPC(0.2);
 969   cuts->SetMaxChi2PerClusterTPC(4.0);
 970   cuts->SetMaxDCAToVertexXY(3.0);
 971   cuts->SetMaxDCAToVertexZ(3.0);
 972   cuts->SetDCAToVertex2D(kTRUE);
 973   cuts->SetAcceptKinkDaughters(kFALSE);
 974   cuts->SetMinimalTPCdedx(10.);
 975
 976   return cuts;
 977 }
 978
 979 //-----------------------------------------------------------------------
 980 AliFlowTrackCuts* AliFlowTrackCuts::GetStandardITSTPCTrackCuts2009(Bool_t selPrimaries)
 981 {
 982   //get standard cuts
 983   AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard global track cuts 2009");
 984   delete cuts->fAliESDtrackCuts;
 985   cuts->fAliESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(selPrimaries);
 986   cuts->SetParamType(kGlobal);
 987   return cuts;
 988 }
 989
 990 //-----------------------------------------------------------------------
 991 Bool_t AliFlowTrackCuts::FillFlowTrackGeneric(AliFlowTrack* flowtrack) const
 992 {
 993   //fill a flow track from tracklet,vzero,pmd,...
 994   TParticle *tmpTParticle=NULL;
 995   AliMCParticle* tmpAliMCParticle=NULL;
 996   switch (fParamMix)
 997   {
 998     case kPure:
 999       flowtrack->SetPhi(fTrackPhi);
1000       flowtrack->SetEta(fTrackEta);
1001       break;
1002     case kTrackWithMCkine:
1003       if (!fMCparticle) return kFALSE;
1004       flowtrack->SetPhi( fMCparticle->Phi() );
1005       flowtrack->SetEta( fMCparticle->Eta() );
1006       flowtrack->SetPt( fMCparticle->Pt() );
1007       break;
1008     case kTrackWithMCpt:
1009       if (!fMCparticle) return kFALSE;
1010       flowtrack->SetPhi(fTrackPhi);
1011       flowtrack->SetEta(fTrackEta);
1012       flowtrack->SetPt(fMCparticle->Pt());
1013       break;
1014     case kTrackWithPtFromFirstMother:
1015       if (!fMCparticle) return kFALSE;
1016       flowtrack->SetPhi(fTrackPhi);
1017       flowtrack->SetEta(fTrackEta);
1018       tmpTParticle = fMCparticle->Particle();
1019       tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
1020       flowtrack->SetPt(tmpAliMCParticle->Pt());
1021       break;
1022     default:
1023       flowtrack->SetPhi(fTrackPhi);
1024       flowtrack->SetEta(fTrackEta);
1025       break;
1026   }
1027   flowtrack->SetSource(AliFlowTrack::kFromTracklet);
1028   return kTRUE;
1029 }
1030
1031 //-----------------------------------------------------------------------
1032 Bool_t AliFlowTrackCuts::FillFlowTrackVParticle(AliFlowTrack* flowtrack) const
1033 {
1034   //fill flow track from AliVParticle (ESD,AOD,MC)
1035   if (!fTrack) return kFALSE;
1036   TParticle *tmpTParticle=NULL;
1037   AliMCParticle* tmpAliMCParticle=NULL;
1038   AliExternalTrackParam* externalParams=NULL;
1039   AliESDtrack* esdtrack=NULL;
1040   switch(fParamMix)
1041   {
1042     case kPure:
1043       flowtrack->Set(fTrack);
1044       break;
1045     case kTrackWithMCkine:
1046       flowtrack->Set(fMCparticle);
1047       break;
1048     case kTrackWithMCPID:
1049       flowtrack->Set(fTrack);
1050       //flowtrack->setPID(...) from mc, when implemented
1051       break;
1052     case kTrackWithMCpt:
1053       if (!fMCparticle) return kFALSE;
1054       flowtrack->Set(fTrack);
1055       flowtrack->SetPt(fMCparticle->Pt());
1056       break;
1057     case kTrackWithPtFromFirstMother:
1058       if (!fMCparticle) return kFALSE;
1059       flowtrack->Set(fTrack);
1060       tmpTParticle = fMCparticle->Particle();
1061       tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
1062       flowtrack->SetPt(tmpAliMCParticle->Pt());
1063       break;
1064     case kTrackWithTPCInnerParams:
1065       esdtrack = dynamic_cast<AliESDtrack*>(fTrack);
1066       if (!esdtrack) return kFALSE;
1067       externalParams = const_cast<AliExternalTrackParam*>(esdtrack->GetTPCInnerParam());
1068       if (!externalParams) return kFALSE;
1069       flowtrack->Set(externalParams);
1070       break;
1071     default:
1072       flowtrack->Set(fTrack);
1073       break;
1074   }
1075   if (fParamType==kMC)
1076   {
1077     flowtrack->SetSource(AliFlowTrack::kFromMC);
1078     flowtrack->SetID(fTrack->GetLabel());
1079   }
1080   else if (dynamic_cast<AliESDtrack*>(fTrack))
1081   {
1082     flowtrack->SetSource(AliFlowTrack::kFromESD);
1083     flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1084   }
1085   else if (dynamic_cast<AliAODTrack*>(fTrack))
1086   {
1087     flowtrack->SetSource(AliFlowTrack::kFromAOD);
1088     flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1089   }
1090   else if (dynamic_cast<AliMCParticle*>(fTrack))
1091   {
1092     flowtrack->SetSource(AliFlowTrack::kFromMC);
1093     flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1094   }
1095   return kTRUE;
1096 }
1097
1098 //-----------------------------------------------------------------------
1099 Bool_t AliFlowTrackCuts::FillFlowTrack(AliFlowTrack* track) const
1100 {
1101   //fill a flow track constructed from whatever we applied cuts on
1102   //return true on success
1103   switch (fParamType)
1104   {
1105     case kSPDtracklet:
1106       return FillFlowTrackGeneric(track);
1107     case kPMD:
1108       return FillFlowTrackGeneric(track);
1109     case kV0:
1110       return FillFlowTrackGeneric(track);
1111     default:
1112       return FillFlowTrackVParticle(track);
1113   }
1114 }
1115
1116 //-----------------------------------------------------------------------
1117 AliFlowTrack* AliFlowTrackCuts::MakeFlowTrackSPDtracklet() const
1118 {
1119   //make a flow track from tracklet
1120   AliFlowTrack* flowtrack=NULL;
1121   TParticle *tmpTParticle=NULL;
1122   AliMCParticle* tmpAliMCParticle=NULL;
1123   switch (fParamMix)
1124   {
1125     case kPure:
1126       flowtrack = new AliFlowTrack();
1127       flowtrack->SetPhi(fTrackPhi);
1128       flowtrack->SetEta(fTrackEta);
1129       break;
1130     case kTrackWithMCkine:
1131       if (!fMCparticle) return NULL;
1132       flowtrack = new AliFlowTrack();
1133       flowtrack->SetPhi( fMCparticle->Phi() );
1134       flowtrack->SetEta( fMCparticle->Eta() );
1135       flowtrack->SetPt( fMCparticle->Pt() );
1136       break;
1137     case kTrackWithMCpt:
1138       if (!fMCparticle) return NULL;
1139       flowtrack = new AliFlowTrack();
1140       flowtrack->SetPhi(fTrackPhi);
1141       flowtrack->SetEta(fTrackEta);
1142       flowtrack->SetPt(fMCparticle->Pt());
1143       break;
1144     case kTrackWithPtFromFirstMother:
1145       if (!fMCparticle) return NULL;
1146       flowtrack = new AliFlowTrack();
1147       flowtrack->SetPhi(fTrackPhi);
1148       flowtrack->SetEta(fTrackEta);
1149       tmpTParticle = fMCparticle->Particle();
1150       tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
1151       flowtrack->SetPt(tmpAliMCParticle->Pt());
1152       break;
1153     default:
1154       flowtrack = new AliFlowTrack();
1155       flowtrack->SetPhi(fTrackPhi);
1156       flowtrack->SetEta(fTrackEta);
1157       break;
1158   }
1159   flowtrack->SetSource(AliFlowTrack::kFromTracklet);
1160   return flowtrack;
1161 }
1162
1163 //-----------------------------------------------------------------------
1164 AliFlowTrack* AliFlowTrackCuts::MakeFlowTrackVParticle() const
1165 {
1166   //make flow track from AliVParticle (ESD,AOD,MC)
1167   if (!fTrack) return NULL;
1168   AliFlowTrack* flowtrack=NULL;
1169   TParticle *tmpTParticle=NULL;
1170   AliMCParticle* tmpAliMCParticle=NULL;
1171   AliExternalTrackParam* externalParams=NULL;
1172   AliESDtrack* esdtrack=NULL;
1173   switch(fParamMix)
1174   {
1175     case kPure:
1176       flowtrack = new AliFlowTrack(fTrack);
1177       break;
1178     case kTrackWithMCkine:
1179       flowtrack = new AliFlowTrack(fMCparticle);
1180       break;
1181     case kTrackWithMCPID:
1182       flowtrack = new AliFlowTrack(fTrack);
1183       //flowtrack->setPID(...) from mc, when implemented
1184       break;
1185     case kTrackWithMCpt:
1186       if (!fMCparticle) return NULL;
1187       flowtrack = new AliFlowTrack(fTrack);
1188       flowtrack->SetPt(fMCparticle->Pt());
1189       break;
1190     case kTrackWithPtFromFirstMother:
1191       if (!fMCparticle) return NULL;
1192       flowtrack = new AliFlowTrack(fTrack);
1193       tmpTParticle = fMCparticle->Particle();
1194       tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
1195       flowtrack->SetPt(tmpAliMCParticle->Pt());
1196       break;
1197     case kTrackWithTPCInnerParams:
1198       esdtrack = dynamic_cast<AliESDtrack*>(fTrack);
1199       if (!esdtrack) return NULL;
1200       externalParams = const_cast<AliExternalTrackParam*>(esdtrack->GetTPCInnerParam());
1201       if (!externalParams) return NULL;
1202       flowtrack = new AliFlowTrack(externalParams);
1203       break;
1204     default:
1205       flowtrack = new AliFlowTrack(fTrack);
1206       break;
1207   }
1208   if (fParamType==kMC)
1209   {
1210     flowtrack->SetSource(AliFlowTrack::kFromMC);
1211     flowtrack->SetID(fTrack->GetLabel());
1212   }
1213   else if (dynamic_cast<AliESDtrack*>(fTrack))
1214   {
1215     flowtrack->SetSource(AliFlowTrack::kFromESD);
1216     flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1217   }
1218   else if (dynamic_cast<AliAODTrack*>(fTrack))
1219   {
1220     flowtrack->SetSource(AliFlowTrack::kFromAOD);
1221     flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1222   }
1223   else if (dynamic_cast<AliMCParticle*>(fTrack))
1224   {
1225     flowtrack->SetSource(AliFlowTrack::kFromMC);
1226     flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
1227   }
1228   return flowtrack;
1229 }
1230
1231 //-----------------------------------------------------------------------
1232 AliFlowTrack* AliFlowTrackCuts::MakeFlowTrackPMDtrack() const
1233 {
1234   //make a flow track from PMD track
1235   AliFlowTrack* flowtrack=NULL;
1236   TParticle *tmpTParticle=NULL;
1237   AliMCParticle* tmpAliMCParticle=NULL;
1238   switch (fParamMix)
1239   {
1240     case kPure:
1241       flowtrack = new AliFlowTrack();
1242       flowtrack->SetPhi(fTrackPhi);
1243       flowtrack->SetEta(fTrackEta);
1244       flowtrack->SetWeight(fTrackWeight);
1245       break;
1246     case kTrackWithMCkine:
1247       if (!fMCparticle) return NULL;
1248       flowtrack = new AliFlowTrack();
1249       flowtrack->SetPhi( fMCparticle->Phi() );
1250       flowtrack->SetEta( fMCparticle->Eta() );
1251       flowtrack->SetWeight(fTrackWeight);
1252       flowtrack->SetPt( fMCparticle->Pt() );
1253       break;
1254     case kTrackWithMCpt:
1255       if (!fMCparticle) return NULL;
1256       flowtrack = new AliFlowTrack();
1257       flowtrack->SetPhi(fTrackPhi);
1258       flowtrack->SetEta(fTrackEta);
1259       flowtrack->SetWeight(fTrackWeight);
1260       flowtrack->SetPt(fMCparticle->Pt());
1261       break;
1262     case kTrackWithPtFromFirstMother:
1263       if (!fMCparticle) return NULL;
1264       flowtrack = new AliFlowTrack();
1265       flowtrack->SetPhi(fTrackPhi);
1266       flowtrack->SetEta(fTrackEta);
1267       flowtrack->SetWeight(fTrackWeight);
1268       tmpTParticle = fMCparticle->Particle();
1269       tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
1270       flowtrack->SetPt(tmpAliMCParticle->Pt());
1271       break;
1272     default:
1273       flowtrack = new AliFlowTrack();
1274       flowtrack->SetPhi(fTrackPhi);
1275       flowtrack->SetEta(fTrackEta);
1276       flowtrack->SetWeight(fTrackWeight);
1277       break;
1278   }
1279
1280   flowtrack->SetSource(AliFlowTrack::kFromPMD);
1281   return flowtrack;
1282 }
1283
1284 //-----------------------------------------------------------------------
1285 AliFlowTrack* AliFlowTrackCuts::MakeFlowTrackV0() const
1286 {
1287   //make a flow track from V0
1288   AliFlowTrack* flowtrack=NULL;
1289   TParticle *tmpTParticle=NULL;
1290   AliMCParticle* tmpAliMCParticle=NULL;
1291   switch (fParamMix)
1292   {
1293     case kPure:
1294       flowtrack = new AliFlowTrack();
1295       flowtrack->SetPhi(fTrackPhi);
1296       flowtrack->SetEta(fTrackEta);
1297       flowtrack->SetWeight(fTrackWeight);
1298       break;
1299     case kTrackWithMCkine:
1300       if (!fMCparticle) return NULL;
1301       flowtrack = new AliFlowTrack();
1302       flowtrack->SetPhi( fMCparticle->Phi() );
1303       flowtrack->SetEta( fMCparticle->Eta() );
1304       flowtrack->SetWeight(fTrackWeight);
1305       flowtrack->SetPt( fMCparticle->Pt() );
1306       break;
1307     case kTrackWithMCpt:
1308       if (!fMCparticle) return NULL;
1309       flowtrack = new AliFlowTrack();
1310       flowtrack->SetPhi(fTrackPhi);
1311       flowtrack->SetEta(fTrackEta);
1312       flowtrack->SetWeight(fTrackWeight);
1313       flowtrack->SetPt(fMCparticle->Pt());
1314       break;
1315     case kTrackWithPtFromFirstMother:
1316       if (!fMCparticle) return NULL;
1317       flowtrack = new AliFlowTrack();
1318       flowtrack->SetPhi(fTrackPhi);
1319       flowtrack->SetEta(fTrackEta);
1320       flowtrack->SetWeight(fTrackWeight);
1321       tmpTParticle = fMCparticle->Particle();
1322       tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
1323       flowtrack->SetPt(tmpAliMCParticle->Pt());
1324       break;
1325     default:
1326       flowtrack = new AliFlowTrack();
1327       flowtrack->SetPhi(fTrackPhi);
1328       flowtrack->SetEta(fTrackEta);
1329       flowtrack->SetWeight(fTrackWeight);
1330       break;
1331   }
1332
1333   flowtrack->SetSource(AliFlowTrack::kFromV0);
1334   return flowtrack;
1335 }
1336
1337 //-----------------------------------------------------------------------
1338 AliFlowTrack* AliFlowTrackCuts::MakeFlowTrack() const
1339 {
1340   //get a flow track constructed from whatever we applied cuts on
1341   //caller is resposible for deletion
1342   //if construction fails return NULL
1343   //TODO: for tracklets, PMD and V0 we probably need just one method,
1344   //something like MakeFlowTrackGeneric(), wait with this until
1345   //requirements quirks are known.
1346   switch (fParamType)
1347   {
1348     case kSPDtracklet:
1349       return MakeFlowTrackSPDtracklet();
1350     case kPMD:
1351       return MakeFlowTrackPMDtrack();
1352     case kV0:
1353       return MakeFlowTrackV0();
1354     default:
1355       return MakeFlowTrackVParticle();
1356   }
1357 }
1358
1359 //-----------------------------------------------------------------------
1360 Bool_t AliFlowTrackCuts::IsPhysicalPrimary() const
1361 {
1362   //check if current particle is a physical primary
1363   if (!fMCevent) return kFALSE;
1364   if (fTrackLabel<0) return kFALSE;
1365   return IsPhysicalPrimary(fMCevent, fTrackLabel, fRequireTransportBitForPrimaries);
1366 }
1367
1368 //-----------------------------------------------------------------------
1369 Bool_t AliFlowTrackCuts::IsPhysicalPrimary(AliMCEvent* mcEvent, Int_t label, Bool_t requiretransported)
1370 {
1371   //check if current particle is a physical primary
1372   Bool_t physprim=mcEvent->IsPhysicalPrimary(label);
1373   AliMCParticle* track = static_cast<AliMCParticle*>(mcEvent->GetTrack(label));
1374   if (!track) return kFALSE;
1375   TParticle* particle = track->Particle();
1376   Bool_t transported = particle->TestBit(kTransportBit);
1377   //printf("label: %i prim: %s, transp: %s, pass: %s\n",label, (physprim)?"YES":"NO ",(transported)?"YES":"NO ",
1378         //(physprim && (transported || !requiretransported))?"YES":"NO"  );
1379   return (physprim && (transported || !requiretransported));
1380 }
1381
1382 //-----------------------------------------------------------------------
1383 void AliFlowTrackCuts::DefineHistograms()
1384 {
1385   //define qa histograms
1386   if (fQA) return;
1387
1388   const Int_t kNbinsP=60;
1389   Double_t binsP[kNbinsP+1];
1390   binsP[0]=0.0;
1391   for(int i=1; i<kNbinsP+1; i++)
1392   {
1393     if(binsP[i-1]+0.05<1.01)
1394       binsP[i]=binsP[i-1]+0.05;
1395     else
1396       binsP[i]=binsP[i-1]+0.1;
1397   }
1398
1399   const Int_t nBinsDCA=1000;
1400   Double_t binsDCA[nBinsDCA+1];
1401   for(int i=0; i<nBinsDCA+1; i++) {binsDCA[i]=0.01*i-5.;}
1402   //for(int i=1; i<41; i++) {binsDCA[i+100]=0.1*i+1.0;}
1403
1404   Bool_t adddirstatus = TH1::AddDirectoryStatus();
1405   TH1::AddDirectory(kFALSE);
1406   fQA=new TList(); fQA->SetOwner();
1407   fQA->SetName(Form("%s QA",GetName()));
1408   TList* before = new TList(); before->SetOwner();
1409   before->SetName("before");
1410   TList* after = new TList(); after->SetOwner();
1411   after->SetName("after");
1412   fQA->Add(before);
1413   fQA->Add(after);
1414   before->Add(new TH2F("TOFbeta",";p [GeV/c];#beta",kNbinsP,binsP,1000,0.4,1.1)); //0
1415   after->Add(new TH2F("TOFbeta",";p [GeV/c];#beta",kNbinsP,binsP,1000,0.4,1.1)); //0
1416   before->Add(new TH2F("TPCdedx",";p [GeV/c];dEdx",kNbinsP,binsP,500,0,500)); //1
1417   after->Add(new TH2F("TPCdedx",";p [GeV/c];dEdx",kNbinsP,binsP,500,0,500)); //1
1418   before->Add(new TH2F("MC pid",";p[GeV/c];species",kNbinsP,binsP,10,-5, 5)); //2
1419   after->Add(new TH2F("MC pid",";p[GeV/c];species",kNbinsP,binsP,10,-5, 5)); //2
1420   before->Add(new TH2F("MC primary",";p[GeV/c];primary",kNbinsP,binsP,2,-1,1)); //3
1421   after->Add(new TH2F("MC primary",";p[GeV/c];primary",kNbinsP,binsP,2,-1,1)); //3
1422   //production process
1423   TH2F* hb = new TH2F("MC production process",";p[GeV/c];",kNbinsP,binsP,kMaxMCProcess,
1424                       -0.5, kMaxMCProcess-0.5);
1425   TH2F* ha = new TH2F("MC production process",";p[GeV/c];",kNbinsP,binsP,kMaxMCProcess,
1426                       -0.5, kMaxMCProcess-0.5);
1427   TAxis* axis = hb->GetYaxis();
1428   for (Int_t i=0; i<kMaxMCProcess; i++)
1429   {
1430     axis->SetBinLabel(i+1,TMCProcessName[i]);
1431   }
1432   axis = ha->GetYaxis();
1433   for (Int_t i=0; i<kMaxMCProcess; i++)
1434   {
1435     axis->SetBinLabel(i+1,TMCProcessName[i]);
1436   }
1437   before->Add(hb); //4
1438   after->Add(ha); //4
1439   //DCA
1440   before->Add(new TH2F("DCAxy",";p_{t}[GeV/c];DCAxy[cm]", 100, 0., 10., nBinsDCA, binsDCA));//5
1441   after->Add(new TH2F("DCAxy",";p_{t}[GeV/c];DCAxy[cm]", 100, 0., 10., nBinsDCA, binsDCA));//5
1442   before->Add(new TH2F("DCAz",";p_{t}[GeV/c];DCAz[cm]", 100, 0., 10., nBinsDCA, binsDCA));//6
1443   after->Add(new TH2F("DCAz",";p_{t}[GeV/c];DCAz[cm]", 100, 0., 10., nBinsDCA, binsDCA));//6
1444   //first mother
1445   hb = new TH2F("MC first mother",";p[GeV/c];",kNbinsP,binsP,4,0.,4.);
1446   ha = new TH2F("MC first mother",";p[GeV/c];",kNbinsP,binsP,4,0.,4.);
1447   hb->GetYaxis()->SetBinLabel(1,"#Lambda");
1448   ha->GetYaxis()->SetBinLabel(1,"#Lambda");
1449   hb->GetYaxis()->SetBinLabel(2,"#Sigma");
1450   ha->GetYaxis()->SetBinLabel(2,"#Sigma");
1451   hb->GetYaxis()->SetBinLabel(3,"#Xi");
1452   ha->GetYaxis()->SetBinLabel(3,"#Xi");
1453   hb->GetYaxis()->SetBinLabel(4,"#Omega");
1454   ha->GetYaxis()->SetBinLabel(4,"#Omega");
1455   before->Add(hb);//7
1456   after->Add(ha);//7
1457
1458   TH1::AddDirectory(adddirstatus);
1459 }
1460
1461 //-----------------------------------------------------------------------
1462 Int_t AliFlowTrackCuts::GetNumberOfInputObjects() const
1463 {
1464   //get the number of tracks in the input event according source
1465   //selection (ESD tracks, tracklets, MC particles etc.)
1466   AliESDEvent* esd=NULL;
1467   switch (fParamType)
1468   {
1469     case kSPDtracklet:
1470       esd = dynamic_cast<AliESDEvent*>(fEvent);
1471       if (!esd) return 0;
1472       return esd->GetMultiplicity()->GetNumberOfTracklets();
1473     case kMC:
1474       if (!fMCevent) return 0;
1475       return fMCevent->GetNumberOfTracks();
1476     case kPMD:
1477       esd = dynamic_cast<AliESDEvent*>(fEvent);
1478       if (!esd) return 0;
1479       return esd->GetNumberOfPmdTracks();
1480     case kV0:
1481       return fgkNumberOfV0tracks;
1482     default:
1483       if (!fEvent) return 0;
1484       return fEvent->GetNumberOfTracks();
1485   }
1486   return 0;
1487 }
1488
1489 //-----------------------------------------------------------------------
1490 TObject* AliFlowTrackCuts::GetInputObject(Int_t i)
1491 {
1492   //get the input object according the data source selection:
1493   //(esd tracks, traclets, mc particles,etc...)
1494   AliESDEvent* esd=NULL;
1495   switch (fParamType)
1496   {
1497     case kSPDtracklet:
1498       esd = dynamic_cast<AliESDEvent*>(fEvent);
1499       if (!esd) return NULL;
1500       return const_cast<AliMultiplicity*>(esd->GetMultiplicity());
1501     case kMC:
1502       if (!fMCevent) return NULL;
1503       return fMCevent->GetTrack(i);
1504     case kPMD:
1505       esd = dynamic_cast<AliESDEvent*>(fEvent);
1506       if (!esd) return NULL;
1507       return esd->GetPmdTrack(i);
1508     case kV0:
1509       esd = dynamic_cast<AliESDEvent*>(fEvent);
1510       if (!esd) //contributed by G.Ortona
1511       {
1512         AliAODEvent* aod = dynamic_cast<AliAODEvent*>(fEvent);
1513         if(!aod)return NULL;
1514         return aod->GetVZEROData();
1515       }
1516       return esd->GetVZEROData();
1517     default:
1518       if (!fEvent) return NULL;
1519       return fEvent->GetTrack(i);
1520   }
1521 }
1522
1523 //-----------------------------------------------------------------------
1524 void AliFlowTrackCuts::Clear(Option_t*)
1525 {
1526   //clean up
1527   fTrack=NULL;
1528   fMCevent=NULL;
1529   fMCparticle=NULL;
1530   fTrackLabel=-1;
1531   fTrackWeight=0.0;
1532   fTrackEta=0.0;
1533   fTrackPhi=0.0;
1534 }
1535
1536 //-----------------------------------------------------------------------
1537 Bool_t AliFlowTrackCuts::PassesTOFbetaSimpleCut(const AliESDtrack* track )
1538 {
1539   //check if passes PID cut using timing in TOF
1540   Bool_t goodtrack = (track->GetStatus() & AliESDtrack::kTOFpid) &&
1541                      (track->GetTOFsignal() > 12000) &&
1542                      (track->GetTOFsignal() < 100000) &&
1543                      (track->GetIntegratedLength() > 365);
1544
1545   if (!fAllowTOFmismatchFlag) {if ((track->GetStatus() & AliESDtrack::kTOFmismatch)) return kFALSE;}
1546
1547   Bool_t statusMatchingHard = TPCTOFagree(track);
1548   if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
1549        return kFALSE;
1550
1551   if (!goodtrack) return kFALSE;
1552
1553   const Float_t c = 2.99792457999999984e-02;
1554   Float_t p = track->GetP();
1555   Float_t l = track->GetIntegratedLength();
1556   Float_t trackT0 = fESDpid.GetTOFResponse().GetStartTime(p);
1557   Float_t timeTOF = track->GetTOFsignal()- trackT0;
1558   Float_t beta = l/timeTOF/c;
1559   Double_t integratedTimes[5] = {-1.0,-1.0,-1.0,-1.0,-1.0};
1560   track->GetIntegratedTimes(integratedTimes);
1561   Float_t betaHypothesis[5] = {0.0,0.0,0.0,0.0,0.0};
1562   Float_t s[5] = {0.0,0.0,0.0,0.0,0.0};
1563   for (Int_t i=0;i<5;i++)
1564   {
1565     betaHypothesis[i] = l/integratedTimes[i]/c;
1566     s[i] = beta-betaHypothesis[i];
1567   }
1568
1569   switch (fParticleID)
1570   {
1571     case AliPID::kPion:
1572       return ( (s[2]<0.015) && (s[2]>-0.015) &&
1573                (s[3]>0.025) &&
1574                (s[4]>0.03) );
1575     case AliPID::kKaon:
1576       return ( (s[3]<0.015) && (s[3]>-0.015) &&
1577                (s[2]<-0.03) &&
1578                (s[4]>0.03) );
1579     case AliPID::kProton:
1580       return ( (s[4]<0.015) && (s[4]>-0.015) &&
1581                (s[3]<-0.025) &&
1582                (s[2]<-0.025) );
1583     default:
1584       return kFALSE;
1585   }
1586   return kFALSE;
1587 }
1588
1589 //-----------------------------------------------------------------------
1590 Float_t AliFlowTrackCuts::GetBeta(const AliESDtrack* track)
1591 {
1592   //get beta
1593   const Float_t c = 2.99792457999999984e-02;
1594   Float_t p = track->GetP();
1595   Float_t l = track->GetIntegratedLength();
1596   Float_t trackT0 = fESDpid.GetTOFResponse().GetStartTime(p);
1597   Float_t timeTOF = track->GetTOFsignal()- trackT0;
1598   return l/timeTOF/c;
1599 }
1600
1601 //-----------------------------------------------------------------------
1602 Bool_t AliFlowTrackCuts::PassesTOFbetaCut(const AliESDtrack* track )
1603 {
1604   if (!fTOFpidCuts)
1605   {
1606     //printf("no TOFpidCuts\n");
1607     return kFALSE;
1608   }
1609
1610   //check if passes PID cut using timing in TOF
1611   Bool_t goodtrack = (track->GetStatus() & AliESDtrack::kTOFpid) &&
1612                      (track->GetTOFsignal() > 12000) &&
1613                      (track->GetTOFsignal() < 100000) &&
1614                      (track->GetIntegratedLength() > 365);
1615
1616   if (!fAllowTOFmismatchFlag) {if ((track->GetStatus() & AliESDtrack::kTOFmismatch)) return kFALSE;}
1617
1618   Bool_t statusMatchingHard = TPCTOFagree(track);
1619   if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
1620        return kFALSE;
1621
1622   if (!goodtrack) return kFALSE;
1623
1624   Float_t beta = GetBeta(track);
1625
1626   Double_t integratedTimes[5] = {-1.0,-1.0,-1.0,-1.0,-1.0};
1627   track->GetIntegratedTimes(integratedTimes);
1628
1629   //construct the pid index because it's not AliPID::EParticleType
1630   Int_t pid = 0;
1631   switch (fParticleID)
1632   {
1633     case AliPID::kPion:
1634       pid=2;
1635       break;
1636     case AliPID::kKaon:
1637       pid=3;
1638       break;
1639     case AliPID::kProton:
1640       pid=4;
1641       break;
1642     default:
1643       return kFALSE;
1644   }
1645
1646   //signal to cut on
1647   const Float_t c = 2.99792457999999984e-02;
1648   Float_t l = track->GetIntegratedLength();
1649   Float_t p = track->GetP();
1650   Float_t betahypothesis = l/integratedTimes[pid]/c;
1651   Float_t betadiff = beta-betahypothesis;
1652
1653   Float_t* arr = fTOFpidCuts->GetMatrixArray();
1654   Int_t col = TMath::BinarySearch(fTOFpidCuts->GetNcols(),arr,static_cast<Float_t>(p));
1655   if (col<0) return kFALSE;
1656   Float_t min = (*fTOFpidCuts)(1,col);
1657   Float_t max = (*fTOFpidCuts)(2,col);
1658
1659   Bool_t pass = (betadiff>min && betadiff<max);
1660
1661   return pass;
1662 }
1663
1664 //-----------------------------------------------------------------------
1665 Bool_t AliFlowTrackCuts::PassesTOFpidCut(const AliESDtrack* track) const
1666 {
1667   //check if passes PID cut using default TOF pid
1668   Double_t pidTOF[AliPID::kSPECIES];
1669   track->GetTOFpid(pidTOF);
1670   if (pidTOF[fParticleID]>=fParticleProbability) return kTRUE;
1671   return kFALSE;
1672 }
1673
1674 //-----------------------------------------------------------------------
1675 Bool_t AliFlowTrackCuts::PassesTPCpidCut(const AliESDtrack* track) const
1676 {
1677   //check if passes PID cut using default TPC pid
1678   Double_t pidTPC[AliPID::kSPECIES];
1679   track->GetTPCpid(pidTPC);
1680   Double_t probablity = 0.;
1681   switch (fParticleID)
1682   {
1683     case AliPID::kPion:
1684       probablity = pidTPC[AliPID::kPion] + pidTPC[AliPID::kMuon];
1685       break;
1686     default:
1687       probablity = pidTPC[fParticleID];
1688   }
1689   if (probablity >= fParticleProbability) return kTRUE;
1690   return kFALSE;
1691 }
1692
1693 //-----------------------------------------------------------------------
1694 Float_t AliFlowTrackCuts::Getdedx(const AliESDtrack* track)
1695 {
1696   //get TPC dedx
1697   return track->GetTPCsignal();
1698 }
1699
1700 //-----------------------------------------------------------------------
1701 Bool_t AliFlowTrackCuts::PassesTPCdedxCut(const AliESDtrack* track)
1702 {
1703   //check if passes PID cut using dedx signal in the TPC
1704   if (!fTPCpidCuts)
1705   {
1706     //printf("no TPCpidCuts\n");
1707     return kFALSE;
1708   }
1709
1710   const AliExternalTrackParam* tpcparam = track->GetInnerParam(); //tpc only params at the inner wall
1711   if (!tpcparam) return kFALSE;
1712   Double_t p = tpcparam->GetP();
1713   Float_t sigExp = fESDpid.GetTPCResponse().GetExpectedSignal(p, fParticleID);
1714   Float_t sigTPC = track->GetTPCsignal();
1715   Float_t s = (sigTPC-sigExp)/sigExp;
1716
1717   Float_t* arr = fTPCpidCuts->GetMatrixArray();
1718   Int_t arrSize = fTPCpidCuts->GetNcols();
1719   Int_t col = TMath::BinarySearch( arrSize, arr, static_cast<Float_t>(p));
1720   if (col<0) return kFALSE;
1721   Float_t min = (*fTPCpidCuts)(1,col);
1722   Float_t max = (*fTPCpidCuts)(2,col);
1723
1724   //printf("------------TPC pid cut %s\n",(s>min && s<max)?"PASS":"FAIL");
1725   return (s>min && s<max);
1726 }
1727
1728 //-----------------------------------------------------------------------
1729 void AliFlowTrackCuts::InitPIDcuts()
1730 {
1731   //init matrices with PID cuts
1732   TMatrixF* t = NULL;
1733   if (!fTPCpidCuts)
1734   {
1735     if (fParticleID==AliPID::kPion)
1736     {
1737       t = new TMatrixF(3,15);
1738       (*t)(0,0)  = 0.20;  (*t)(1,0)  = -0.4;  (*t)(2,0)  =   0.0;
1739       (*t)(0,1)  = 0.25;  (*t)(1,1)  = -0.4;  (*t)(2,1)  =   0.1;
1740       (*t)(0,2)  = 0.30;  (*t)(1,2)  = -0.4;  (*t)(2,2)  =  0.2;
1741       (*t)(0,3)  = 0.35;  (*t)(1,3)  = -0.4;  (*t)(2,3)  =  0.2;
1742       (*t)(0,4)  = 0.40;  (*t)(1,4)  = -0.4;  (*t)(2,4)  =   0.3;
1743       (*t)(0,5)  = 0.45;  (*t)(1,5)  = -0.4;  (*t)(2,5)  =   0.3;
1744       (*t)(0,6)  = 0.50;  (*t)(1,6)  = -0.4;  (*t)(2,6)  =  0.25;
1745       (*t)(0,7)  = 0.55;  (*t)(1,7)  = -0.4;  (*t)(2,7)  =  0.15;
1746       (*t)(0,8)  = 0.60;  (*t)(1,8)  = -0.4;  (*t)(2,8)  =   0.1;
1747       (*t)(0,9)  = 0.65;  (*t)(1,9)  = -0.4;  (*t)(2,9)  =  0.05;
1748       (*t)(0,10)  = 0.70;  (*t)(1,10)  = -0.4;  (*t)(2,10)  =     0;
1749       (*t)(0,11)  = 0.75;  (*t)(1,11)  = -0.4;  (*t)(2,11)  =     0;
1750       (*t)(0,12)  = 0.80;  (*t)(1,12)  = -0.4;  (*t)(2,12)  = -0.05;
1751       (*t)(0,13)  = 0.85;  (*t)(1,13)  = -0.4;  (*t)(2,13)  = -0.1;
1752       (*t)(0,14)  = 0.90;  (*t)(1,14)  = 0;     (*t)(2,14)  =     0;
1753     }
1754     else
1755     if (fParticleID==AliPID::kKaon)
1756     {
1757       t = new TMatrixF(3,12);
1758       (*t)(0,0)  = 0.20;  (*t)(1,0)  = -0.2;  (*t)(2,0)  = 0.2;
1759       (*t)(0,1)  = 0.25;  (*t)(1,1)  = -0.2;  (*t)(2,1)  = 0.2;
1760       (*t)(0,2)  = 0.30;  (*t)(1,2)  = -0.2;  (*t)(2,2)  = 0.2;
1761       (*t)(0,3)  = 0.35;  (*t)(1,3)  = -0.2;  (*t)(2,3)  = 0.2;
1762       (*t)(0,4)  = 0.40;  (*t)(1,4)  = -0.1;  (*t)(2,4)  = 0.2;
1763       (*t)(0,5)  = 0.45;  (*t)(1,5)  = -0.1;  (*t)(2,5)  = 0.2;
1764       (*t)(0,6)  = 0.50;  (*t)(1,6)  =-0.05;  (*t)(2,6)  = 0.2;
1765       (*t)(0,7)  = 0.55;  (*t)(1,7)  = -0.1;  (*t)(2,7)  = 0.1;
1766       (*t)(0,8)  = 0.60;  (*t)(1,8)  =-0.05;  (*t)(2,8)  = 0.1;
1767       (*t)(0,9)  = 0.65;  (*t)(1,9)  =    0;  (*t)(2,9)  = 0.15;
1768       (*t)(0,10)  = 0.70;  (*t)(1,10)  = 0.05;  (*t)(2,10)  = 0.2;
1769       (*t)(0,11)  = 0.75;  (*t)(1,11)  =    0;  (*t)(2,11)  = 0;
1770     }
1771     else
1772     if (fParticleID==AliPID::kProton)
1773     {
1774       t = new TMatrixF(3,9);
1775       (*t)(0,0)  = 0.20;  (*t)(1,0)  = -0.1;  (*t)(2,0)  =  0.1;
1776       (*t)(0,1)  = 0.25;  (*t)(1,1)  = -0.2;  (*t)(2,1)  =  0.2;
1777       (*t)(0,2)  = 0.80;  (*t)(1,2)  = -0.1;  (*t)(2,2)  =  0.2;
1778       (*t)(0,3)  = 0.85;  (*t)(1,3)  =-0.05;  (*t)(2,3)  =  0.2;
1779       (*t)(0,4)  = 0.90;  (*t)(1,4)  =-0.05;  (*t)(2,4)  = 0.25;
1780       (*t)(0,5)  = 0.95;  (*t)(1,5)  =-0.05;  (*t)(2,5)  = 0.25;
1781       (*t)(0,6)  = 1.00;  (*t)(1,6)  = -0.1;  (*t)(2,6)  = 0.25;
1782       (*t)(0,7)  = 1.10;  (*t)(1,7)  =-0.05;  (*t)(2,7)  =  0.3;
1783       (*t)(0,8) = 1.20;   (*t)(1,8)  =    0;  (*t)(2,8) =    0;
1784     }
1785     delete fTPCpidCuts;
1786     fTPCpidCuts=t;
1787   }
1788   t = NULL;
1789   if (!fTOFpidCuts)
1790   {
1791     if (fParticleID==AliPID::kPion)
1792     {
1793       //TOF pions, 0.9 purity
1794       t = new TMatrixF(3,61);
1795       (*t)(0,0)  = 0.000;  (*t)(1,0)  = 0.000;  (*t)(2,0)  =   0.000;
1796       (*t)(0,1)  = 0.050;  (*t)(1,1)  = 0.000;  (*t)(2,1)  =   0.000;
1797       (*t)(0,2)  = 0.100;  (*t)(1,2)  = 0.000;  (*t)(2,2)  =   0.000;
1798       (*t)(0,3)  = 0.150;  (*t)(1,3)  = 0.000;  (*t)(2,3)  =   0.000;
1799       (*t)(0,4)  = 0.200;  (*t)(1,4)  = -0.030;  (*t)(2,4)  =   0.030;
1800       (*t)(0,5)  = 0.250;  (*t)(1,5)  = -0.036;  (*t)(2,5)  =   0.032;
1801       (*t)(0,6)  = 0.300;  (*t)(1,6)  = -0.038;  (*t)(2,6)  =   0.032;
1802       (*t)(0,7)  = 0.350;  (*t)(1,7)  = -0.034;  (*t)(2,7)  =   0.032;
1803       (*t)(0,8)  = 0.400;  (*t)(1,8)  = -0.032;  (*t)(2,8)  =   0.020;
1804       (*t)(0,9)  = 0.450;  (*t)(1,9)  = -0.030;  (*t)(2,9)  =   0.020;
1805       (*t)(0,10)  = 0.500;  (*t)(1,10)  = -0.030;  (*t)(2,10)  =   0.020;
1806       (*t)(0,11)  = 0.550;  (*t)(1,11)  = -0.030;  (*t)(2,11)  =   0.020;
1807       (*t)(0,12)  = 0.600;  (*t)(1,12)  = -0.030;  (*t)(2,12)  =   0.020;
1808       (*t)(0,13)  = 0.650;  (*t)(1,13)  = -0.030;  (*t)(2,13)  =   0.020;
1809       (*t)(0,14)  = 0.700;  (*t)(1,14)  = -0.030;  (*t)(2,14)  =   0.020;
1810       (*t)(0,15)  = 0.750;  (*t)(1,15)  = -0.030;  (*t)(2,15)  =   0.020;
1811       (*t)(0,16)  = 0.800;  (*t)(1,16)  = -0.030;  (*t)(2,16)  =   0.020;
1812       (*t)(0,17)  = 0.850;  (*t)(1,17)  = -0.030;  (*t)(2,17)  =   0.020;
1813       (*t)(0,18)  = 0.900;  (*t)(1,18)  = -0.030;  (*t)(2,18)  =   0.020;
1814       (*t)(0,19)  = 0.950;  (*t)(1,19)  = -0.028;  (*t)(2,19)  =   0.028;
1815       (*t)(0,20)  = 1.000;  (*t)(1,20)  = -0.028;  (*t)(2,20)  =   0.028;
1816       (*t)(0,21)  = 1.100;  (*t)(1,21)  = -0.028;  (*t)(2,21)  =   0.028;
1817       (*t)(0,22)  = 1.200;  (*t)(1,22)  = -0.026;  (*t)(2,22)  =   0.028;
1818       (*t)(0,23)  = 1.300;  (*t)(1,23)  = -0.024;  (*t)(2,23)  =   0.028;
1819       (*t)(0,24)  = 1.400;  (*t)(1,24)  = -0.020;  (*t)(2,24)  =   0.028;
1820       (*t)(0,25)  = 1.500;  (*t)(1,25)  = -0.018;  (*t)(2,25)  =   0.028;
1821       (*t)(0,26)  = 1.600;  (*t)(1,26)  = -0.016;  (*t)(2,26)  =   0.028;
1822       (*t)(0,27)  = 1.700;  (*t)(1,27)  = -0.014;  (*t)(2,27)  =   0.028;
1823       (*t)(0,28)  = 1.800;  (*t)(1,28)  = -0.012;  (*t)(2,28)  =   0.026;
1824       (*t)(0,29)  = 1.900;  (*t)(1,29)  = -0.010;  (*t)(2,29)  =   0.026;
1825       (*t)(0,30)  = 2.000;  (*t)(1,30)  = -0.008;  (*t)(2,30)  =   0.026;
1826       (*t)(0,31)  = 2.100;  (*t)(1,31)  = -0.008;  (*t)(2,31)  =   0.024;
1827       (*t)(0,32)  = 2.200;  (*t)(1,32)  = -0.006;  (*t)(2,32)  =   0.024;
1828       (*t)(0,33)  = 2.300;  (*t)(1,33)  = -0.004;  (*t)(2,33)  =   0.024;
1829       (*t)(0,34)  = 2.400;  (*t)(1,34)  = -0.004;  (*t)(2,34)  =   0.024;
1830       (*t)(0,35)  = 2.500;  (*t)(1,35)  = -0.002;  (*t)(2,35)  =   0.024;
1831       (*t)(0,36)  = 2.600;  (*t)(1,36)  = -0.002;  (*t)(2,36)  =   0.024;
1832       (*t)(0,37)  = 2.700;  (*t)(1,37)  = 0.000;  (*t)(2,37)  =   0.024;
1833       (*t)(0,38)  = 2.800;  (*t)(1,38)  = 0.000;  (*t)(2,38)  =   0.026;
1834       (*t)(0,39)  = 2.900;  (*t)(1,39)  = 0.000;  (*t)(2,39)  =   0.024;
1835       (*t)(0,40)  = 3.000;  (*t)(1,40)  = 0.002;  (*t)(2,40)  =   0.026;
1836       (*t)(0,41)  = 3.100;  (*t)(1,41)  = 0.002;  (*t)(2,41)  =   0.026;
1837       (*t)(0,42)  = 3.200;  (*t)(1,42)  = 0.002;  (*t)(2,42)  =   0.026;
1838       (*t)(0,43)  = 3.300;  (*t)(1,43)  = 0.002;  (*t)(2,43)  =   0.026;
1839       (*t)(0,44)  = 3.400;  (*t)(1,44)  = 0.002;  (*t)(2,44)  =   0.026;
1840       (*t)(0,45)  = 3.500;  (*t)(1,45)  = 0.002;  (*t)(2,45)  =   0.026;
1841       (*t)(0,46)  = 3.600;  (*t)(1,46)  = 0.002;  (*t)(2,46)  =   0.026;
1842       (*t)(0,47)  = 3.700;  (*t)(1,47)  = 0.002;  (*t)(2,47)  =   0.026;
1843       (*t)(0,48)  = 3.800;  (*t)(1,48)  = 0.002;  (*t)(2,48)  =   0.026;
1844       (*t)(0,49)  = 3.900;  (*t)(1,49)  = 0.004;  (*t)(2,49)  =   0.024;
1845       (*t)(0,50)  = 4.000;  (*t)(1,50)  = 0.004;  (*t)(2,50)  =   0.026;
1846       (*t)(0,51)  = 4.100;  (*t)(1,51)  = 0.004;  (*t)(2,51)  =   0.026;
1847       (*t)(0,52)  = 4.200;  (*t)(1,52)  = 0.004;  (*t)(2,52)  =   0.024;
1848       (*t)(0,53)  = 4.300;  (*t)(1,53)  = 0.006;  (*t)(2,53)  =   0.024;
1849       (*t)(0,54)  = 4.400;  (*t)(1,54)  = 0.000;  (*t)(2,54)  =   0.000;
1850       (*t)(0,55)  = 4.500;  (*t)(1,55)  = 0.000;  (*t)(2,55)  =   0.000;
1851       (*t)(0,56)  = 4.600;  (*t)(1,56)  = 0.000;  (*t)(2,56)  =   0.000;
1852       (*t)(0,57)  = 4.700;  (*t)(1,57)  = 0.000;  (*t)(2,57)  =   0.000;
1853       (*t)(0,58)  = 4.800;  (*t)(1,58)  = 0.000;  (*t)(2,58)  =   0.000;
1854       (*t)(0,59)  = 4.900;  (*t)(1,59)  = 0.000;  (*t)(2,59)  =   0.000;
1855       (*t)(0,60)  = 5.900;  (*t)(1,60)  = 0.000;  (*t)(2,60)  =   0.000;
1856     }
1857     else
1858     if (fParticleID==AliPID::kProton)
1859     {
1860       //TOF protons, 0.9 purity
1861       t = new TMatrixF(3,61);
1862       (*t)(0,0)  = 0.000;  (*t)(1,0)  = 0.000;  (*t)(2,0)  =   0.000;
1863       (*t)(0,1)  = 0.050;  (*t)(1,1)  = 0.000;  (*t)(2,1)  =   0.000;
1864       (*t)(0,2)  = 0.100;  (*t)(1,2)  = 0.000;  (*t)(2,2)  =   0.000;
1865       (*t)(0,3)  = 0.150;  (*t)(1,3)  = 0.000;  (*t)(2,3)  =   0.000;
1866       (*t)(0,4)  = 0.200;  (*t)(1,4)  = -0.07;  (*t)(2,4)  =   0.07;
1867       (*t)(0,5)  = 0.200;  (*t)(1,5)  = -0.07;  (*t)(2,5)  =   0.07;
1868       (*t)(0,6)  = 0.200;  (*t)(1,6)  = -0.07;  (*t)(2,6)  =   0.07;
1869       (*t)(0,7)  = 0.200;  (*t)(1,7)  = -0.07;  (*t)(2,7)  =   0.07;
1870       (*t)(0,8)  = 0.200;  (*t)(1,8)  = -0.07;  (*t)(2,8)  =   0.07;
1871       (*t)(0,9)  = 0.200;  (*t)(1,9)  = -0.07;  (*t)(2,9)  =   0.07;
1872       (*t)(0,10)  = 0.200;  (*t)(1,10)  = -0.07;  (*t)(2,10)  =   0.07;
1873       (*t)(0,11)  = 0.200;  (*t)(1,11)  = -0.07;  (*t)(2,11)  =   0.07;
1874       (*t)(0,12)  = 0.200;  (*t)(1,12)  = -0.07;  (*t)(2,12)  =   0.07;
1875       (*t)(0,13)  = 0.200;  (*t)(1,13)  = -0.07;  (*t)(2,13)  =   0.07;
1876       (*t)(0,14)  = 0.200;  (*t)(1,14)  = -0.07;  (*t)(2,14)  =   0.07;
1877       (*t)(0,15)  = 0.200;  (*t)(1,15)  = -0.07;  (*t)(2,15)  =   0.07;
1878       (*t)(0,16)  = 0.200;  (*t)(1,16)  = -0.07;  (*t)(2,16)  =   0.07;
1879       (*t)(0,17)  = 0.850;  (*t)(1,17)  = -0.070;  (*t)(2,17)  =   0.070;
1880       (*t)(0,18)  = 0.900;  (*t)(1,18)  = -0.072;  (*t)(2,18)  =   0.072;
1881       (*t)(0,19)  = 0.950;  (*t)(1,19)  = -0.072;  (*t)(2,19)  =   0.072;
1882       (*t)(0,20)  = 1.000;  (*t)(1,20)  = -0.074;  (*t)(2,20)  =   0.074;
1883       (*t)(0,21)  = 1.100;  (*t)(1,21)  = -0.032;  (*t)(2,21)  =   0.032;
1884       (*t)(0,22)  = 1.200;  (*t)(1,22)  = -0.026;  (*t)(2,22)  =   0.026;
1885       (*t)(0,23)  = 1.300;  (*t)(1,23)  = -0.026;  (*t)(2,23)  =   0.026;
1886       (*t)(0,24)  = 1.400;  (*t)(1,24)  = -0.024;  (*t)(2,24)  =   0.024;
1887       (*t)(0,25)  = 1.500;  (*t)(1,25)  = -0.024;  (*t)(2,25)  =   0.024;
1888       (*t)(0,26)  = 1.600;  (*t)(1,26)  = -0.026;  (*t)(2,26)  =   0.026;
1889       (*t)(0,27)  = 1.700;  (*t)(1,27)  = -0.026;  (*t)(2,27)  =   0.026;
1890       (*t)(0,28)  = 1.800;  (*t)(1,28)  = -0.026;  (*t)(2,28)  =   0.026;
1891       (*t)(0,29)  = 1.900;  (*t)(1,29)  = -0.026;  (*t)(2,29)  =   0.026;
1892       (*t)(0,30)  = 2.000;  (*t)(1,30)  = -0.026;  (*t)(2,30)  =   0.026;
1893       (*t)(0,31)  = 2.100;  (*t)(1,31)  = -0.026;  (*t)(2,31)  =   0.026;
1894       (*t)(0,32)  = 2.200;  (*t)(1,32)  = -0.026;  (*t)(2,32)  =   0.024;
1895       (*t)(0,33)  = 2.300;  (*t)(1,33)  = -0.028;  (*t)(2,33)  =   0.022;
1896       (*t)(0,34)  = 2.400;  (*t)(1,34)  = -0.028;  (*t)(2,34)  =   0.020;
1897       (*t)(0,35)  = 2.500;  (*t)(1,35)  = -0.028;  (*t)(2,35)  =   0.018;
1898       (*t)(0,36)  = 2.600;  (*t)(1,36)  = -0.028;  (*t)(2,36)  =   0.016;
1899       (*t)(0,37)  = 2.700;  (*t)(1,37)  = -0.028;  (*t)(2,37)  =   0.016;
1900       (*t)(0,38)  = 2.800;  (*t)(1,38)  = -0.030;  (*t)(2,38)  =   0.014;
1901       (*t)(0,39)  = 2.900;  (*t)(1,39)  = -0.030;  (*t)(2,39)  =   0.012;
1902       (*t)(0,40)  = 3.000;  (*t)(1,40)  = -0.030;  (*t)(2,40)  =   0.012;
1903       (*t)(0,41)  = 3.100;  (*t)(1,41)  = -0.030;  (*t)(2,41)  =   0.010;
1904       (*t)(0,42)  = 3.200;  (*t)(1,42)  = -0.030;  (*t)(2,42)  =   0.010;
1905       (*t)(0,43)  = 3.300;  (*t)(1,43)  = -0.030;  (*t)(2,43)  =   0.010;
1906       (*t)(0,44)  = 3.400;  (*t)(1,44)  = -0.030;  (*t)(2,44)  =   0.008;
1907       (*t)(0,45)  = 3.500;  (*t)(1,45)  = -0.030;  (*t)(2,45)  =   0.008;
1908       (*t)(0,46)  = 3.600;  (*t)(1,46)  = -0.030;  (*t)(2,46)  =   0.008;
1909       (*t)(0,47)  = 3.700;  (*t)(1,47)  = -0.030;  (*t)(2,47)  =   0.006;
1910       (*t)(0,48)  = 3.800;  (*t)(1,48)  = -0.030;  (*t)(2,48)  =   0.006;
1911       (*t)(0,49)  = 3.900;  (*t)(1,49)  = -0.030;  (*t)(2,49)  =   0.006;
1912       (*t)(0,50)  = 4.000;  (*t)(1,50)  = -0.028;  (*t)(2,50)  =   0.004;
1913       (*t)(0,51)  = 4.100;  (*t)(1,51)  = -0.030;  (*t)(2,51)  =   0.004;
1914       (*t)(0,52)  = 4.200;  (*t)(1,52)  = -0.030;  (*t)(2,52)  =   0.004;
1915       (*t)(0,53)  = 4.300;  (*t)(1,53)  = -0.028;  (*t)(2,53)  =   0.002;
1916       (*t)(0,54)  = 4.400;  (*t)(1,54)  = -0.030;  (*t)(2,54)  =   0.002;
1917       (*t)(0,55)  = 4.500;  (*t)(1,55)  = -0.028;  (*t)(2,55)  =   0.002;
1918       (*t)(0,56)  = 4.600;  (*t)(1,56)  = -0.028;  (*t)(2,56)  =   0.002;
1919       (*t)(0,57)  = 4.700;  (*t)(1,57)  = -0.028;  (*t)(2,57)  =   0.000;
1920       (*t)(0,58)  = 4.800;  (*t)(1,58)  = -0.028;  (*t)(2,58)  =   0.002;
1921       (*t)(0,59)  = 4.900;  (*t)(1,59)  = 0.000;  (*t)(2,59)  =   0.000;
1922       (*t)(0,60)  = 5.900;  (*t)(1,60)  = 0.000;  (*t)(2,60)  =   0.000;
1923     }
1924     else
1925     if (fParticleID==AliPID::kKaon)
1926     {
1927       //TOF kaons, 0.9 purity
1928       t = new TMatrixF(3,61);
1929       (*t)(0,0)  = 0.000;  (*t)(1,0)  = 0.000;  (*t)(2,0)  =   0.000;
1930       (*t)(0,1)  = 0.050;  (*t)(1,1)  = 0.000;  (*t)(2,1)  =   0.000;
1931       (*t)(0,2)  = 0.100;  (*t)(1,2)  = 0.000;  (*t)(2,2)  =   0.000;
1932       (*t)(0,3)  = 0.150;  (*t)(1,3)  = 0.000;  (*t)(2,3)  =   0.000;
1933       (*t)(0,4)  = 0.200;  (*t)(1,4)  = -0.05;  (*t)(2,4)  =   0.05;
1934       (*t)(0,5)  = 0.200;  (*t)(1,5)  = -0.05;  (*t)(2,5)  =   0.05;
1935       (*t)(0,6)  = 0.200;  (*t)(1,6)  = -0.05;  (*t)(2,6)  =   0.05;
1936       (*t)(0,7)  = 0.200;  (*t)(1,7)  = -0.05;  (*t)(2,7)  =   0.05;
1937       (*t)(0,8)  = 0.200;  (*t)(1,8)  = -0.05;  (*t)(2,8)  =   0.05;
1938       (*t)(0,9)  = 0.200;  (*t)(1,9)  = -0.05;  (*t)(2,9)  =   0.05;
1939       (*t)(0,10)  = 0.200;  (*t)(1,10)  = -0.05;  (*t)(2,10)  =   0.05;
1940       (*t)(0,11)  = 0.550;  (*t)(1,11)  = -0.026;  (*t)(2,11)  =   0.026;
1941       (*t)(0,12)  = 0.600;  (*t)(1,12)  = -0.026;  (*t)(2,12)  =   0.026;
1942       (*t)(0,13)  = 0.650;  (*t)(1,13)  = -0.026;  (*t)(2,13)  =   0.026;
1943       (*t)(0,14)  = 0.700;  (*t)(1,14)  = -0.026;  (*t)(2,14)  =   0.026;
1944       (*t)(0,15)  = 0.750;  (*t)(1,15)  = -0.026;  (*t)(2,15)  =   0.026;
1945       (*t)(0,16)  = 0.800;  (*t)(1,16)  = -0.026;  (*t)(2,16)  =   0.026;
1946       (*t)(0,17)  = 0.850;  (*t)(1,17)  = -0.024;  (*t)(2,17)  =   0.024;
1947       (*t)(0,18)  = 0.900;  (*t)(1,18)  = -0.024;  (*t)(2,18)  =   0.024;
1948       (*t)(0,19)  = 0.950;  (*t)(1,19)  = -0.024;  (*t)(2,19)  =   0.024;
1949       (*t)(0,20)  = 1.000;  (*t)(1,20)  = -0.024;  (*t)(2,20)  =   0.024;
1950       (*t)(0,21)  = 1.100;  (*t)(1,21)  = -0.024;  (*t)(2,21)  =   0.024;
1951       (*t)(0,22)  = 1.200;  (*t)(1,22)  = -0.024;  (*t)(2,22)  =   0.022;
1952       (*t)(0,23)  = 1.300;  (*t)(1,23)  = -0.024;  (*t)(2,23)  =   0.020;
1953       (*t)(0,24)  = 1.400;  (*t)(1,24)  = -0.026;  (*t)(2,24)  =   0.016;
1954       (*t)(0,25)  = 1.500;  (*t)(1,25)  = -0.028;  (*t)(2,25)  =   0.014;
1955       (*t)(0,26)  = 1.600;  (*t)(1,26)  = -0.028;  (*t)(2,26)  =   0.012;
1956       (*t)(0,27)  = 1.700;  (*t)(1,27)  = -0.028;  (*t)(2,27)  =   0.010;
1957       (*t)(0,28)  = 1.800;  (*t)(1,28)  = -0.028;  (*t)(2,28)  =   0.010;
1958       (*t)(0,29)  = 1.900;  (*t)(1,29)  = -0.028;  (*t)(2,29)  =   0.008;
1959       (*t)(0,30)  = 2.000;  (*t)(1,30)  = -0.028;  (*t)(2,30)  =   0.006;
1960       (*t)(0,31)  = 2.100;  (*t)(1,31)  = -0.026;  (*t)(2,31)  =   0.006;
1961       (*t)(0,32)  = 2.200;  (*t)(1,32)  = -0.024;  (*t)(2,32)  =   0.004;
1962       (*t)(0,33)  = 2.300;  (*t)(1,33)  = -0.020;  (*t)(2,33)  =   0.002;
1963       (*t)(0,34)  = 2.400;  (*t)(1,34)  = -0.020;  (*t)(2,34)  =   0.002;
1964       (*t)(0,35)  = 2.500;  (*t)(1,35)  = -0.018;  (*t)(2,35)  =   0.000;
1965       (*t)(0,36)  = 2.600;  (*t)(1,36)  = -0.016;  (*t)(2,36)  =   0.000;
1966       (*t)(0,37)  = 2.700;  (*t)(1,37)  = -0.014;  (*t)(2,37)  =   -0.002;
1967       (*t)(0,38)  = 2.800;  (*t)(1,38)  = -0.014;  (*t)(2,38)  =   -0.004;
1968       (*t)(0,39)  = 2.900;  (*t)(1,39)  = -0.012;  (*t)(2,39)  =   -0.004;
1969       (*t)(0,40)  = 3.000;  (*t)(1,40)  = -0.010;  (*t)(2,40)  =   -0.006;
1970       (*t)(0,41)  = 3.100;  (*t)(1,41)  = 0.000;  (*t)(2,41)  =   0.000;
1971       (*t)(0,42)  = 3.200;  (*t)(1,42)  = 0.000;  (*t)(2,42)  =   0.000;
1972       (*t)(0,43)  = 3.300;  (*t)(1,43)  = 0.000;  (*t)(2,43)  =   0.000;
1973       (*t)(0,44)  = 3.400;  (*t)(1,44)  = 0.000;  (*t)(2,44)  =   0.000;
1974       (*t)(0,45)  = 3.500;  (*t)(1,45)  = 0.000;  (*t)(2,45)  =   0.000;
1975       (*t)(0,46)  = 3.600;  (*t)(1,46)  = 0.000;  (*t)(2,46)  =   0.000;
1976       (*t)(0,47)  = 3.700;  (*t)(1,47)  = 0.000;  (*t)(2,47)  =   0.000;
1977       (*t)(0,48)  = 3.800;  (*t)(1,48)  = 0.000;  (*t)(2,48)  =   0.000;
1978       (*t)(0,49)  = 3.900;  (*t)(1,49)  = 0.000;  (*t)(2,49)  =   0.000;
1979       (*t)(0,50)  = 4.000;  (*t)(1,50)  = 0.000;  (*t)(2,50)  =   0.000;
1980       (*t)(0,51)  = 4.100;  (*t)(1,51)  = 0.000;  (*t)(2,51)  =   0.000;
1981       (*t)(0,52)  = 4.200;  (*t)(1,52)  = 0.000;  (*t)(2,52)  =   0.000;
1982       (*t)(0,53)  = 4.300;  (*t)(1,53)  = 0.000;  (*t)(2,53)  =   0.000;
1983       (*t)(0,54)  = 4.400;  (*t)(1,54)  = 0.000;  (*t)(2,54)  =   0.000;
1984       (*t)(0,55)  = 4.500;  (*t)(1,55)  = 0.000;  (*t)(2,55)  =   0.000;
1985       (*t)(0,56)  = 4.600;  (*t)(1,56)  = 0.000;  (*t)(2,56)  =   0.000;
1986       (*t)(0,57)  = 4.700;  (*t)(1,57)  = 0.000;  (*t)(2,57)  =   0.000;
1987       (*t)(0,58)  = 4.800;  (*t)(1,58)  = 0.000;  (*t)(2,58)  =   0.000;
1988       (*t)(0,59)  = 4.900;  (*t)(1,59)  = 0.000;  (*t)(2,59)  =   0.000;
1989       (*t)(0,60)  = 5.900;  (*t)(1,60)  = 0.000;  (*t)(2,60)  =   0.000;
1990     }
1991     delete fTOFpidCuts;
1992     fTOFpidCuts=t;
1993   }
1994 }
1995
1996 //-----------------------------------------------------------------------
1997 Bool_t AliFlowTrackCuts::PassesTPCbayesianCut(const AliESDtrack* track)
1998 {
1999   //cut on TPC bayesian pid
2000   //TODO: maybe join all bayesian methods, make GetESDPdg aware of pid mode selected
2001   //if (! fAllowTOFmismatchFlag) {if (track->GetStatus() & AliESDtrack::kTOFmismatch) return kFALSE;}
2002
2003   //Bool_t statusMatchingHard = TPCTOFagree(track);
2004   //if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
2005   //     return kFALSE;
2006
2007   Int_t pdg = GetESDPdg(track,"bayesianTPC");
2008   //  printf("pdg set to %i\n",pdg);
2009
2010   Int_t pid = 0;
2011   Float_t prob = 0;
2012   switch (fParticleID)
2013   {
2014     case AliPID::kPion:
2015       pid=211;
2016       prob = fProbBayes[2];
2017       break;
2018     case AliPID::kKaon:
2019       pid=321;
2020       prob = fProbBayes[3];
2021      break;
2022     case AliPID::kProton:
2023       pid=2212;
2024       prob = fProbBayes[4];
2025       break;
2026     case AliPID::kElectron:
2027       pid=-11;
2028        prob = fProbBayes[0];
2029      break;
2030     default:
2031       return kFALSE;
2032   }
2033
2034   if(TMath::Abs(pdg) == TMath::Abs(pid) && prob > fParticleProbability)
2035   {
2036     if(!fCutCharge)
2037       return kTRUE;
2038     else if (fCutCharge && fCharge * track->GetSign() > 0)
2039       return kTRUE;
2040   }
2041   return kFALSE;
2042 }
2043
2044 //-----------------------------------------------------------------------
2045 // part added by F. Noferini (some methods)
2046 Bool_t AliFlowTrackCuts::PassesTOFbayesianCut(const AliESDtrack* track)
2047 {
2048
2049   Bool_t goodtrack = (track->GetStatus() & AliESDtrack::kTOFpid) &&
2050                      (track->GetTOFsignal() > 12000) &&
2051                      (track->GetTOFsignal() < 100000) &&
2052                      (track->GetIntegratedLength() > 365);
2053
2054   if (! goodtrack)
2055        return kFALSE;
2056
2057   if (! fAllowTOFmismatchFlag) {if (track->GetStatus() & AliESDtrack::kTOFmismatch) return kFALSE;}
2058
2059   Bool_t statusMatchingHard = TPCTOFagree(track);
2060   if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
2061        return kFALSE;
2062
2063   Int_t pdg = GetESDPdg(track,"bayesianALL");
2064   //  printf("pdg set to %i\n",pdg);
2065
2066   Int_t pid = 0;
2067   Float_t prob = 0;
2068   switch (fParticleID)
2069   {
2070     case AliPID::kPion:
2071       pid=211;
2072       prob = fProbBayes[2];
2073       break;
2074     case AliPID::kKaon:
2075       pid=321;
2076       prob = fProbBayes[3];
2077      break;
2078     case AliPID::kProton:
2079       pid=2212;
2080       prob = fProbBayes[4];
2081       break;
2082     case AliPID::kElectron:
2083       pid=-11;
2084        prob = fProbBayes[0];
2085      break;
2086     default:
2087       return kFALSE;
2088   }
2089
2090   //  printf("pt = %f -- all prob = [%4.2f,%4.2f,%4.2f,%4.2f,%4.2f] -- prob = %f\n",track->Pt(),fProbBayes[0],fProbBayes[1],fProbBayes[2],fProbBayes[3],fProbBayes[4],prob);
2091   if(TMath::Abs(pdg) == TMath::Abs(pid) && prob > fParticleProbability){
2092     if(!fCutCharge)
2093       return kTRUE;
2094     else if (fCutCharge && fCharge * track->GetSign() > 0)
2095       return kTRUE;
2096   }
2097   return kFALSE;
2098 }
2099
2100
2101 //-----------------------------------------------------------------------
2102  // part added by Natasha
2103 Bool_t AliFlowTrackCuts::PassesNucleiSelection(const AliESDtrack* track)
2104 {
2105   Bool_t select=kFALSE;
2106
2107   //if (!track) continue;
2108
2109   if (!track->GetInnerParam())
2110     return kFALSE;    //break;
2111
2112   const AliExternalTrackParam* tpcTrack = track->GetInnerParam();
2113
2114   Double_t ptotTPC = tpcTrack->GetP();
2115   Double_t sigTPC = track->GetTPCsignal();
2116   Double_t dEdxBBA = 0.;
2117   Double_t dSigma = 0.;
2118
2119   switch (fParticleID)
2120     {
2121     case AliPID::kDeuteron:
2122       //pid=10;
2123       dEdxBBA = AliExternalTrackParam::BetheBlochAleph(ptotTPC/1.8756,
2124                                                        4.60e+00,
2125                                                        8.9684e+00,
2126                                                        1.640e-05,
2127                                                        2.35e+00,
2128                                                        2.35e+00);
2129       dSigma = (sigTPC -  dEdxBBA)/dEdxBBA;
2130
2131       if( ptotTPC<=1.1 && (dSigma < (0.5 - (0.1818*ptotTPC)) ) && (dSigma > ( (0.218*ptotTPC - 0.4) ) )  )
2132         {select=kTRUE;}
2133       break;
2134
2135     case AliPID::kTriton:
2136       //pid=11;
2137       select=kFALSE;
2138       break;
2139
2140     case AliPID::kHe3:
2141       //pid=12;
2142       select=kFALSE;
2143       break;
2144
2145     case AliPID::kAlpha:
2146       //pid=13;
2147       select=kFALSE;
2148       break;
2149
2150     default:
2151       return kFALSE;
2152     }
2153
2154   return select;
2155 }
2156 // end part added by Natasha
2157
2158
2159
2160 //-----------------------------------------------------------------------
2161 Int_t AliFlowTrackCuts::GetESDPdg(const AliESDtrack *track,Option_t *option,Int_t ipart,Float_t cPi,Float_t cKa,Float_t cPr)
2162 {
2163   //Get ESD Pdg
2164   Int_t pdg = 0;
2165   Int_t pdgvalues[5] = {-11,-13,211,321,2212};
2166   Float_t mass[5] = {5.10998909999999971e-04,1.05658000000000002e-01,1.39570000000000000e-01,4.93676999999999977e-01,9.38271999999999995e-01};
2167
2168   if(strstr(option,"bayesianTOF")){ // Bayesian TOF PID
2169     Double_t c[5]={0.01, 0.01, 0.85, 0.1, 0.05};
2170     Double_t rcc=0.;
2171
2172     Float_t pt = track->Pt();
2173
2174     Int_t iptesd = 0;
2175     while(pt > fBinLimitPID[iptesd] && iptesd < fgkPIDptBin-1) iptesd++;
2176
2177     if(cPi < 0){
2178       c[0] = fC[iptesd][0];
2179       c[1] = fC[iptesd][1];
2180       c[2] = fC[iptesd][2];
2181       c[3] = fC[iptesd][3];
2182       c[4] = fC[iptesd][4];
2183     }
2184     else{
2185       c[0] = 0.0;
2186       c[1] = 0.0;
2187       c[2] = cPi;
2188       c[3] = cKa;
2189       c[4] = cPr;
2190     }
2191
2192     Double_t r1[10]; track->GetTOFpid(r1);
2193
2194     Int_t i;
2195     for (i=0; i<5; i++) rcc+=(c[i]*r1[i]);
2196
2197     Double_t w[10];
2198     for (i=0; i<5; i++){
2199         w[i]=c[i]*r1[i]/rcc;
2200         fProbBayes[i] = w[i];
2201     }
2202     if (w[2]>=w[3] && w[2]>=w[4] && w[2]>=w[1] && w[2]>=w[0]) {//pion
2203       pdg = 211*Int_t(track->GetSign());
2204     }
2205     else if (w[4]>=w[3] && w[4]>=w[1] && w[4]>=w[0]) {//proton
2206       pdg = 2212*Int_t(track->GetSign());
2207     }
2208     else if (w[3]>=w[1] && w[3]>=w[0]){//kaon
2209       pdg = 321*Int_t(track->GetSign());
2210     }
2211     else if (w[0]>=w[1]) { //electrons
2212       pdg = -11*Int_t(track->GetSign());
2213     }
2214     else{ // muon
2215       pdg = -13*Int_t(track->GetSign());
2216     }
2217   }
2218
2219   else if(strstr(option,"bayesianTPC")){ // Bayesian TPC PID
2220     Double_t c[5]={0.01, 0.01, 0.85, 0.1, 0.05};
2221     Double_t rcc=0.;
2222
2223     Float_t pt = track->Pt();
2224
2225     Int_t iptesd = 0;
2226     while(pt > fBinLimitPID[iptesd] && iptesd < fgkPIDptBin-1) iptesd++;
2227
2228     if(cPi < 0){
2229       c[0] = fC[iptesd][0];
2230       c[1] = fC[iptesd][1];
2231       c[2] = fC[iptesd][2];
2232       c[3] = fC[iptesd][3];
2233       c[4] = fC[iptesd][4];
2234     }
2235     else{
2236       c[0] = 0.0;
2237       c[1] = 0.0;
2238       c[2] = cPi;
2239       c[3] = cKa;
2240       c[4] = cPr;
2241     }
2242
2243     Double_t r1[10]; track->GetTPCpid(r1);
2244
2245     Int_t i;
2246     for (i=0; i<5; i++) rcc+=(c[i]*r1[i]);
2247
2248     Double_t w[10];
2249     for (i=0; i<5; i++){
2250         w[i]=c[i]*r1[i]/rcc;
2251         fProbBayes[i] = w[i];
2252     }
2253     if (w[2]>=w[3] && w[2]>=w[4] && w[2]>=w[1] && w[2]>=w[0]) {//pion
2254       pdg = 211*Int_t(track->GetSign());
2255     }
2256     else if (w[4]>=w[3] && w[4]>=w[1] && w[4]>=w[0]) {//proton
2257       pdg = 2212*Int_t(track->GetSign());
2258     }
2259     else if (w[3]>=w[1] && w[3]>=w[0]){//kaon
2260       pdg = 321*Int_t(track->GetSign());
2261     }
2262     else if (w[0]>=w[1]) { //electrons
2263       pdg = -11*Int_t(track->GetSign());
2264     }
2265     else{ // muon
2266       pdg = -13*Int_t(track->GetSign());
2267     }
2268   }
2269
2270   else if(strstr(option,"bayesianALL")){
2271     Double_t c[5]={0.01, 0.01, 0.85, 0.1, 0.05};
2272     Double_t rcc=0.;
2273
2274     Float_t pt = track->Pt();
2275
2276     Int_t iptesd = 0;
2277     while(pt > fBinLimitPID[iptesd] && iptesd < fgkPIDptBin-1) iptesd++;
2278
2279     if(cPi < 0){
2280       c[0] = fC[iptesd][0];
2281       c[1] = fC[iptesd][1];
2282       c[2] = fC[iptesd][2];
2283       c[3] = fC[iptesd][3];
2284       c[4] = fC[iptesd][4];
2285     }
2286     else{
2287       c[0] = 0.0;
2288       c[1] = 0.0;
2289       c[2] = cPi;
2290       c[3] = cKa;
2291       c[4] = cPr;
2292     }
2293
2294     Double_t r1[10]; track->GetTOFpid(r1);
2295     Double_t r2[10]; track->GetTPCpid(r2);
2296
2297     r1[0] = TMath::Min(r1[2],r1[0]);
2298     r1[1] = TMath::Min(r1[2],r1[1]);
2299
2300     Int_t i;
2301     for (i=0; i<5; i++) rcc+=(c[i]*r1[i]*r2[i]);
2302
2303
2304     Double_t w[10];
2305     for (i=0; i<5; i++){
2306         w[i]=c[i]*r1[i]*r2[i]/rcc;
2307         fProbBayes[i] = w[i];
2308     }
2309
2310     if (w[2]>=w[3] && w[2]>=w[4] && w[2]>=w[1] && w[2]>=w[0]) {//pion
2311       pdg = 211*Int_t(track->GetSign());
2312     }
2313     else if (w[4]>=w[3] && w[4]>=w[1] && w[4]>=w[0]) {//proton
2314       pdg = 2212*Int_t(track->GetSign());
2315     }
2316     else if (w[3]>=w[1] && w[3]>=w[0]){//kaon
2317       pdg = 321*Int_t(track->GetSign());
2318     }
2319     else if (w[0]>=w[1]) { //electrons
2320       pdg = -11*Int_t(track->GetSign());
2321     }
2322     else{ // muon
2323       pdg = -13*Int_t(track->GetSign());
2324     }
2325   }
2326
2327   else if(strstr(option,"sigmacutTOF")){
2328     printf("PID not implemented yet: %s\nNO PID!!!!\n",option);
2329     Float_t p = track->P();
2330
2331     // Take expected times
2332     Double_t exptimes[5];
2333     track->GetIntegratedTimes(exptimes);
2334
2335     // Take resolution for TOF response
2336     // like fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[ipart], mass[ipart]);
2337     Float_t resolution = fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[ipart], mass[ipart]);
2338
2339     if(TMath::Abs(exptimes[ipart] - track->GetTOFsignal()) < 3 * resolution){
2340       pdg = pdgvalues[ipart] * Int_t(track->GetSign());
2341     }
2342   }
2343
2344   else{
2345     printf("Invalid PID option: %s\nNO PID!!!!\n",option);
2346   }
2347
2348   return pdg;
2349 }
2350
2351 //-----------------------------------------------------------------------
2352 void AliFlowTrackCuts::SetPriors(Float_t centrCur){
2353   fBinLimitPID[0] = 0.300000;
2354   fBinLimitPID[1] = 0.400000;
2355   fBinLimitPID[2] = 0.500000;
2356   fBinLimitPID[3] = 0.600000;
2357   fBinLimitPID[4] = 0.700000;
2358   fBinLimitPID[5] = 0.800000;
2359   fBinLimitPID[6] = 0.900000;
2360   fBinLimitPID[7] = 1.000000;
2361   fBinLimitPID[8] = 1.200000;
2362   fBinLimitPID[9] = 1.400000;
2363   fBinLimitPID[10] = 1.600000;
2364   fBinLimitPID[11] = 1.800000;
2365   fBinLimitPID[12] = 2.000000;
2366   fBinLimitPID[13] = 2.200000;
2367   fBinLimitPID[14] = 2.400000;
2368   fBinLimitPID[15] = 2.600000;
2369   fBinLimitPID[16] = 2.800000;
2370   fBinLimitPID[17] = 3.000000;
2371
2372   // 0-10%
2373   if(centrCur < 10){
2374       fC[0][0] = 0.005;
2375       fC[0][1] = 0.005;
2376       fC[0][2] = 1.0000;
2377       fC[0][3] = 0.0010;
2378       fC[0][4] = 0.0010;
2379
2380       fC[1][0] = 0.005;
2381       fC[1][1] = 0.005;
2382       fC[1][2] = 1.0000;
2383       fC[1][3] = 0.0168;
2384       fC[1][4] = 0.0122;
2385
2386       fC[2][0] = 0.005;
2387       fC[2][1] = 0.005;
2388       fC[2][2] = 1.0000;
2389       fC[2][3] = 0.0272;
2390       fC[2][4] = 0.0070;
2391
2392       fC[3][0] = 0.005;
2393       fC[3][1] = 0.005;
2394       fC[3][2] = 1.0000;
2395       fC[3][3] = 0.0562;
2396       fC[3][4] = 0.0258;
2397
2398       fC[4][0] = 0.005;
2399       fC[4][1] = 0.005;
2400       fC[4][2] = 1.0000;
2401       fC[4][3] = 0.0861;
2402       fC[4][4] = 0.0496;
2403
2404       fC[5][0] = 0.005;
2405       fC[5][1] = 0.005;
2406       fC[5][2] = 1.0000;
2407       fC[5][3] = 0.1168;
2408       fC[5][4] = 0.0740;
2409
2410       fC[6][0] = 0.005;
2411       fC[6][1] = 0.005;
2412       fC[6][2] = 1.0000;
2413       fC[6][3] = 0.1476;
2414       fC[6][4] = 0.0998;
2415
2416       fC[7][0] = 0.005;
2417       fC[7][1] = 0.005;
2418       fC[7][2] = 1.0000;
2419       fC[7][3] = 0.1810;
2420       fC[7][4] = 0.1296;
2421
2422       fC[8][0] = 0.005;
2423       fC[8][1] = 0.005;
2424       fC[8][2] = 1.0000;
2425       fC[8][3] = 0.2240;
2426       fC[8][4] = 0.1827;
2427
2428       fC[9][0] = 0.005;
2429       fC[9][1] = 0.005;
2430       fC[9][2] = 1.0000;
2431       fC[9][3] = 0.2812;
2432       fC[9][4] = 0.2699;
2433
2434       fC[10][0] = 0.005;
2435       fC[10][1] = 0.005;
2436       fC[10][2] = 1.0000;
2437       fC[10][3] = 0.3328;
2438       fC[10][4] = 0.3714;
2439
2440       fC[11][0] = 0.005;
2441       fC[11][1] = 0.005;
2442       fC[11][2] = 1.0000;
2443       fC[11][3] = 0.3780;
2444       fC[11][4] = 0.4810;
2445
2446       fC[12][0] = 0.005;
2447       fC[12][1] = 0.005;
2448       fC[12][2] = 1.0000;
2449       fC[12][3] = 0.4125;
2450       fC[12][4] = 0.5771;
2451
2452       fC[13][0] = 0.005;
2453       fC[13][1] = 0.005;
2454       fC[13][2] = 1.0000;
2455       fC[13][3] = 0.4486;
2456       fC[13][4] = 0.6799;
2457
2458       fC[14][0] = 0.005;
2459       fC[14][1] = 0.005;
2460       fC[14][2] = 1.0000;
2461       fC[14][3] = 0.4840;
2462       fC[14][4] = 0.7668;
2463
2464       fC[15][0] = 0.005;
2465       fC[15][1] = 0.005;
2466       fC[15][2] = 1.0000;
2467       fC[15][3] = 0.4971;
2468       fC[15][4] = 0.8288;
2469
2470       fC[16][0] = 0.005;
2471       fC[16][1] = 0.005;
2472       fC[16][2] = 1.0000;
2473       fC[16][3] = 0.4956;
2474       fC[16][4] = 0.8653;
2475
2476       fC[17][0] = 0.005;
2477       fC[17][1] = 0.005;
2478       fC[17][2] = 1.0000;
2479       fC[17][3] = 0.5173;
2480       fC[17][4] = 0.9059;
2481   }
2482   // 10-20%
2483   else if(centrCur < 20){
2484      fC[0][0] = 0.005;
2485       fC[0][1] = 0.005;
2486       fC[0][2] = 1.0000;
2487       fC[0][3] = 0.0010;
2488       fC[0][4] = 0.0010;
2489
2490       fC[1][0] = 0.005;
2491       fC[1][1] = 0.005;
2492       fC[1][2] = 1.0000;
2493       fC[1][3] = 0.0132;
2494       fC[1][4] = 0.0088;
2495
2496       fC[2][0] = 0.005;
2497       fC[2][1] = 0.005;
2498       fC[2][2] = 1.0000;
2499       fC[2][3] = 0.0283;
2500       fC[2][4] = 0.0068;
2501
2502       fC[3][0] = 0.005;
2503       fC[3][1] = 0.005;
2504       fC[3][2] = 1.0000;
2505       fC[3][3] = 0.0577;
2506       fC[3][4] = 0.0279;
2507
2508       fC[4][0] = 0.005;
2509       fC[4][1] = 0.005;
2510       fC[4][2] = 1.0000;
2511       fC[4][3] = 0.0884;
2512       fC[4][4] = 0.0534;
2513
2514       fC[5][0] = 0.005;
2515       fC[5][1] = 0.005;
2516       fC[5][2] = 1.0000;
2517       fC[5][3] = 0.1179;
2518       fC[5][4] = 0.0794;
2519
2520       fC[6][0] = 0.005;
2521       fC[6][1] = 0.005;
2522       fC[6][2] = 1.0000;
2523       fC[6][3] = 0.1480;
2524       fC[6][4] = 0.1058;
2525
2526       fC[7][0] = 0.005;
2527       fC[7][1] = 0.005;
2528       fC[7][2] = 1.0000;
2529       fC[7][3] = 0.1807;
2530       fC[7][4] = 0.1366;
2531
2532       fC[8][0] = 0.005;
2533       fC[8][1] = 0.005;
2534       fC[8][2] = 1.0000;
2535       fC[8][3] = 0.2219;
2536       fC[8][4] = 0.1891;
2537
2538       fC[9][0] = 0.005;
2539       fC[9][1] = 0.005;
2540       fC[9][2] = 1.0000;
2541       fC[9][3] = 0.2804;
2542       fC[9][4] = 0.2730;
2543
2544       fC[10][0] = 0.005;
2545       fC[10][1] = 0.005;
2546       fC[10][2] = 1.0000;
2547       fC[10][3] = 0.3283;
2548       fC[10][4] = 0.3660;
2549
2550       fC[11][0] = 0.005;
2551       fC[11][1] = 0.005;
2552       fC[11][2] = 1.0000;
2553       fC[11][3] = 0.3710;
2554       fC[11][4] = 0.4647;
2555
2556       fC[12][0] = 0.005;
2557       fC[12][1] = 0.005;
2558       fC[12][2] = 1.0000;
2559       fC[12][3] = 0.4093;
2560       fC[12][4] = 0.5566;
2561
2562       fC[13][0] = 0.005;
2563       fC[13][1] = 0.005;
2564       fC[13][2] = 1.0000;
2565       fC[13][3] = 0.4302;
2566       fC[13][4] = 0.6410;
2567
2568       fC[14][0] = 0.005;
2569       fC[14][1] = 0.005;
2570       fC[14][2] = 1.0000;
2571       fC[14][3] = 0.4649;
2572       fC[14][4] = 0.7055;
2573
2574       fC[15][0] = 0.005;
2575       fC[15][1] = 0.005;
2576       fC[15][2] = 1.0000;
2577       fC[15][3] = 0.4523;
2578       fC[15][4] = 0.7440;
2579
2580       fC[16][0] = 0.005;
2581       fC[16][1] = 0.005;
2582       fC[16][2] = 1.0000;
2583       fC[16][3] = 0.4591;
2584       fC[16][4] = 0.7799;
2585
2586       fC[17][0] = 0.005;
2587       fC[17][1] = 0.005;
2588       fC[17][2] = 1.0000;
2589       fC[17][3] = 0.4804;
2590       fC[17][4] = 0.8218;
2591   }
2592   // 20-30%
2593   else if(centrCur < 30){
2594      fC[0][0] = 0.005;
2595       fC[0][1] = 0.005;
2596       fC[0][2] = 1.0000;
2597       fC[0][3] = 0.0010;
2598       fC[0][4] = 0.0010;
2599
2600       fC[1][0] = 0.005;
2601       fC[1][1] = 0.005;
2602       fC[1][2] = 1.0000;
2603       fC[1][3] = 0.0102;
2604       fC[1][4] = 0.0064;
2605
2606       fC[2][0] = 0.005;
2607       fC[2][1] = 0.005;
2608       fC[2][2] = 1.0000;
2609       fC[2][3] = 0.0292;
2610       fC[2][4] = 0.0066;
2611
2612       fC[3][0] = 0.005;
2613       fC[3][1] = 0.005;
2614       fC[3][2] = 1.0000;
2615       fC[3][3] = 0.0597;
2616       fC[3][4] = 0.0296;
2617
2618       fC[4][0] = 0.005;
2619       fC[4][1] = 0.005;
2620       fC[4][2] = 1.0000;
2621       fC[4][3] = 0.0900;
2622       fC[4][4] = 0.0589;
2623
2624       fC[5][0] = 0.005;
2625       fC[5][1] = 0.005;
2626       fC[5][2] = 1.0000;
2627       fC[5][3] = 0.1199;
2628       fC[5][4] = 0.0859;
2629
2630       fC[6][0] = 0.005;
2631       fC[6][1] = 0.005;
2632       fC[6][2] = 1.0000;
2633       fC[6][3] = 0.1505;
2634       fC[6][4] = 0.1141;
2635
2636       fC[7][0] = 0.005;
2637       fC[7][1] = 0.005;
2638       fC[7][2] = 1.0000;
2639       fC[7][3] = 0.1805;
2640       fC[7][4] = 0.1454;
2641
2642       fC[8][0] = 0.005;
2643       fC[8][1] = 0.005;
2644       fC[8][2] = 1.0000;
2645       fC[8][3] = 0.2221;
2646       fC[8][4] = 0.2004;
2647
2648       fC[9][0] = 0.005;
2649       fC[9][1] = 0.005;
2650       fC[9][2] = 1.0000;
2651       fC[9][3] = 0.2796;
2652       fC[9][4] = 0.2838;
2653
2654       fC[10][0] = 0.005;
2655       fC[10][1] = 0.005;
2656       fC[10][2] = 1.0000;
2657       fC[10][3] = 0.3271;
2658       fC[10][4] = 0.3682;
2659
2660       fC[11][0] = 0.005;
2661       fC[11][1] = 0.005;
2662       fC[11][2] = 1.0000;
2663       fC[11][3] = 0.3648;
2664       fC[11][4] = 0.4509;
2665
2666       fC[12][0] = 0.005;
2667       fC[12][1] = 0.005;
2668       fC[12][2] = 1.0000;
2669       fC[12][3] = 0.3988;
2670       fC[12][4] = 0.5339;
2671
2672       fC[13][0] = 0.005;
2673       fC[13][1] = 0.005;
2674       fC[13][2] = 1.0000;
2675       fC[13][3] = 0.4315;
2676       fC[13][4] = 0.5995;
2677
2678       fC[14][0] = 0.005;
2679       fC[14][1] = 0.005;
2680       fC[14][2] = 1.0000;
2681       fC[14][3] = 0.4548;
2682       fC[14][4] = 0.6612;
2683
2684       fC[15][0] = 0.005;
2685       fC[15][1] = 0.005;
2686       fC[15][2] = 1.0000;
2687       fC[15][3] = 0.4744;
2688       fC[15][4] = 0.7060;
2689
2690       fC[16][0] = 0.005;
2691       fC[16][1] = 0.005;
2692       fC[16][2] = 1.0000;
2693       fC[16][3] = 0.4899;
2694       fC[16][4] = 0.7388;
2695
2696       fC[17][0] = 0.005;
2697       fC[17][1] = 0.005;
2698       fC[17][2] = 1.0000;
2699       fC[17][3] = 0.4411;
2700       fC[17][4] = 0.7293;
2701   }
2702   // 30-40%
2703   else if(centrCur < 40){
2704       fC[0][0] = 0.005;
2705       fC[0][1] = 0.005;
2706       fC[0][2] = 1.0000;
2707       fC[0][3] = 0.0010;
2708       fC[0][4] = 0.0010;
2709
2710       fC[1][0] = 0.005;
2711       fC[1][1] = 0.005;
2712       fC[1][2] = 1.0000;
2713       fC[1][3] = 0.0102;
2714       fC[1][4] = 0.0048;
2715
2716       fC[2][0] = 0.005;
2717       fC[2][1] = 0.005;
2718       fC[2][2] = 1.0000;
2719       fC[2][3] = 0.0306;
2720       fC[2][4] = 0.0079;
2721
2722       fC[3][0] = 0.005;
2723       fC[3][1] = 0.005;
2724       fC[3][2] = 1.0000;
2725       fC[3][3] = 0.0617;
2726       fC[3][4] = 0.0338;
2727
2728       fC[4][0] = 0.005;
2729       fC[4][1] = 0.005;
2730       fC[4][2] = 1.0000;
2731       fC[4][3] = 0.0920;
2732       fC[4][4] = 0.0652;
2733
2734       fC[5][0] = 0.005;
2735       fC[5][1] = 0.005;
2736       fC[5][2] = 1.0000;
2737       fC[5][3] = 0.1211;
2738       fC[5][4] = 0.0955;
2739
2740       fC[6][0] = 0.005;
2741       fC[6][1] = 0.005;
2742       fC[6][2] = 1.0000;
2743       fC[6][3] = 0.1496;
2744       fC[6][4] = 0.1242;
2745
2746       fC[7][0] = 0.005;
2747       fC[7][1] = 0.005;
2748       fC[7][2] = 1.0000;
2749       fC[7][3] = 0.1807;
2750       fC[7][4] = 0.1576;
2751
2752       fC[8][0] = 0.005;
2753       fC[8][1] = 0.005;
2754       fC[8][2] = 1.0000;
2755       fC[8][3] = 0.2195;
2756       fC[8][4] = 0.2097;
2757
2758       fC[9][0] = 0.005;
2759       fC[9][1] = 0.005;
2760       fC[9][2] = 1.0000;
2761       fC[9][3] = 0.2732;
2762       fC[9][4] = 0.2884;
2763
2764       fC[10][0] = 0.005;
2765       fC[10][1] = 0.005;
2766       fC[10][2] = 1.0000;
2767       fC[10][3] = 0.3204;
2768       fC[10][4] = 0.3679;
2769
2770       fC[11][0] = 0.005;
2771       fC[11][1] = 0.005;
2772       fC[11][2] = 1.0000;
2773       fC[11][3] = 0.3564;
2774       fC[11][4] = 0.4449;
2775
2776       fC[12][0] = 0.005;
2777       fC[12][1] = 0.005;
2778       fC[12][2] = 1.0000;
2779       fC[12][3] = 0.3791;
2780       fC[12][4] = 0.5052;
2781
2782       fC[13][0] = 0.005;
2783       fC[13][1] = 0.005;
2784       fC[13][2] = 1.0000;
2785       fC[13][3] = 0.4062;
2786       fC[13][4] = 0.5647;
2787
2788       fC[14][0] = 0.005;
2789       fC[14][1] = 0.005;
2790       fC[14][2] = 1.0000;
2791       fC[14][3] = 0.4234;
2792       fC[14][4] = 0.6203;
2793
2794       fC[15][0] = 0.005;
2795       fC[15][1] = 0.005;
2796       fC[15][2] = 1.0000;
2797       fC[15][3] = 0.4441;
2798       fC[15][4] = 0.6381;
2799
2800       fC[16][0] = 0.005;
2801       fC[16][1] = 0.005;
2802       fC[16][2] = 1.0000;
2803       fC[16][3] = 0.4629;
2804       fC[16][4] = 0.6496;
2805
2806       fC[17][0] = 0.005;
2807       fC[17][1] = 0.005;
2808       fC[17][2] = 1.0000;
2809       fC[17][3] = 0.4293;
2810       fC[17][4] = 0.6491;
2811   }
2812   // 40-50%
2813   else if(centrCur < 50){
2814       fC[0][0] = 0.005;
2815       fC[0][1] = 0.005;
2816       fC[0][2] = 1.0000;
2817       fC[0][3] = 0.0010;
2818       fC[0][4] = 0.0010;
2819
2820       fC[1][0] = 0.005;
2821       fC[1][1] = 0.005;
2822       fC[1][2] = 1.0000;
2823       fC[1][3] = 0.0093;
2824       fC[1][4] = 0.0057;
2825
2826       fC[2][0] = 0.005;
2827       fC[2][1] = 0.005;
2828       fC[2][2] = 1.0000;
2829       fC[2][3] = 0.0319;
2830       fC[2][4] = 0.0075;
2831
2832       fC[3][0] = 0.005;
2833       fC[3][1] = 0.005;
2834       fC[3][2] = 1.0000;
2835       fC[3][3] = 0.0639;
2836       fC[3][4] = 0.0371;
2837
2838       fC[4][0] = 0.005;
2839       fC[4][1] = 0.005;
2840       fC[4][2] = 1.0000;
2841       fC[4][3] = 0.0939;
2842       fC[4][4] = 0.0725;
2843
2844       fC[5][0] = 0.005;
2845       fC[5][1] = 0.005;
2846       fC[5][2] = 1.0000;
2847       fC[5][3] = 0.1224;
2848       fC[5][4] = 0.1045;
2849
2850       fC[6][0] = 0.005;
2851       fC[6][1] = 0.005;
2852       fC[6][2] = 1.0000;
2853       fC[6][3] = 0.1520;
2854       fC[6][4] = 0.1387;
2855
2856       fC[7][0] = 0.005;
2857       fC[7][1] = 0.005;
2858       fC[7][2] = 1.0000;
2859       fC[7][3] = 0.1783;
2860       fC[7][4] = 0.1711;
2861
2862       fC[8][0] = 0.005;
2863       fC[8][1] = 0.005;
2864       fC[8][2] = 1.0000;
2865       fC[8][3] = 0.2202;
2866       fC[8][4] = 0.2269;
2867
2868       fC[9][0] = 0.005;
2869       fC[9][1] = 0.005;
2870       fC[9][2] = 1.0000;
2871       fC[9][3] = 0.2672;
2872       fC[9][4] = 0.2955;
2873
2874       fC[10][0] = 0.005;
2875       fC[10][1] = 0.005;
2876       fC[10][2] = 1.0000;
2877       fC[10][3] = 0.3191;
2878       fC[10][4] = 0.3676;
2879
2880       fC[11][0] = 0.005;
2881       fC[11][1] = 0.005;
2882       fC[11][2] = 1.0000;
2883       fC[11][3] = 0.3434;
2884       fC[11][4] = 0.4321;
2885
2886       fC[12][0] = 0.005;
2887       fC[12][1] = 0.005;
2888       fC[12][2] = 1.0000;
2889       fC[12][3] = 0.3692;
2890       fC[12][4] = 0.4879;
2891
2892       fC[13][0] = 0.005;
2893       fC[13][1] = 0.005;
2894       fC[13][2] = 1.0000;
2895       fC[13][3] = 0.3993;
2896       fC[13][4] = 0.5377;
2897
2898       fC[14][0] = 0.005;
2899       fC[14][1] = 0.005;
2900       fC[14][2] = 1.0000;
2901       fC[14][3] = 0.3818;
2902       fC[14][4] = 0.5547;
2903
2904       fC[15][0] = 0.005;
2905       fC[15][1] = 0.005;
2906       fC[15][2] = 1.0000;
2907       fC[15][3] = 0.4003;
2908       fC[15][4] = 0.5484;
2909
2910       fC[16][0] = 0.005;
2911       fC[16][1] = 0.005;
2912       fC[16][2] = 1.0000;
2913       fC[16][3] = 0.4281;
2914       fC[16][4] = 0.5383;
2915
2916       fC[17][0] = 0.005;
2917       fC[17][1] = 0.005;
2918       fC[17][2] = 1.0000;
2919       fC[17][3] = 0.3960;
2920       fC[17][4] = 0.5374;
2921   }
2922   // 50-60%
2923   else if(centrCur < 60){
2924       fC[0][0] = 0.005;
2925       fC[0][1] = 0.005;
2926       fC[0][2] = 1.0000;
2927       fC[0][3] = 0.0010;
2928       fC[0][4] = 0.0010;
2929
2930       fC[1][0] = 0.005;
2931       fC[1][1] = 0.005;
2932       fC[1][2] = 1.0000;
2933       fC[1][3] = 0.0076;
2934       fC[1][4] = 0.0032;
2935
2936       fC[2][0] = 0.005;
2937       fC[2][1] = 0.005;
2938       fC[2][2] = 1.0000;
2939       fC[2][3] = 0.0329;
2940       fC[2][4] = 0.0085;
2941
2942       fC[3][0] = 0.005;
2943       fC[3][1] = 0.005;
2944       fC[3][2] = 1.0000;
2945       fC[3][3] = 0.0653;
2946       fC[3][4] = 0.0423;
2947
2948       fC[4][0] = 0.005;
2949       fC[4][1] = 0.005;
2950       fC[4][2] = 1.0000;
2951       fC[4][3] = 0.0923;
2952       fC[4][4] = 0.0813;
2953
2954       fC[5][0] = 0.005;
2955       fC[5][1] = 0.005;
2956       fC[5][2] = 1.0000;
2957       fC[5][3] = 0.1219;
2958       fC[5][4] = 0.1161;
2959
2960       fC[6][0] = 0.005;
2961       fC[6][1] = 0.005;
2962       fC[6][2] = 1.0000;
2963       fC[6][3] = 0.1519;
2964       fC[6][4] = 0.1520;
2965
2966       fC[7][0] = 0.005;
2967       fC[7][1] = 0.005;
2968       fC[7][2] = 1.0000;
2969       fC[7][3] = 0.1763;
2970       fC[7][4] = 0.1858;
2971
2972       fC[8][0] = 0.005;
2973       fC[8][1] = 0.005;
2974       fC[8][2] = 1.0000;
2975       fC[8][3] = 0.2178;
2976       fC[8][4] = 0.2385;
2977
2978       fC[9][0] = 0.005;
2979       fC[9][1] = 0.005;
2980       fC[9][2] = 1.0000;
2981       fC[9][3] = 0.2618;
2982       fC[9][4] = 0.3070;
2983
2984       fC[10][0] = 0.005;
2985       fC[10][1] = 0.005;
2986       fC[10][2] = 1.0000;
2987       fC[10][3] = 0.3067;
2988       fC[10][4] = 0.3625;
2989
2990       fC[11][0] = 0.005;
2991       fC[11][1] = 0.005;
2992       fC[11][2] = 1.0000;
2993       fC[11][3] = 0.3336;
2994       fC[11][4] = 0.4188;
2995
2996       fC[12][0] = 0.005;
2997       fC[12][1] = 0.005;
2998       fC[12][2] = 1.0000;
2999       fC[12][3] = 0.3706;
3000       fC[12][4] = 0.4511;
3001
3002       fC[13][0] = 0.005;
3003       fC[13][1] = 0.005;
3004       fC[13][2] = 1.0000;
3005       fC[13][3] = 0.3765;
3006       fC[13][4] = 0.4729;
3007
3008       fC[14][0] = 0.005;
3009       fC[14][1] = 0.005;
3010       fC[14][2] = 1.0000;
3011       fC[14][3] = 0.3942;
3012       fC[14][4] = 0.4855;
3013
3014       fC[15][0] = 0.005;
3015       fC[15][1] = 0.005;
3016       fC[15][2] = 1.0000;
3017       fC[15][3] = 0.4051;
3018       fC[15][4] = 0.4762;
3019
3020       fC[16][0] = 0.005;
3021       fC[16][1] = 0.005;
3022       fC[16][2] = 1.0000;
3023       fC[16][3] = 0.3843;
3024       fC[16][4] = 0.4763;
3025
3026       fC[17][0] = 0.005;
3027       fC[17][1] = 0.005;
3028       fC[17][2] = 1.0000;
3029       fC[17][3] = 0.4237;
3030       fC[17][4] = 0.4773;
3031   }
3032   // 60-70%
3033   else if(centrCur < 70){
3034          fC[0][0] = 0.005;
3035       fC[0][1] = 0.005;
3036       fC[0][2] = 1.0000;
3037       fC[0][3] = 0.0010;
3038       fC[0][4] = 0.0010;
3039
3040       fC[1][0] = 0.005;
3041       fC[1][1] = 0.005;
3042       fC[1][2] = 1.0000;
3043       fC[1][3] = 0.0071;
3044       fC[1][4] = 0.0012;
3045
3046       fC[2][0] = 0.005;
3047       fC[2][1] = 0.005;
3048       fC[2][2] = 1.0000;
3049       fC[2][3] = 0.0336;
3050       fC[2][4] = 0.0097;
3051
3052       fC[3][0] = 0.005;
3053       fC[3][1] = 0.005;
3054       fC[3][2] = 1.0000;
3055       fC[3][3] = 0.0662;
3056       fC[3][4] = 0.0460;
3057
3058       fC[4][0] = 0.005;
3059       fC[4][1] = 0.005;
3060       fC[4][2] = 1.0000;
3061       fC[4][3] = 0.0954;
3062       fC[4][4] = 0.0902;
3063
3064       fC[5][0] = 0.005;
3065       fC[5][1] = 0.005;
3066       fC[5][2] = 1.0000;
3067       fC[5][3] = 0.1181;
3068       fC[5][4] = 0.1306;
3069
3070       fC[6][0] = 0.005;
3071       fC[6][1] = 0.005;
3072       fC[6][2] = 1.0000;
3073       fC[6][3] = 0.1481;
3074       fC[6][4] = 0.1662;
3075
3076       fC[7][0] = 0.005;
3077       fC[7][1] = 0.005;
3078       fC[7][2] = 1.0000;
3079       fC[7][3] = 0.1765;
3080       fC[7][4] = 0.1963;
3081
3082       fC[8][0] = 0.005;
3083       fC[8][1] = 0.005;
3084       fC[8][2] = 1.0000;
3085       fC[8][3] = 0.2155;
3086       fC[8][4] = 0.2433;
3087
3088       fC[9][0] = 0.005;
3089       fC[9][1] = 0.005;
3090       fC[9][2] = 1.0000;
3091       fC[9][3] = 0.2580;
3092       fC[9][4] = 0.3022;
3093
3094       fC[10][0] = 0.005;
3095       fC[10][1] = 0.005;
3096       fC[10][2] = 1.0000;
3097       fC[10][3] = 0.2872;
3098       fC[10][4] = 0.3481;
3099
3100       fC[11][0] = 0.005;
3101       fC[11][1] = 0.005;
3102       fC[11][2] = 1.0000;
3103       fC[11][3] = 0.3170;
3104       fC[11][4] = 0.3847;
3105
3106       fC[12][0] = 0.005;
3107       fC[12][1] = 0.005;
3108       fC[12][2] = 1.0000;
3109       fC[12][3] = 0.3454;
3110       fC[12][4] = 0.4258;
3111
3112       fC[13][0] = 0.005;
3113       fC[13][1] = 0.005;
3114       fC[13][2] = 1.0000;
3115       fC[13][3] = 0.3580;
3116       fC[13][4] = 0.4299;
3117
3118       fC[14][0] = 0.005;
3119       fC[14][1] = 0.005;
3120       fC[14][2] = 1.0000;
3121       fC[14][3] = 0.3903;
3122       fC[14][4] = 0.4326;
3123
3124       fC[15][0] = 0.005;
3125       fC[15][1] = 0.005;
3126       fC[15][2] = 1.0000;
3127       fC[15][3] = 0.3690;
3128       fC[15][4] = 0.4491;
3129
3130       fC[16][0] = 0.005;
3131       fC[16][1] = 0.005;
3132       fC[16][2] = 1.0000;
3133       fC[16][3] = 0.4716;
3134       fC[16][4] = 0.4298;
3135
3136       fC[17][0] = 0.005;
3137       fC[17][1] = 0.005;
3138       fC[17][2] = 1.0000;
3139       fC[17][3] = 0.3875;
3140       fC[17][4] = 0.4083;
3141   }
3142   // 70-80%
3143   else if(centrCur < 80){
3144       fC[0][0] = 0.005;
3145       fC[0][1] = 0.005;
3146       fC[0][2] = 1.0000;
3147       fC[0][3] = 0.0010;
3148       fC[0][4] = 0.0010;
3149
3150       fC[1][0] = 0.005;
3151       fC[1][1] = 0.005;
3152       fC[1][2] = 1.0000;
3153       fC[1][3] = 0.0075;
3154       fC[1][4] = 0.0007;
3155
3156       fC[2][0] = 0.005;
3157       fC[2][1] = 0.005;
3158       fC[2][2] = 1.0000;
3159       fC[2][3] = 0.0313;
3160       fC[2][4] = 0.0124;
3161
3162       fC[3][0] = 0.005;
3163       fC[3][1] = 0.005;
3164       fC[3][2] = 1.0000;
3165       fC[3][3] = 0.0640;
3166       fC[3][4] = 0.0539;
3167
3168       fC[4][0] = 0.005;
3169       fC[4][1] = 0.005;
3170       fC[4][2] = 1.0000;
3171       fC[4][3] = 0.0923;
3172       fC[4][4] = 0.0992;
3173
3174       fC[5][0] = 0.005;
3175       fC[5][1] = 0.005;
3176       fC[5][2] = 1.0000;
3177       fC[5][3] = 0.1202;
3178       fC[5][4] = 0.1417;
3179
3180       fC[6][0] = 0.005;
3181       fC[6][1] = 0.005;
3182       fC[6][2] = 1.0000;
3183       fC[6][3] = 0.1413;
3184       fC[6][4] = 0.1729;
3185
3186       fC[7][0] = 0.005;
3187       fC[7][1] = 0.005;
3188       fC[7][2] = 1.0000;
3189       fC[7][3] = 0.1705;
3190       fC[7][4] = 0.1999;
3191
3192       fC[8][0] = 0.005;
3193       fC[8][1] = 0.005;
3194       fC[8][2] = 1.0000;
3195       fC[8][3] = 0.2103;
3196       fC[8][4] = 0.2472;
3197
3198       fC[9][0] = 0.005;
3199       fC[9][1] = 0.005;
3200       fC[9][2] = 1.0000;
3201       fC[9][3] = 0.2373;
3202       fC[9][4] = 0.2916;
3203
3204       fC[10][0] = 0.005;
3205       fC[10][1] = 0.005;
3206       fC[10][2] = 1.0000;
3207       fC[10][3] = 0.2824;
3208       fC[10][4] = 0.3323;
3209
3210       fC[11][0] = 0.005;
3211       fC[11][1] = 0.005;
3212       fC[11][2] = 1.0000;
3213       fC[11][3] = 0.3046;
3214       fC[11][4] = 0.3576;
3215
3216       fC[12][0] = 0.005;
3217       fC[12][1] = 0.005;
3218       fC[12][2] = 1.0000;
3219       fC[12][3] = 0.3585;
3220       fC[12][4] = 0.4003;
3221
3222       fC[13][0] = 0.005;
3223       fC[13][1] = 0.005;
3224       fC[13][2] = 1.0000;
3225       fC[13][3] = 0.3461;
3226       fC[13][4] = 0.3982;
3227
3228       fC[14][0] = 0.005;
3229       fC[14][1] = 0.005;
3230       fC[14][2] = 1.0000;
3231       fC[14][3] = 0.3362;
3232       fC[14][4] = 0.3776;
3233
3234       fC[15][0] = 0.005;
3235       fC[15][1] = 0.005;
3236       fC[15][2] = 1.0000;
3237       fC[15][3] = 0.3071;
3238       fC[15][4] = 0.3500;
3239
3240       fC[16][0] = 0.005;
3241       fC[16][1] = 0.005;
3242       fC[16][2] = 1.0000;
3243       fC[16][3] = 0.2914;
3244       fC[16][4] = 0.3937;
3245
3246       fC[17][0] = 0.005;
3247       fC[17][1] = 0.005;
3248       fC[17][2] = 1.0000;
3249       fC[17][3] = 0.3727;
3250       fC[17][4] = 0.3877;
3251   }
3252   // 80-100%
3253   else{
3254       fC[0][0] = 0.005;
3255       fC[0][1] = 0.005;
3256       fC[0][2] = 1.0000;
3257       fC[0][3] = 0.0010;
3258       fC[0][4] = 0.0010;
3259
3260       fC[1][0] = 0.005;
3261       fC[1][1] = 0.005;
3262       fC[1][2] = 1.0000;
3263       fC[1][3] = 0.0060;
3264       fC[1][4] = 0.0035;
3265
3266       fC[2][0] = 0.005;
3267       fC[2][1] = 0.005;
3268       fC[2][2] = 1.0000;
3269       fC[2][3] = 0.0323;
3270       fC[2][4] = 0.0113;
3271
3272       fC[3][0] = 0.005;
3273       fC[3][1] = 0.005;
3274       fC[3][2] = 1.0000;
3275       fC[3][3] = 0.0609;
3276       fC[3][4] = 0.0653;
3277
3278       fC[4][0] = 0.005;
3279       fC[4][1] = 0.005;
3280       fC[4][2] = 1.0000;
3281       fC[4][3] = 0.0922;
3282       fC[4][4] = 0.1076;
3283
3284       fC[5][0] = 0.005;
3285       fC[5][1] = 0.005;
3286       fC[5][2] = 1.0000;
3287       fC[5][3] = 0.1096;
3288       fC[5][4] = 0.1328;
3289
3290       fC[6][0] = 0.005;
3291       fC[6][1] = 0.005;
3292       fC[6][2] = 1.0000;
3293       fC[6][3] = 0.1495;
3294       fC[6][4] = 0.1779;
3295
3296       fC[7][0] = 0.005;
3297       fC[7][1] = 0.005;
3298       fC[7][2] = 1.0000;
3299       fC[7][3] = 0.1519;
3300       fC[7][4] = 0.1989;
3301
3302       fC[8][0] = 0.005;
3303       fC[8][1] = 0.005;
3304       fC[8][2] = 1.0000;
3305       fC[8][3] = 0.1817;
3306       fC[8][4] = 0.2472;
3307
3308       fC[9][0] = 0.005;
3309       fC[9][1] = 0.005;
3310       fC[9][2] = 1.0000;
3311       fC[9][3] = 0.2429;
3312       fC[9][4] = 0.2684;
3313
3314       fC[10][0] = 0.005;
3315       fC[10][1] = 0.005;
3316       fC[10][2] = 1.0000;
3317       fC[10][3] = 0.2760;
3318       fC[10][4] = 0.3098;
3319
3320       fC[11][0] = 0.005;
3321       fC[11][1] = 0.005;
3322       fC[11][2] = 1.0000;
3323       fC[11][3] = 0.2673;
3324       fC[11][4] = 0.3198;
3325
3326       fC[12][0] = 0.005;
3327       fC[12][1] = 0.005;
3328       fC[12][2] = 1.0000;
3329       fC[12][3] = 0.3165;
3330       fC[12][4] = 0.3564;
3331
3332       fC[13][0] = 0.005;
3333       fC[13][1] = 0.005;
3334       fC[13][2] = 1.0000;
3335       fC[13][3] = 0.3526;
3336       fC[13][4] = 0.3011;
3337
3338       fC[14][0] = 0.005;
3339       fC[14][1] = 0.005;
3340       fC[14][2] = 1.0000;
3341       fC[14][3] = 0.3788;
3342       fC[14][4] = 0.3011;
3343
3344       fC[15][0] = 0.005;
3345       fC[15][1] = 0.005;
3346       fC[15][2] = 1.0000;
3347       fC[15][3] = 0.3788;
3348       fC[15][4] = 0.3011;
3349
3350       fC[16][0] = 0.005;
3351       fC[16][1] = 0.005;
3352       fC[16][2] = 1.0000;
3353       fC[16][3] = 0.3788;
3354       fC[16][4] = 0.3011;
3355
3356       fC[17][0] = 0.005;
3357       fC[17][1] = 0.005;
3358       fC[17][2] = 1.0000;
3359       fC[17][3] = 0.3788;
3360       fC[17][4] = 0.3011;
3361   }
3362
3363   for(Int_t i=18;i<fgkPIDptBin;i++){
3364     fBinLimitPID[i] = 3.0 + 0.2 * (i-18);
3365     fC[i][0] = fC[17][0];
3366     fC[i][1] = fC[17][1];
3367     fC[i][2] = fC[17][2];
3368     fC[i][3] = fC[17][3];
3369     fC[i][4] = fC[17][4];
3370   }
3371 }
3372
3373 //---------------------------------------------------------------//
3374 Bool_t AliFlowTrackCuts::TPCTOFagree(const AliESDtrack *track)
3375 {
3376   Bool_t status = kFALSE;
3377
3378   Float_t mass[5] = {5.10998909999999971e-04,1.05658000000000002e-01,1.39570000000000000e-01,4.93676999999999977e-01,9.38271999999999995e-01};
3379
3380
3381   Double_t exptimes[5];
3382   track->GetIntegratedTimes(exptimes);
3383
3384   Float_t dedx = track->GetTPCsignal();
3385
3386   Float_t p = track->P();
3387   Float_t time = track->GetTOFsignal()- fESDpid.GetTOFResponse().GetStartTime(p);
3388   Float_t tl = track->GetIntegratedLength();
3389
3390   Float_t betagammares =  fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[4], mass[4]);
3391
3392   Float_t betagamma1 = tl/(time-5 *betagammares) * 33.3564095198152043;
3393
3394 //  printf("betagamma1 = %f\n",betagamma1);
3395
3396   if(betagamma1 < 0.1) betagamma1 = 0.1;
3397
3398   if(betagamma1 < 0.99999) betagamma1 /= TMath::Sqrt(1-betagamma1*betagamma1);
3399   else betagamma1 = 100;
3400
3401   Float_t betagamma2 = tl/(time+5 *betagammares) * 33.3564095198152043;
3402 //  printf("betagamma2 = %f\n",betagamma2);
3403
3404   if(betagamma2 < 0.1) betagamma2 = 0.1;
3405
3406   if(betagamma2 < 0.99999) betagamma2 /= TMath::Sqrt(1-betagamma2*betagamma2);
3407   else betagamma2 = 100;
3408
3409
3410   Double_t ptpc[3];
3411   track->GetInnerPxPyPz(ptpc);
3412   Float_t momtpc=TMath::Sqrt(ptpc[0]*ptpc[0] + ptpc[1]*ptpc[1] + ptpc[2]*ptpc[2]);
3413
3414   for(Int_t i=0;i < 5;i++){
3415     Float_t resolutionTOF =  fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[i], mass[i]);
3416     if(TMath::Abs(exptimes[i] - time) < 5 * resolutionTOF){
3417       Float_t dedxExp = 0;
3418       if(i==0) dedxExp =  fESDpid.GetTPCResponse().GetExpectedSignal(momtpc,AliPID::kElectron);
3419       else if(i==1) dedxExp =  fESDpid.GetTPCResponse().GetExpectedSignal(momtpc,AliPID::kMuon);
3420       else if(i==2) dedxExp =  fESDpid.GetTPCResponse().GetExpectedSignal(momtpc,AliPID::kPion);
3421       else if(i==3) dedxExp =  fESDpid.GetTPCResponse().GetExpectedSignal(momtpc,AliPID::kKaon);
3422       else if(i==4) dedxExp =  fESDpid.GetTPCResponse().GetExpectedSignal(momtpc,AliPID::kProton);
3423
3424       Float_t resolutionTPC = 2;
3425       if(i==0) resolutionTPC =   fESDpid.GetTPCResponse().GetExpectedSigma(momtpc,track->GetTPCsignalN(),AliPID::kElectron);
3426       else if(i==1) resolutionTPC =   fESDpid.GetTPCResponse().GetExpectedSigma(momtpc,track->GetTPCsignalN(),AliPID::kMuon);
3427       else if(i==2) resolutionTPC =   fESDpid.GetTPCResponse().GetExpectedSigma(momtpc,track->GetTPCsignalN(),AliPID::kPion);
3428       else if(i==3) resolutionTPC =   fESDpid.GetTPCResponse().GetExpectedSigma(momtpc,track->GetTPCsignalN(),AliPID::kKaon);
3429       else if(i==4) resolutionTPC =   fESDpid.GetTPCResponse().GetExpectedSigma(momtpc,track->GetTPCsignalN(),AliPID::kProton);
3430
3431       if(TMath::Abs(dedx - dedxExp) < 3 * resolutionTPC){
3432         status = kTRUE;
3433       }
3434     }
3435   }
3436
3437   Float_t bb1 =  fESDpid.GetTPCResponse().Bethe(betagamma1);
3438   Float_t bb2 =  fESDpid.GetTPCResponse().Bethe(betagamma2);
3439   Float_t bbM =  fESDpid.GetTPCResponse().Bethe((betagamma1+betagamma2)*0.5);
3440
3441
3442   //  status = kFALSE;
3443   // for nuclei
3444   Float_t resolutionTOFpr =   fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[4], mass[4]);
3445   Float_t resolutionTPCpr =   fESDpid.GetTPCResponse().GetExpectedSigma(momtpc,track->GetTPCsignalN(),AliPID::kProton);
3446   if(TMath::Abs(dedx-bb1) < resolutionTPCpr*3 && exptimes[4] < time-7*resolutionTOFpr){
3447      status = kTRUE;
3448   }
3449   else if(TMath::Abs(dedx-bb2) < resolutionTPCpr*3 && exptimes[4] < time-7*resolutionTOFpr){
3450      status = kTRUE;
3451   }
3452   else if(TMath::Abs(dedx-bbM) < resolutionTPCpr*3 && exptimes[4] < time-7*resolutionTOFpr){
3453      status = kTRUE;
3454   }
3455
3456   return status;
3457 }
3458 // end part added by F. Noferini
3459 //-----------------------------------------------------------------------
3460
3461 //-----------------------------------------------------------------------
3462 const char* AliFlowTrackCuts::PIDsourceName(PIDsource s)
3463 {
3464   //get the name of the particle id source
3465   switch (s)
3466   {
3467     case kTPCdedx:
3468       return "TPCdedx";
3469     case kTPCbayesian:
3470       return "TPCbayesian";
3471     case kTOFbeta:
3472       return "TOFbeta";
3473     case kTPCpid:
3474       return "TPCpid";
3475     case kTOFpid:
3476       return "TOFpid";
3477     case kTOFbayesian:
3478       return "TOFbayesianPID";
3479     case kTOFbetaSimple:
3480       return "TOFbetaSimple";
3481     case kTPCNuclei:
3482       return "TPCnuclei";
3483     default:
3484       return "NOPID";
3485   }
3486 }
3487
3488 //-----------------------------------------------------------------------
3489 const char* AliFlowTrackCuts::GetParamTypeName(trackParameterType type)
3490 {
3491   //return the name of the selected parameter type
3492   switch (type)
3493   {
3494     case kMC:
3495       return "MC";
3496     case kGlobal:
3497       return "Global";
3498     case kTPCstandalone:
3499       return "TPCstandalone";
3500     case kSPDtracklet:
3501       return "SPDtracklets";
3502     case kPMD:
3503       return "PMD";
3504     case kV0:
3505       return "V0";
3506     default:
3507       return "unknown";
3508   }
3509 }
3510
3511 //-----------------------------------------------------------------------
3512 Bool_t AliFlowTrackCuts::PassesPMDcuts(AliESDPmdTrack* track )
3513 {
3514   //check PMD specific cuts
3515   //clean up from last iteration, and init label
3516   Int_t   det   = track->GetDetector();
3517   //Int_t   smn   = track->GetSmn();
3518   Float_t clsX  = track->GetClusterX();
3519   Float_t clsY  = track->GetClusterY();
3520   Float_t clsZ  = track->GetClusterZ();
3521   Float_t ncell = track->GetClusterCells();
3522   Float_t adc   = track->GetClusterADC();
3523
3524   fTrack = NULL;
3525   fMCparticle=NULL;
3526   fTrackLabel=-1;
3527
3528   fTrackEta = GetPmdEta(clsX,clsY,clsZ);
3529   fTrackPhi = GetPmdPhi(clsX,clsY);
3530   fTrackWeight = 1.0;
3531
3532   Bool_t pass=kTRUE;
3533   if (fCutEta) {if (  fTrackEta < fEtaMin || fTrackEta >= fEtaMax ) pass = kFALSE;}
3534   if (fCutPhi) {if ( fTrackPhi < fPhiMin || fTrackPhi >= fPhiMax ) pass = kFALSE;}
3535   if (fCutPmdDet) {if(det != fPmdDet) pass = kFALSE;}
3536   if (fCutPmdAdc) {if(adc < fPmdAdc) pass = kFALSE;}
3537   if (fCutPmdNcell) {if(ncell < fPmdNcell) pass = kFALSE;}
3538
3539   return pass;
3540 }
3541
3542 //-----------------------------------------------------------------------
3543 Bool_t AliFlowTrackCuts::PassesV0cuts(AliVVZERO* vzero, Int_t id)
3544 {
3545   //check V0 cuts
3546
3547   //clean up from last iter
3548   fTrack = NULL;
3549   fMCparticle=NULL;
3550   fTrackLabel=-1;
3551
3552   fTrackPhi = TMath::PiOver4()*(0.5+id%8);
3553   if(id<32)
3554     fTrackEta = -3.45+0.5*(id/8); // taken from PPR
3555   else
3556     fTrackEta = +4.8-0.5*((id/8)-4); // taken from PPR
3557   fTrackWeight = vzero->GetMultiplicity(id); // not corrected yet: we should use AliESDUtils
3558
3559   Bool_t pass=kTRUE;
3560   if (fCutEta) {if (  fTrackEta < fEtaMin || fTrackEta >= fEtaMax ) pass = kFALSE;}
3561   if (fCutPhi) {if ( fTrackPhi < fPhiMin || fTrackPhi >= fPhiMax ) pass = kFALSE;}
3562
3563   return pass;
3564 }
3565
3566 //----------------------------------------------------------------------------//
3567 Double_t AliFlowTrackCuts::GetPmdEta(Float_t xPos, Float_t yPos, Float_t zPos)
3568 {
3569   Float_t rpxpy, theta, eta;
3570   rpxpy  = TMath::Sqrt(xPos*xPos + yPos*yPos);
3571   theta  = TMath::ATan2(rpxpy,zPos);
3572   eta    = -TMath::Log(TMath::Tan(0.5*theta));
3573   return eta;
3574 }
3575
3576 //--------------------------------------------------------------------------//
3577 Double_t AliFlowTrackCuts::GetPmdPhi(Float_t xPos, Float_t yPos)
3578 {
3579   Float_t pybypx, phi = 0., phi1;
3580   if(xPos==0)
3581     {
3582       if(yPos>0) phi = 90.;
3583       if(yPos<0) phi = 270.;
3584     }
3585   if(xPos != 0)
3586     {
3587       pybypx = yPos/xPos;
3588       if(pybypx < 0) pybypx = - pybypx;
3589       phi1 = TMath::ATan(pybypx)*180./3.14159;
3590
3591       if(xPos > 0 && yPos > 0) phi = phi1;        // 1st Quadrant
3592       if(xPos < 0 && yPos > 0) phi = 180 - phi1;  // 2nd Quadrant
3593       if(xPos < 0 && yPos < 0) phi = 180 + phi1;  // 3rd Quadrant
3594       if(xPos > 0 && yPos < 0) phi = 360 - phi1;  // 4th Quadrant
3595
3596     }
3597   phi = phi*3.14159/180.;
3598   return   phi;
3599 }
3600
3601 //---------------------------------------------------------------//
3602 void AliFlowTrackCuts::Browse(TBrowser* b)
3603 {
3604   //some browsing capabilities
3605   if (fQA) b->Add(fQA);
3606 }
3607
3608 //---------------------------------------------------------------//
3609 Long64_t AliFlowTrackCuts::Merge(TCollection* list)
3610 {
3611   //merge
3612   if (!fQA || !list) return 0;
3613   if (list->IsEmpty()) return 0;
3614   AliFlowTrackCuts* obj=NULL;
3615   TList tmplist;
3616   TIter next(list);
3617   while ( (obj = dynamic_cast<AliFlowTrackCuts*>(next())) )
3618   {
3619     if (obj==this) continue;
3620     tmplist.Add(obj->GetQA());
3621   }
3622   return fQA->Merge(&tmplist);
3623 }
3624