PWGCF/FEMTOSCOPY/AliFemtoUser/AliFemtoMJTrackCut.cxx

   1 /*
   2 ***************************************************************************
   3 *
   4 * $Id$
   5 *
   6 *
   7 ***************************************************************************
   8 *
   9 *
  10 *
  11 *
  12 ***************************************************************************
  13 *
  14 * $Log$
  15 * Revision 1.3  2007/05/22 09:01:42  akisiel
  16 * Add the possibiloity to save cut settings in the ROOT file
  17 *
  18 * Revision 1.2  2007/05/21 10:38:25  akisiel
  19 * More coding rule conformance
  20 *
  21 * Revision 1.1  2007/05/16 10:25:06  akisiel
  22 * Making the directory structure of AliFemtoUser flat. All files go into one common directory
  23 *
  24 * Revision 1.4  2007/05/03 09:46:10  akisiel
  25 * Fixing Effective C++ warnings
  26 *
  27 * Revision 1.3  2007/04/27 07:25:59  akisiel
  28 * Make revisions needed for compilation from the main AliRoot tree
  29 *
  30 * Revision 1.1.1.1  2007/04/25 15:38:41  panos
  31 * Importing the HBT code dir
  32 *
  33 * Revision 1.4  2007-04-03 16:00:08  mchojnacki
  34 * Changes to iprove memory managing
  35 *
  36 * Revision 1.3  2007/03/13 15:30:03  mchojnacki
  37 * adding reader for simulated data
  38 *
  39 * Revision 1.2  2007/03/08 14:58:03  mchojnacki
  40 * adding some alice stuff
  41 *
  42 * Revision 1.1.1.1  2007/03/07 10:14:49  mchojnacki
  43 * First version on CVS
  44 *
  45 **************************************************************************/
  46
  47 #include "AliFemtoMJTrackCut.h"
  48 #include <cstdio>
  49
  50 #ifdef __ROOT__
  51 ClassImp(AliFemtoMJTrackCut)
  52 #endif
  53
  54
  55 // electron
  56 // 0.13 - 1.8
  57 // 0       7.594129e-02    8.256141e-03
  58 // 1       -5.535827e-01   8.170825e-02
  59 // 2       1.728591e+00    3.104210e-01
  60 // 3       -2.827893e+00   5.827802e-01
  61 // 4       2.503553e+00    5.736207e-01
  62 // 5       -1.125965e+00   2.821170e-01
  63 // 6       2.009036e-01    5.438876e-02
  64
  65 // pion
  66 // 0.13 - 2.0
  67 // 0       1.063457e+00    8.872043e-03
  68 // 1       -4.222208e-01   2.534402e-02
  69 // 2       1.042004e-01    1.503945e-02
  70
  71 // kaon
  72 // 0.18 - 2.0
  73 // 0       -7.289406e-02   1.686074e-03
  74 // 1       4.415666e-01    1.143939e-02
  75 // 2       -2.996790e-01   1.840964e-02
  76 // 3       6.704652e-02    7.783990e-03
  77
  78 // proton
  79 // 0.26 - 2.0
  80 // 0       -3.730200e-02   2.347311e-03
  81 // 1       1.163684e-01    1.319316e-02
  82 // 2       8.354116e-02    1.997948e-02
  83 // 3       -4.608098e-02   8.336400e-03
  84
  85
  86   AliFemtoMJTrackCut::AliFemtoMJTrackCut() :
  87     fCharge(0),
  88     fLabel(0),
  89     fStatus(0),
  90     fPIDMethod(knSigma),
  91   fNsigmaTPCTOF(kFALSE),
  92   fNsigmaTPConly(kFALSE),
  93   fNsigma(3.),
  94     fminTPCclsF(0),
  95     fminTPCncls(0),
  96     fminITScls(0),
  97     fMaxITSchiNdof(1000.0),
  98     fMaxTPCchiNdof(1000.0),
  99     fMaxSigmaToVertex(1000.0),
 100     fNTracksPassed(0),
 101     fNTracksFailed(0),
 102     fRemoveKinks(kFALSE),
 103     fRemoveITSFake(kFALSE),
 104     fMostProbable(0),
 105     fMaxImpactXY(1000.0),
 106   fMinImpactXY(-1000.0),
 107     fMaxImpactZ(1000.0),
 108     fMaxImpactXYPtOff(1000.0),
 109     fMaxImpactXYPtNrm(1000.0),
 110     fMaxImpactXYPtPow(1000.0),
 111     fMinPforTOFpid(0.0),
 112     fMaxPforTOFpid(10000.0),
 113     fMinPforTPCpid(0.0),
 114     fMaxPforTPCpid(10000.0),
 115     fMinPforITSpid(0.0),
 116    fMaxPforITSpid(10000.0),
 117    fElectronRejection(0)
 118 {
 119   // Default constructor
 120   fNTracksPassed = fNTracksFailed = 0;
 121   fCharge = 0;  // takes both charges 0
 122   fPt[0]=0.0;              fPt[1] = 100.0;//100
 123   fRapidity[0]=-2;       fRapidity[1]=2;//-2 2
 124   fEta[0]=-2;       fEta[1]=2;//-2 2
 125   fPidProbElectron[0]=-1;fPidProbElectron[1]=2;
 126   fPidProbPion[0]=-1;    fPidProbPion[1]=2;
 127   fPidProbKaon[0]=-1;fPidProbKaon[1]=2;
 128   fPidProbProton[0]=-1;fPidProbProton[1]=2;
 129   fPidProbMuon[0]=-1;fPidProbMuon[1]=2;
 130   for (Int_t i = 0; i < 3; i++)
 131     fCutClusterRequirementITS[i] = AliESDtrackCuts::kOff;
 132   fLabel=false;
 133   fStatus=0;
 134   fminTPCclsF=0;
 135   fminITScls=0;
 136   fPIDMethod=knSigma;
 137   fNsigmaTPCTOF=kFALSE;
 138   fNsigmaTPConly=kFALSE;
 139   fNsigma=3.;
 140 }
 141 //------------------------------
 142 AliFemtoMJTrackCut::~AliFemtoMJTrackCut(){
 143   /* noop */
 144 }
 145 //------------------------------
 146 bool AliFemtoMJTrackCut::Pass(const AliFemtoTrack* track)
 147 {
 148   //cout<<"AliFemtoMJTrackCut::Pass"<<endl;
 149
 150   // test the particle and return
 151   // true if it meets all the criteria
 152   // false if it doesn't meet at least one of the criteria
 153   float tMost[5];
 154
 155   //cout<<"AliFemtoESD  cut"<<endl;
 156   //cout<<fPidProbPion[0]<<" < pi ="<<track->PidProbPion()<<" <"<<fPidProbPion[1]<<endl;
 157   if (fStatus!=0)
 158     {
 159       //cout<<" status "<<track->Label()<<" "<<track->Flags()<<" "<<track->TPCnclsF()<<" "<<track->ITSncls()<<endl;
 160       if ((track->Flags()&fStatus)!=fStatus)
 161         {
 162           //      cout<<track->Flags()<<" "<<fStatus<<" no go through status"<<endl;
 163           return false;
 164         }
 165
 166     }
 167   if (fRemoveKinks) {
 168     if ((track->KinkIndex(0)) || (track->KinkIndex(1)) || (track->KinkIndex(2)))
 169       return false;
 170   }
 171   if (fRemoveITSFake) {
 172     if (track->ITSncls() < 0)
 173       return false;
 174   }
 175   if (fminTPCclsF>track->TPCnclsF())
 176     {
 177       //cout<<" No go because TPC Number of ClsF"<<fminTPCclsF<< " "<<track->TPCnclsF()<<endl;
 178       return false;
 179     }
 180   if (fminTPCncls>track->TPCncls())
 181     {
 182       //cout<<" No go because TPC Number of ClsF"<<fminTPCclsF<< " "<<track->TPCnclsF()<<endl;
 183       return false;
 184     }
 185   if (fminITScls>track->ITSncls())
 186     {
 187       //cout<<" No go because ITS Number of Cls"<<fminITScls<< " "<<track->ITSncls()<<endl;
 188       return false;
 189     }
 190
 191   if (fMaxImpactXY < TMath::Abs(track->ImpactD()))
 192     return false;
 193
 194   if (fMinImpactXY > TMath::Abs(track->ImpactD()))
 195     return false;
 196
 197   if (fMaxImpactZ < TMath::Abs(track->ImpactZ()))
 198     return false;
 199
 200   if (fMaxSigmaToVertex < track->SigmaToVertex()) {
 201     return false;
 202   }
 203
 204   if (track->ITSncls() > 0)
 205     if ((track->ITSchi2()/track->ITSncls()) > fMaxITSchiNdof) {
 206       return false;
 207     }
 208
 209   if (track->TPCncls() > 0)
 210     if ((track->TPCchi2()/track->TPCncls()) > fMaxTPCchiNdof) {
 211       return false;
 212     }
 213   //ITS cluster requirenments
 214   for (Int_t i = 0; i < 3; i++)
 215     if(!CheckITSClusterRequirement(fCutClusterRequirementITS[i], track->HasPointOnITSLayer(i*2), track->HasPointOnITSLayer(i*2+1)))
 216       return false;
 217
 218   if (fLabel)
 219     {
 220       //cout<<"labels"<<endl;
 221       if(track->Label()<0)
 222         {
 223           fNTracksFailed++;
 224           //   cout<<"No Go Through the cut"<<endl;
 225           //  cout<<fLabel<<" Label="<<track->Label()<<endl;
 226           return false;
 227         }
 228     }
 229   if (fCharge!=0)
 230     {
 231       //cout<<"AliFemtoESD  cut ch "<<endl;
 232       //cout<<fCharge<<" Charge="<<track->Charge()<<endl;
 233       if (track->Charge()!= fCharge)
 234         {
 235           fNTracksFailed++;
 236           //  cout<<"No Go Through the cut"<<endl;
 237           // cout<<fCharge<<" Charge="<<track->Charge()<<endl;
 238           return false;
 239         }
 240     }
 241
 242
 243
 244
 245   Bool_t tTPCPidIn = (track->Flags()&AliFemtoTrack::kTPCpid)>0;
 246   Bool_t tITSPidIn = (track->Flags()&AliFemtoTrack::kITSpid)>0;
 247   Bool_t tTOFPidIn = (track->Flags()&AliFemtoTrack::kTOFpid)>0;
 248
 249   if(fMinPforTOFpid > 0 && track->P().Mag() > fMinPforTOFpid &&
 250      track->P().Mag() < fMaxPforTOFpid && !tTOFPidIn)
 251     {
 252       fNTracksFailed++;
 253       return false;
 254     }
 255
 256   if(fMinPforTPCpid > 0 && track->P().Mag() > fMinPforTPCpid &&
 257      track->P().Mag() < fMaxPforTPCpid && !tTPCPidIn)
 258     {
 259       fNTracksFailed++;
 260       return false;
 261     }
 262
 263   if(fMinPforITSpid > 0 && track->P().Mag() > fMinPforITSpid &&
 264      track->P().Mag() < fMaxPforITSpid && !tITSPidIn)
 265     {
 266       fNTracksFailed++;
 267       return false;
 268     }
 269
 270
 271   float tEnergy = ::sqrt(track->P().Mag2()+fMass*fMass);
 272   float tRapidity = 0;
 273   if(tEnergy-track->P().z()!=0 && (tEnergy+track->P().z())/(tEnergy-track->P().z())>0)
 274     tRapidity = 0.5*::log((tEnergy+track->P().z())/(tEnergy-track->P().z()));
 275   float tPt = ::sqrt((track->P().x())*(track->P().x())+(track->P().y())*(track->P().y()));
 276   float tEta = track->P().PseudoRapidity();
 277
 278   if (fMaxImpactXYPtOff < 999.0) {
 279     if ((fMaxImpactXYPtOff + fMaxImpactXYPtNrm*TMath::Power(tPt, fMaxImpactXYPtPow)) < TMath::Abs(track->ImpactD())) {
 280       fNTracksFailed++;
 281       return false;
 282     }
 283   }
 284
 285   if ((tRapidity<fRapidity[0])||(tRapidity>fRapidity[1]))
 286     {
 287       fNTracksFailed++;
 288       //cout<<"No Go Through the cut"<<endl;
 289       //cout<<fRapidity[0]<<" < Rapidity ="<<tRapidity<<" <"<<fRapidity[1]<<endl;
 290       return false;
 291     }
 292   if ((tEta<fEta[0])||(tEta>fEta[1]))
 293     {
 294       fNTracksFailed++;
 295       //cout<<"No Go Through the cut"<<endl;
 296       //cout<<fEta[0]<<" < Eta ="<<tEta<<" <"<<fEta[1]<<endl;
 297       return false;
 298     }
 299   if ((tPt<fPt[0])||(tPt>fPt[1]))
 300     {
 301       fNTracksFailed++;
 302       //cout<<"No Go Through the cut"<<endl;
 303       //cout<<fPt[0]<<" < Pt ="<<Pt<<" <"<<fPt[1]<<endl;
 304       return false;
 305     }
 306
 307
 308
 309
 310   //   cout << "Track has pids: "
 311   //        << track->PidProbElectron() << " "
 312   //        << track->PidProbMuon() << " "
 313   //        << track->PidProbPion() << " "
 314   //        << track->PidProbKaon() << " "
 315   //        << track->PidProbProton() << " "
 316   //        << track->PidProbElectron()+track->PidProbMuon()+track->PidProbPion()+track->PidProbKaon()+track->PidProbProton() << endl;
 317
 318
 319   if ((track->PidProbElectron()<fPidProbElectron[0])||(track->PidProbElectron()>fPidProbElectron[1]))
 320     {
 321       fNTracksFailed++;
 322       //cout<<"No Go Through the cut"<<endl;
 323       //cout<<fPidProbElectron[0]<<" < e ="<<track->PidProbElectron()<<" <"<<fPidProbElectron[1]<<endl;
 324       return false;
 325     }
 326   if ((track->PidProbPion()<fPidProbPion[0])||(track->PidProbPion()>fPidProbPion[1]))
 327     {
 328       fNTracksFailed++;
 329       //cout<<"No Go Through the cut"<<endl;
 330       //cout<<fPidProbPion[0]<<" < pi ="<<track->PidProbPion()<<" <"<<fPidProbPion[1]<<endl;
 331       return false;
 332     }
 333   if ((track->PidProbKaon()<fPidProbKaon[0])||(track->PidProbKaon()>fPidProbKaon[1]))
 334     {
 335       fNTracksFailed++;
 336       //cout<<"No Go Through the cut"<<endl;
 337       //cout<<fPidProbKaon[0]<<" < k ="<<track->PidProbKaon()<<" <"<<fPidProbKaon[1]<<endl;
 338       return false;
 339     }
 340   if ((track->PidProbProton()<fPidProbProton[0])||(track->PidProbProton()>fPidProbProton[1]))
 341     {
 342       fNTracksFailed++;
 343       //cout<<"No Go Through the cut"<<endl;
 344       //cout<<fPidProbProton[0]<<" < p  ="<<track->PidProbProton()<<" <"<<fPidProbProton[1]<<endl;
 345       return false;
 346     }
 347   if ((track->PidProbMuon()<fPidProbMuon[0])||(track->PidProbMuon()>fPidProbMuon[1]))
 348     {
 349       fNTracksFailed++;
 350       //cout<<"No Go Through the cut"<<endl;
 351       //cout<<fPidProbMuon[0]<<" <  mi="<<track->PidProbMuon()<<" <"<<fPidProbMuon[1]<<endl;
 352       return false;
 353     }
 354
 355   //****N Sigma Method -- electron rejection****
 356   if(fElectronRejection)
 357     if(!IsElectron(track->NSigmaTPCE(),track->NSigmaTPCPi(),track->NSigmaTPCK(), track->NSigmaTPCP()))
 358       return false;
 359
 360
 361   if (fMostProbable) {
 362
 363     int imost=0;
 364     tMost[0] = track->PidProbElectron()*PidFractionElectron(track->P().Mag());
 365     tMost[1] = 0.0;
 366     tMost[2] = track->PidProbPion()*PidFractionPion(track->P().Mag());
 367     tMost[3] = track->PidProbKaon()*PidFractionKaon(track->P().Mag());
 368     tMost[4] = track->PidProbProton()*PidFractionProton(track->P().Mag());
 369     float ipidmax = 0.0;
 370
 371     //****N Sigma Method****
 372         if(fPIDMethod==0){
 373           // Looking for pions
 374           if (fMostProbable == 2) {
 375             if (IsPionNSigma(track->P().Mag(), track->NSigmaTPCPi(), track->NSigmaTOFPi()))
 376               imost = 2;
 377
 378           }
 379           else if (fMostProbable == 3) {
 380
 381
 382             if (IsKaonNSigma(track->P().Mag(), track->NSigmaTPCK(), track->NSigmaTOFK())){
 383
 384               imost = 3;
 385             }
 386           }
 387           else if (fMostProbable == 4) { // proton nsigma-PID required contour adjusting (in LHC10h)
 388             if ( IsProtonNSigma(track->P().Mag(), track->NSigmaTPCP(), track->NSigmaTOFP()) // && (TMath::Abs(track->NSigmaTPCP()) < TMath::Abs(track->NSigmaTPCPi())) && (TMath::Abs(track->NSigmaTPCP()) < TMath::Abs(track->NSigmaTPCK())) && (TMath::Abs(track->NSigmaTOFP()) < TMath::Abs(track->NSigmaTOFPi())) && (TMath::Abs(track->NSigmaTOFP()) < TMath::Abs(track->NSigmaTOFK()))
 389                  // && IsProtonTPCdEdx(track->P().Mag(), track->TPCsignal())
 390                  )
 391               imost = 4;
 392           }
 393           else if (fMostProbable == 5) { // no-protons
 394             if ( !IsProtonNSigma(track->P().Mag(), track->NSigmaTPCP(), track->NSigmaTOFP()) )
 395               imost = 5;
 396           }
 397           else if (fMostProbable == 6) { //pions OR kaons OR protons
 398             if (IsPionNSigma(track->P().Mag(), track->NSigmaTPCPi(), track->NSigmaTOFPi()))
 399               imost = 6;
 400             else if (IsKaonNSigma(track->P().Mag(), track->NSigmaTPCK(), track->NSigmaTOFK()))
 401               imost = 6;
 402             else if (IsProtonNSigma(track->P().Mag(), track->NSigmaTPCP(), track->NSigmaTOFP()) )
 403               imost = 6;
 404           }
 405           else if (fMostProbable == 7) { // pions OR kaons OR protons OR electrons
 406             if (IsPionNSigma(track->P().Mag(), track->NSigmaTPCPi(), track->NSigmaTOFPi()))
 407               imost = 7;
 408             else if (IsKaonNSigma(track->P().Mag(), track->NSigmaTPCK(), track->NSigmaTOFK()))
 409               imost = 7;
 410             else if (IsProtonNSigma(track->P().Mag(), track->NSigmaTPCP(), track->NSigmaTOFP()) )
 411               imost = 7;
 412             else if (TMath::Abs(track->NSigmaTPCE())<3)
 413               imost = 7;
 414
 415           }
 416           else if (fMostProbable == 8) { // TOF matching
 417             if(track->NSigmaTOFPi() != -1000 || track->Pt()<0.5){
 418               imost = 8;
 419             }
 420           }
 421           else if (fMostProbable == 9) { // Other: no kaons, no pions, no protons
 422             if (IsPionNSigma(track->P().Mag(), track->NSigmaTPCPi(), track->NSigmaTOFPi()))
 423               imost = -1;
 424             else if (IsKaonNSigma(track->P().Mag(), track->NSigmaTPCK(), track->NSigmaTOFK()))
 425               imost = -1;
 426             else if (IsProtonNSigma(track->P().Mag(), track->NSigmaTPCP(), track->NSigmaTOFP()) )
 427               imost = -1;
 428             else if(track->NSigmaTOFPi() != -1000 || track->Pt()<0.5){
 429               imost = 9;
 430             }
 431           }
 432           if (fMostProbable == 10) {//cut on Nsigma in pT not p
 433             if (IsPionNSigma(track->Pt(), track->NSigmaTPCPi(), track->NSigmaTOFPi()))
 434               imost = 10;
 435           }
 436           else if (fMostProbable == 11) {//cut on Nsigma in pT not p
 437             if (IsKaonNSigma(track->Pt(), track->NSigmaTPCK(), track->NSigmaTOFK())){
 438               imost = 11;
 439             }
 440           }
 441           else if (fMostProbable == 12) { //cut on Nsigma in pT not p
 442             if ( IsProtonNSigma(track->Pt(), track->NSigmaTPCP(), track->NSigmaTOFP()) )
 443               imost = 12;
 444           }
 445           else if (fMostProbable == 13) { //cut on Nsigma in pT not p, EXCLUSIVE PID
 446             if ((IsPionNSigma(track->Pt(),track->NSigmaTPCPi(), track->NSigmaTOFPi()) && !IsKaonNSigma(track->Pt(),track->NSigmaTPCK(), track->NSigmaTOFK()) && !IsProtonNSigma(track->Pt(),track->NSigmaTPCP(), track->NSigmaTOFP())))
 447               imost = 13;
 448           }
 449           else if (fMostProbable == 14) { //cut on Nsigma in pT not p, EXCLUSIVE PID
 450             bool isPionNsigma = (IsPionNSigma(track->Pt(),track->NSigmaTPCPi(), track->NSigmaTOFPi()) && !IsKaonNSigma(track->Pt(),track->NSigmaTPCK(), track->NSigmaTOFK()) && !IsProtonNSigma(track->Pt(),track->NSigmaTPCP(), track->NSigmaTOFP()));
 451             if ( !isPionNsigma && IsKaonNSigma(track->Pt(),track->NSigmaTPCK(), track->NSigmaTOFK()) && !IsProtonNSigma(track->Pt(),track->NSigmaTPCP(), track->NSigmaTOFP()))
 452               imost = 14;
 453           }
 454           else if (fMostProbable == 15) { //cut on Nsigma in pT not p, EXCLUSIVE PID
 455             bool isPionNsigma = (IsPionNSigma(track->Pt(),track->NSigmaTPCPi(), track->NSigmaTOFPi()) && !IsKaonNSigma(track->Pt(),track->NSigmaTPCK(), track->NSigmaTOFK()) && !IsProtonNSigma(track->Pt(),track->NSigmaTPCP(), track->NSigmaTOFP()));
 456             bool isKaonNsigma = (!isPionNsigma && IsKaonNSigma(track->Pt(),track->NSigmaTPCK(), track->NSigmaTOFK()) && !IsProtonNSigma(track->Pt(),track->NSigmaTPCP(), track->NSigmaTOFP()));
 457             if (!isPionNsigma && !isKaonNsigma && IsProtonNSigma(track->Pt(),track->NSigmaTPCP(), track->NSigmaTOFP()))
 458               imost = 15;
 459           }
 460           //*************** Without double counting, pions ************************
 461           else if (fMostProbable == 16) { //cut on Nsigma in pT not p,  without double counting
 462             double nSigmaPIDPi = 0, nSigmaPIDK = 0, nSigmaPIDP = 0;
 463             if(track->Pt()<0.5){
 464               nSigmaPIDPi = abs(track->NSigmaTPCPi());
 465               nSigmaPIDK  = abs(track->NSigmaTPCK());
 466               nSigmaPIDP  = abs(track->NSigmaTPCP());
 467             }
 468             else{
 469               nSigmaPIDPi = TMath::Hypot(track->NSigmaTOFPi(), track->NSigmaTPCPi());
 470               nSigmaPIDK= TMath::Hypot(track->NSigmaTOFK(), track->NSigmaTPCK());
 471               nSigmaPIDP= TMath::Hypot(track->NSigmaTOFP(), track->NSigmaTPCP());
 472             }
 473             bool isPionNsigma = 0;
 474
 475             if(nSigmaPIDPi<nSigmaPIDK && nSigmaPIDPi<nSigmaPIDP){
 476               isPionNsigma = (IsPionNSigma(track->Pt(),track->NSigmaTPCPi(), track->NSigmaTOFPi()));
 477               if(isPionNsigma) imost=16;
 478             }
 479           }
 480
 481           else if (fMostProbable == 17) { //cut on Nsigma in pT not p,  without double counting
 482             double nSigmaPIDPi = 0, nSigmaPIDK = 0, nSigmaPIDP = 0;
 483             if(track->Pt()<0.5){
 484               nSigmaPIDPi = abs(track->NSigmaTPCPi());
 485               nSigmaPIDK  = abs(track->NSigmaTPCK());
 486               nSigmaPIDP  = abs(track->NSigmaTPCP());
 487             }
 488             else{
 489               nSigmaPIDPi = TMath::Hypot(track->NSigmaTOFPi(), track->NSigmaTPCPi());
 490               nSigmaPIDK= TMath::Hypot(track->NSigmaTOFK(), track->NSigmaTPCK());
 491               nSigmaPIDP= TMath::Hypot(track->NSigmaTOFP(), track->NSigmaTPCP());
 492             }
 493
 494             bool isKaonNsigma = 0;
 495             if(nSigmaPIDK<nSigmaPIDPi && nSigmaPIDK<nSigmaPIDP){
 496               isKaonNsigma = (IsKaonNSigma(track->Pt(),track->NSigmaTPCK(), track->NSigmaTOFK()));
 497               if(isKaonNsigma) imost=17;
 498             }
 499           }
 500
 501           else if (fMostProbable == 18) { //cut on Nsigma in pT not p,  without double counting
 502             double nSigmaPIDPi = 0, nSigmaPIDK = 0, nSigmaPIDP = 0;
 503             if(track->Pt()<0.5){
 504               nSigmaPIDPi = abs(track->NSigmaTPCPi());
 505               nSigmaPIDK  = abs(track->NSigmaTPCK());
 506               nSigmaPIDP  = abs(track->NSigmaTPCP());
 507             }
 508             else{
 509               nSigmaPIDPi = TMath::Hypot(track->NSigmaTOFPi(), track->NSigmaTPCPi());
 510               nSigmaPIDK= TMath::Hypot(track->NSigmaTOFK(), track->NSigmaTPCK());
 511               nSigmaPIDP= TMath::Hypot(track->NSigmaTOFP(), track->NSigmaTPCP());
 512             }
 513             bool isProtonNsigma  = 0;
 514
 515             if(nSigmaPIDP<nSigmaPIDPi && nSigmaPIDP<nSigmaPIDK){
 516               isProtonNsigma = (IsProtonNSigma(track->Pt(),track->NSigmaTPCP(), track->NSigmaTOFP()));
 517             }
 518             if(isProtonNsigma) imost=18;
 519           }
 520         }
 521
 522
 523
 524
 525     //****Contour Method****
 526         if(fPIDMethod==1){
 527           for (int ip=0; ip<5; ip++)
 528             if (tMost[ip] > ipidmax) { ipidmax = tMost[ip]; imost = ip; };
 529
 530           // Looking for pions
 531           if (fMostProbable == 2) {
 532             if (imost == 2) {
 533               // Using the TPC to reject non-pions
 534               if (!(IsPionTPCdEdx(track->P().Mag(), track->TPCsignal())))
 535                 imost = 0;
 536               if (0) {
 537                 // Using the TOF to reject non-pions
 538                 if (track->P().Mag() < 0.6) {
 539                   if (tTOFPidIn)
 540                     if (!IsPionTOFTime(track->P().Mag(), track->TOFpionTime()))
 541                       imost = 0;
 542                 }
 543                 else {
 544                   if (tTOFPidIn) {
 545                     if (!IsPionTOFTime(track->P().Mag(), track->TOFpionTime()))
 546                       imost = 0;
 547                   }
 548                   else {
 549                     imost = 0;
 550                   }
 551                 }
 552               }
 553             }
 554           }
 555
 556           // Looking for kaons
 557           else if (fMostProbable == 3) {
 558             //       if (imost == 3) {
 559             // Using the TPC to reject non-kaons
 560             if (track->P().Mag() < 0.6) {
 561               if (!(IsKaonTPCdEdx(track->P().Mag(), track->TPCsignal())))
 562                 imost = 0;
 563               else imost = 3;
 564               if (1) {
 565                 // Using the TOF to reject non-kaons
 566                 if (tTOFPidIn)
 567                   if (!IsKaonTOFTime(track->P().Mag(), track->TOFkaonTime()))
 568                     imost = 0;
 569               }
 570             }
 571             else {
 572               if (1) {
 573                 if (tTOFPidIn) {
 574                   if (!IsKaonTOFTime(track->P().Mag(), track->TOFkaonTime()))
 575                     imost = 0;
 576                   else
 577                     imost = 3;
 578                 }
 579                 else {
 580                   if (!(IsKaonTPCdEdx(track->P().Mag(), track->TPCsignal())))
 581                     imost = 0;
 582                   else
 583                     imost = 3;
 584                 }
 585               }
 586             }
 587             //       }
 588           }
 589
 590           // Looking for protons
 591           else if (fMostProbable == 4) {
 592             //       if (imost == 3) {
 593             // Using the TPC to reject non-kaons
 594             if (track->P().Mag() < 0.8) {
 595               if (!(IsProtonTPCdEdx(track->P().Mag(), track->TPCsignal())))
 596                 imost = 0;
 597               else imost = 4;
 598               if (0) {
 599                 // Using the TOF to reject non-kaons
 600                 if (tTOFPidIn)
 601                   if (!IsKaonTOFTime(track->P().Mag(), track->TOFkaonTime()))
 602                     imost = 0;
 603               }
 604             }
 605             else {
 606               if (0) {
 607                 if (tTOFPidIn) {
 608                   if (!IsKaonTOFTime(track->P().Mag(), track->TOFkaonTime()))
 609                     imost = 0;
 610                   else
 611                     imost = 3;
 612                 }
 613                 else {
 614                   if (!(IsKaonTPCdEdx(track->P().Mag(), track->TPCsignal())))
 615                     imost = 0;
 616                   else
 617                     imost = 3;
 618                 }
 619               }
 620             }
 621             //       }
 622           }
 623         }
 624     if (imost != fMostProbable) return false;
 625   }
 626
 627   //fan
 628   //cout<<"****** Go Through the cut ******"<<endl;
 629   // cout<<fLabel<<" Label="<<track->Label()<<endl;
 630   // cout<<fCharge<<" Charge="<<track->Charge()<<endl;
 631   // cout<<fPt[0]<<" < Pt ="<<Pt<<" <"<<fPt[1]<<endl;
 632   //cout<<fRapidity[0]<<" < Rapidity ="<<tRapidity<<" <"<<fRapidity[1]<<endl;
 633   //cout<<fPidProbElectron[0]<<" <  e="<<track->PidProbElectron()<<" <"<<fPidProbElectron[1]<<endl;
 634   //cout<<fPidProbPion[0]<<" <  pi="<<track->PidProbPion()<<" <"<<fPidProbPion[1]<<endl;
 635   //cout<<fPidProbKaon[0]<<" <  k="<<track->PidProbKaon()<<" <"<<fPidProbKaon[1]<<endl;
 636   //cout<<fPidProbProton[0]<<" <  p="<<track->PidProbProton()<<" <"<<fPidProbProton[1]<<endl;
 637   //cout<<fPidProbMuon[0]<<" <  mi="<<track->PidProbMuon()<<" <"<<fPidProbMuon[1]<<endl;
 638   fNTracksPassed++ ;
 639   return true;
 640
 641
 642 }
 643 //------------------------------
 644 AliFemtoString AliFemtoMJTrackCut::Report()
 645 {
 646   // Prepare report from the execution
 647   string tStemp;
 648   char tCtemp[100];
 649   snprintf(tCtemp , 100, "Particle mass:\t%E\n",this->Mass());
 650   tStemp=tCtemp;
 651   snprintf(tCtemp , 100, "Particle charge:\t%d\n",fCharge);
 652   tStemp+=tCtemp;
 653   snprintf(tCtemp , 100, "Particle pT:\t%E - %E\n",fPt[0],fPt[1]);
 654   tStemp+=tCtemp;
 655   snprintf(tCtemp , 100, "Particle rapidity:\t%E - %E\n",fRapidity[0],fRapidity[1]);
 656   tStemp+=tCtemp;
 657   snprintf(tCtemp , 100, "Particle eta:\t%E - %E\n",fEta[0],fEta[1]);
 658   tStemp+=tCtemp;
 659   snprintf(tCtemp , 100, "Number of tracks which passed:\t%ld  Number which failed:\t%ld\n",fNTracksPassed,fNTracksFailed);
 660   tStemp += tCtemp;
 661   AliFemtoString returnThis = tStemp;
 662   return returnThis;
 663 }
 664 TList *AliFemtoMJTrackCut::ListSettings()
 665 {
 666   // return a list of settings in a writable form
 667   TList *tListSetttings = new TList();
 668   char buf[200];
 669   snprintf(buf, 200, "AliFemtoMJTrackCut.mass=%f", this->Mass());
 670   tListSetttings->AddLast(new TObjString(buf));
 671
 672   snprintf(buf, 200, "AliFemtoMJTrackCut.charge=%i", fCharge);
 673   tListSetttings->AddLast(new TObjString(buf));
 674   snprintf(buf, 200, "AliFemtoMJTrackCut.pidprobpion.minimum=%f", fPidProbPion[0]);
 675   tListSetttings->AddLast(new TObjString(buf));
 676   snprintf(buf, 200, "AliFemtoMJTrackCut.pidprobpion.maximum=%f", fPidProbPion[1]);
 677   tListSetttings->AddLast(new TObjString(buf));
 678   snprintf(buf, 200, "AliFemtoMJTrackCut.pidprobkaon.minimum=%f", fPidProbKaon[0]);
 679   tListSetttings->AddLast(new TObjString(buf));
 680   snprintf(buf, 200, "AliFemtoMJTrackCut.pidprobkaon.maximum=%f", fPidProbKaon[1]);
 681   tListSetttings->AddLast(new TObjString(buf));
 682   snprintf(buf, 200, "AliFemtoMJTrackCut.pidprobproton.minimum=%f", fPidProbProton[0]);
 683   tListSetttings->AddLast(new TObjString(buf));
 684   snprintf(buf, 200, "AliFemtoMJTrackCut.pidprobproton.maximum=%f", fPidProbProton[1]);
 685   tListSetttings->AddLast(new TObjString(buf));
 686   snprintf(buf, 200, "AliFemtoMJTrackCut.pidprobelectron.minimum=%f", fPidProbElectron[0]);
 687   tListSetttings->AddLast(new TObjString(buf));
 688   snprintf(buf, 200, "AliFemtoMJTrackCut.pidprobelectron.maximum=%f", fPidProbElectron[1]);
 689   tListSetttings->AddLast(new TObjString(buf));
 690   snprintf(buf, 200, "AliFemtoMJTrackCut.pidprobMuon.minimum=%f", fPidProbMuon[0]);
 691   tListSetttings->AddLast(new TObjString(buf));
 692   snprintf(buf, 200, "AliFemtoMJTrackCut.pidprobMuon.maximum=%f", fPidProbMuon[1]);
 693   tListSetttings->AddLast(new TObjString(buf));
 694   snprintf(buf, 200, "AliFemtoMJTrackCut.minimumtpcclusters=%i", fminTPCclsF);
 695   tListSetttings->AddLast(new TObjString(buf));
 696   snprintf(buf, 200, "AliFemtoMJTrackCut.minimumitsclusters=%i", fminTPCclsF);
 697   tListSetttings->AddLast(new TObjString(buf));
 698   snprintf(buf, 200, "AliFemtoMJTrackCut.pt.minimum=%f", fPt[0]);
 699   tListSetttings->AddLast(new TObjString(buf));
 700   snprintf(buf, 200, "AliFemtoMJTrackCut.pt.maximum=%f", fPt[1]);
 701   tListSetttings->AddLast(new TObjString(buf));
 702   snprintf(buf, 200, "AliFemtoMJTrackCut.rapidity.minimum=%f", fRapidity[0]);
 703   tListSetttings->AddLast(new TObjString(buf));
 704   snprintf(buf, 200, "AliFemtoMJTrackCut.rapidity.maximum=%f", fRapidity[1]);
 705   tListSetttings->AddLast(new TObjString(buf));
 706   snprintf(buf, 200, "AliFemtoMJTrackCut.removekinks=%i", fRemoveKinks);
 707   tListSetttings->AddLast(new TObjString(buf));
 708   snprintf(buf, 200, "AliFemtoMJTrackCut.maxitschindof=%f", fMaxITSchiNdof);
 709   tListSetttings->AddLast(new TObjString(buf));
 710   snprintf(buf, 200, "AliFemtoMJTrackCut.maxtpcchindof=%f", fMaxTPCchiNdof);
 711   tListSetttings->AddLast(new TObjString(buf));
 712   snprintf(buf, 200, "AliFemtoMJTrackCut.maxsigmatovertex=%f", fMaxSigmaToVertex);
 713   tListSetttings->AddLast(new TObjString(buf));
 714   snprintf(buf, 200, "AliFemtoMJTrackCut.maximpactxy=%f", fMaxImpactXY);
 715   tListSetttings->AddLast(new TObjString(buf));
 716   snprintf(buf, 200, "AliFemtoMJTrackCut.maximpactz=%f", fMaxImpactZ);
 717   tListSetttings->AddLast(new TObjString(buf));
 718   if (fMostProbable) {
 719     if (fMostProbable == 2)
 720       snprintf(buf, 200, "AliFemtoMJTrackCut.mostprobable=%s", "Pion");
 721     if (fMostProbable == 3)
 722       snprintf(buf, 200, "AliFemtoMJTrackCut.mostprobable=%s", "Kaon");
 723     if (fMostProbable == 4)
 724       snprintf(buf, 200, "AliFemtoMJTrackCut.mostprobable=%s", "Proton");
 725     tListSetttings->AddLast(new TObjString(buf));
 726   }
 727   return tListSetttings;
 728 }
 729 void AliFemtoMJTrackCut::SetRemoveKinks(const bool& flag)
 730 {
 731   fRemoveKinks = flag;
 732 }
 733
 734 void AliFemtoMJTrackCut::SetRemoveITSFake(const bool& flag)
 735 {
 736   fRemoveITSFake = flag;
 737 }
 738
 739 // electron
 740 // 0.13 - 1.8
 741 // 0       7.594129e-02    8.256141e-03
 742 // 1       -5.535827e-01   8.170825e-02
 743 // 2       1.728591e+00    3.104210e-01
 744 // 3       -2.827893e+00   5.827802e-01
 745 // 4       2.503553e+00    5.736207e-01
 746 // 5       -1.125965e+00   2.821170e-01
 747 // 6       2.009036e-01    5.438876e-02
 748 float AliFemtoMJTrackCut::PidFractionElectron(float mom) const
 749 {
 750   // Provide a parameterized fraction of electrons dependent on momentum
 751   if (mom<0.13)
 752     return (7.594129e-02
 753             -5.535827e-01*0.13
 754             +1.728591e+00*0.13*0.13
 755             -2.827893e+00*0.13*0.13*0.13
 756             +2.503553e+00*0.13*0.13*0.13*0.13
 757             -1.125965e+00*0.13*0.13*0.13*0.13*0.13
 758             +2.009036e-01*0.13*0.13*0.13*0.13*0.13*0.13);
 759
 760   if (mom>1.8)
 761     return (7.594129e-02
 762             -5.535827e-01*1.8
 763             +1.728591e+00*1.8*1.8
 764             -2.827893e+00*1.8*1.8*1.8
 765             +2.503553e+00*1.8*1.8*1.8*1.8
 766             -1.125965e+00*1.8*1.8*1.8*1.8*1.8
 767             +2.009036e-01*1.8*1.8*1.8*1.8*1.8*1.8);
 768   return (7.594129e-02
 769           -5.535827e-01*mom
 770           +1.728591e+00*mom*mom
 771           -2.827893e+00*mom*mom*mom
 772           +2.503553e+00*mom*mom*mom*mom
 773           -1.125965e+00*mom*mom*mom*mom*mom
 774           +2.009036e-01*mom*mom*mom*mom*mom*mom);
 775 }
 776
 777 // pion
 778 // 0.13 - 2.0
 779 // 0       1.063457e+00    8.872043e-03
 780 // 1       -4.222208e-01   2.534402e-02
 781 // 2       1.042004e-01    1.503945e-02
 782 float AliFemtoMJTrackCut::PidFractionPion(float mom) const
 783 {
 784   // Provide a parameterized fraction of pions dependent on momentum
 785   if (mom<0.13)
 786     return ( 1.063457e+00
 787              -4.222208e-01*0.13
 788              +1.042004e-01*0.0169);
 789   if (mom>2.0)
 790     return ( 1.063457e+00
 791              -4.222208e-01*2.0
 792              +1.042004e-01*4.0);
 793   return ( 1.063457e+00
 794            -4.222208e-01*mom
 795            +1.042004e-01*mom*mom);
 796 }
 797
 798 // kaon
 799 // 0.18 - 2.0
 800 // 0       -7.289406e-02   1.686074e-03
 801 // 1       4.415666e-01    1.143939e-02
 802 // 2       -2.996790e-01   1.840964e-02
 803 // 3       6.704652e-02    7.783990e-03
 804 float AliFemtoMJTrackCut::PidFractionKaon(float mom) const
 805 {
 806   // Provide a parameterized fraction of kaons dependent on momentum
 807   if (mom<0.18)
 808     return (-7.289406e-02
 809             +4.415666e-01*0.18
 810             -2.996790e-01*0.18*0.18
 811             +6.704652e-02*0.18*0.18*0.18);
 812   if (mom>2.0)
 813     return (-7.289406e-02
 814             +4.415666e-01*2.0
 815             -2.996790e-01*2.0*2.0
 816             +6.704652e-02*2.0*2.0*2.0);
 817   return (-7.289406e-02
 818           +4.415666e-01*mom
 819           -2.996790e-01*mom*mom
 820           +6.704652e-02*mom*mom*mom);
 821 }
 822
 823 // proton
 824 // 0.26 - 2.0
 825 // 0       -3.730200e-02   2.347311e-03
 826 // 1       1.163684e-01    1.319316e-02
 827 // 2       8.354116e-02    1.997948e-02
 828 // 3       -4.608098e-02   8.336400e-03
 829 float AliFemtoMJTrackCut::PidFractionProton(float mom) const
 830 {
 831   // Provide a parameterized fraction of protons dependent on momentum
 832   if (mom<0.26) return  0.0;
 833   if (mom>2.0)
 834     return (-3.730200e-02
 835             +1.163684e-01*2.0
 836             +8.354116e-02*2.0*2.0
 837             -4.608098e-02*2.0*2.0*2.0);
 838   return (-3.730200e-02
 839           +1.163684e-01*mom
 840           +8.354116e-02*mom*mom
 841           -4.608098e-02*mom*mom*mom);
 842 }
 843
 844 void AliFemtoMJTrackCut::SetMomRangeTOFpidIs(const float& minp, const float& maxp)
 845 {
 846   fMinPforTOFpid = minp;
 847   fMaxPforTOFpid = maxp;
 848 }
 849
 850 void AliFemtoMJTrackCut::SetMomRangeTPCpidIs(const float& minp, const float& maxp)
 851 {
 852   fMinPforTPCpid = minp;
 853   fMaxPforTPCpid = maxp;
 854 }
 855
 856 void AliFemtoMJTrackCut::SetMomRangeITSpidIs(const float& minp, const float& maxp)
 857 {
 858   fMinPforITSpid = minp;
 859   fMaxPforITSpid = maxp;
 860 }
 861
 862 bool AliFemtoMJTrackCut::IsPionTPCdEdx(float mom, float dEdx)
 863 {
 864   //   double a1 = -95.4545, b1 = 86.5455;
 865   //   double a2 = 0.0,      b2 = 56.0;
 866   double a1 = -343.75,  b1 = 168.125;
 867   double a2 = 0.0,      b2 = 65.0;
 868
 869   if (mom < 0.32) {
 870     if (dEdx < a1*mom+b1) return true;
 871   }
 872   if (dEdx < a2*mom+b2) return true;
 873
 874   return false;
 875 }
 876
 877 bool AliFemtoMJTrackCut::IsKaonTPCdEdx(float mom, float dEdx)
 878 {
 879
 880 //   double a1 = -547.0; double b1 =  297.0;
 881 //   double a2 = -125.0; double b2 =  145.0;
 882 //   double a3 = -420.0; double b3 =  357.0;
 883 //   double a4 = -110.0; double b4 =  171.0;
 884 //   double b5 =   72.0;
 885
 886 //   if (mom<0.2) return false;
 887
 888 //   if (mom<0.36) {
 889 //     if (dEdx < a1*mom+b1) return false;
 890 //     if (dEdx > a3*mom+b3) return false;
 891 //   }
 892 //   else if (mom<0.6) {
 893 //     if (dEdx < a2*mom+b2) return false;
 894 //     if (dEdx > a3*mom+b3) return false;
 895 //   }
 896 //   else if (mom<0.9) {
 897 //     if (dEdx > a4*mom+b4) return false;
 898 //     if (dEdx <        b5) return false;
 899 //   }
 900 //   else
 901 //     return false;
 902 //   //   else {
 903 //   //     if (dEdx > b5) return false;
 904 //   //   }
 905
 906 //   return true;
 907
 908   double a1 = -268.896; double b1 =  198.669;
 909   double a2 = -49.0012;  double b2 =  88.7214;
 910
 911   if (mom<0.2) return false;
 912
 913   if (mom>0.3 && mom<0.5) {
 914     if (dEdx < a1*mom+b1) return false;
 915   }
 916   else  if (mom<1.2) {
 917     if (dEdx < a2*mom+b2) return false;
 918   }
 919
 920   return true;
 921
 922 }
 923
 924 bool AliFemtoMJTrackCut::IsProtonTPCdEdx(float mom, float dEdx)
 925 {
 926   double a1 = -1800.0; double b1 =  940.0;
 927   double a2 = -500.0;  double b2 =  420.0;
 928   double a3 = -216.7;  double b3 =  250.0;
 929
 930   if (mom<0.2) return false;
 931
 932   if (mom>0.3 && mom<0.4) {
 933     if (dEdx < a1*mom+b1) return false;
 934   }
 935   else  if (mom<0.6) {
 936     if (dEdx < a2*mom+b2) return false;
 937   }
 938   else  if (mom<0.9) {
 939     if (dEdx < a3*mom+b3) return false;
 940   }
 941
 942   return true;
 943
 944 }
 945
 946 bool AliFemtoMJTrackCut::IsPionTOFTime(float mom, float ttof)
 947 {
 948   double a1 = -427.0; double b1 =  916.0;
 949   double a2 =  327.0; double b2 = -888.0;
 950   if (mom<0.3) return kFALSE;
 951   if (mom>2.0) return kFALSE;
 952   if (ttof > a1*mom+b1) return kFALSE;
 953   if (ttof < a2*mom+b2) return kFALSE;
 954
 955   return kTRUE;
 956 }
 957
 958 bool AliFemtoMJTrackCut::IsKaonTOFTime(float mom, float ttof)
 959 {
 960   double a1 =   000.0; double b1 =  -500.0;
 961   double a2 =   000.0; double b2 =   500.0;
 962   double a3 =   850.0; double b3 = -1503.0;
 963   double a4 = -1637.0; double b4 =  3621.0;
 964
 965   if (mom<0.3) return kFALSE;
 966   if (mom>2.06) return kFALSE;
 967   if (mom<1.2) {
 968     if (ttof > a2*mom+b2) return kFALSE;
 969     if (ttof < a1*mom+b1) return kFALSE;
 970   }
 971   if (mom<1.9) {
 972     if (ttof > a2*mom+b2) return kFALSE;
 973     if (ttof < a3*mom+b3) return kFALSE;
 974   }
 975   if (mom<2.06) {
 976     if (ttof > a4*mom+b4) return kFALSE;
 977     if (ttof < a3*mom+b3) return kFALSE;
 978   }
 979   return kTRUE;
 980 }
 981
 982 bool AliFemtoMJTrackCut::IsProtonTOFTime(float mom, float ttof)
 983 {
 984   double a1 =   000.0; double b1 =  -915.0;
 985   double a2 =   000.0; double b2 =   600.0;
 986   double a3 =   572.0; double b3 = -1715.0;
 987
 988   if (mom<0.3) return kFALSE;
 989   if (mom>3.0) return kFALSE;
 990   if (mom<1.4) {
 991     if (ttof > a2*mom+b2) return kFALSE;
 992     if (ttof < a1*mom+b1) return kFALSE;
 993   }
 994   if (mom<3.0) {
 995     if (ttof > a2*mom+b2) return kFALSE;
 996     if (ttof < a3*mom+b3) return kFALSE;
 997   }
 998   return kTRUE;
 999 }
