]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG2/FLOW/AliFlowTasks/AliFlowTrackCuts.cxx
hooks for PMD flow analysis
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowTasks / AliFlowTrackCuts.cxx
index 014ffa57152b1197038674d8f8f232fc7010cc27..ccfad35cc124af51392a8f722a1e7104d0965ba6 100644 (file)
 #include "AliFlowTrackCuts.h"
 #include "AliLog.h"
 #include "AliESDpid.h"
+#include "AliESDPmdTrack.h"
 
 ClassImp(AliFlowTrackCuts)
 
 //-----------------------------------------------------------------------
 AliFlowTrackCuts::AliFlowTrackCuts():
+  AliFlowTrackSimpleCuts(),
+  fAliESDtrackCuts(NULL),
+  fQA(NULL),
+  fCutMC(kFALSE),
+  fCutMCprocessType(kFALSE),
+  fMCprocessType(kPNoProcess),
+  fCutMCPID(kFALSE),
+  fMCPID(0),
+  fIgnoreSignInPID(kFALSE),
+  fCutMCisPrimary(kFALSE),
+  fRequireTransportBitForPrimaries(kTRUE),
+  fMCisPrimary(kFALSE),
+  fRequireCharge(kFALSE),
+  fFakesAreOK(kTRUE),
+  fCutSPDtrackletDeltaPhi(kFALSE),
+  fSPDtrackletDeltaPhiMax(FLT_MAX),
+  fSPDtrackletDeltaPhiMin(-FLT_MAX),
+  fIgnoreTPCzRange(kFALSE),
+  fIgnoreTPCzRangeMax(FLT_MAX),
+  fIgnoreTPCzRangeMin(-FLT_MAX),
+  fCutChi2PerClusterTPC(kFALSE),
+  fMaxChi2PerClusterTPC(FLT_MAX),
+  fMinChi2PerClusterTPC(-FLT_MAX),
+  fCutNClustersTPC(kFALSE),
+  fNClustersTPCMax(INT_MAX),
+  fNClustersTPCMin(INT_MIN),  
+  fCutNClustersITS(kFALSE),
+  fNClustersITSMax(INT_MAX),
+  fNClustersITSMin(INT_MIN),  
+  fUseAODFilterBit(kFALSE),
+  fAODFilterBit(0),
+  fCutDCAToVertexXY(kFALSE),
+  fCutDCAToVertexZ(kFALSE),
+  fCutMinimalTPCdedx(kFALSE),
+  fMinimalTPCdedx(0.),
+  fParamType(kGlobal),
+  fParamMix(kPure),
+  fTrack(NULL),
+  fTrackPhi(0.),
+  fTrackEta(0.),
+  fTrackWeight(0.),
+  fTrackLabel(INT_MIN),
+  fMCevent(NULL),
+  fMCparticle(NULL),
+  fEvent(NULL),
+  fTPCtrack(),
+  fESDpid(),
+  fPIDsource(kTOFpid),
+  fTPCpidCuts(NULL),
+  fTOFpidCuts(NULL),
+  fParticleID(AliPID::kPion),
+  fParticleProbability(.9)
+{
+  //io constructor 
+  for ( Int_t i=0; i<5; i++ ) { fProbBayes[i]=0.0; }
+  SetPriors(); //init arrays
+}
+
+//-----------------------------------------------------------------------
+AliFlowTrackCuts::AliFlowTrackCuts(const char* name):
   AliFlowTrackSimpleCuts(),
   fAliESDtrackCuts(new AliESDtrackCuts()),
   fQA(NULL),
@@ -66,6 +127,7 @@ AliFlowTrackCuts::AliFlowTrackCuts():
   fMCPID(0),
   fIgnoreSignInPID(kFALSE),
   fCutMCisPrimary(kFALSE),
+  fRequireTransportBitForPrimaries(kTRUE),
   fMCisPrimary(kFALSE),
   fRequireCharge(kFALSE),
   fFakesAreOK(kTRUE),
@@ -81,6 +143,15 @@ AliFlowTrackCuts::AliFlowTrackCuts():
   fCutNClustersTPC(kFALSE),
   fNClustersTPCMax(INT_MAX),
   fNClustersTPCMin(INT_MIN),  
+  fCutNClustersITS(kFALSE),
+  fNClustersITSMax(INT_MAX),
+  fNClustersITSMin(INT_MIN),
+  fUseAODFilterBit(kFALSE),
+  fAODFilterBit(0),
+  fCutDCAToVertexXY(kFALSE),
+  fCutDCAToVertexZ(kFALSE),
+  fCutMinimalTPCdedx(kFALSE),
+  fMinimalTPCdedx(0.),
   fParamType(kGlobal),
   fParamMix(kPure),
   fTrack(NULL),
@@ -92,20 +163,29 @@ AliFlowTrackCuts::AliFlowTrackCuts():
   fMCparticle(NULL),
   fEvent(NULL),
   fTPCtrack(),
-  fESDpid(NULL),
-  fPIDsource(kTPCTOFpid),
+  fESDpid(),
+  fPIDsource(kTOFpid),
   fTPCpidCuts(NULL),
   fTOFpidCuts(NULL),
-  fTPCTOFpidCrossOverPt(0.4),
-  fAliPID(AliPID::kPion)
+  fParticleID(AliPID::kPion),
+  fParticleProbability(.9)
 {
   //constructor 
+  SetName(name);
+  SetTitle("AliFlowTrackCuts");
+  fESDpid.GetTPCResponse().SetBetheBlochParameters( 0.0283086,
+                                                    2.63394e+01,
+                                                    5.04114e-11,
+                                                    2.12543e+00,
+                                                    4.88663e+00 );
+  for ( Int_t i=0; i<5; i++ ) { fProbBayes[i]=0.0; }
+  SetPriors(); //init arrays
 }
 
 //-----------------------------------------------------------------------
 AliFlowTrackCuts::AliFlowTrackCuts(const AliFlowTrackCuts& that):
   AliFlowTrackSimpleCuts(that),
-  fAliESDtrackCuts(new AliESDtrackCuts(*(that.fAliESDtrackCuts))),
+  fAliESDtrackCuts(NULL),
   fQA(NULL),
   fCutMC(that.fCutMC),
   fCutMCprocessType(that.fCutMCprocessType),
@@ -114,6 +194,7 @@ AliFlowTrackCuts::AliFlowTrackCuts(const AliFlowTrackCuts& that):
   fMCPID(that.fMCPID),
   fIgnoreSignInPID(that.fIgnoreSignInPID),
   fCutMCisPrimary(that.fCutMCisPrimary),
+  fRequireTransportBitForPrimaries(that.fRequireTransportBitForPrimaries),
   fMCisPrimary(that.fMCisPrimary),
   fRequireCharge(that.fRequireCharge),
   fFakesAreOK(that.fFakesAreOK),
@@ -129,6 +210,15 @@ AliFlowTrackCuts::AliFlowTrackCuts(const AliFlowTrackCuts& that):
   fCutNClustersTPC(that.fCutNClustersTPC),
   fNClustersTPCMax(that.fNClustersTPCMax),
   fNClustersTPCMin(that.fNClustersTPCMin),
+  fCutNClustersITS(that.fCutNClustersITS),
+  fNClustersITSMax(that.fNClustersITSMax),
+  fNClustersITSMin(that.fNClustersITSMin),
+  fUseAODFilterBit(that.fUseAODFilterBit),
+  fAODFilterBit(that.fAODFilterBit),
+  fCutDCAToVertexXY(that.fCutDCAToVertexXY),
+  fCutDCAToVertexZ(that.fCutDCAToVertexZ),
+  fCutMinimalTPCdedx(that.fCutMinimalTPCdedx),
+  fMinimalTPCdedx(that.fMinimalTPCdedx),
   fParamType(that.fParamType),
   fParamMix(that.fParamMix),
   fTrack(NULL),
@@ -144,20 +234,25 @@ AliFlowTrackCuts::AliFlowTrackCuts(const AliFlowTrackCuts& that):
   fPIDsource(that.fPIDsource),
   fTPCpidCuts(NULL),
   fTOFpidCuts(NULL),
-  fTPCTOFpidCrossOverPt(that.fTPCTOFpidCrossOverPt),
-  fAliPID(that.fAliPID)
+  fParticleID(that.fParticleID),
+  fParticleProbability(that.fParticleProbability)
 {
   //copy constructor
   if (that.fTPCpidCuts) fTPCpidCuts = new TMatrixF(*(that.fTPCpidCuts));
   if (that.fTOFpidCuts) fTOFpidCuts = new TMatrixF(*(that.fTOFpidCuts));
+  if (that.fAliESDtrackCuts) fAliESDtrackCuts = new AliESDtrackCuts(*(that.fAliESDtrackCuts));
+  memcpy(fProbBayes,that.fProbBayes,sizeof(fProbBayes));
+  SetPriors(); //init arrays
 }
 
 //-----------------------------------------------------------------------
 AliFlowTrackCuts& AliFlowTrackCuts::operator=(const AliFlowTrackCuts& that)
 {
   //assignment
+  if (this==&that) return *this;
+
   AliFlowTrackSimpleCuts::operator=(that);
-  *fAliESDtrackCuts=*(that.fAliESDtrackCuts);
+  if (that.fAliESDtrackCuts) *fAliESDtrackCuts=*(that.fAliESDtrackCuts);
   fQA=NULL;
   fCutMC=that.fCutMC;
   fCutMCprocessType=that.fCutMCprocessType;
@@ -166,6 +261,7 @@ AliFlowTrackCuts& AliFlowTrackCuts::operator=(const AliFlowTrackCuts& that)
   fMCPID=that.fMCPID;
   fIgnoreSignInPID=that.fIgnoreSignInPID,
   fCutMCisPrimary=that.fCutMCisPrimary;
+  fRequireTransportBitForPrimaries=that.fRequireTransportBitForPrimaries;
   fMCisPrimary=that.fMCisPrimary;
   fRequireCharge=that.fRequireCharge;
   fFakesAreOK=that.fFakesAreOK;
@@ -181,6 +277,15 @@ AliFlowTrackCuts& AliFlowTrackCuts::operator=(const AliFlowTrackCuts& that)
   fCutNClustersTPC=that.fCutNClustersTPC;
   fNClustersTPCMax=that.fNClustersTPCMax;
   fNClustersTPCMin=that.fNClustersTPCMin;  
+  fCutNClustersITS=that.fCutNClustersITS;
+  fNClustersITSMax=that.fNClustersITSMax;
+  fNClustersITSMin=that.fNClustersITSMin;  
+  fUseAODFilterBit=that.fUseAODFilterBit;
+  fAODFilterBit=that.fAODFilterBit;
+  fCutDCAToVertexXY=that.fCutDCAToVertexXY;
+  fCutDCAToVertexZ=that.fCutDCAToVertexZ;
+  fCutMinimalTPCdedx=that.fCutMinimalTPCdedx;
+  fMinimalTPCdedx=that.fMinimalTPCdedx;
   fParamType=that.fParamType;
   fParamMix=that.fParamMix;
 
@@ -196,11 +301,14 @@ AliFlowTrackCuts& AliFlowTrackCuts::operator=(const AliFlowTrackCuts& that)
   fESDpid = that.fESDpid;
   fPIDsource = that.fPIDsource;
 
+  delete fTPCpidCuts;
+  delete fTOFpidCuts;
   if (that.fTPCpidCuts) fTPCpidCuts = new TMatrixF(*(that.fTPCpidCuts));
   if (that.fTOFpidCuts) fTOFpidCuts = new TMatrixF(*(that.fTOFpidCuts));
-  fTPCTOFpidCrossOverPt=that.fTPCTOFpidCrossOverPt;
 
-  fAliPID=that.fAliPID;
+  fParticleID=that.fParticleID;
+  fParticleProbability=that.fParticleProbability;
+  memcpy(fProbBayes,that.fProbBayes,sizeof(fProbBayes));
 
   return *this;
 }
@@ -214,6 +322,45 @@ AliFlowTrackCuts::~AliFlowTrackCuts()
   delete fTOFpidCuts;
 }
 
