PWG2/FEMTOSCOPY/AliFemto/AliFemtoAODTrackCut.cxx

   1 /***************************************************************************
   2  *
   3  * $Id: AliFemtoAODTrackCut.cxx 23733 2008-02-05 16:51:41Z akisiel $
   4  *
   5  *
   6  ***************************************************************************
   7  *
   8  *
   9  *
  10  *
  11  ***************************************************************************
  12  *
  13  * $Log$
  14  * Revision 1.3  2007/05/22 09:01:42  akisiel
  15  * Add the possibiloity to save cut settings in the ROOT file
  16  *
  17  * Revision 1.2  2007/05/21 10:38:25  akisiel
  18  * More coding rule conformance
  19  *
  20  * Revision 1.1  2007/05/16 10:25:06  akisiel
  21  * Making the directory structure of AliFemtoUser flat. All files go into one common directory
  22  *
  23  * Revision 1.4  2007/05/03 09:46:10  akisiel
  24  * Fixing Effective C++ warnings
  25  *
  26  * Revision 1.3  2007/04/27 07:25:59  akisiel
  27  * Make revisions needed for compilation from the main AliRoot tree
  28  *
  29  * Revision 1.1.1.1  2007/04/25 15:38:41  panos
  30  * Importing the HBT code dir
  31  *
  32  * Revision 1.4  2007-04-03 16:00:08  mchojnacki
  33  * Changes to iprove memory managing
  34  *
  35  * Revision 1.3  2007/03/13 15:30:03  mchojnacki
  36  * adding reader for simulated data
  37  *
  38  * Revision 1.2  2007/03/08 14:58:03  mchojnacki
  39  * adding some alice stuff
  40  *
  41  * Revision 1.1.1.1  2007/03/07 10:14:49  mchojnacki
  42  * First version on CVS
  43  *
  44  **************************************************************************/
  45
  46 #include "AliFemtoAODTrackCut.h"
  47 #include <cstdio>
  48
  49 #ifdef __ROOT__
  50 ClassImp(AliFemtoAODTrackCut)
  51 #endif
  52
  53
  54 // electron
  55 // 0.13 - 1.8
  56 // 0       7.594129e-02    8.256141e-03
  57 // 1       -5.535827e-01   8.170825e-02
  58 // 2       1.728591e+00    3.104210e-01
  59 // 3       -2.827893e+00   5.827802e-01
  60 // 4       2.503553e+00    5.736207e-01
  61 // 5       -1.125965e+00   2.821170e-01
  62 // 6       2.009036e-01    5.438876e-02
  63
  64 // pion
  65 // 0.13 - 2.0
  66 // 0       1.063457e+00    8.872043e-03
  67 // 1       -4.222208e-01   2.534402e-02
  68 // 2       1.042004e-01    1.503945e-02
  69
  70 // kaon
  71 // 0.18 - 2.0
  72 // 0       -7.289406e-02   1.686074e-03
  73 // 1       4.415666e-01    1.143939e-02
  74 // 2       -2.996790e-01   1.840964e-02
  75 // 3       6.704652e-02    7.783990e-03
  76
  77 // proton
  78 // 0.26 - 2.0
  79 // 0       -3.730200e-02   2.347311e-03
  80 // 1       1.163684e-01    1.319316e-02
  81 // 2       8.354116e-02    1.997948e-02
  82 // 3       -4.608098e-02   8.336400e-03
  83
  84
  85 AliFemtoAODTrackCut::AliFemtoAODTrackCut() :
  86     fCharge(0),
  87     fLabel(0),
  88     fMaxchiNdof(1000.0),
  89     fMaxSigmaToVertex(1000.0),
  90     fNTracksPassed(0),
  91     fNTracksFailed(0),
  92     fMostProbable(0)
  93 {
  94   // Default constructor
  95   fNTracksPassed = fNTracksFailed = 0;
  96   fCharge = 0;  // takes both charges 0
  97   fPt[0]=0.0;              fPt[1] = 100.0;//100
  98   fRapidity[0]=-2;       fRapidity[1]=2;//-2 2
  99   fPidProbElectron[0]=-1;fPidProbElectron[1]=2;
 100   fPidProbPion[0]=-1;    fPidProbPion[1]=2;
 101   fPidProbKaon[0]=-1;fPidProbKaon[1]=2;
 102   fPidProbProton[0]=-1;fPidProbProton[1]=2;
 103   fPidProbMuon[0]=-1;fPidProbMuon[1]=2;
 104   fLabel=false;
 105 }
 106 //------------------------------
 107 AliFemtoAODTrackCut::~AliFemtoAODTrackCut(){
 108   /* noop */
 109 }
 110 //------------------------------
 111 bool AliFemtoAODTrackCut::Pass(const AliFemtoTrack* track)
 112 {
 113   // test the particle and return
 114   // true if it meets all the criteria
 115   // false if it doesn't meet at least one of the criteria
 116   float tMost[5];
 117
 118   if (((track->ITSchi2() + track->TPCchi2())/(track->ITSncls() + track->TPCncls())) > fMaxchiNdof) {
 119     return false;
 120   }
 121
 122   if (fMaxSigmaToVertex < track->SigmaToVertex()) {
 123     return false;
 124   }
 125
 126   if (fLabel)
 127     {
 128       //cout<<"labels"<<endl;
 129       if(track->Label()<0)
 130         {
 131           fNTracksFailed++;
 132           //   cout<<"No Go Through the cut"<<endl;
 133           //  cout<<fLabel<<" Label="<<track->Label()<<endl;
 134           return false;
 135         }
 136     }
 137   if (fCharge!=0)
 138     {
 139       //cout<<"AliFemtoESD  cut ch "<<endl;
 140       //cout<<fCharge<<" Charge="<<track->Charge()<<endl;
 141       if (track->Charge()!= fCharge)
 142         {
 143           fNTracksFailed++;
 144           //  cout<<"No Go Through the cut"<<endl;
 145           // cout<<fCharge<<" Charge="<<track->Charge()<<endl;
 146           return false;
 147         }
 148     }
 149   float tEnergy = ::sqrt(track->P().mag2()+fMass*fMass);
 150   float tRapidity = 0.5*::log((tEnergy+track->P().z())/(tEnergy-track->P().z()));
 151   float tPt = ::sqrt((track->P().x())*(track->P().x())+(track->P().y())*(track->P().y()));
 152   if ((tRapidity<fRapidity[0])||(tRapidity>fRapidity[1]))
 153     {
 154       fNTracksFailed++;
 155       //cout<<"No Go Through the cut"<<endl;
 156       //cout<<fRapidity[0]<<" < Rapidity ="<<tRapidity<<" <"<<fRapidity[1]<<endl;
 157       return false;
 158     }
 159   if ((tPt<fPt[0])||(tPt>fPt[1]))
 160     {
 161       fNTracksFailed++;
 162       //cout<<"No Go Through the cut"<<endl;
 163       //cout<<fPt[0]<<" < Pt ="<<Pt<<" <"<<fPt[1]<<endl;
 164       return false;
 165     }
 166
 167
 168   if ((track->PidProbElectron()<fPidProbElectron[0])||(track->PidProbElectron()>fPidProbElectron[1]))
 169     {
 170       fNTracksFailed++;
 171       //cout<<"No Go Through the cut"<<endl;
 172       //cout<<fPidProbElectron[0]<<" < e ="<<track->PidProbElectron()<<" <"<<fPidProbElectron[1]<<endl;
 173       return false;
 174     }
 175   if ((track->PidProbPion()<fPidProbPion[0])||(track->PidProbPion()>fPidProbPion[1]))
 176     {
 177       fNTracksFailed++;
 178       //cout<<"No Go Through the cut"<<endl;
 179       //cout<<fPidProbPion[0]<<" < pi ="<<track->PidProbPion()<<" <"<<fPidProbPion[1]<<endl;
 180       return false;
 181     }
 182   if ((track->PidProbKaon()<fPidProbKaon[0])||(track->PidProbKaon()>fPidProbKaon[1]))
 183     {
 184       fNTracksFailed++;
 185       //cout<<"No Go Through the cut"<<endl;
 186       //cout<<fPidProbKaon[0]<<" < k ="<<track->PidProbKaon()<<" <"<<fPidProbKaon[1]<<endl;
 187       return false;
 188     }
 189   if ((track->PidProbProton()<fPidProbProton[0])||(track->PidProbProton()>fPidProbProton[1]))
 190     {
 191       fNTracksFailed++;
 192       //cout<<"No Go Through the cut"<<endl;
 193       //cout<<fPidProbProton[0]<<" < p  ="<<track->PidProbProton()<<" <"<<fPidProbProton[1]<<endl;
 194       return false;
 195     }
 196   if ((track->PidProbMuon()<fPidProbMuon[0])||(track->PidProbMuon()>fPidProbMuon[1]))
 197     {
 198       fNTracksFailed++;
 199       //cout<<"No Go Through the cut"<<endl;
 200       //cout<<fPidProbMuon[0]<<" <  mi="<<track->PidProbMuon()<<" <"<<fPidProbMuon[1]<<endl;
 201       return false;
 202     }
 203
 204   if (fMostProbable) {
 205     tMost[0] = track->PidProbElectron()*PidFractionElectron(track->P().mag());
 206     tMost[1] = 0.0;
 207     tMost[2] = track->PidProbPion()*PidFractionPion(track->P().mag());
 208     tMost[3] = track->PidProbKaon()*PidFractionKaon(track->P().mag());
 209     tMost[4] = track->PidProbProton()*PidFractionProton(track->P().mag());
 210     int imost=0;
 211     float ipidmax = 0.0;
 212     for (int ip=0; ip<5; ip++)
 213       if (tMost[ip] > ipidmax) { ipidmax = tMost[ip]; imost = ip; };
 214     if (imost != fMostProbable) return false;
 215   }
 216
 217   // cout<<"Go Through the cut"<<endl;
 218   // cout<<fLabel<<" Label="<<track->Label()<<endl;
 219   // cout<<fCharge<<" Charge="<<track->Charge()<<endl;
 220   // cout<<fPt[0]<<" < Pt ="<<Pt<<" <"<<fPt[1]<<endl;
 221   //cout<<fRapidity[0]<<" < Rapidity ="<<tRapidity<<" <"<<fRapidity[1]<<endl;
 222   //cout<<fPidProbElectron[0]<<" <  e="<<track->PidProbElectron()<<" <"<<fPidProbElectron[1]<<endl;
 223   //cout<<fPidProbPion[0]<<" <  pi="<<track->PidProbPion()<<" <"<<fPidProbPion[1]<<endl;
 224   //cout<<fPidProbKaon[0]<<" <  k="<<track->PidProbKaon()<<" <"<<fPidProbKaon[1]<<endl;
 225   //cout<<fPidProbProton[0]<<" <  p="<<track->PidProbProton()<<" <"<<fPidProbProton[1]<<endl;
 226   //cout<<fPidProbMuon[0]<<" <  mi="<<track->PidProbMuon()<<" <"<<fPidProbMuon[1]<<endl;
 227   fNTracksPassed++ ;
 228   return true;
 229
 230
 231 }
 232 //------------------------------
 233 AliFemtoString AliFemtoAODTrackCut::Report()
 234 {
 235   // Prepare report from the execution
 236   string tStemp;
 237   char tCtemp[100];
 238   sprintf(tCtemp,"Particle mass:\t%E\n",this->Mass());
 239   tStemp=tCtemp;
 240   sprintf(tCtemp,"Particle charge:\t%d\n",fCharge);
 241   tStemp+=tCtemp;
 242   sprintf(tCtemp,"Particle pT:\t%E - %E\n",fPt[0],fPt[1]);
 243   tStemp+=tCtemp;
 244   sprintf(tCtemp,"Particle rapidity:\t%E - %E\n",fRapidity[0],fRapidity[1]);
 245   tStemp+=tCtemp;
 246   sprintf(tCtemp,"Number of tracks which passed:\t%ld  Number which failed:\t%ld\n",fNTracksPassed,fNTracksFailed);
 247   tStemp += tCtemp;
 248   AliFemtoString returnThis = tStemp;
 249   return returnThis;
 250 }
 251 TList *AliFemtoAODTrackCut::ListSettings()
 252 {
 253   // return a list of settings in a writable form
 254   TList *tListSetttings = new TList();
 255   char buf[200];
 256   snprintf(buf, 200, "AliFemtoAODTrackCut.mass=%lf", this->Mass());
 257   tListSetttings->AddLast(new TObjString(buf));
 258
 259   snprintf(buf, 200, "AliFemtoAODTrackCut.charge=%i", fCharge);
 260   tListSetttings->AddLast(new TObjString(buf));
 261   snprintf(buf, 200, "AliFemtoAODTrackCut.pidprobpion.minimum=%lf", fPidProbPion[0]);
 262   tListSetttings->AddLast(new TObjString(buf));
 263   snprintf(buf, 200, "AliFemtoAODTrackCut.pidprobpion.maximum=%lf", fPidProbPion[1]);
 264   tListSetttings->AddLast(new TObjString(buf));
 265   snprintf(buf, 200, "AliFemtoAODTrackCut.pidprobkaon.minimum=%lf", fPidProbKaon[0]);
 266   tListSetttings->AddLast(new TObjString(buf));
 267   snprintf(buf, 200, "AliFemtoAODTrackCut.pidprobkaon.maximum=%lf", fPidProbKaon[1]);
 268   tListSetttings->AddLast(new TObjString(buf));
 269   snprintf(buf, 200, "AliFemtoAODTrackCut.pidprobproton.minimum=%lf", fPidProbProton[0]);
 270   tListSetttings->AddLast(new TObjString(buf));
 271   snprintf(buf, 200, "AliFemtoAODTrackCut.pidprobproton.maximum=%lf", fPidProbProton[1]);
 272   tListSetttings->AddLast(new TObjString(buf));
 273   snprintf(buf, 200, "AliFemtoAODTrackCut.pidprobelectron.minimum=%lf", fPidProbElectron[0]);
 274   tListSetttings->AddLast(new TObjString(buf));
 275   snprintf(buf, 200, "AliFemtoAODTrackCut.pidprobelectron.maximum=%lf", fPidProbElectron[1]);
 276   tListSetttings->AddLast(new TObjString(buf));
 277   snprintf(buf, 200, "AliFemtoAODTrackCut.pidprobMuon.minimum=%lf", fPidProbMuon[0]);
 278   tListSetttings->AddLast(new TObjString(buf));
 279   snprintf(buf, 200, "AliFemtoAODTrackCut.pidprobMuon.maximum=%lf", fPidProbMuon[1]);
 280   tListSetttings->AddLast(new TObjString(buf));
 281   snprintf(buf, 200, "AliFemtoAODTrackCut.pt.minimum=%lf", fPt[0]);
 282   tListSetttings->AddLast(new TObjString(buf));
 283   snprintf(buf, 200, "AliFemtoAODTrackCut.pt.maximum=%lf", fPt[1]);
 284   tListSetttings->AddLast(new TObjString(buf));
 285   snprintf(buf, 200, "AliFemtoAODTrackCut.rapidity.minimum=%lf", fRapidity[0]);
 286   tListSetttings->AddLast(new TObjString(buf));
 287   snprintf(buf, 200, "AliFemtoAODTrackCut.rapidity.maximum=%lf", fRapidity[1]);
 288   tListSetttings->AddLast(new TObjString(buf));
 289   snprintf(buf, 200, "AliFemtoAODTrackCut.maxchindof=%lf", fMaxchiNdof);
 290   tListSetttings->AddLast(new TObjString(buf));
 291   snprintf(buf, 200, "AliFemtoAODTrackCut.maxsigmatovertex=%lf", fMaxSigmaToVertex);
 292   tListSetttings->AddLast(new TObjString(buf));
 293   if (fMostProbable) {
 294     if (fMostProbable == 2)
 295       snprintf(buf, 200, "AliFemtoAODTrackCut.mostprobable=%s", "Pion");
 296     if (fMostProbable == 3)
 297       snprintf(buf, 200, "AliFemtoAODTrackCut.mostprobable=%s", "Kaon");
 298     if (fMostProbable == 4)
 299       snprintf(buf, 200, "AliFemtoAODTrackCut.mostprobable=%s", "Proton");
 300     tListSetttings->AddLast(new TObjString(buf));
 301   }
 302   return tListSetttings;
 303 }
 304                             // electron
 305 // 0.13 - 1.8
 306 // 0       7.594129e-02    8.256141e-03
 307 // 1       -5.535827e-01   8.170825e-02
 308 // 2       1.728591e+00    3.104210e-01
 309 // 3       -2.827893e+00   5.827802e-01
 310 // 4       2.503553e+00    5.736207e-01
 311 // 5       -1.125965e+00   2.821170e-01
 312 // 6       2.009036e-01    5.438876e-02
 313 float AliFemtoAODTrackCut::PidFractionElectron(float mom) const
 314 {
 315   // Provide a parameterized fraction of electrons dependent on momentum
 316   if (mom<0.13) return 0.0;
 317   if (mom>1.8) return 0.0;
 318   return (7.594129e-02
 319           -5.535827e-01*mom
 320           +1.728591e+00*mom*mom
 321           -2.827893e+00*mom*mom*mom
 322           +2.503553e+00*mom*mom*mom*mom
 323           -1.125965e+00*mom*mom*mom*mom*mom
 324           +2.009036e-01*mom*mom*mom*mom*mom*mom);
 325 }
 326
 327 // pion
 328 // 0.13 - 2.0
 329 // 0       1.063457e+00    8.872043e-03
 330 // 1       -4.222208e-01   2.534402e-02
 331 // 2       1.042004e-01    1.503945e-02
 332 float AliFemtoAODTrackCut::PidFractionPion(float mom) const
 333 {
 334   // Provide a parameterized fraction of pions dependent on momentum
 335   if (mom<0.13) return 0.0;
 336   if (mom>2.0) return 0.0;
 337   return ( 1.063457e+00
 338            -4.222208e-01*mom
 339            +1.042004e-01*mom*mom);
 340 }
 341
 342 // kaon
 343 // 0.18 - 2.0
 344 // 0       -7.289406e-02   1.686074e-03
 345 // 1       4.415666e-01    1.143939e-02
 346 // 2       -2.996790e-01   1.840964e-02
 347 // 3       6.704652e-02    7.783990e-03
 348 float AliFemtoAODTrackCut::PidFractionKaon(float mom) const
 349 {
 350   // Provide a parameterized fraction of kaons dependent on momentum
 351   if (mom<0.18) return 0.0;
 352   if (mom>2.0) return 0.0;
 353   return (-7.289406e-02
 354           +4.415666e-01*mom
 355           -2.996790e-01*mom*mom
 356           +6.704652e-02*mom*mom*mom);
 357 }
 358
 359 // proton
 360 // 0.26 - 2.0
 361 // 0       -3.730200e-02   2.347311e-03
 362 // 1       1.163684e-01    1.319316e-02
 363 // 2       8.354116e-02    1.997948e-02
 364 // 3       -4.608098e-02   8.336400e-03
 365 float AliFemtoAODTrackCut::PidFractionProton(float mom) const
 366 {
 367   // Provide a parameterized fraction of protons dependent on momentum
 368   if (mom<0.26) return  0.0;
 369   if (mom>2.0) return 0.0;
 370   return (-3.730200e-02
 371           +1.163684e-01*mom
 372           +8.354116e-02*mom*mom
 373           -4.608098e-02*mom*mom*mom);
 374 }