1000
1001
1002
1003
1004 bool AliFemtoMJTrackCut::IsKaonTPCdEdxNSigma(float mom, float nsigmaK)
1005 {
1006 //  cout<<" AliFemtoMJTrackCut::IsKaonTPCdEdxNSigma "<<mom<<" "<<nsigmaK<<endl;
1007
1008
1009   if(mom<0.35 && TMath::Abs(nsigmaK)<5.0)return true;
1010   if(mom>=0.35 && mom<0.5 && TMath::Abs(nsigmaK)<3.0)return true;
1011   if(mom>=0.5 && mom<0.7 && TMath::Abs(nsigmaK)<2.0)return true;
1012
1013   return false;
1014 }
1015
1016 bool AliFemtoMJTrackCut::IsKaonTOFNSigma(float mom, float nsigmaK)
1017 {
1018 //  cout<<" AliFemtoMJTrackCut::IsKaonTPCdEdxNSigma "<<mom<<" "<<nsigmaK<<endl;
1019   //fan
1020   //  if(mom<1.5 && TMath::Abs(nsigmaK)<3.0)return true;
1021   if(mom>=1.5 && TMath::Abs(nsigmaK)<2.0)return true;
1022   return false;
1023 }
1024
1025 /*
1026 bool AliFemtoMJTrackCut::IsKaonNSigma(float mom, float nsigmaTPCK, float nsigmaTOFK)
1027 {
1028
1029
1030   if(mom<0.5)
1031     {
1032           if(TMath::Abs(nsigmaTPCK)<2.0)
1033            {
1034            return true;
1035            }
1036            else
1037            {
1038            return false;
1039            }
1040     }
1041
1042
1043    if(mom>=0.5)
1044     {
1045          if(TMath::Abs(nsigmaTOFK)<3.0 && TMath::Abs(nsigmaTPCK)<3.0)
1046          {
1047          return true;
1048          }
1049          else
1050          {
1051          return false;
1052          }
1053     }
1054
1055 //   if(mom>1.5 || mom<0.15)return false;
1056
1057
1058 }
1059
1060 */
1061
1062 //old
1063 bool AliFemtoMJTrackCut::IsKaonNSigma(float mom, float nsigmaTPCK, float nsigmaTOFK)
1064 {
1065   if (fNsigmaTPCTOF) {
1066     if (mom > 0.5) {
1067       //        if (TMath::Hypot( nsigmaTOFP, nsigmaTPCP )/TMath::Sqrt(2) < 3.0)
1068       if (TMath::Hypot( nsigmaTOFK, nsigmaTPCK ) < fNsigma)
1069         return true;
1070     }
1071     else {
1072       if (TMath::Abs(nsigmaTPCK) < fNsigma)
1073         return true;
1074     }
1075   }
1076   else {
1077
1078     if(mom<0.4)
1079       {
1080         if(nsigmaTOFK<-999.)
1081           {
1082             if(TMath::Abs(nsigmaTPCK)<2.0) return true;
1083           }
1084         else if(TMath::Abs(nsigmaTOFK)<3.0 && TMath::Abs(nsigmaTPCK)<3.0) return true;
1085       }
1086     else if(mom>=0.4 && mom<=0.6)
1087       {
1088         if(nsigmaTOFK<-999.)
1089           {
1090             if(TMath::Abs(nsigmaTPCK)<2.0) return true;
1091           }
1092         else if(TMath::Abs(nsigmaTOFK)<3.0 && TMath::Abs(nsigmaTPCK)<3.0) return true;
1093       }
1094     else if(nsigmaTOFK<-999.)
1095       {
1096         return false;
1097       }
1098     else if(TMath::Abs(nsigmaTOFK)<3.0 && TMath::Abs(nsigmaTPCK)<3.0) return true;
1099   }
1100   return false;
1101 }
1102
1103
1104
1105 bool AliFemtoMJTrackCut::IsPionNSigma(float mom, float nsigmaTPCPi, float nsigmaTOFPi)
1106 {
1107   if (fNsigmaTPCTOF) {
1108     if (mom > 0.5) {
1109       //        if (TMath::Hypot( nsigmaTOFP, nsigmaTPCP )/TMath::Sqrt(2) < 3.0)
1110       if (TMath::Hypot( nsigmaTOFPi, nsigmaTPCPi ) < fNsigma)
1111         return true;
1112     }
1113     else {
1114       if (TMath::Abs(nsigmaTPCPi) < fNsigma)
1115         return true;
1116     }
1117   }
1118   else {
1119     if(mom<0.65)
1120       {
1121         if(nsigmaTOFPi<-999.)
1122           {
1123             if(mom<0.35 && TMath::Abs(nsigmaTPCPi)<3.0) return true;
1124             else if(mom<0.5 && mom>=0.35 && TMath::Abs(nsigmaTPCPi)<3.0) return true;
1125             else if(mom>=0.5 && TMath::Abs(nsigmaTPCPi)<2.0) return true;
1126             else return false;
1127           }
1128         else if(TMath::Abs(nsigmaTOFPi)<3.0 && TMath::Abs(nsigmaTPCPi)<3.0) return true;
1129       }
1130     else if(nsigmaTOFPi<-999.)
1131       {
1132         return false;
1133       }
1134     else if(mom<1.5 && TMath::Abs(nsigmaTOFPi)<3.0 && TMath::Abs(nsigmaTPCPi)<5.0) return true;
1135     else if(mom>=1.5 && TMath::Abs(nsigmaTOFPi)<2.0 && TMath::Abs(nsigmaTPCPi)<5.0) return true;
1136   }
1137   return false;
1138 }
1139
1140
1141 bool AliFemtoMJTrackCut::IsProtonNSigma(float mom, float nsigmaTPCP, float nsigmaTOFP)
1142 {
1143
1144   if (fNsigmaTPCTOF) {
1145     if (mom > 0.5) {
1146 //        if (TMath::Hypot( nsigmaTOFP, nsigmaTPCP )/TMath::Sqrt(2) < 3.0)
1147         if (TMath::Hypot( nsigmaTOFP, nsigmaTPCP ) < fNsigma)
1148             return true;
1149         }
1150     else {
1151         if (TMath::Abs(nsigmaTPCP) < fNsigma)
1152             return true;
1153     }
1154   }
1155   else if (fNsigmaTPConly) {
1156     if (TMath::Abs(nsigmaTPCP) < fNsigma)
1157       return true;
1158   }
1159   else {
1160     if (mom > 0.8 && mom < 2.5) {
1161       if ( TMath::Abs(nsigmaTPCP) < 3.0 && TMath::Abs(nsigmaTOFP) < 3.0)
1162         return true;
1163     }
1164     else if (mom > 2.5) {
1165       if ( TMath::Abs(nsigmaTPCP) < 3.0 && TMath::Abs(nsigmaTOFP) < 2.0)
1166         return true;
1167     }
1168     else {
1169       if (TMath::Abs(nsigmaTPCP) < 3.0)
1170         return true;
1171     }
1172   }
1173
1174   return false;
1175 }
1176
1177
1178 void AliFemtoMJTrackCut::SetPIDMethod(ReadPIDMethodType newMethod)
1179 {
1180   fPIDMethod = newMethod;
1181 }
1182
1183
1184 void AliFemtoMJTrackCut::SetNsigmaTPCTOF(Bool_t nsigma)
1185 {
1186   fNsigmaTPCTOF = nsigma;
1187 }
1188
1189 void AliFemtoMJTrackCut::SetNsigmaTPConly(Bool_t nsigma)
1190 {
1191   fNsigmaTPConly = nsigma;
1192 }
1193
1194 void AliFemtoMJTrackCut::SetNsigma(Double_t nsigma)
1195 {
1196   fNsigma = nsigma;
1197 }
1198
1199
1200 void AliFemtoMJTrackCut::SetClusterRequirementITS(AliESDtrackCuts::Detector det, AliESDtrackCuts::ITSClusterRequirement req)
1201 {
1202   fCutClusterRequirementITS[det] = req;
1203 }
1204
1205 Bool_t AliFemtoMJTrackCut::CheckITSClusterRequirement(AliESDtrackCuts::ITSClusterRequirement req, Bool_t clusterL1, Bool_t clusterL2)
1206 {
1207   // checks if the cluster requirement is fullfilled (in this case: return kTRUE)
1208
1209   switch (req)
1210     {
1211     case AliESDtrackCuts::kOff:        return kTRUE;
1212     case AliESDtrackCuts::kNone:       return !clusterL1 && !clusterL2;
1213     case AliESDtrackCuts::kAny:        return clusterL1 || clusterL2;
1214     case AliESDtrackCuts::kFirst:      return clusterL1;
1215     case AliESDtrackCuts::kOnlyFirst:  return clusterL1 && !clusterL2;
1216     case AliESDtrackCuts::kSecond:     return clusterL2;
1217     case AliESDtrackCuts::kOnlySecond: return clusterL2 && !clusterL1;
1218     case AliESDtrackCuts::kBoth:       return clusterL1 && clusterL2;
1219   }
1220
1221   return kFALSE;
1222 }
1223
1224 bool AliFemtoMJTrackCut::IsElectron(float nsigmaTPCE, float nsigmaTPCPi,float nsigmaTPCK, float nsigmaTPCP)
1225 {
1226   if(TMath::Abs(nsigmaTPCE)<3 && TMath::Abs(nsigmaTPCPi)>3 && TMath::Abs(nsigmaTPCK)>3 && TMath::Abs(nsigmaTPCP)>3)
1227       return false;
1228    else
1229      return true;
1230 }