+//-----------------------------------------------------------------------
+void AliFlowTrackCuts::SetEvent(AliVEvent* event, AliMCEvent* mcEvent)
+{
+  //set the event
+  Clear();
+  fEvent=event;
+  fMCevent=mcEvent;
+
+  //do the magic for ESD
+  AliESDEvent* myESD = dynamic_cast<AliESDEvent*>(event);
+  if (fCutPID && myESD)
+  {
+    //TODO: maybe call it only for the TOF options?
+    // Added by F. Noferini for TOF PID
+    fESDpid.SetTOFResponse(myESD,AliESDpid::kTOF_T0);
+    fESDpid.MakePID(myESD,kFALSE);
+    // End F. Noferini added part
+  }
+
+  //TODO: AOD
+}
+
+//-----------------------------------------------------------------------
+void AliFlowTrackCuts::SetCutMC( Bool_t b )
+{
+  //will we be cutting on MC information?
+  fCutMC=b;
+
+  //if we cut on MC info then also the Bethe Bloch should be the one tuned for MC
+  if (fCutMC)
+  {
+    fESDpid.GetTPCResponse().SetBetheBlochParameters( 2.15898e+00/50.,
+                                                       1.75295e+01,
+                                                       3.40030e-09,
+                                                       1.96178e+00,
+                                                       3.91720e+00);
+  }
+}
+
 //-----------------------------------------------------------------------
 Bool_t AliFlowTrackCuts::IsSelected(TObject* obj, Int_t id)
 {
@@ -224,6 +371,8 @@ Bool_t AliFlowTrackCuts::IsSelected(TObject* obj, Int_t id)
   if (flowtrack) return PassesCuts(flowtrack);
   AliMultiplicity* tracklets = dynamic_cast<AliMultiplicity*>(obj);
   if (tracklets) return PassesCuts(tracklets,id);
+  AliESDPmdTrack* pmdtrack = dynamic_cast<AliESDPmdTrack*>(obj);
+  if (pmdtrack) return PassesPMDcuts(pmdtrack);
   return kFALSE;  //default when passed wrong type of object
 }
 
@@ -298,7 +447,7 @@ Bool_t AliFlowTrackCuts::PassesMCcuts(AliMCEvent* mcEvent, Int_t label)
 
   if (fCutMCisPrimary)
   {
-    if (IsPhysicalPrimary(mcEvent,label) != fMCisPrimary) return kFALSE;
+    if (IsPhysicalPrimary(mcEvent,label,fRequireTransportBitForPrimaries) != fMCisPrimary) return kFALSE;
   }
   if (fCutMCPID)
   {
@@ -347,18 +496,24 @@ Bool_t AliFlowTrackCuts::PassesCuts(AliVParticle* vparticle)
 
   Bool_t isMCparticle = kFALSE; //some things are different for MC particles, check!
   AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*>(vparticle);
+  AliAODTrack* aodTrack = NULL;
   if (esdTrack)
+    //for an ESD track we do some magic sometimes like constructing TPC only parameters
+    //or doing some hybrid, handle that here
     HandleESDtrack(esdTrack);
   else
   {
     HandleVParticle(vparticle);
     //now check if produced particle is MC
     isMCparticle = (dynamic_cast<AliMCParticle*>(fTrack))!=NULL;
+    aodTrack = dynamic_cast<AliAODTrack*>(vparticle); //keep the additional dynamic cast out of the way for ESDs
   }
   ////////////////////////////////////////////////////////////////
   ////////////////////////////////////////////////////////////////
 
   if (!fTrack) return kFALSE;
+  //because it may be different from global, not needed for aodTrack because we dont do anything funky there
+  if (esdTrack) esdTrack = static_cast<AliESDtrack*>(fTrack);
   
   Bool_t pass=kTRUE;
   //check the common cuts for the current particle fTrack (MC,AOD,ESD)
@@ -380,70 +535,128 @@ Bool_t AliFlowTrackCuts::PassesCuts(AliVParticle* vparticle)
   //when additionally MC info is required
   if (fCutMC && !PassesMCcuts()) pass=kFALSE;
 
-  //check all else for ESDs using aliesdtrackcuts
-  if (esdTrack && (fParamType!=kMC) ) 
+  //the case of ESD or AOD
+  if (esdTrack) { if (!PassesESDcuts(esdTrack)) { pass=kFALSE; } }
+  if (aodTrack) { if (!PassesAODcuts(aodTrack)) { pass=kFALSE; } }
+
+  //true by default, if we didn't set any cuts
+  return pass;
+}
+
+//_______________________________________________________________________
+Bool_t AliFlowTrackCuts::PassesAODcuts(AliAODTrack* track)
+{
+  Bool_t pass = kTRUE;
+
+  if (fCutNClustersTPC)
+  {
+    Int_t ntpccls = track->GetTPCNcls();
+    if (ntpccls < fNClustersTPCMin || ntpccls > fNClustersTPCMax) pass=kFALSE;
+  }
+
+  if (fCutNClustersITS)
+  {
+    Int_t nitscls = track->GetITSNcls();
+    if (nitscls < fNClustersITSMin || nitscls > fNClustersITSMax) pass=kFALSE;
+  }
+  
+   if (fCutChi2PerClusterTPC)
+  {
+    Double_t chi2tpc = track->Chi2perNDF();
+    if (chi2tpc < fMinChi2PerClusterTPC || chi2tpc > fMaxChi2PerClusterTPC) pass=kFALSE;
+  }
+  
+  if (GetRequireTPCRefit() && !(track->GetStatus() & AliESDtrack::kTPCrefit) ) pass=kFALSE;
+  if (GetRequireITSRefit() && !(track->GetStatus() & AliESDtrack::kITSrefit) ) pass=kFALSE;
+  
+  if (fUseAODFilterBit && !track->TestFilterBit(fAODFilterBit)) pass=kFALSE;
+  
+  if (fCutDCAToVertexXY && track->DCA()>GetMaxDCAToVertexXY()) pass=kFALSE;
+    
+
+  return pass;
+}
+
+//_______________________________________________________________________
+Bool_t AliFlowTrackCuts::PassesESDcuts(AliESDtrack* track)
+{
+  Bool_t pass=kTRUE;
+  if (fIgnoreTPCzRange)
   {
-    if (fIgnoreTPCzRange)
+    const AliExternalTrackParam* pin = track->GetOuterParam();
+    const AliExternalTrackParam* pout = track->GetInnerParam();
+    if (pin&&pout)
     {
-      const AliExternalTrackParam* pin = esdTrack->GetOuterParam();
-      const AliExternalTrackParam* pout = esdTrack->GetInnerParam();
-      if (pin&&pout)
-      {
-        Double_t zin = pin->GetZ();
-        Double_t zout = pout->GetZ();
-        if (zin*zout<0) pass=kFALSE;   //reject if cross the membrane
-        if (zin < fIgnoreTPCzRangeMin || zin > fIgnoreTPCzRangeMax) pass=kFALSE;
-        if (zout < fIgnoreTPCzRangeMin || zout > fIgnoreTPCzRangeMax) pass=kFALSE;
-      }
+      Double_t zin = pin->GetZ();
+      Double_t zout = pout->GetZ();
+      if (zin*zout<0) pass=kFALSE;   //reject if cross the membrane
+      if (zin < fIgnoreTPCzRangeMin || zin > fIgnoreTPCzRangeMax) pass=kFALSE;
+      if (zout < fIgnoreTPCzRangeMin || zout > fIgnoreTPCzRangeMax) pass=kFALSE;
     }
+  }
  
-    if (!fAliESDtrackCuts->IsSelected(static_cast<AliESDtrack*>(fTrack))) pass=kFALSE;
-    Int_t ntpccls = ( fParamType==kESD_TPConly )?
-                      esdTrack->GetTPCNclsIter1():esdTrack->GetTPCNcls();    
-    if (fCutChi2PerClusterTPC)
-    {
-      Float_t tpcchi2 = (fParamType==kESD_TPConly)?
-                         esdTrack->GetTPCchi2Iter1():esdTrack->GetTPCchi2();
-      tpcchi2 = (ntpccls>0)?tpcchi2/ntpccls:-FLT_MAX;
-      if (tpcchi2<fMinChi2PerClusterTPC || tpcchi2 >=fMaxChi2PerClusterTPC)
-        pass=kFALSE;
-    }
+  Int_t ntpccls = ( fParamType==kESD_TPConly )?
+                    track->GetTPCNclsIter1():track->GetTPCNcls();    
+  if (fCutChi2PerClusterTPC)
+  {
+    Float_t tpcchi2 = (fParamType==kESD_TPConly)?
+                       track->GetTPCchi2Iter1():track->GetTPCchi2();
+    tpcchi2 = (ntpccls>0)?tpcchi2/ntpccls:-FLT_MAX;
+    if (tpcchi2<fMinChi2PerClusterTPC || tpcchi2 >=fMaxChi2PerClusterTPC)
+      pass=kFALSE;
+  }
+
+  if (fCutMinimalTPCdedx) 
+  {
+    if (track->GetTPCsignal() < fMinimalTPCdedx) pass=kFALSE;
+  }
+
+  if (fCutNClustersTPC)
+  {
+    if (ntpccls < fNClustersTPCMin || ntpccls > fNClustersTPCMax) pass=kFALSE;
+  }
+
+  Int_t nitscls = track->GetNcls(0);
+  if (fCutNClustersITS)
+  {
+    if (nitscls < fNClustersITSMin || nitscls > fNClustersITSMax) pass=kFALSE;
+  }
 
-    if (fCutNClustersTPC)
+  if (fCutPID)
+  {
+    switch (fPIDsource)    
     {
-      if (ntpccls < fNClustersTPCMin || ntpccls > fNClustersTPCMax) pass=kFALSE;
+      case kTPCpid:
+        if (!PassesTPCpidCut(track)) pass=kFALSE;
+        break;
+      case kTPCdedx:
+        if (!PassesTPCdedxCut(track)) pass=kFALSE;
+        break;
+      case kTOFpid:
+        if (!PassesTOFpidCut(track)) pass=kFALSE;
+        break;
+      case kTOFbeta:
+        if (!PassesTOFbetaCut(track)) pass=kFALSE;
+        break;
+           // part added by F. Noferini
+      case kTOFbayesian:
+             if (!PassesTOFbayesianCut(track)) pass=kFALSE;
+             break;
+           // end part added by F. Noferini
+      default:
+        printf("AliFlowTrackCuts::PassesCuts() this should never be called!\n");
+        pass=kFALSE;
+        break;
     }
+  }    
 
-    if (fCutPID)
-    {
-      switch (fPIDsource)    
-      {
-        case kTPCpid:
-          if (!PassesTPCpidCut(esdTrack)) pass=kFALSE;
-          break;
-        case kTOFpid:
-          if (!PassesTOFpidCut(esdTrack)) pass=kFALSE;
-          break;
-        case kTPCTOFpid:
-          if (pt< fTPCTOFpidCrossOverPt)
-          {
-            if (!PassesTPCpidCut(esdTrack)) pass=kFALSE;
-          }
-          else //if (pt>=fTPCTOFpidCrossOverPt)
-          {
-            if (!PassesTOFpidCut(esdTrack)) pass=kFALSE;
-          }
-          break;
-        default:
-          printf("AliFlowTrackCuts::PassesCuts() this should never be called!\n");
-          pass=kFALSE;
-          break;
-      }
-    }    
+  //some stuff is still handled by AliESDtrackCuts class - delegate
+  if (fAliESDtrackCuts)
+  {
+    if (!fAliESDtrackCuts->IsSelected(track)) pass=kFALSE;
   }
-
-  return pass; //true by default, if we didn't set any cuts
+  return pass;
 }
 
 //-----------------------------------------------------------------------
@@ -468,7 +681,7 @@ void AliFlowTrackCuts::HandleESDtrack(AliESDtrack* track)
       fTrack = track;
       break;
     case kESD_TPConly:
-      if (!track->FillTPCOnlyTrack(fTPCtrack))
+      if (!track->FillTPCOnlyTrack(fTPCtrack)) 
       {
         fTrack=NULL;
         fMCparticle=NULL;
@@ -487,15 +700,89 @@ void AliFlowTrackCuts::HandleESDtrack(AliESDtrack* track)
   }
 }
 
+//-----------------------------------------------------------------------
+Int_t AliFlowTrackCuts::Count(AliVEvent* event)
+{
+  //calculate the number of track in given event.
+  //if argument is NULL(default) take the event attached
+  //by SetEvent()
+  Int_t multiplicity = 0;
+  if (!event)
+  {
+    for (Int_t i=0; i<GetNumberOfInputObjects(); i++)
+    {
+      if (IsSelected(GetInputObject(i))) multiplicity++;
+    }
+  }
+  else
+  {
+    for (Int_t i=0; i<event->GetNumberOfTracks(); i++)
+    {
+      if (IsSelected(event->GetTrack(i))) multiplicity++;
+    }
+  }
+  return multiplicity;
+}
+
+//-----------------------------------------------------------------------
+AliFlowTrackCuts* AliFlowTrackCuts::GetStandardGlobalTrackCuts2010()
+{
+  //get standard cuts
+  AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard Global flow cuts");
+  cuts->SetParamType(kGlobal);
+  cuts->SetPtRange(0.2,5.);
+  cuts->SetEtaRange(-0.8,0.8);
+  cuts->SetMinNClustersTPC(70);
+  cuts->SetMinChi2PerClusterTPC(0.1);
+  cuts->SetMaxChi2PerClusterTPC(4.0);
+  cuts->SetMinNClustersITS(2);
+  cuts->SetRequireITSRefit(kTRUE);
+  cuts->SetRequireTPCRefit(kTRUE);
+  cuts->SetMaxDCAToVertexXY(0.3);
+  cuts->SetMaxDCAToVertexZ(0.3);
+  cuts->SetAcceptKinkDaughters(kFALSE);
+  cuts->SetMinimalTPCdedx(10.);
+
+  return cuts;
+}
+
+//-----------------------------------------------------------------------
+AliFlowTrackCuts* AliFlowTrackCuts::GetStandardTPCOnlyTrackCuts2010()
+{
+  //get standard cuts
+  AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard TPConly flow cuts");
+  cuts->SetParamType(kESD_TPConly);
+  cuts->SetPtRange(0.2,5.);
+  cuts->SetEtaRange(-0.8,0.8);
+  cuts->SetMinNClustersTPC(70);
+  cuts->SetMinChi2PerClusterTPC(0.2);
+  cuts->SetMaxChi2PerClusterTPC(4.0);
+  cuts->SetMaxDCAToVertexXY(3.0);
+  cuts->SetMaxDCAToVertexZ(3.0);
+  cuts->SetDCAToVertex2D(kTRUE);
+  cuts->SetAcceptKinkDaughters(kFALSE);
+  cuts->SetMinimalTPCdedx(10.);
+
+  return cuts;
+}
+
 //-----------------------------------------------------------------------
 AliFlowTrackCuts* AliFlowTrackCuts::GetStandardTPCOnlyTrackCuts()
 {
   //get standard cuts
-  AliFlowTrackCuts* cuts = new AliFlowTrackCuts();
-  cuts->SetName("standard TPConly cuts");
-  delete cuts->fAliESDtrackCuts;
-  cuts->fAliESDtrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
+  AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard TPConly flow cuts");
   cuts->SetParamType(kESD_TPConly);
+  cuts->SetPtRange(0.2,5.);
+  cuts->SetEtaRange(-0.8,0.8);
+  cuts->SetMinNClustersTPC(70);
+  cuts->SetMinChi2PerClusterTPC(0.2);
+  cuts->SetMaxChi2PerClusterTPC(4.0);
+  cuts->SetMaxDCAToVertexXY(3.0);
+  cuts->SetMaxDCAToVertexZ(3.0);
+  cuts->SetDCAToVertex2D(kTRUE);
+  cuts->SetAcceptKinkDaughters(kFALSE);
+  cuts->SetMinimalTPCdedx(10.);
+
   return cuts;
 }
 
@@ -503,8 +790,7 @@ AliFlowTrackCuts* AliFlowTrackCuts::GetStandardTPCOnlyTrackCuts()
 AliFlowTrackCuts* AliFlowTrackCuts::GetStandardITSTPCTrackCuts2009(Bool_t selPrimaries)
 {
   //get standard cuts
-  AliFlowTrackCuts* cuts = new AliFlowTrackCuts();
-  cuts->SetName("standard global track cuts 2009");
+  AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard global track cuts 2009");
   delete cuts->fAliESDtrackCuts;
   cuts->fAliESDtrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2009(selPrimaries);
   cuts->SetParamType(kGlobal);
@@ -512,139 +798,192 @@ AliFlowTrackCuts* AliFlowTrackCuts::GetStandardITSTPCTrackCuts2009(Bool_t selPri
 }
 
 //-----------------------------------------------------------------------
-AliFlowTrack* AliFlowTrackCuts::MakeFlowTrack() const
+AliFlowTrack* AliFlowTrackCuts::MakeFlowTrackSPDtracklet() const
 {
-  //get a flow track constructed from whatever we applied cuts on
-  //caller is resposible for deletion
-  //if construction fails return NULL
+  //make a flow track from tracklet
   AliFlowTrack* flowtrack=NULL;
   TParticle *tmpTParticle=NULL;
   AliMCParticle* tmpAliMCParticle=NULL;
-  if (fParamType==kESD_SPDtracklet)
+  switch (fParamMix)
   {
-    switch (fParamMix)
-    {
-      case kPure:
-        flowtrack = new AliFlowTrack();
-        flowtrack->SetPhi(fTrackPhi);
-        flowtrack->SetEta(fTrackEta);
-        break;
-      case kTrackWithMCkine:
-        if (!fMCparticle) return NULL;
-        flowtrack = new AliFlowTrack();
-        flowtrack->SetPhi( fMCparticle->Phi() );
-        flowtrack->SetEta( fMCparticle->Eta() );
-        flowtrack->SetPt( fMCparticle->Pt() );
-        break;
-      case kTrackWithMCpt:
-        if (!fMCparticle) return NULL;
-        flowtrack = new AliFlowTrack();
-        flowtrack->SetPhi(fTrackPhi);
-        flowtrack->SetEta(fTrackEta);
-        flowtrack->SetPt(fMCparticle->Pt());
-        break;
-      case kTrackWithPtFromFirstMother:
-        if (!fMCparticle) return NULL;
-        flowtrack = new AliFlowTrack();
-        flowtrack->SetPhi(fTrackPhi);
-        flowtrack->SetEta(fTrackEta);
-        tmpTParticle = fMCparticle->Particle();
-        tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
-        flowtrack->SetPt(tmpAliMCParticle->Pt());
-        break;
-      default:
-        flowtrack = new AliFlowTrack();
-        flowtrack->SetPhi(fTrackPhi);
-        flowtrack->SetEta(fTrackEta);
-        break;
-    }
-    flowtrack->SetSource(AliFlowTrack::kFromTracklet);
+    case kPure:
+      flowtrack = new AliFlowTrack();
+      flowtrack->SetPhi(fTrackPhi);
+      flowtrack->SetEta(fTrackEta);
+      break;
+    case kTrackWithMCkine:
+      if (!fMCparticle) return NULL;
+      flowtrack = new AliFlowTrack();
+      flowtrack->SetPhi( fMCparticle->Phi() );
+      flowtrack->SetEta( fMCparticle->Eta() );
+      flowtrack->SetPt( fMCparticle->Pt() );
+      break;
+    case kTrackWithMCpt:
+      if (!fMCparticle) return NULL;
+      flowtrack = new AliFlowTrack();
+      flowtrack->SetPhi(fTrackPhi);
+      flowtrack->SetEta(fTrackEta);
+      flowtrack->SetPt(fMCparticle->Pt());
+      break;
+    case kTrackWithPtFromFirstMother:
+      if (!fMCparticle) return NULL;
+      flowtrack = new AliFlowTrack();
+      flowtrack->SetPhi(fTrackPhi);
+      flowtrack->SetEta(fTrackEta);
+      tmpTParticle = fMCparticle->Particle();
+      tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
+      flowtrack->SetPt(tmpAliMCParticle->Pt());
+      break;
+    default:
+      flowtrack = new AliFlowTrack();
+      flowtrack->SetPhi(fTrackPhi);
+      flowtrack->SetEta(fTrackEta);
+      break;
   }
-  else
+  flowtrack->SetSource(AliFlowTrack::kFromTracklet);
+  return flowtrack;
+}
+
+//-----------------------------------------------------------------------
+AliFlowTrack* AliFlowTrackCuts::MakeFlowTrackVParticle() const
+{
+  //make flow track from AliVParticle (ESD,AOD,MC)
+  if (!fTrack) return NULL;
+  AliFlowTrack* flowtrack=NULL;
+  TParticle *tmpTParticle=NULL;
+  AliMCParticle* tmpAliMCParticle=NULL;
+  switch(fParamMix)
   {
-    if (!fTrack) return NULL;
-    switch(fParamMix)
-    {
-      case kPure:
-        flowtrack = new AliFlowTrack(fTrack);
-        break;
-      case kTrackWithMCkine:
-        flowtrack = new AliFlowTrack(fMCparticle);
-        break;
-      case kTrackWithMCPID:
-        flowtrack = new AliFlowTrack(fTrack);
-        //flowtrack->setPID(...) from mc, when implemented
-        break;
-      case kTrackWithMCpt:
-        if (!fMCparticle) return NULL;
-        flowtrack = new AliFlowTrack(fTrack);
-        flowtrack->SetPt(fMCparticle->Pt());
-        break;
-      case kTrackWithPtFromFirstMother:
-        if (!fMCparticle) return NULL;
-        flowtrack = new AliFlowTrack(fTrack);
-        tmpTParticle = fMCparticle->Particle();
-        tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
-        flowtrack->SetPt(tmpAliMCParticle->Pt());
-        break;
-      default:
-        flowtrack = new AliFlowTrack(fTrack);
-        break;
-    }
-    if (fParamType==kMC) flowtrack->SetSource(AliFlowTrack::kFromMC);
-    else if (dynamic_cast<AliESDtrack*>(fTrack)) flowtrack->SetSource(AliFlowTrack::kFromESD);
-    else if (dynamic_cast<AliAODTrack*>(fTrack)) flowtrack->SetSource(AliFlowTrack::kFromAOD);
-    else if (dynamic_cast<AliMCParticle*>(fTrack)) flowtrack->SetSource(AliFlowTrack::kFromMC);
+    case kPure:
+      flowtrack = new AliFlowTrack(fTrack);
+      break;
+    case kTrackWithMCkine:
+      flowtrack = new AliFlowTrack(fMCparticle);
+      break;
+    case kTrackWithMCPID:
+      flowtrack = new AliFlowTrack(fTrack);
+      //flowtrack->setPID(...) from mc, when implemented
+      break;
+    case kTrackWithMCpt:
+      if (!fMCparticle) return NULL;
+      flowtrack = new AliFlowTrack(fTrack);
+      flowtrack->SetPt(fMCparticle->Pt());
+      break;
+    case kTrackWithPtFromFirstMother:
+      if (!fMCparticle) return NULL;
+      flowtrack = new AliFlowTrack(fTrack);
+      tmpTParticle = fMCparticle->Particle();
+      tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
+      flowtrack->SetPt(tmpAliMCParticle->Pt());
+      break;
+    default:
+      flowtrack = new AliFlowTrack(fTrack);
+      break;
+  }
+  if (fParamType==kMC) flowtrack->SetSource(AliFlowTrack::kFromMC);
+  else if (dynamic_cast<AliESDtrack*>(fTrack)) flowtrack->SetSource(AliFlowTrack::kFromESD);
+  else if (dynamic_cast<AliAODTrack*>(fTrack)) flowtrack->SetSource(AliFlowTrack::kFromAOD);
+  else if (dynamic_cast<AliMCParticle*>(fTrack)) flowtrack->SetSource(AliFlowTrack::kFromMC);
+  return flowtrack;
+}
+
+//-----------------------------------------------------------------------
+AliFlowTrack* AliFlowTrackCuts::MakeFlowTrackPMDtrack() const
+{
+  //make a flow track from PMD track
+  AliFlowTrack* flowtrack=NULL;
+  TParticle *tmpTParticle=NULL;
+  AliMCParticle* tmpAliMCParticle=NULL;
+  switch (fParamMix)
+  {
+    case kPure:
+      flowtrack = new AliFlowTrack();
+      flowtrack->SetPhi(fTrackPhi);
+      flowtrack->SetEta(fTrackEta);
+      flowtrack->SetWeight(fTrackWeight);
+      break;
+    case kTrackWithMCkine:
+      if (!fMCparticle) return NULL;
+      flowtrack = new AliFlowTrack();
+      flowtrack->SetPhi( fMCparticle->Phi() );
+      flowtrack->SetEta( fMCparticle->Eta() );
+      flowtrack->SetWeight(fTrackWeight);
+      flowtrack->SetPt( fMCparticle->Pt() );
+      break;
+    case kTrackWithMCpt:
+      if (!fMCparticle) return NULL;
+      flowtrack = new AliFlowTrack();
+      flowtrack->SetPhi(fTrackPhi);
+      flowtrack->SetEta(fTrackEta);
+      flowtrack->SetWeight(fTrackWeight);
+      flowtrack->SetPt(fMCparticle->Pt());
+      break;
+    case kTrackWithPtFromFirstMother:
+      if (!fMCparticle) return NULL;
+      flowtrack = new AliFlowTrack();
+      flowtrack->SetPhi(fTrackPhi);
+      flowtrack->SetEta(fTrackEta);
+      flowtrack->SetWeight(fTrackWeight);
+      tmpTParticle = fMCparticle->Particle();
+      tmpAliMCParticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(tmpTParticle->GetFirstMother()));
+      flowtrack->SetPt(tmpAliMCParticle->Pt());
+      break;
+    default:
+      flowtrack = new AliFlowTrack();
+      flowtrack->SetPhi(fTrackPhi);
+      flowtrack->SetEta(fTrackEta);
+      flowtrack->SetWeight(fTrackWeight);
+      break;
   }
+
+  flowtrack->SetSource(AliFlowTrack::kFromPMD);
   return flowtrack;
 }
 
+//-----------------------------------------------------------------------
+AliFlowTrack* AliFlowTrackCuts::MakeFlowTrack() const
+{
+  //get a flow track constructed from whatever we applied cuts on
+  //caller is resposible for deletion
+  //if construction fails return NULL
+  switch (fParamType)
+  {
+    case kESD_SPDtracklet:
+      return MakeFlowTrackSPDtracklet();
+    case kPMD:
+      return MakeFlowTrackPMDtrack();
+    default:
+      return MakeFlowTrackVParticle();
+  }
+}
+
 //-----------------------------------------------------------------------
 Bool_t AliFlowTrackCuts::IsPhysicalPrimary() const
 {
   //check if current particle is a physical primary
   if (!fMCevent) return kFALSE;
   if (fTrackLabel<0) return kFALSE;
-  return IsPhysicalPrimary(fMCevent, fTrackLabel);
+  return IsPhysicalPrimary(fMCevent, fTrackLabel, fRequireTransportBitForPrimaries);
 }
 
 //-----------------------------------------------------------------------
-Bool_t AliFlowTrackCuts::IsPhysicalPrimary(AliMCEvent* mcEvent, Int_t label)
+Bool_t AliFlowTrackCuts::IsPhysicalPrimary(AliMCEvent* mcEvent, Int_t label, Bool_t requiretransported)
 {
   //check if current particle is a physical primary
   Bool_t physprim=mcEvent->IsPhysicalPrimary(label);
-  if (!physprim) return kFALSE;
   AliMCParticle* track = static_cast<AliMCParticle*>(mcEvent->GetTrack(label));
   if (!track) return kFALSE;
   TParticle* particle = track->Particle();
   Bool_t transported = particle->TestBit(kTransportBit);
-  //printf("prim: %s, transp: %s\n",(physprim)?"YES":"NO ",(transported)?"YES":"NO ");
-  return (physprim && transported);
-}
-
-//-----------------------------------------------------------------------
-const char* AliFlowTrackCuts::GetParamTypeName(trackParameterType type) 
-{
-  //return the name of the selected parameter type
-  switch (type)
-  {
-    case kMC:
-      return "MC";
-    case kGlobal:
-      return "ESD global";
-    case kESD_TPConly:
-      return "TPC only";
-    case kESD_SPDtracklet:
-        return "SPD tracklet";
-    default:
-        return "unknown";
-  }
-  return "unknown";
+  //printf("label: %i prim: %s, transp: %s, pass: %s\n",label, (physprim)?"YES":"NO ",(transported)?"YES":"NO ",
+        //(physprim && (transported || !requiretransported))?"YES":"NO"  );
+  return (physprim && (transported || !requiretransported));
 }
 
 //-----------------------------------------------------------------------
 void AliFlowTrackCuts::DefineHistograms()
 {
+  //define qa histograms
 }
 
 //-----------------------------------------------------------------------
@@ -662,6 +1001,10 @@ Int_t AliFlowTrackCuts::GetNumberOfInputObjects() const
     case kMC:
       if (!fMCevent) return 0;
       return fMCevent->GetNumberOfTracks();
+    case kPMD:
+      esd = dynamic_cast<AliESDEvent*>(fEvent);
+      if (!esd) return 0;
+      return esd->GetNumberOfPmdTracks();
     default:
       if (!fEvent) return 0;
       return fEvent->GetNumberOfTracks();
@@ -684,6 +1027,10 @@ TObject* AliFlowTrackCuts::GetInputObject(Int_t i)
     case kMC:
       if (!fMCevent) return NULL;
       return fMCevent->GetTrack(i);
+    case kPMD:
+      esd = dynamic_cast<AliESDEvent*>(fEvent);
+      if (!esd) return NULL;
+      return esd->GetPmdTrack(i);
     default:
       if (!fEvent) return NULL;
       return fEvent->GetTrack(i);
@@ -704,23 +1051,30 @@ void AliFlowTrackCuts::Clear(Option_t*)
 }
 
 //-----------------------------------------------------------------------
-Bool_t AliFlowTrackCuts::PassesTOFpidCut(AliESDtrack* t )
+Bool_t AliFlowTrackCuts::PassesTOFbetaCut(AliESDtrack* track )
 {
   //check if passes PID cut using timing in TOF
-  if (!fESDpid) return kFALSE;
-  if (!(t && (t->GetStatus() & AliESDtrack::kTOFout) && (t->GetStatus() & AliESDtrack::kTIME)
-       && (t->GetTOFsignal() > 12000) && (t->GetTOFsignal() < 100000) && (t->GetIntegratedLength() > 365)))
-       return kFALSE;
-  Float_t pt = t->Pt();
-  Float_t p = t->GetP();
-  Float_t trackT0 = fESDpid->GetTOFResponse().GetStartTime(p);
-  Float_t timeTOF = t->GetTOFsignal()- trackT0; 
-  //2=pion 3=kaon 4=protons
-  Double_t inttimes[5] = {-1.0,-1.0,-1.0,-1.0,-1.0};
-  t->GetIntegratedTimes(inttimes);
-  //construct the pid index because it's screwed up in TOF
+  Bool_t goodtrack = (track) && 
+                     (track->GetStatus() & AliESDtrack::kTOFpid) && 
+                     (track->GetTOFsignal() > 12000) && 
+                     (track->GetTOFsignal() < 100000) && 
+                     (track->GetIntegratedLength() > 365) && 
+                    !(track->GetStatus() & AliESDtrack::kTOFmismatch);
+
+  if (!goodtrack) return kFALSE;
+  
+  const Float_t c = 2.99792457999999984e-02;  
+  Float_t p = track->GetP();
+  Float_t L = track->GetIntegratedLength();  
+  Float_t trackT0 = fESDpid.GetTOFResponse().GetStartTime(p);
+  Float_t timeTOF = track->GetTOFsignal()- trackT0; 
+  Float_t beta = L/timeTOF/c;
+  Double_t integratedTimes[5] = {-1.0,-1.0,-1.0,-1.0,-1.0};
+  track->GetIntegratedTimes(integratedTimes);
+
+  //construct the pid index because it's not AliPID::EParticleType
   Int_t pid = 0;
-  switch (fAliPID)
+  switch (fParticleID)
   {
     case AliPID::kPion:
       pid=2;
@@ -734,35 +1088,62 @@ Bool_t AliFlowTrackCuts::PassesTOFpidCut(AliESDtrack* t )
     default:
       return kFALSE;
   }
-  Float_t s = timeTOF-inttimes[pid];
+
+  //signal to cut on
+  Float_t s = beta-L/integratedTimes[pid]/c;
 
   Float_t* arr = fTOFpidCuts->GetMatrixArray();
-  Int_t col = TMath::BinarySearch(fTOFpidCuts->GetNcols(),arr,static_cast<Float_t>(pt));
+  Int_t col = TMath::BinarySearch(fTOFpidCuts->GetNcols(),arr,static_cast<Float_t>(p));
   if (col<0) return kFALSE;
   Float_t min = (*fTOFpidCuts)(1,col);
   Float_t max = (*fTOFpidCuts)(2,col);
 
-  //printf("TOF pid cut %s\n",(s>min && s<max)?"PASS":"FAIL");
+  //printf("--------------TOF beta cut %s\n",(s>min && s<max)?"PASS":"FAIL");
   return (s>min && s<max);
 }
 
+//-----------------------------------------------------------------------
+Bool_t AliFlowTrackCuts::PassesTOFpidCut(AliESDtrack* track)
+{
+  //check if passes PID cut using default TOF pid
+  Double_t pidTOF[AliPID::kSPECIES];
+  track->GetTOFpid(pidTOF);
+  if (pidTOF[fParticleID]>=fParticleProbability) return kTRUE;
+  return kFALSE;
+}
+
 //-----------------------------------------------------------------------
 Bool_t AliFlowTrackCuts::PassesTPCpidCut(AliESDtrack* track)
 {
-  //check if passes PID cut using dedx signal in the TPC
-  if (!fESDpid) 
+  //check if passes PID cut using default TPC pid
+  Double_t pidTPC[AliPID::kSPECIES];
+  track->GetTPCpid(pidTPC);
+  Double_t probablity = 0.;
+  switch (fParticleID)
   {
-    return kFALSE;
+    case AliPID::kPion:
+      probablity = pidTPC[AliPID::kPion] + pidTPC[AliPID::kMuon];
+      break;
+    default:
+      probablity = pidTPC[fParticleID];
   }
+  if (probablity >= fParticleProbability) return kTRUE;
+  return kFALSE;
+}
+
+//-----------------------------------------------------------------------
+Bool_t AliFlowTrackCuts::PassesTPCdedxCut(AliESDtrack* track)
+{
+  //check if passes PID cut using dedx signal in the TPC
   if (!fTPCpidCuts)
   {
     printf("no TPCpidCuts\n");
     return kFALSE;
   }
 
-  const AliExternalTrackParam* tpcparam = track->GetInnerParam();
+  const AliExternalTrackParam* tpcparam = track->GetInnerParam(); //tpc only params at the inner wall
   if (!tpcparam) return kFALSE;
-  Float_t sigExp = fESDpid->GetTPCResponse().GetExpectedSignal(tpcparam->GetP(), fAliPID);
+  Float_t sigExp = fESDpid.GetTPCResponse().GetExpectedSignal(tpcparam->GetP(), fParticleID);
   Float_t sigTPC = track->GetTPCsignal();
   Float_t s = (sigTPC-sigExp)/sigExp;
   Double_t pt = track->Pt();
@@ -784,7 +1165,7 @@ void AliFlowTrackCuts::InitPIDcuts()
   TMatrixF* t = NULL;
   if (!fTPCpidCuts)
   {
-    if (fAliPID==AliPID::kPion)
+    if (fParticleID==AliPID::kPion)
     {
       t = new TMatrixF(3,10);
       (*t)(0,0)  = 0.20;  (*t)(1,0)  = -0.4;  (*t)(2,0)  =   0.2;
@@ -799,17 +1180,29 @@ void AliFlowTrackCuts::InitPIDcuts()
       (*t)(0,9)  = 0.65;  (*t)(1,9)  =    0;  (*t)(2,9)  =     0;
     }
     else
-    if (fAliPID==AliPID::kProton)
+    if (fParticleID==AliPID::kKaon)
+    {
+      t = new TMatrixF(3,7);
+      (*t)(0,0)  = 0.20;  (*t)(1,0)  = -0.2;  (*t)(2,0)  = 0.4; 
+      (*t)(0,1)  = 0.25;  (*t)(1,1)  =-0.15;  (*t)(2,1)  = 0.4;
+      (*t)(0,2)  = 0.30;  (*t)(1,2)  = -0.1;  (*t)(2,2)  = 0.4;
+      (*t)(0,3)  = 0.35;  (*t)(1,3)  = -0.1;  (*t)(2,3)  = 0.4;
+      (*t)(0,4)  = 0.40;  (*t)(1,4)  = -0.1;  (*t)(2,4)  = 0.6;
+      (*t)(0,5)  = 0.45;  (*t)(1,5)  = -0.1;  (*t)(2,5)  = 0.6;
+      (*t)(0,6)  = 0.50;  (*t)(1,6)  =    0;  (*t)(2,6)  =   0;
+    }
+    else
+    if (fParticleID==AliPID::kProton)
     {
       t = new TMatrixF(3,16);
-      (*t)(0,0)  = 0.20;  (*t)(1,0)  =     0;  (*t)(2,0)  =  0.3
-      (*t)(0,1)  = 0.25;  (*t)(1,1)  =  -0.2;  (*t)(2,1)  =  0.6
+      (*t)(0,0)  = 0.20;  (*t)(1,0)  =     0;  (*t)(2,0)  =    0
+      (*t)(0,1)  = 0.25;  (*t)(1,1)  =  -0.2;  (*t)(2,1)  =  0.3
       (*t)(0,2)  = 0.30;  (*t)(1,2)  =  -0.2;  (*t)(2,2)  =  0.6; 
       (*t)(0,3)  = 0.35;  (*t)(1,3)  =  -0.2;  (*t)(2,3)  =  0.6; 
       (*t)(0,4)  = 0.40;  (*t)(1,4)  =  -0.2;  (*t)(2,4)  =  0.6; 
       (*t)(0,5)  = 0.45;  (*t)(1,5)  = -0.15;  (*t)(2,5)  =  0.6; 
       (*t)(0,6)  = 0.50;  (*t)(1,6)  =  -0.1;  (*t)(2,6)  =  0.6; 
-      (*t)(0,7)  = 0.55;  (*t)(1,7)  = -0.05;  (*t)(2,7)  = 0.45
+      (*t)(0,7)  = 0.55;  (*t)(1,7)  = -0.05;  (*t)(2,7)  =  0.6
       (*t)(0,8)  = 0.60;  (*t)(1,8)  = -0.05;  (*t)(2,8)  = 0.45; 
       (*t)(0,9)  = 0.65;  (*t)(1,9)  = -0.05;  (*t)(2,9)  = 0.45; 
       (*t)(0,10) = 0.70;  (*t)(1,10) = -0.05;  (*t)(2,10) = 0.45; 
@@ -819,124 +1212,623 @@ void AliFlowTrackCuts::InitPIDcuts()
       (*t)(0,14) = 0.90;  (*t)(1,14) =     0;  (*t)(2,14) = 0.45;
       (*t)(0,15) = 0.95;  (*t)(1,15) =     0;  (*t)(2,15) =    0;
     }
-    else
-    if (fAliPID==AliPID::kKaon)
-    {
-      t = new TMatrixF(3,7);
-      (*t)(0,0)  = 0.20;  (*t)(1,0)  = -0.2;  (*t)(2,0)  = 0.4; 
-      (*t)(0,1)  = 0.25;  (*t)(1,1)  =-0.15;  (*t)(2,1)  = 0.4;
-      (*t)(0,2)  = 0.30;  (*t)(1,2)  = -0.1;  (*t)(2,2)  = 0.4;
-      (*t)(0,3)  = 0.35;  (*t)(1,3)  = -0.1;  (*t)(2,3)  = 0.4;
-      (*t)(0,4)  = 0.40;  (*t)(1,4)  = -0.1;  (*t)(2,4)  = 0.6;
-      (*t)(0,5)  = 0.45;  (*t)(1,5)  = -0.1;  (*t)(2,5)  = 0.6;
-      (*t)(0,6)  = 0.50;  (*t)(1,6)  =    0;  (*t)(2,6)  =0;
-    }
     fTPCpidCuts=t;
   }
+  t = NULL;
   if (!fTOFpidCuts)
   {
-    if (fAliPID==AliPID::kPion)
+    if (fParticleID==AliPID::kPion)
     {
-      t = new TMatrixF(3,26);
-      (*t)(0,0)  = 0.3;   (*t)(1,0)  = -700;  (*t)(2,0)  = 700;
-      (*t)(0,1)  = 0.35;  (*t)(1,1)  = -800;  (*t)(2,1)  = 800;
-      (*t)(0,2)  = 0.40;  (*t)(1,2)  = -600;  (*t)(2,2)  = 800;
-      (*t)(0,3)  = 0.45;  (*t)(1,3)  = -500;  (*t)(2,3)  = 700;
-      (*t)(0,4)  = 0.50;  (*t)(1,4)  = -400;  (*t)(2,4)  = 700;
-      (*t)(0,5)  = 0.55;  (*t)(1,5)  = -400;  (*t)(2,5)  = 700;
-      (*t)(0,6)  = 0.60;  (*t)(1,6)  = -400;  (*t)(2,6)  = 700;
-      (*t)(0,7)  = 0.65;  (*t)(1,7)  = -400;  (*t)(2,7)  = 700;
-      (*t)(0,8)  = 0.70;  (*t)(1,8)  = -400;  (*t)(2,8)  = 700;
-      (*t)(0,9)  = 0.75;  (*t)(1,9)  = -400;  (*t)(2,9)  = 700;
-      (*t)(0,10) = 0.80;  (*t)(1,10) = -400;  (*t)(2,10) = 600;
-      (*t)(0,11) = 0.85;  (*t)(1,11) = -400;  (*t)(2,11) = 600;
-      (*t)(0,12) = 0.90;  (*t)(1,12) = -400;  (*t)(2,12) = 600;
-      (*t)(0,13) = 0.95;  (*t)(1,13) = -400;  (*t)(2,13) = 600;
-      (*t)(0,14) = 1.00;  (*t)(1,14) = -400;  (*t)(2,14) = 550;
-      (*t)(0,15) = 1.10;  (*t)(1,15) = -400;  (*t)(2,15) = 450;
-      (*t)(0,16) = 1.20;  (*t)(1,16) = -400;  (*t)(2,16) = 400;
-      (*t)(0,17) = 1.30;  (*t)(1,17) = -400;  (*t)(2,17) = 300;
-      (*t)(0,18) = 1.40;  (*t)(1,18) = -400;  (*t)(2,18) = 300;
-      (*t)(0,19) = 1.50;  (*t)(1,19) = -400;  (*t)(2,19) = 250;
-      (*t)(0,20) = 1.60;  (*t)(1,20) = -400;  (*t)(2,20) = 200;
-      (*t)(0,21) = 1.70;  (*t)(1,21) = -400;  (*t)(2,21) = 150;
-      (*t)(0,22) = 1.80;  (*t)(1,22) = -400;  (*t)(2,22) = 100;
-      (*t)(0,23) = 1.90;  (*t)(1,23) = -400;  (*t)(2,23) =  70;
-      (*t)(0,24) = 2.00;  (*t)(1,24) = -400;  (*t)(2,24) =  50;
-      (*t)(0,25) = 2.10;  (*t)(1,25) =    0;  (*t)(2,25) =   0;
+      //TOF pions, 0.9 purity
+      t = new TMatrixF(3,61);
+      (*t)(0,0)  = 0.000;  (*t)(2,0)  = 0.000;  (*t)(2,0)  =   0.000;
+      (*t)(0,1)  = 0.050;  (*t)(2,1)  = 0.000;  (*t)(2,1)  =   0.000;
+      (*t)(0,2)  = 0.100;  (*t)(2,2)  = 0.000;  (*t)(2,2)  =   0.000;
+      (*t)(0,3)  = 0.150;  (*t)(2,3)  = 0.000;  (*t)(2,3)  =   0.000;
+      (*t)(0,4)  = 0.200;  (*t)(2,4)  = 0.000;  (*t)(2,4)  =   0.000;
+      (*t)(0,5)  = 0.250;  (*t)(2,5)  = -0.046;  (*t)(2,5)  =   0.046;
+      (*t)(0,6)  = 0.300;  (*t)(2,6)  = -0.038;  (*t)(2,6)  =   0.038;
+      (*t)(0,7)  = 0.350;  (*t)(2,7)  = -0.034;  (*t)(2,7)  =   0.034;
+      (*t)(0,8)  = 0.400;  (*t)(2,8)  = -0.032;  (*t)(2,8)  =   0.032;
+      (*t)(0,9)  = 0.450;  (*t)(2,9)  = -0.030;  (*t)(2,9)  =   0.030;
+      (*t)(0,10)  = 0.500;  (*t)(2,10)  = -0.030;  (*t)(2,10)  =   0.030;
+      (*t)(0,11)  = 0.550;  (*t)(2,11)  = -0.030;  (*t)(2,11)  =   0.030;
+      (*t)(0,12)  = 0.600;  (*t)(2,12)  = -0.030;  (*t)(2,12)  =   0.030;
+      (*t)(0,13)  = 0.650;  (*t)(2,13)  = -0.030;  (*t)(2,13)  =   0.030;
+      (*t)(0,14)  = 0.700;  (*t)(2,14)  = -0.030;  (*t)(2,14)  =   0.030;
+      (*t)(0,15)  = 0.750;  (*t)(2,15)  = -0.030;  (*t)(2,15)  =   0.030;
+      (*t)(0,16)  = 0.800;  (*t)(2,16)  = -0.030;  (*t)(2,16)  =   0.030;
+      (*t)(0,17)  = 0.850;  (*t)(2,17)  = -0.030;  (*t)(2,17)  =   0.030;
+      (*t)(0,18)  = 0.900;  (*t)(2,18)  = -0.030;  (*t)(2,18)  =   0.030;
+      (*t)(0,19)  = 0.950;  (*t)(2,19)  = -0.028;  (*t)(2,19)  =   0.028;
+      (*t)(0,20)  = 1.000;  (*t)(2,20)  = -0.028;  (*t)(2,20)  =   0.028;
+      (*t)(0,21)  = 1.100;  (*t)(2,21)  = -0.028;  (*t)(2,21)  =   0.028;
+      (*t)(0,22)  = 1.200;  (*t)(2,22)  = -0.026;  (*t)(2,22)  =   0.028;
+      (*t)(0,23)  = 1.300;  (*t)(2,23)  = -0.024;  (*t)(2,23)  =   0.028;
+      (*t)(0,24)  = 1.400;  (*t)(2,24)  = -0.020;  (*t)(2,24)  =   0.028;
+      (*t)(0,25)  = 1.500;  (*t)(2,25)  = -0.018;  (*t)(2,25)  =   0.028;
+      (*t)(0,26)  = 1.600;  (*t)(2,26)  = -0.016;  (*t)(2,26)  =   0.028;
+      (*t)(0,27)  = 1.700;  (*t)(2,27)  = -0.014;  (*t)(2,27)  =   0.028;
+      (*t)(0,28)  = 1.800;  (*t)(2,28)  = -0.012;  (*t)(2,28)  =   0.026;
+      (*t)(0,29)  = 1.900;  (*t)(2,29)  = -0.010;  (*t)(2,29)  =   0.026;
+      (*t)(0,30)  = 2.000;  (*t)(2,30)  = -0.008;  (*t)(2,30)  =   0.026;
+      (*t)(0,31)  = 2.100;  (*t)(2,31)  = -0.008;  (*t)(2,31)  =   0.024;
+      (*t)(0,32)  = 2.200;  (*t)(2,32)  = -0.006;  (*t)(2,32)  =   0.024;
+      (*t)(0,33)  = 2.300;  (*t)(2,33)  = -0.004;  (*t)(2,33)  =   0.024;
+      (*t)(0,34)  = 2.400;  (*t)(2,34)  = -0.004;  (*t)(2,34)  =   0.024;
+      (*t)(0,35)  = 2.500;  (*t)(2,35)  = -0.002;  (*t)(2,35)  =   0.024;
+      (*t)(0,36)  = 2.600;  (*t)(2,36)  = -0.002;  (*t)(2,36)  =   0.024;
+      (*t)(0,37)  = 2.700;  (*t)(2,37)  = 0.000;  (*t)(2,37)  =   0.024;
+      (*t)(0,38)  = 2.800;  (*t)(2,38)  = 0.000;  (*t)(2,38)  =   0.026;
+      (*t)(0,39)  = 2.900;  (*t)(2,39)  = 0.000;  (*t)(2,39)  =   0.024;
+      (*t)(0,40)  = 3.000;  (*t)(2,40)  = 0.002;  (*t)(2,40)  =   0.026;
+      (*t)(0,41)  = 3.100;  (*t)(2,41)  = 0.002;  (*t)(2,41)  =   0.026;
+      (*t)(0,42)  = 3.200;  (*t)(2,42)  = 0.002;  (*t)(2,42)  =   0.026;
+      (*t)(0,43)  = 3.300;  (*t)(2,43)  = 0.002;  (*t)(2,43)  =   0.026;
+      (*t)(0,44)  = 3.400;  (*t)(2,44)  = 0.002;  (*t)(2,44)  =   0.026;
+      (*t)(0,45)  = 3.500;  (*t)(2,45)  = 0.002;  (*t)(2,45)  =   0.026;
+      (*t)(0,46)  = 3.600;  (*t)(2,46)  = 0.002;  (*t)(2,46)  =   0.026;
+      (*t)(0,47)  = 3.700;  (*t)(2,47)  = 0.002;  (*t)(2,47)  =   0.026;
+      (*t)(0,48)  = 3.800;  (*t)(2,48)  = 0.002;  (*t)(2,48)  =   0.026;
+      (*t)(0,49)  = 3.900;  (*t)(2,49)  = 0.004;  (*t)(2,49)  =   0.024;
+      (*t)(0,50)  = 4.000;  (*t)(2,50)  = 0.004;  (*t)(2,50)  =   0.026;
+      (*t)(0,51)  = 4.100;  (*t)(2,51)  = 0.004;  (*t)(2,51)  =   0.026;
+      (*t)(0,52)  = 4.200;  (*t)(2,52)  = 0.004;  (*t)(2,52)  =   0.024;
+      (*t)(0,53)  = 4.300;  (*t)(2,53)  = 0.006;  (*t)(2,53)  =   0.024;
+      (*t)(0,54)  = 4.400;  (*t)(2,54)  = 0.000;  (*t)(2,54)  =   0.000;
+      (*t)(0,55)  = 4.500;  (*t)(2,55)  = 0.000;  (*t)(2,55)  =   0.000;
+      (*t)(0,56)  = 4.600;  (*t)(2,56)  = 0.000;  (*t)(2,56)  =   0.000;
+      (*t)(0,57)  = 4.700;  (*t)(2,57)  = 0.000;  (*t)(2,57)  =   0.000;
+      (*t)(0,58)  = 4.800;  (*t)(2,58)  = 0.000;  (*t)(2,58)  =   0.000;
+      (*t)(0,59)  = 4.900;  (*t)(2,59)  = 0.000;  (*t)(2,59)  =   0.000;
+      (*t)(0,60)  = 5.900;  (*t)(2,60)  = 0.000;  (*t)(2,60)  =   0.000;
     }
     else
-    if (fAliPID==AliPID::kProton)
+    if (fParticleID==AliPID::kProton)
     {
-      t = new TMatrixF(3,39);
-      (*t)(0,0)  = 0.3;  (*t)(1,0)   = 0; (*t)(2,0) = 0;
-      (*t)(0,1)  = 0.35;  (*t)(1,1)  = 0; (*t)(2,1) = 0;
-      (*t)(0,2)  = 0.40;  (*t)(1,2)  = 0; (*t)(2,2) = 0;
-      (*t)(0,3)  = 0.45;  (*t)(1,3)  = 0; (*t)(2,3) = 0;
-      (*t)(0,4)  = 0.50;  (*t)(1,4)  = 0; (*t)(2,4) = 0;
-      (*t)(0,5)  = 0.55;  (*t)(1,5)  = -900;  (*t)(2,5)  = 600;
-      (*t)(0,6)  = 0.60;  (*t)(1,6)  = -800;  (*t)(2,6)  = 600;
-      (*t)(0,7)  = 0.65;  (*t)(1,7)  = -800;  (*t)(2,7)  = 600;
-      (*t)(0,8)  = 0.70;  (*t)(1,8)  = -800;  (*t)(2,8)  = 600;
-      (*t)(0,9)  = 0.75;  (*t)(1,9)  = -700;  (*t)(2,9)  = 500;
-      (*t)(0,10) = 0.80;  (*t)(1,10) = -700;  (*t)(2,10) = 500;
-      (*t)(0,11) = 0.85;  (*t)(1,11) = -700;  (*t)(2,11) = 500;
-      (*t)(0,12) = 0.90;  (*t)(1,12) = -600;  (*t)(2,12) = 500;
-      (*t)(0,13) = 0.95;  (*t)(1,13) = -600;  (*t)(2,13) = 500;
-      (*t)(0,14) = 1.00;  (*t)(1,14) = -600;  (*t)(2,14) = 500;
-      (*t)(0,15) = 1.10;  (*t)(1,15) = -600;  (*t)(2,15) = 500;
-      (*t)(0,16) = 1.20;  (*t)(1,16) = -500;  (*t)(2,16) = 500;
-      (*t)(0,17) = 1.30;  (*t)(1,17) = -500;  (*t)(2,17) = 500;
-      (*t)(0,18) = 1.40;  (*t)(1,18) = -500;  (*t)(2,18) = 500;
-      (*t)(0,19) = 1.50;  (*t)(1,19) = -500;  (*t)(2,19) = 500;
-      (*t)(0,20) = 1.60;  (*t)(1,20) = -400;  (*t)(2,20) = 500;
-      (*t)(0,21) = 1.70;  (*t)(1,21) = -400;  (*t)(2,21) = 500;
-      (*t)(0,22) = 1.80;  (*t)(1,22) = -400;  (*t)(2,22) = 500;
-      (*t)(0,23) = 1.90;  (*t)(1,23) = -400;  (*t)(2,23) = 500;
-      (*t)(0,24) = 2.00;  (*t)(1,24) = -400;  (*t)(2,24) = 500;
-      (*t)(0,25) = 2.10;  (*t)(1,25) = -350;  (*t)(2,25) = 500;
-      (*t)(0,26) = 2.20;  (*t)(1,26) = -350;  (*t)(2,26) = 500;
-      (*t)(0,27) = 2.30;  (*t)(1,27) = -300;  (*t)(2,27) = 500;
-      (*t)(0,28) = 2.40;  (*t)(1,28) = -300;  (*t)(2,28) = 500;
-      (*t)(0,29) = 3.30;  (*t)(1,29) = -300;  (*t)(2,29) = 500;
-      (*t)(0,30) = 3.30;  (*t)(1,30) = -250;  (*t)(2,30) = 500;
-      (*t)(0,31) = 3.30;  (*t)(1,31) = -200;  (*t)(2,31) = 500;
-      (*t)(0,32) = 3.30;  (*t)(1,32) = -150;  (*t)(2,32) = 500;
-      (*t)(0,33) = 3.30;  (*t)(1,33) = -150;  (*t)(2,33) = 500;
-      (*t)(0,34) = 3.30;  (*t)(1,34) = -100;  (*t)(2,34) = 400;
-      (*t)(0,35) = 3.30;  (*t)(1,35) = -100;  (*t)(2,35) = 400;
-      (*t)(0,36) = 3.30;  (*t)(1,36) =    0;  (*t)(2,36) = 0;
-      (*t)(0,37) = 3.30;  (*t)(1,37) =    0;  (*t)(2,37) = 0;
-      (*t)(0,38) = 3.30;  (*t)(1,38) =    0;  (*t)(2,38) = 0;
+      //TOF protons, 0.9 purity
+      t = new TMatrixF(3,61);
+      (*t)(0,0)  = 0.000;  (*t)(2,0)  = 0.000;  (*t)(2,0)  =   0.000;
+      (*t)(0,1)  = 0.050;  (*t)(2,1)  = 0.000;  (*t)(2,1)  =   0.000;
+      (*t)(0,2)  = 0.100;  (*t)(2,2)  = 0.000;  (*t)(2,2)  =   0.000;
+      (*t)(0,3)  = 0.150;  (*t)(2,3)  = 0.000;  (*t)(2,3)  =   0.000;
+      (*t)(0,4)  = 0.200;  (*t)(2,4)  = 0.000;  (*t)(2,4)  =   0.000;
+      (*t)(0,5)  = 0.250;  (*t)(2,5)  = 0.000;  (*t)(2,5)  =   0.000;
+      (*t)(0,6)  = 0.300;  (*t)(2,6)  = 0.000;  (*t)(2,6)  =   0.000;
+      (*t)(0,7)  = 0.350;  (*t)(2,7)  = 0.000;  (*t)(2,7)  =   0.000;
+      (*t)(0,8)  = 0.400;  (*t)(2,8)  = 0.000;  (*t)(2,8)  =   0.000;
+      (*t)(0,9)  = 0.450;  (*t)(2,9)  = 0.000;  (*t)(2,9)  =   0.000;
+      (*t)(0,10)  = 0.500;  (*t)(2,10)  = 0.000;  (*t)(2,10)  =   0.000;
+      (*t)(0,11)  = 0.550;  (*t)(2,11)  = 0.000;  (*t)(2,11)  =   0.000;
+      (*t)(0,12)  = 0.600;  (*t)(2,12)  = 0.000;  (*t)(2,12)  =   0.000;
+      (*t)(0,13)  = 0.650;  (*t)(2,13)  = 0.000;  (*t)(2,13)  =   0.000;
+      (*t)(0,14)  = 0.700;  (*t)(2,14)  = 0.000;  (*t)(2,14)  =   0.000;
+      (*t)(0,15)  = 0.750;  (*t)(2,15)  = 0.000;  (*t)(2,15)  =   0.000;
+      (*t)(0,16)  = 0.800;  (*t)(2,16)  = 0.000;  (*t)(2,16)  =   0.000;
+      (*t)(0,17)  = 0.850;  (*t)(2,17)  = -0.070;  (*t)(2,17)  =   0.070;
+      (*t)(0,18)  = 0.900;  (*t)(2,18)  = -0.072;  (*t)(2,18)  =   0.072;
+      (*t)(0,19)  = 0.950;  (*t)(2,19)  = -0.072;  (*t)(2,19)  =   0.072;
+      (*t)(0,20)  = 1.000;  (*t)(2,20)  = -0.074;  (*t)(2,20)  =   0.074;
+      (*t)(0,21)  = 1.100;  (*t)(2,21)  = -0.032;  (*t)(2,21)  =   0.032;
+      (*t)(0,22)  = 1.200;  (*t)(2,22)  = -0.026;  (*t)(2,22)  =   0.026;
+      (*t)(0,23)  = 1.300;  (*t)(2,23)  = -0.026;  (*t)(2,23)  =   0.026;
+      (*t)(0,24)  = 1.400;  (*t)(2,24)  = -0.024;  (*t)(2,24)  =   0.024;
+      (*t)(0,25)  = 1.500;  (*t)(2,25)  = -0.024;  (*t)(2,25)  =   0.024;
+      (*t)(0,26)  = 1.600;  (*t)(2,26)  = -0.026;  (*t)(2,26)  =   0.026;
+      (*t)(0,27)  = 1.700;  (*t)(2,27)  = -0.026;  (*t)(2,27)  =   0.026;
+      (*t)(0,28)  = 1.800;  (*t)(2,28)  = -0.026;  (*t)(2,28)  =   0.026;
+      (*t)(0,29)  = 1.900;  (*t)(2,29)  = -0.026;  (*t)(2,29)  =   0.026;
+      (*t)(0,30)  = 2.000;  (*t)(2,30)  = -0.026;  (*t)(2,30)  =   0.026;
+      (*t)(0,31)  = 2.100;  (*t)(2,31)  = -0.026;  (*t)(2,31)  =   0.026;
+      (*t)(0,32)  = 2.200;  (*t)(2,32)  = -0.026;  (*t)(2,32)  =   0.024;
+      (*t)(0,33)  = 2.300;  (*t)(2,33)  = -0.028;  (*t)(2,33)  =   0.022;
+      (*t)(0,34)  = 2.400;  (*t)(2,34)  = -0.028;  (*t)(2,34)  =   0.020;
+      (*t)(0,35)  = 2.500;  (*t)(2,35)  = -0.028;  (*t)(2,35)  =   0.018;
+      (*t)(0,36)  = 2.600;  (*t)(2,36)  = -0.028;  (*t)(2,36)  =   0.016;
+      (*t)(0,37)  = 2.700;  (*t)(2,37)  = -0.028;  (*t)(2,37)  =   0.016;
+      (*t)(0,38)  = 2.800;  (*t)(2,38)  = -0.030;  (*t)(2,38)  =   0.014;
+      (*t)(0,39)  = 2.900;  (*t)(2,39)  = -0.030;  (*t)(2,39)  =   0.012;
+      (*t)(0,40)  = 3.000;  (*t)(2,40)  = -0.030;  (*t)(2,40)  =   0.012;
+      (*t)(0,41)  = 3.100;  (*t)(2,41)  = -0.030;  (*t)(2,41)  =   0.010;
+      (*t)(0,42)  = 3.200;  (*t)(2,42)  = -0.030;  (*t)(2,42)  =   0.010;
+      (*t)(0,43)  = 3.300;  (*t)(2,43)  = -0.030;  (*t)(2,43)  =   0.010;
+      (*t)(0,44)  = 3.400;  (*t)(2,44)  = -0.030;  (*t)(2,44)  =   0.008;
+      (*t)(0,45)  = 3.500;  (*t)(2,45)  = -0.030;  (*t)(2,45)  =   0.008;
+      (*t)(0,46)  = 3.600;  (*t)(2,46)  = -0.030;  (*t)(2,46)  =   0.008;
+      (*t)(0,47)  = 3.700;  (*t)(2,47)  = -0.030;  (*t)(2,47)  =   0.006;
+      (*t)(0,48)  = 3.800;  (*t)(2,48)  = -0.030;  (*t)(2,48)  =   0.006;
+      (*t)(0,49)  = 3.900;  (*t)(2,49)  = -0.030;  (*t)(2,49)  =   0.006;
+      (*t)(0,50)  = 4.000;  (*t)(2,50)  = -0.028;  (*t)(2,50)  =   0.004;
+      (*t)(0,51)  = 4.100;  (*t)(2,51)  = -0.030;  (*t)(2,51)  =   0.004;
+      (*t)(0,52)  = 4.200;  (*t)(2,52)  = -0.030;  (*t)(2,52)  =   0.004;
+      (*t)(0,53)  = 4.300;  (*t)(2,53)  = -0.028;  (*t)(2,53)  =   0.002;
+      (*t)(0,54)  = 4.400;  (*t)(2,54)  = -0.030;  (*t)(2,54)  =   0.002;
+      (*t)(0,55)  = 4.500;  (*t)(2,55)  = -0.028;  (*t)(2,55)  =   0.002;
+      (*t)(0,56)  = 4.600;  (*t)(2,56)  = -0.028;  (*t)(2,56)  =   0.002;
+      (*t)(0,57)  = 4.700;  (*t)(2,57)  = -0.028;  (*t)(2,57)  =   0.000;
+      (*t)(0,58)  = 4.800;  (*t)(2,58)  = -0.028;  (*t)(2,58)  =   0.002;
+      (*t)(0,59)  = 4.900;  (*t)(2,59)  = 0.000;  (*t)(2,59)  =   0.000;
+      (*t)(0,60)  = 5.900;  (*t)(2,60)  = 0.000;  (*t)(2,60)  =   0.000; 
     }
     else
-    if (fAliPID==AliPID::kKaon)
+    if (fParticleID==AliPID::kKaon)
     {
-      t = new TMatrixF(3,23);
-      (*t)(0,0)  = 0.3;   (*t)(1,0)  =    0;  (*t)(2,0)  =    0;
-      (*t)(0,1)  = 0.35;  (*t)(1,1)  =    0;  (*t)(2,1)  =    0;
-      (*t)(0,2)  = 0.40;  (*t)(1,2)  = -800;  (*t)(2,2)  =  600;
-      (*t)(0,3)  = 0.45;  (*t)(1,3)  = -800;  (*t)(2,3)  =  600;
-      (*t)(0,4)  = 0.50;  (*t)(1,4)  = -800;  (*t)(2,4)  =  600;
-      (*t)(0,5)  = 0.55;  (*t)(1,5)  = -800;  (*t)(2,5)  =  600;
-      (*t)(0,6)  = 0.60;  (*t)(1,6)  = -800;  (*t)(2,6)  =  600;
-      (*t)(0,7)  = 0.65;  (*t)(1,7)  = -700;  (*t)(2,7)  =  600;
-      (*t)(0,8)  = 0.70;  (*t)(1,8)  = -600;  (*t)(2,8)  =  600;
-      (*t)(0,9)  = 0.75;  (*t)(1,9)  = -600;  (*t)(2,9)  =  500;
-      (*t)(0,10) = 0.80;  (*t)(1,10) = -500;  (*t)(2,10) =  500;
-      (*t)(0,11) = 0.85;  (*t)(1,11) = -500;  (*t)(2,11) =  500;
-      (*t)(0,12) = 0.90;  (*t)(1,12) = -400;  (*t)(2,12) =  500;
-      (*t)(0,13) = 0.95;  (*t)(1,13) = -400;  (*t)(2,13) =  500;
-      (*t)(0,14) = 1.00;  (*t)(1,14) = -400;  (*t)(2,14) =  500;
-      (*t)(0,15) = 1.10;  (*t)(1,15) = -350;  (*t)(2,15) =  450;
-      (*t)(0,16) = 1.20;  (*t)(1,16) = -300;  (*t)(2,16) =  400;
-      (*t)(0,17) = 1.30;  (*t)(1,17) = -300;  (*t)(2,17) =  400;
-      (*t)(0,18) = 1.40;  (*t)(1,18) = -250;  (*t)(2,18) =  400;
-      (*t)(0,19) = 1.50;  (*t)(1,19) = -200;  (*t)(2,19) =  400;
-      (*t)(0,20) = 1.60;  (*t)(1,20) = -150;  (*t)(2,20) =  400;
-      (*t)(0,21) = 1.70;  (*t)(1,21) = -100;  (*t)(2,21) =  400;
-      (*t)(0,22) = 1.80;  (*t)(1,22) =    0;  (*t)(2,22) =    0;
+      //TOF kaons, 0.9 purity
+      t = new TMatrixF(3,61);
+      (*t)(0,0)  = 0.000;  (*t)(2,0)  = 0.000;  (*t)(2,0)  =   0.000;
+      (*t)(0,1)  = 0.050;  (*t)(2,1)  = 0.000;  (*t)(2,1)  =   0.000;
+      (*t)(0,2)  = 0.100;  (*t)(2,2)  = 0.000;  (*t)(2,2)  =   0.000;
+      (*t)(0,3)  = 0.150;  (*t)(2,3)  = 0.000;  (*t)(2,3)  =   0.000;
+      (*t)(0,4)  = 0.200;  (*t)(2,4)  = 0.000;  (*t)(2,4)  =   0.000;
+      (*t)(0,5)  = 0.250;  (*t)(2,5)  = 0.000;  (*t)(2,5)  =   0.000;
+      (*t)(0,6)  = 0.300;  (*t)(2,6)  = 0.000;  (*t)(2,6)  =   0.000;
+      (*t)(0,7)  = 0.350;  (*t)(2,7)  = 0.000;  (*t)(2,7)  =   0.000;
+      (*t)(0,8)  = 0.400;  (*t)(2,8)  = 0.000;  (*t)(2,8)  =   0.000;
+      (*t)(0,9)  = 0.450;  (*t)(2,9)  = 0.000;  (*t)(2,9)  =   0.000;
+      (*t)(0,10)  = 0.500;  (*t)(2,10)  = 0.000;  (*t)(2,10)  =   0.000;
+      (*t)(0,11)  = 0.550;  (*t)(2,11)  = -0.026;  (*t)(2,11)  =   0.026;
+      (*t)(0,12)  = 0.600;  (*t)(2,12)  = -0.026;  (*t)(2,12)  =   0.026;
+      (*t)(0,13)  = 0.650;  (*t)(2,13)  = -0.026;  (*t)(2,13)  =   0.026;
+      (*t)(0,14)  = 0.700;  (*t)(2,14)  = -0.026;  (*t)(2,14)  =   0.026;
+      (*t)(0,15)  = 0.750;  (*t)(2,15)  = -0.026;  (*t)(2,15)  =   0.026;
+      (*t)(0,16)  = 0.800;  (*t)(2,16)  = -0.026;  (*t)(2,16)  =   0.026;
+      (*t)(0,17)  = 0.850;  (*t)(2,17)  = -0.024;  (*t)(2,17)  =   0.024;
+      (*t)(0,18)  = 0.900;  (*t)(2,18)  = -0.024;  (*t)(2,18)  =   0.024;
+      (*t)(0,19)  = 0.950;  (*t)(2,19)  = -0.024;  (*t)(2,19)  =   0.024;
+      (*t)(0,20)  = 1.000;  (*t)(2,20)  = -0.024;  (*t)(2,20)  =   0.024;
+      (*t)(0,21)  = 1.100;  (*t)(2,21)  = -0.024;  (*t)(2,21)  =   0.024;
+      (*t)(0,22)  = 1.200;  (*t)(2,22)  = -0.024;  (*t)(2,22)  =   0.022;
+      (*t)(0,23)  = 1.300;  (*t)(2,23)  = -0.024;  (*t)(2,23)  =   0.020;
+      (*t)(0,24)  = 1.400;  (*t)(2,24)  = -0.026;  (*t)(2,24)  =   0.016;
+      (*t)(0,25)  = 1.500;  (*t)(2,25)  = -0.028;  (*t)(2,25)  =   0.014;
+      (*t)(0,26)  = 1.600;  (*t)(2,26)  = -0.028;  (*t)(2,26)  =   0.012;
+      (*t)(0,27)  = 1.700;  (*t)(2,27)  = -0.028;  (*t)(2,27)  =   0.010;
+      (*t)(0,28)  = 1.800;  (*t)(2,28)  = -0.028;  (*t)(2,28)  =   0.010;
+      (*t)(0,29)  = 1.900;  (*t)(2,29)  = -0.028;  (*t)(2,29)  =   0.008;
+      (*t)(0,30)  = 2.000;  (*t)(2,30)  = -0.028;  (*t)(2,30)  =   0.006;
+      (*t)(0,31)  = 2.100;  (*t)(2,31)  = -0.026;  (*t)(2,31)  =   0.006;
+      (*t)(0,32)  = 2.200;  (*t)(2,32)  = -0.024;  (*t)(2,32)  =   0.004;
+      (*t)(0,33)  = 2.300;  (*t)(2,33)  = -0.020;  (*t)(2,33)  =   0.002;
+      (*t)(0,34)  = 2.400;  (*t)(2,34)  = -0.020;  (*t)(2,34)  =   0.002;
+      (*t)(0,35)  = 2.500;  (*t)(2,35)  = -0.018;  (*t)(2,35)  =   0.000;
+      (*t)(0,36)  = 2.600;  (*t)(2,36)  = -0.016;  (*t)(2,36)  =   0.000;
+      (*t)(0,37)  = 2.700;  (*t)(2,37)  = -0.014;  (*t)(2,37)  =   -0.002;
+      (*t)(0,38)  = 2.800;  (*t)(2,38)  = -0.014;  (*t)(2,38)  =   -0.004;
+      (*t)(0,39)  = 2.900;  (*t)(2,39)  = -0.012;  (*t)(2,39)  =   -0.004;
+      (*t)(0,40)  = 3.000;  (*t)(2,40)  = -0.010;  (*t)(2,40)  =   -0.006;
+      (*t)(0,41)  = 3.100;  (*t)(2,41)  = 0.000;  (*t)(2,41)  =   0.000;
+      (*t)(0,42)  = 3.200;  (*t)(2,42)  = 0.000;  (*t)(2,42)  =   0.000;
+      (*t)(0,43)  = 3.300;  (*t)(2,43)  = 0.000;  (*t)(2,43)  =   0.000;
+      (*t)(0,44)  = 3.400;  (*t)(2,44)  = 0.000;  (*t)(2,44)  =   0.000;
+      (*t)(0,45)  = 3.500;  (*t)(2,45)  = 0.000;  (*t)(2,45)  =   0.000;
+      (*t)(0,46)  = 3.600;  (*t)(2,46)  = 0.000;  (*t)(2,46)  =   0.000;
+      (*t)(0,47)  = 3.700;  (*t)(2,47)  = 0.000;  (*t)(2,47)  =   0.000;
+      (*t)(0,48)  = 3.800;  (*t)(2,48)  = 0.000;  (*t)(2,48)  =   0.000;
+      (*t)(0,49)  = 3.900;  (*t)(2,49)  = 0.000;  (*t)(2,49)  =   0.000;
+      (*t)(0,50)  = 4.000;  (*t)(2,50)  = 0.000;  (*t)(2,50)  =   0.000;
+      (*t)(0,51)  = 4.100;  (*t)(2,51)  = 0.000;  (*t)(2,51)  =   0.000;
+      (*t)(0,52)  = 4.200;  (*t)(2,52)  = 0.000;  (*t)(2,52)  =   0.000;
+      (*t)(0,53)  = 4.300;  (*t)(2,53)  = 0.000;  (*t)(2,53)  =   0.000;
+      (*t)(0,54)  = 4.400;  (*t)(2,54)  = 0.000;  (*t)(2,54)  =   0.000;
+      (*t)(0,55)  = 4.500;  (*t)(2,55)  = 0.000;  (*t)(2,55)  =   0.000;
+      (*t)(0,56)  = 4.600;  (*t)(2,56)  = 0.000;  (*t)(2,56)  =   0.000;
+      (*t)(0,57)  = 4.700;  (*t)(2,57)  = 0.000;  (*t)(2,57)  =   0.000;
+      (*t)(0,58)  = 4.800;  (*t)(2,58)  = 0.000;  (*t)(2,58)  =   0.000;
+      (*t)(0,59)  = 4.900;  (*t)(2,59)  = 0.000;  (*t)(2,59)  =   0.000;
+      (*t)(0,60)  = 5.900;  (*t)(2,60)  = 0.000;  (*t)(2,60)  =   0.000;
     }
     fTOFpidCuts=t;
   }
 }
+
+//-----------------------------------------------------------------------
+// part added by F. Noferini (some methods)
+Bool_t AliFlowTrackCuts::PassesTOFbayesianCut(AliESDtrack* track){
+
+  Bool_t goodtrack = track && (track->GetStatus() & AliESDtrack::kTOFpid) && (track->GetTOFsignal() > 12000) && (track->GetTOFsignal() < 100000) && (track->GetIntegratedLength() > 365) && !(track->GetStatus() & AliESDtrack::kTOFmismatch);
+
+  if (! goodtrack)
+       return kFALSE;
+
+  Int_t pdg = GetESDPdg(track,"bayesianALL");
+  //  printf("pdg set to %i\n",pdg);
+
+  Int_t pid = 0;
+  Float_t prob = 0;
+  switch (fParticleID)
+  {
+    case AliPID::kPion:
+      pid=211;
+      prob = fProbBayes[2];
+      break;
+    case AliPID::kKaon:
+      pid=321;
+      prob = fProbBayes[3];
+     break;
+    case AliPID::kProton:
+      pid=2212;
+      prob = fProbBayes[4];
+      break;
+    case AliPID::kElectron:
+      pid=-11;
+       prob = fProbBayes[0];
+     break;
+    default:
+      return kFALSE;
+  }
+
+  //  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);
+  if(TMath::Abs(pdg) == TMath::Abs(pid) && prob > 0.8){
+    if(!fCutCharge)
+      return kTRUE;
+    else if (fCutCharge && fCharge * track->GetSign() > 0)
+      return kTRUE;
+  }
+  return kFALSE;
+}
+//-----------------------------------------------------------------------
+Int_t AliFlowTrackCuts::GetESDPdg(AliESDtrack *track,Option_t *option,Int_t ipart,Float_t cPi,Float_t cKa,Float_t cPr){
+  Int_t pdg = 0;
+  Int_t pdgvalues[5] = {-11,-13,211,321,2212};
+  Float_t mass[5] = {5.10998909999999971e-04,1.05658000000000002e-01,1.39570000000000000e-01,4.93676999999999977e-01,9.38271999999999995e-01};
+
+  if(strstr(option,"bayesianTOF")){ // Bayesian TOF PID
+    Double_t c[5]={0.01, 0.01, 0.85, 0.1, 0.05};
+    Double_t rcc=0.;
+    
+    Float_t pt = track->Pt();
+    
+    Int_t iptesd = 0;
+    while(pt > fBinLimitPID[iptesd] && iptesd < fnPIDptBin-1) iptesd++;
+  
+    if(cPi < 0){
+      c[0] = fC[iptesd][0];
+      c[1] = fC[iptesd][1];
+      c[2] = fC[iptesd][2];
+      c[3] = fC[iptesd][3];
+      c[4] = fC[iptesd][4];
+    }
+    else{
+      c[0] = 0.0;
+      c[1] = 0.0;
+      c[2] = cPi;
+      c[3] = cKa;
+      c[4] = cPr;      
+    }
+
+    Double_t r1[10]; track->GetTOFpid(r1);
+    
+    Int_t i;
+    for (i=0; i<5; i++) rcc+=(c[i]*r1[i]);
+    
+    Double_t w[10];
+    for (i=0; i<5; i++){
+       w[i]=c[i]*r1[i]/rcc;
+       fProbBayes[i] = w[i];
+    }
+    if (w[2]>=w[3] && w[2]>=w[4] && w[2]>=w[1] && w[2]>=w[0]) {//pion
+      pdg = 211*Int_t(track->GetSign());
+    }
+    else if (w[4]>=w[3] && w[4]>=w[1] && w[4]>=w[0]) {//proton
+      pdg = 2212*Int_t(track->GetSign());
+    }
+    else if (w[3]>=w[1] && w[3]>=w[0]){//kaon
+      pdg = 321*Int_t(track->GetSign());
+    }
+    else if (w[0]>=w[1]) { //electrons
+      pdg = -11*Int_t(track->GetSign());
+    }
+    else{ // muon
+      pdg = -13*Int_t(track->GetSign());
+    }
+  }
+
+  else if(strstr(option,"bayesianTPC")){ // Bayesian TPC PID
+    Double_t c[5]={0.01, 0.01, 0.85, 0.1, 0.05};
+    Double_t rcc=0.;
+    
+    Float_t pt = track->Pt();
+    
+    Int_t iptesd = 0;
+    while(pt > fBinLimitPID[iptesd] && iptesd < fnPIDptBin-1) iptesd++;
+  
+    if(cPi < 0){
+      c[0] = fC[iptesd][0];
+      c[1] = fC[iptesd][1];
+      c[2] = fC[iptesd][2];
+      c[3] = fC[iptesd][3];
+      c[4] = fC[iptesd][4];
+    }
+    else{
+      c[0] = 0.0;
+      c[1] = 0.0;
+      c[2] = cPi;
+      c[3] = cKa;
+      c[4] = cPr;      
+    }
+
+    Double_t r1[10]; track->GetTPCpid(r1);
+    
+    Int_t i;
+    for (i=0; i<5; i++) rcc+=(c[i]*r1[i]);
+    
+    Double_t w[10];
+    for (i=0; i<5; i++){
+       w[i]=c[i]*r1[i]/rcc;
+       fProbBayes[i] = w[i];
+    }
+    if (w[2]>=w[3] && w[2]>=w[4] && w[2]>=w[1] && w[2]>=w[0]) {//pion
+      pdg = 211*Int_t(track->GetSign());
+    }
+    else if (w[4]>=w[3] && w[4]>=w[1] && w[4]>=w[0]) {//proton
+      pdg = 2212*Int_t(track->GetSign());
+    }
+    else if (w[3]>=w[1] && w[3]>=w[0]){//kaon
+      pdg = 321*Int_t(track->GetSign());
+    }
+    else if (w[0]>=w[1]) { //electrons
+      pdg = -11*Int_t(track->GetSign());
+    }
+    else{ // muon
+      pdg = -13*Int_t(track->GetSign());
+    }
+  }
+  
+  else if(strstr(option,"bayesianALL")){
+    Double_t c[5]={0.01, 0.01, 0.85, 0.1, 0.05};
+    Double_t rcc=0.;
+    
+    Float_t pt = track->Pt();
+    
+    Int_t iptesd = 0;
+    while(pt > fBinLimitPID[iptesd] && iptesd < fnPIDptBin-1) iptesd++;
+
+    if(cPi < 0){
+      c[0] = fC[iptesd][0];
+      c[1] = fC[iptesd][1];
+      c[2] = fC[iptesd][2];
+      c[3] = fC[iptesd][3];
+      c[4] = fC[iptesd][4];
+    }
+    else{
+      c[0] = 0.0;
+      c[1] = 0.0;
+      c[2] = cPi;
+      c[3] = cKa;
+      c[4] = cPr;      
+    }
+
+    Double_t r1[10]; track->GetTOFpid(r1);
+    Double_t r2[10]; track->GetTPCpid(r2);
+
+    Int_t i;
+    for (i=0; i<5; i++) rcc+=(c[i]*r1[i]*r2[i]);
+    
+
+    Double_t w[10];
+    for (i=0; i<5; i++){
+       w[i]=c[i]*r1[i]*r2[i]/rcc;
+       fProbBayes[i] = w[i];
+    }
+
+    if (w[2]>=w[3] && w[2]>=w[4] && w[2]>=w[1] && w[2]>=w[0]) {//pion
+      pdg = 211*Int_t(track->GetSign());
+    }
+    else if (w[4]>=w[3] && w[4]>=w[1] && w[4]>=w[0]) {//proton
+      pdg = 2212*Int_t(track->GetSign());
+    }
+    else if (w[3]>=w[1] && w[3]>=w[0]){//kaon
+      pdg = 321*Int_t(track->GetSign());
+    }
+    else if (w[0]>=w[1]) { //electrons
+      pdg = -11*Int_t(track->GetSign());
+    }
+    else{ // muon
+      pdg = -13*Int_t(track->GetSign());
+    }
+  }
+
+  else if(strstr(option,"sigmacutTOF")){
+    printf("PID not implemented yet: %s\nNO PID!!!!\n",option);
+    Float_t p = track->P();
+
+    // Take expected times
+    Double_t exptimes[5];
+    track->GetIntegratedTimes(exptimes);
+
+    // Take resolution for TOF response
+    // like fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[ipart], mass[ipart]);
+    Float_t resolution = fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[ipart], mass[ipart]);
+
+    if(TMath::Abs(exptimes[ipart] - track->GetTOFsignal()) < 3 * resolution){
+      pdg = pdgvalues[ipart] * Int_t(track->GetSign());
+    }
+  }
+
+  else{
+    printf("Invalid PID option: %s\nNO PID!!!!\n",option);
+  }
+
+  return pdg;
+}
+//-----------------------------------------------------------------------
+void AliFlowTrackCuts::SetPriors(){
+  // set abbundancies
+    fBinLimitPID[0] = 0.30;
+    fC[0][0] = 0.015;
+    fC[0][1] = 0.015;
+    fC[0][2] = 1;
+    fC[0][3] = 0.0025;
+    fC[0][4] = 0.000015;
+    fBinLimitPID[1] = 0.35;
+    fC[1][0] = 0.015;
+    fC[1][1] = 0.015;
+    fC[1][2] = 1;
+    fC[1][3] = 0.01;
+    fC[1][4] = 0.001;
+    fBinLimitPID[2] = 0.40;
+    fC[2][0] = 0.015;
+    fC[2][1] = 0.015;
+    fC[2][2] = 1;
+    fC[2][3] = 0.026;
+    fC[2][4] = 0.004;
+    fBinLimitPID[3] = 0.45;
+    fC[3][0] = 0.015;
+    fC[3][1] = 0.015;
+    fC[3][2] = 1;
+    fC[3][3] = 0.026;
+    fC[3][4] = 0.004;
+    fBinLimitPID[4] = 0.50;
+    fC[4][0] = 0.015;
+    fC[4][1] = 0.015;
+    fC[4][2] = 1.000000;
+    fC[4][3] = 0.05;
+    fC[4][4] = 0.01;
+    fBinLimitPID[5] = 0.60;
+    fC[5][0] = 0.012;
+    fC[5][1] = 0.012;
+    fC[5][2] = 1;
+    fC[5][3] = 0.085;
+    fC[5][4] = 0.022;
+    fBinLimitPID[6] = 0.70;
+    fC[6][0] = 0.01;
+    fC[6][1] = 0.01;
+    fC[6][2] = 1;
+    fC[6][3] = 0.12;
+    fC[6][4] = 0.036;
+    fBinLimitPID[7] = 0.80;
+    fC[7][0] = 0.0095;
+    fC[7][1] = 0.0095;
+    fC[7][2] = 1;
+    fC[7][3] = 0.15;
+    fC[7][4] = 0.05;
+    fBinLimitPID[8] = 0.90;
+    fC[8][0] = 0.0085;
+    fC[8][1] = 0.0085;
+    fC[8][2] = 1;
+    fC[8][3] = 0.18;
+    fC[8][4] = 0.074;
+    fBinLimitPID[9] = 1;
+    fC[9][0] = 0.008;
+    fC[9][1] = 0.008;
+    fC[9][2] = 1;
+    fC[9][3] = 0.22;
+    fC[9][4] = 0.1;
+    fBinLimitPID[10] = 1.20;
+    fC[10][0] = 0.007;
+    fC[10][1] = 0.007;
+    fC[10][2] = 1;
+    fC[10][3] = 0.28;
+    fC[10][4] = 0.16;
+    fBinLimitPID[11] = 1.40;
+    fC[11][0] = 0.0066;
+    fC[11][1] = 0.0066;
+    fC[11][2] = 1;
+    fC[11][3] = 0.35;
+    fC[11][4] = 0.23;
+    fBinLimitPID[12] = 1.60;
+    fC[12][0] = 0.0075;
+    fC[12][1] = 0.0075;
+    fC[12][2] = 1;
+    fC[12][3] = 0.40;
+    fC[12][4] = 0.31;
+    fBinLimitPID[13] = 1.80;
+    fC[13][0] = 0.0062;
+    fC[13][1] = 0.0062;
+    fC[13][2] = 1;
+    fC[13][3] = 0.45;
+    fC[13][4] = 0.39;
+    fBinLimitPID[14] = 2.00;
+    fC[14][0] = 0.005;
+    fC[14][1] = 0.005;
+    fC[14][2] = 1;
+    fC[14][3] = 0.46;
+    fC[14][4] = 0.47;
+    fBinLimitPID[15] = 2.20;
+    fC[15][0] = 0.0042;
+    fC[15][1] = 0.0042;
+    fC[15][2] = 1;
+    fC[15][3] = 0.5;
+    fC[15][4] = 0.55;
+    fBinLimitPID[16] = 2.40;
+    fC[16][0] = 0.007;
+    fC[16][1] = 0.007;
+    fC[16][2] = 1;
+    fC[16][3] = 0.5;
+    fC[16][4] = 0.6;
+    
+    for(Int_t i=17;i<fnPIDptBin;i++){
+       fBinLimitPID[i] = 2.0 + 0.2 * (i-14);
+       fC[i][0] = fC[13][0];
+       fC[i][1] = fC[13][1];
+       fC[i][2] = fC[13][2];
+       fC[i][3] = fC[13][3];
+       fC[i][4] = fC[13][4];
+    }  
+}
+// end part added by F. Noferini
+//-----------------------------------------------------------------------
+
+
+//-----------------------------------------------------------------------
+const char* AliFlowTrackCuts::PIDsourceName(PIDsource s)
+{
+  //get the name of the particle id source
+  switch (s)
+  {
+    case kTPCdedx:
+      return "TPCdedx";
+    case kTOFbeta:
+      return "TOFbeta";
+    case kTPCpid:
+      return "TPCpid";
+    case kTOFpid:
+      return "TOFpid";
+    case kTOFbayesian:
+      return "TOFbayesianPID";
+    default:
+      return "NOPID";
+  }
+}
+
+//-----------------------------------------------------------------------
+const char* AliFlowTrackCuts::GetParamTypeName(trackParameterType type) 
+{
+  //return the name of the selected parameter type
+  switch (type)
+  {
+    case kMC:
+      return "MC";
+    case kGlobal:
+      return "ESD global";
+    case kESD_TPConly:
+      return "TPC only";
+    case kESD_SPDtracklet:
+      return "SPD tracklet";
+    case kPMD:
+      return "PMD";
+    default:
+      return "unknown";
+  }
+}
+
+//-----------------------------------------------------------------------
+Bool_t AliFlowTrackCuts::PassesPMDcuts(AliESDPmdTrack* /*track*/ )
+{
+  //check PMD specific cuts
+  //clean up from last iteration, and init label
+
+  fTrack = NULL;
+  fMCparticle=NULL;
+  fTrackLabel=-1;
+
+  fTrackPhi = 0.;
+  fTrackEta = 0.;
+  fTrackWeight = 1.0;
+
+  Bool_t pass=kTRUE;
+  if (fCutEta) {if (  fTrackEta < fEtaMin || fTrackEta >= fEtaMax ) pass = kFALSE;}
+  if (fCutPhi) {if ( fTrackPhi < fPhiMin || fTrackPhi >= fPhiMax ) pass = kFALSE;}
+
+  return kFALSE;
+}