[u/mrichter/AliRoot.git] / ITS / AliITSPid.cxx

//------------------------------------------------------------//
// Class for identification of pions,kaons and protons in ITS //
// Prior particles population (probabilities) are taken from  //
// Hijing event generator.                                   //
//------------------------------------------------------------//
// #include <stdlib.h>
#include "AliITSPid.h"
//#include "TMath.h"
#include <Riostream.h>
#include <TClonesArray.h>
#include <TVector.h>
#include "AliKalmanTrack.h"
#include "AliITSIOTrack.h"
#include "AliITStrackV2.h"
#include <TF1.h>

ClassImp(AliITSPid)


AliITSPid::AliITSPid(const AliITSPid &source) : TObject(source),
fMxtrs(source.fMxtrs),
fTrs(source.fTrs),
fWpi(source.fWpi),
fWk(source.fWk),
fWp(source.fWp),
fRpik(source.fRpik),
fRppi(source.fRppi),
fRpka(source.fRpka),
fRp(source.fRp),
fPcode(source.fPcode),
fSigmin(source.fSigmin),
fSilent(source.fSilent),
fCutKa(source.fCutKa),
fCutPr(source.fCutPr),
fggpi(source.fggpi),
fggka(source.fggka),
fggpr(source.fggpr){
    // Copy constructor. This is a function which is not allowed to be

  
}
  
//______________________________________________________________________
AliITSPid& AliITSPid::operator=(const AliITSPid& source){
  // Assignment operator. This is a function which is not allowed to be
  this->~AliITSPid();
  new(this) AliITSPid(source);
  return *this;
}

//
  Float_t AliITSPid::Qtrm(Int_t track) {
//
// This calculates truncated mean signal for given track.
//
    TVector q(*( this->GetVec(track)  ));
    Int_t ml=(Int_t)q(0);
    if(ml<1)return 0.;
    if(ml>6)ml=6;
    float vf[6];
    Int_t nl=0; for(Int_t i=0;i<ml;i++){if(q(i)>fSigmin){vf[nl]=q(i+1);nl++;}}
    if(nl==0)return 0.;
    switch(nl){
      case  1:q(6)=q(1); break;
      case  2:q(6)=(q(1)+q(2))/2.; break;
      default:
	for(int fi=0;fi<2;fi++){
	  Int_t swap;
	  do{ swap=0; float qmin=vf[fi];
	  for(int j=fi+1;j<nl;j++)
	    if(qmin>vf[j]){qmin=vf[j];vf[j]=vf[fi];vf[fi]=qmin;swap=1;};
	  }while(swap==1);
	}
        q(6)= (vf[0]+vf[1])/2.;
        break;
    }
    for(Int_t i=0;i<nl;i++){q(i+1)=vf[i];} this->SetVec(track,q);
    return (q(6));
}

Float_t AliITSPid::Qtrm(Float_t qarr[6],Int_t narr) const{
//
//This calculates truncated mean signal for given signal set.
//
  Float_t q[6],qm,qmin;
  Int_t nl,ml;
  if(narr>0&&narr<7){ml=narr;}else{return 0;};
  nl=0; for(Int_t i=0;i<ml;i++){if(qarr[i]>fSigmin){q[nl]=qarr[i];nl++;}}
  if(nl==0)return 0.;
    switch(nl){
      case  1:qm=q[0]; break;
      case  2:qm=(q[0]+q[1])/2.; break;
      default:
	Int_t swap;
	for(int fi=0;fi<2;fi++){
	  do{ swap=0; qmin=q[fi];
	  for(int j=fi+1;j<nl;j++)
	    if(qmin>q[j]){qmin=q[j];q[j]=q[fi];q[fi]=qmin;swap=1;};
	  }while(swap==1);
	}
        qm= (q[0]+q[1])/2.;
        break;
    }
    return qm;
}

Int_t	AliITSPid::Wpik(Float_t pm,Float_t q){
  //Calcutates probabilityes of pions and kaons
  //Returns particle code for dominant probability.
  Double_t par[6];
  for(int i=0;i<6;i++){par[i]=fGGpi[i]->Eval(pm);}
  fggpi->SetParameters(par);

  for(int i=0;i<3;i++){par[i]=fGGka[i]->Eval(pm);}
  fggka->SetParameters(par);

  Float_t ppi=fggpi->Eval(q);
  Float_t pka=fggka->Eval(q);
  Float_t p=ppi+pka;
  /*
  if(!fSilent){
    fggka->Print();
    fggpi->Print();
    if(p>0)cout<<" ppi,pka="<<ppi/p<<"  "<<pka/p<<endl;
  }
  */

  if(p>0){
    ppi=ppi/p; 
    pka=pka/p;
    fWp=0.; fWpi=ppi; fWk=pka;
    if( pka>ppi){return fPcode=321;}else{return fPcode=211;}
  }else{return 0;}
}
//-----------------------------------------------------------
Int_t	AliITSPid::Wpikp(Float_t pm,Float_t q){
  //
  //Calcutates probabilityes of pions,kaons and protons.
  //Returns particle code for dominant probability.
  Double_t par[6];
  for(int i=0;i<6;i++){par[i]=fGGpi[i]->Eval(pm);}
  fggpi->SetParameters(par);

  for(int i=0;i<3;i++){par[i]=fGGka[i]->Eval(pm);}
  fggka->SetParameters(par);

  for(int i=0;i<3;i++){par[i]=fGGpr[i]->Eval(pm);}
  fggpr->SetParameters(par);

  Float_t p,ppi,pka,ppr;
  if( q>(fggpr->GetParameter(1)+fggpr->GetParameter(2)) )
      { p=1.0; ppr=1.0; ppi=pka=0.0;
    }else{ 
    ppi=fggpi->Eval(q);
    pka=fggka->Eval(q);
    ppr=fggpr->Eval(q);
    p=ppi+pka+ppr;
  }
  if(p>0){
    ppi=ppi/p; 
    pka=pka/p;
    ppr=ppr/p;
    fWp=ppr; fWpi=ppi; fWk=pka;
    //if(!fSilent)cout<<" ppi,pka,ppr="<<ppi<<"  "<<pka<<" "<<ppr<<endl;

   if( ppi>pka&&ppi>ppr )
           {return fPcode=211;}
   else{ if(pka>ppr){return fPcode=321;}else{return fPcode=2212;}
   }

  }else{return 0;}
}
//-----------------------------------------------------------
Int_t	AliITSPid::GetPcode(TClonesArray* rps,Float_t pm)
{
  //Returns particle code
  Info("GetPcode","method not implemented - Inputs TClonesArray *%x , Float_t %f",rps,pm); 
    return 0;    
}
//-----------------------------------------------------------
Int_t   AliITSPid::GetPcode(AliKalmanTrack *track)
{
  //Returns particle code for given track.
      Double_t xk,par[5]; track->GetExternalParameters(xk,par);
      Float_t phi=TMath::ASin(par[2]) + track->GetAlpha();
      if (phi<-TMath::Pi()) phi+=2*TMath::Pi();
      if (phi>=TMath::Pi()) phi-=2*TMath::Pi();
      Float_t lam=TMath::ATan(par[3]); 
      Float_t pt1=TMath::Abs(par[4]);
      Float_t mom=1./(pt1*TMath::Cos(lam));
      Float_t dedx=track->GetPIDsignal();
    Int_t pcode=GetPcode(dedx/40.,mom);
//    cout<<"TPCtrack dedx,mom,pcode="<<dedx<<","<<mom<<","<<pcode<<endl;
    return pcode?pcode:211;
    }
//------------------------------------------------------------
Int_t   AliITSPid::GetPcode(AliITSIOTrack *track)
{
  //Returns particle code for given track(V1).
    Double_t px,py,pz;
    px=track->GetPx();
    py=track->GetPy();
    pz=track->GetPz();
    Float_t mom=TMath::Sqrt(px*px+py*py+pz*pz);
//???????????????????
    // Float_t dedx=1.0;
    Float_t dedx=track->GetdEdx();
//???????????????????    
    Int_t pcode=GetPcode(dedx,mom);
//    cout<<"ITSV1 dedx,mom,pcode="<<dedx<<","<<mom<<","<<pcode<<endl;
return pcode?pcode:211;
}
//-----------------------------------------------------------
Int_t   AliITSPid::GetPcode(AliITStrackV2 *track)
{
//Returns particle code for given track(V2).
  if(track==0)return 0;
  //      track->Propagate(track->GetAlpha(),3.,0.1/65.19*1.848,0.1*1.848);
      track->PropagateTo(3.,0.0028,65.19);
      //track->PropagateToVertex();          Not needed. (I.B.)
    Double_t xk,par[5]; track->GetExternalParameters(xk,par);
    Float_t lam=TMath::ATan(par[3]);
    Float_t pt1=TMath::Abs(par[4]);
    Float_t mom=0.;
    if( (pt1*TMath::Cos(lam))!=0. ){ mom=1./(pt1*TMath::Cos(lam)); }else{mom=0.;};
    Float_t dedx=track->GetdEdx();
//    cout<<"lam,pt1,mom,dedx="<<lam<<","<<pt1<<","<<mom<<","<<dedx<<endl;
    Int_t pcode=GetPcode(dedx,mom);
//    cout<<"ITS V2 dedx,mom,pcode="<<dedx<<","<<mom<<","<<pcode<<endl;
return pcode?pcode:211;
}
//-----------------------------------------------------------
Int_t	AliITSPid::GetPcode(Float_t q,Float_t pm)
{
  //Returns particle code for given signal and momentum.
    fWpi=fWk=fWp=0.;     fPcode=0;

    if ( pm<=0.400 )
	{ if( q<fCutKa->Eval(pm) )
	    {return Pion();}
	else{ if( q<fCutPr->Eval(pm) )
		{return Kaon();}
	    else{return Proton();}
	    } 
	}
    if ( pm<=0.750 ){
	if ( q>fCutPr->Eval(pm)  ){
	  return Proton();
	}
	else {
	  return Wpik(pm,q);
	}
    }
    if( pm<=1.10 ){ return Wpikp(pm,q); }
    return fPcode;    
}
//-----------------------------------------------------------
void	AliITSPid::SetCut(Int_t n,Float_t pm,Float_t pilo,Float_t pihi,
			Float_t klo,Float_t khi,Float_t plo,Float_t phi)
{
  // The cut-table initializer method.
    fCut[n][0]=pm;
    fCut[n][1]=pilo;
    fCut[n][2]=pihi;
    fCut[n][3]=klo;
    fCut[n][4]=khi;
    fCut[n][5]=plo;
    fCut[n][6]=phi;
    return ;    
}
//------------------------------------------------------------
void AliITSPid::SetVec(Int_t ntrack,const TVector& info) const
{
  //Store track info in tracls table
TClonesArray& arr=*fTrs;
    new( arr[ntrack] ) TVector(info);
}
//-----------------------------------------------------------
TVector* AliITSPid::GetVec(Int_t ntrack) const
{
  //Get given track from track table 
TClonesArray& arr=*fTrs;
    return (TVector*)arr[ntrack];
}
//-----------------------------------------------------------
void AliITSPid::SetEdep(Int_t track,Float_t Edep)
{
  //Set dEdx for given track
    TVector xx(0,11);
    if( ((TVector*)fTrs->At(track))->IsValid() )
	{TVector yy( *((TVector*)fTrs->At(track)) );xx=yy; }
    Int_t j=(Int_t)xx(0); if(j>4)return;
    xx(++j)=Edep;xx(0)=j;
    TClonesArray &arr=*fTrs;
    new(arr[track])TVector(xx);
}
//-----------------------------------------------------------
void AliITSPid::SetPmom(Int_t track,Float_t Pmom)
{
  //Set momentum for given track
    TVector xx(0,11);
    if( ((TVector*)fTrs->At(track))->IsValid() )
	{TVector yy( *((TVector*)fTrs->At(track)) );xx=yy; }
    xx(10)=Pmom;
    TClonesArray &arr=*fTrs;
    new(arr[track])TVector(xx);
}
//-----------------------------------------------------------
void AliITSPid::SetPcod(Int_t track,Int_t partcode)
{
  //Set particle code for given track
    TVector xx(0,11);
    if( ((TVector*)fTrs->At(track))->IsValid() )
	{TVector yy( *((TVector*)fTrs->At(track)) );xx=yy; }
    if(xx(11)==0)
	{xx(11)=partcode; fMxtrs++;
	TClonesArray &arr=*fTrs;
	new(arr[track])TVector(xx);
	}
}
//-----------------------------------------------------------
void AliITSPid::Print(Int_t track)
{
  //Prints information for given track
cout<<fMxtrs<<" tracks in AliITSPid obj."<<endl;
    if( ((TVector*)fTrs->At(track))->IsValid() )
	{TVector xx( *((TVector*)fTrs->At(track)) );
	 xx.Print();
	 }
    else 
	{cout<<"No data for track "<<track<<endl;return;}
}
//-----------------------------------------------------------
void AliITSPid::Tab(void)
{
  //Make PID for tracks stored in tracks table
  if(fTrs->GetEntries()==0){cout<<"No entries in TAB"<<endl;return;}
  cout<<"------------------------------------------------------------------------"<<endl;
  cout<<"Nq"<<"   q1  "<<"   q2  "<<"   q3  "<<"   q4  "<<"   q5   "<<
    " Qtrm    "    <<"  Wpi  "<<"  Wk   "<<"  Wp  "<<"Pmom  "<<endl;
  cout<<"------------------------------------------------------------------------"<<endl;
  for(Int_t i=0;i<fTrs->GetEntries();i++)
    {
      TVector xxx( *((TVector*)fTrs->At(i)) );     
      if( xxx.IsValid() && xxx(0)>0 )
	{
	  TVector xx( *((TVector*)fTrs->At(i)) );
	  if(xx(0)>=2)
	    {
	      //       1)Calculate Qtrm	
	      xx(6)=(this->Qtrm(i));

	    }else{
	      xx(6)=xx(1);
	    }
	  //	 2)Calculate Wpi,Wk,Wp
	  this->GetPcode(xx(6),xx(10)/1000.);
	  xx(7)=GetWpi();
	  xx(8)=GetWk();
	  xx(9)=GetWp();
	  //       3)Print table
	  if(xx(0)>0){
	    //		    cout<<xx(0)<<" ";
	    for(Int_t j=1;j<11;j++){
	      if(i<7){ cout.width(7);cout.precision(4);cout<<xx(j);}
	      if(i>7){ cout.width(7);cout.precision(5);cout<<xx(j);}
	    }
	    cout<<endl;
	  }
	  //	  4)Update data in TVector
	  TClonesArray &arr=*fTrs;
	  new(arr[i])TVector(xx);	 
	}
      else 
	{/*cout<<"No data for track "<<i<<endl;*/}
    }// End loop for tracks
}
void AliITSPid::Reset(void)
{
  //Reset tracks table
  for(Int_t i=0;i<fTrs->GetEntries();i++){
    TVector xx(0,11);
    TClonesArray &arr=*fTrs;
    new(arr[i])TVector(xx);
  }
}
//-----------------------------------------------------------
AliITSPid::AliITSPid(Int_t ntrack):
fMxtrs(0),
fTrs(0),
fWpi(0),
fWk(0),
fWp(0),
fRpik(0),
fRppi(0),
fRpka(0),
fRp(0),
fPcode(0),
fSigmin(0.01),
fSilent(0),
fCutKa(0),
fCutPr(0),
fggpi(0),
fggka(0),
fggpr(0){
  //Constructor for AliITSPid class
 
    fTrs = new TClonesArray("TVector",ntrack);
    TClonesArray &arr=*fTrs;
    for(Int_t i=0;i<ntrack;i++)new(arr[i])TVector(0,11);
     //   
    fCutKa=new TF1("fcutka","pol4",0.05,0.4);
    Double_t ka[5]={25.616, -161.59, 408.97, -462.17, 192.86};
    fCutKa->SetParameters(ka);
    //
    fCutPr=new TF1("fcutpr","[0]/x/x+[1]",0.05,1.1);
    Double_t pr[2]={0.70675,0.4455};
    fCutPr->SetParameters(pr);
    //
    //---------- signal fit ----------
{//Pions
fGGpi[0]=new TF1("fp1pi","pol4",0.34,1.2);
  Double_t parpi0[10]={ -1.9096471071e+03, 4.5354331545e+04, -1.1860738840e+05,
   1.1405329025e+05, -3.8289694496e+04  };
  fGGpi[0]->SetParameters(parpi0);
fGGpi[1]=new TF1("fp2pi","[0]/x/x+[1]",0.34,1.2);
  Double_t parpi1[10]={ 1.0791668283e-02, 9.7347716496e-01  };
  fGGpi[1]->SetParameters(parpi1);
fGGpi[2]=new TF1("fp3pi","[0]/x/x+[1]",0.34,1.2);
  Double_t parpi2[10]={ 5.8191602279e-04, 9.7285601334e-02  };
  fGGpi[2]->SetParameters(parpi2);
fGGpi[3]=new TF1("fp4pi","pol4",0.34,1.2);
  Double_t parpi3[10]={ 6.6267353195e+02, 7.1595101104e+02, -5.3095111914e+03,
   6.2900977606e+03, -2.2935862292e+03  };
  fGGpi[3]->SetParameters(parpi3);
fGGpi[4]=new TF1("fp5pi","[0]/x/x+[1]",0.34,1.2);
  Double_t parpi4[10]={ 9.0419011783e-03, 1.1628922525e+00  };
  fGGpi[4]->SetParameters(parpi4);
fGGpi[5]=new TF1("fp6pi","[0]/x/x+[1]",0.34,1.2);
  Double_t parpi5[10]={ 1.8324872519e-03, 2.1503968838e-01  };
  fGGpi[5]->SetParameters(parpi5);
}//End Pions
{//Kaons
fGGka[0]=new TF1("fp1ka","pol4",0.24,1.2);
  Double_t parka0[20]={
  -1.1204243395e+02,4.6716191428e+01,2.2584059281e+03,
  -3.7123338009e+03,1.6003647641e+03  };
  fGGka[0]->SetParameters(parka0);
fGGka[1]=new TF1("fp2ka","[0]/x/x+[1]",0.24,1.2);
  Double_t parka1[20]={
  2.5181172905e-01,8.7566001814e-01  };
  fGGka[1]->SetParameters(parka1);
fGGka[2]=new TF1("fp3ka","pol6",0.24,1.2);
  Double_t parka2[20]={
  8.6236021573e+00,-7.0970427531e+01,2.4846827669e+02,
  -4.6094401290e+02,4.7546751408e+02,-2.5807112462e+02,
  5.7545491696e+01  };
  fGGka[2]->SetParameters(parka2);
}//End Kaons
{//Protons
fGGpr[0]=new TF1("fp1pr","pol4",0.4,1.2);
  Double_t parpr0[10]={
  6.0150106543e+01,-8.8176206410e+02,3.1222644604e+03,
  -3.5269200901e+03,1.2859128345e+03  };
  fGGpr[0]->SetParameters(parpr0);
fGGpr[1]=new TF1("fp2pr","[0]/x/x+[1]",0.4,1.2);
  Double_t parpr1[10]={
  9.4970837607e-01,7.3573504201e-01  };
  fGGpr[1]->SetParameters(parpr1);
fGGpr[2]=new TF1("fp3pr","[0]/x/x+[1]",0.4,1.2);
  Double_t parpr2[10]={
  1.2498403757e-01,2.7845072306e-02  };
  fGGpr[2]->SetParameters(parpr2);
}//End Protons
    //----------- end fit -----------

    fggpr=new TF1("ggpr","gaus",0.4,1.2);
    fggpi=new TF1("ggpi","gaus+gaus(3)",0.4,1.2);
    fggka=new TF1("ggka","gaus",0.4,1.2);

    //-------------------------------------------------
const int kInf=10;
//         Ncut Pmom   pilo  pihi    klo    khi     plo    phi
//       cut[j] [0]    [1]    [2]    [3]    [4]     [5]    [6]
//----------------------------------------------------------------
    SetCut(  1, 0.12 ,  0.  ,  0.  , kInf  , kInf   , kInf  , kInf  );
    SetCut(  2, 0.20 ,  0.  ,  6.0 , 6.0  , kInf   , kInf  , kInf  );
    SetCut(  3, 0.30 ,  0.  ,  3.5 , 3.5  , 9.0   , 9.0  , kInf  );
    SetCut(  4, 0.41 ,  0.  ,  1.9 , 1.9  , 4.0   , 4.0  , kInf  );
//----------------------------------------------------------------
    SetCut(  5, 0.47 , 0.935, 0.139, 1.738 , 0.498  , 3.5  , kInf  );  //410-470
    SetCut(  6, 0.53 , 0.914, 0.136, 1.493 , 0.436  , 3.0  , kInf  );  //470-530
//----------------------------------------------------------------    
    SetCut(  7, 0.59 , 0.895, 0.131, 1.384 , 0.290 , 2.7  , kInf  );    //530-590
    SetCut(  8, 0.65 , 0.887, 0.121, 1.167 , 0.287 , 2.5  , kInf  );     //590-650
    SetCut(  9, 0.73 , 0.879, 0.120, 1.153 , 0.257 , 2.0  , kInf  );     //650-730
//----------------------------------------------------------------    
    SetCut( 10, 0.83 , 0.880, 0.126, 1.164 , 0.204 , 2.308 , 0.297 );       //730-830
    SetCut( 11, 0.93 , 0.918, 0.145, 1.164 , 0.204 , 2.00 , 0.168 );        //830-930
    SetCut( 12, 1.03 , 0.899, 0.128, 1.164 , 0.204  ,1.80 , 0.168);
    //------------------------ pi,K ---------------------
    fAprob[0][0]=1212;     fAprob[1][0]=33.;   // fAprob[0][i] - const for pions,cut[i+5] 
    fAprob[0][1]=1022;     fAprob[1][1]=46.2 ; // fAprob[1][i] -           kaons
    //---------------------------------------------------
    fAprob[0][2]= 889.7;   fAprob[1][2]=66.58; fAprob[2][2]=14.53;
    fAprob[0][3]= 686.;    fAprob[1][3]=88.8;  fAprob[2][3]=19.27;   
    fAprob[0][4]= 697.;    fAprob[1][4]=125.6; fAprob[2][4]=28.67;
    //------------------------ pi,K,p -------------------
    fAprob[0][5]= 633.7;   fAprob[1][5]=100.1;   fAprob[2][5]=37.99;   // fAprob[2][i] -  protons
    fAprob[0][6]= 469.5;   fAprob[1][6]=20.74;   fAprob[2][6]=25.43;
    fAprob[0][7]= 355.;    fAprob[1][7]=
                          355.*(20.74/469.5);  fAprob[2][7]=34.08;
}
//End AliITSPid.cxx
Commit	Line	Data
8057c644	1	//------------------------------------------------------------//
	2	// Class for identification of pions,kaons and protons in ITS //
	3	// Prior particles population (probabilities) are taken from //
	4	// Hijing event generator. //
	5	//------------------------------------------------------------//
	6	// #include <stdlib.h>
23efe5f1	7	#include "AliITSPid.h"
7d62fb64	8	//#include "TMath.h"
4ae5bbc4	9	#include <Riostream.h>
e75c69cd	10	#include <TClonesArray.h>
e75c69cd	11	#include <TVector.h>
5aa865de	12	#include "AliKalmanTrack.h"
8057c644	13	#include "AliITSIOTrack.h"
e75c69cd	14	#include "AliITStrackV2.h"
	15	#include <TF1.h>
	16
23efe5f1	17	ClassImp(AliITSPid)
8057c644	18
8057c644	19
8221b41b	20	AliITSPid::AliITSPid(const AliITSPid &source) : TObject(source),
	21	fMxtrs(source.fMxtrs),
	22	fTrs(source.fTrs),
	23	fWpi(source.fWpi),
	24	fWk(source.fWk),
	25	fWp(source.fWp),
	26	fRpik(source.fRpik),
	27	fRppi(source.fRppi),
	28	fRpka(source.fRpka),
	29	fRp(source.fRp),
	30	fPcode(source.fPcode),
	31	fSigmin(source.fSigmin),
	32	fSilent(source.fSilent),
	33	fCutKa(source.fCutKa),
	34	fCutPr(source.fCutPr),
	35	fggpi(source.fggpi),
	36	fggka(source.fggka),
	37	fggpr(source.fggpr){
8057c644	38	// Copy constructor. This is a function which is not allowed to be
8057c644	39
8221b41b	40
8057c644	41	}
8221b41b	42
8057c644	43	//______________________________________________________________________
8221b41b	44	AliITSPid& AliITSPid::operator=(const AliITSPid& source){
	45	// Assignment operator. This is a function which is not allowed to be
	46	this->~AliITSPid();
	47	new(this) AliITSPid(source);
8057c644	48	return *this;
	49	}
	50
	51	//
	52	Float_t AliITSPid::Qtrm(Int_t track) {
	53	//
	54	// This calculates truncated mean signal for given track.
11b3ae13	55	//
23efe5f1	56	TVector q(*( this->GetVec(track) ));
38f7c306	57	Int_t ml=(Int_t)q(0);
	58	if(ml<1)return 0.;
	59	if(ml>6)ml=6;
	60	float vf[6];
	61	Int_t nl=0; for(Int_t i=0;i<ml;i++){if(q(i)>fSigmin){vf[nl]=q(i+1);nl++;}}
	62	if(nl==0)return 0.;
23efe5f1	63	switch(nl){
	64	case 1:q(6)=q(1); break;
	65	case 2:q(6)=(q(1)+q(2))/2.; break;
	66	default:
38f7c306	67	for(int fi=0;fi<2;fi++){
	68	Int_t swap;
	69	do{ swap=0; float qmin=vf[fi];
	70	for(int j=fi+1;j<nl;j++)
	71	if(qmin>vf[j]){qmin=vf[j];vf[j]=vf[fi];vf[fi]=qmin;swap=1;};
	72	}while(swap==1);
	73	}
23efe5f1	74	q(6)= (vf[0]+vf[1])/2.;
	75	break;
	76	}
	77	for(Int_t i=0;i<nl;i++){q(i+1)=vf[i];} this->SetVec(track,q);
	78	return (q(6));
	79	}
38f7c306	80
8057c644	81	Float_t AliITSPid::Qtrm(Float_t qarr[6],Int_t narr) const{
	82	//
	83	//This calculates truncated mean signal for given signal set.
	84	//
38f7c306	85	Float_t q[6],qm,qmin;
79b921e1	86	Int_t nl,ml;
38f7c306	87	if(narr>0&&narr<7){ml=narr;}else{return 0;};
38f7c306	88	nl=0; for(Int_t i=0;i<ml;i++){if(qarr[i]>fSigmin){q[nl]=qarr[i];nl++;}}
79b921e1	89	if(nl==0)return 0.;
79b921e1	90	switch(nl){
	91	case 1:qm=q[0]; break;
	92	case 2:qm=(q[0]+q[1])/2.; break;
	93	default:
38f7c306	94	Int_t swap;
	95	for(int fi=0;fi<2;fi++){
	96	do{ swap=0; qmin=q[fi];
	97	for(int j=fi+1;j<nl;j++)
	98	if(qmin>q[j]){qmin=q[j];q[j]=q[fi];q[fi]=qmin;swap=1;};
	99	}while(swap==1);
	100	}
	101	qm= (q[0]+q[1])/2.;
79b921e1	102	break;
	103	}
	104	return qm;
	105	}
11b3ae13	106
8057c644	107	Int_t AliITSPid::Wpik(Float_t pm,Float_t q){
	108	//Calcutates probabilityes of pions and kaons
	109	//Returns particle code for dominant probability.
11b3ae13	110	Double_t par[6];
	111	for(int i=0;i<6;i++){par[i]=fGGpi[i]->Eval(pm);}
	112	fggpi->SetParameters(par);
23efe5f1	113
11b3ae13	114	for(int i=0;i<3;i++){par[i]=fGGka[i]->Eval(pm);}
11b3ae13	115	fggka->SetParameters(par);
38f7c306	116
11b3ae13	117	Float_t ppi=fggpi->Eval(q);
	118	Float_t pka=fggka->Eval(q);
	119	Float_t p=ppi+pka;
	120	/*
	121	if(!fSilent){
	122	fggka->Print();
	123	fggpi->Print();
	124	if(p>0)cout<<" ppi,pka="<<ppi/p<<" "<<pka/p<<endl;
	125	}
	126	*/
38f7c306	127
11b3ae13	128	if(p>0){
	129	ppi=ppi/p;
	130	pka=pka/p;
	131	fWp=0.; fWpi=ppi; fWk=pka;
	132	if( pka>ppi){return fPcode=321;}else{return fPcode=211;}
	133	}else{return 0;}
23efe5f1	134	}
23efe5f1	135	//-----------------------------------------------------------
8057c644	136	Int_t AliITSPid::Wpikp(Float_t pm,Float_t q){
	137	//
	138	//Calcutates probabilityes of pions,kaons and protons.
	139	//Returns particle code for dominant probability.
11b3ae13	140	Double_t par[6];
	141	for(int i=0;i<6;i++){par[i]=fGGpi[i]->Eval(pm);}
	142	fggpi->SetParameters(par);
23efe5f1	143
11b3ae13	144	for(int i=0;i<3;i++){par[i]=fGGka[i]->Eval(pm);}
11b3ae13	145	fggka->SetParameters(par);
23efe5f1	146
11b3ae13	147	for(int i=0;i<3;i++){par[i]=fGGpr[i]->Eval(pm);}
11b3ae13	148	fggpr->SetParameters(par);
38f7c306	149
11b3ae13	150	Float_t p,ppi,pka,ppr;
	151	if( q>(fggpr->GetParameter(1)+fggpr->GetParameter(2)) )
	152	{ p=1.0; ppr=1.0; ppi=pka=0.0;
	153	}else{
	154	ppi=fggpi->Eval(q);
	155	pka=fggka->Eval(q);
	156	ppr=fggpr->Eval(q);
	157	p=ppi+pka+ppr;
	158	}
	159	if(p>0){
	160	ppi=ppi/p;
	161	pka=pka/p;
	162	ppr=ppr/p;
	163	fWp=ppr; fWpi=ppi; fWk=pka;
	164	//if(!fSilent)cout<<" ppi,pka,ppr="<<ppi<<" "<<pka<<" "<<ppr<<endl;
38f7c306	165
11b3ae13	166	if( ppi>pka&&ppi>ppr )
	167	{return fPcode=211;}
	168	else{ if(pka>ppr){return fPcode=321;}else{return fPcode=2212;}
	169	}
38f7c306	170
11b3ae13	171	}else{return 0;}
23efe5f1	172	}
	173	//-----------------------------------------------------------
	174	Int_t AliITSPid::GetPcode(TClonesArray* rps,Float_t pm)
	175	{
8057c644	176	//Returns particle code
088e0b8d	177	Info("GetPcode","method not implemented - Inputs TClonesArray *%x , Float_t %f",rps,pm);
23efe5f1	178	return 0;
	179	}
	180	//-----------------------------------------------------------
5aa865de	181	Int_t AliITSPid::GetPcode(AliKalmanTrack *track)
79b921e1	182	{
8057c644	183	//Returns particle code for given track.
79b921e1	184	Double_t xk,par[5]; track->GetExternalParameters(xk,par);
	185	Float_t phi=TMath::ASin(par[2]) + track->GetAlpha();
	186	if (phi<-TMath::Pi()) phi+=2*TMath::Pi();
	187	if (phi>=TMath::Pi()) phi-=2*TMath::Pi();
	188	Float_t lam=TMath::ATan(par[3]);
e75c69cd	189	Float_t pt1=TMath::Abs(par[4]);
e75c69cd	190	Float_t mom=1./(pt1*TMath::Cos(lam));
bbc6cd2c	191	Float_t dedx=track->GetPIDsignal();
79b921e1	192	Int_t pcode=GetPcode(dedx/40.,mom);
def162f4	193	// cout<<"TPCtrack dedx,mom,pcode="<<dedx<<","<<mom<<","<<pcode<<endl;
79b921e1	194	return pcode?pcode:211;
	195	}
	196	//------------------------------------------------------------
79b921e1	197	Int_t AliITSPid::GetPcode(AliITSIOTrack *track)
79b921e1	198	{
8057c644	199	//Returns particle code for given track(V1).
79b921e1	200	Double_t px,py,pz;
	201	px=track->GetPx();
	202	py=track->GetPy();
	203	pz=track->GetPz();
	204	Float_t mom=TMath::Sqrt(pxpx+pypy+pz*pz);
	205	//???????????????????
def162f4	206	// Float_t dedx=1.0;
def162f4	207	Float_t dedx=track->GetdEdx();
79b921e1	208	//???????????????????
79b921e1	209	Int_t pcode=GetPcode(dedx,mom);
def162f4	210	// cout<<"ITSV1 dedx,mom,pcode="<<dedx<<","<<mom<<","<<pcode<<endl;
79b921e1	211	return pcode?pcode:211;
79b921e1	212	}
79b921e1	213	//-----------------------------------------------------------
	214	Int_t AliITSPid::GetPcode(AliITStrackV2 *track)
	215	{
8057c644	216	//Returns particle code for given track(V2).
79b921e1	217	if(track==0)return 0;
38f7c306	218	// track->Propagate(track->GetAlpha(),3.,0.1/65.191.848,0.11.848);
79b921e1	219	track->PropagateTo(3.,0.0028,65.19);
6c94f330	220	//track->PropagateToVertex(); Not needed. (I.B.)
79b921e1	221	Double_t xk,par[5]; track->GetExternalParameters(xk,par);
79b921e1	222	Float_t lam=TMath::ATan(par[3]);
e75c69cd	223	Float_t pt1=TMath::Abs(par[4]);
79b921e1	224	Float_t mom=0.;
e75c69cd	225	if( (pt1TMath::Cos(lam))!=0. ){ mom=1./(pt1TMath::Cos(lam)); }else{mom=0.;};
79b921e1	226	Float_t dedx=track->GetdEdx();
e75c69cd	227	// cout<<"lam,pt1,mom,dedx="<<lam<<","<<pt1<<","<<mom<<","<<dedx<<endl;
79b921e1	228	Int_t pcode=GetPcode(dedx,mom);
def162f4	229	// cout<<"ITS V2 dedx,mom,pcode="<<dedx<<","<<mom<<","<<pcode<<endl;
79b921e1	230	return pcode?pcode:211;
	231	}
	232	//-----------------------------------------------------------
23efe5f1	233	Int_t AliITSPid::GetPcode(Float_t q,Float_t pm)
23efe5f1	234	{
8057c644	235	//Returns particle code for given signal and momentum.
38f7c306	236	fWpi=fWk=fWp=0.; fPcode=0;
23efe5f1	237
11b3ae13	238	if ( pm<=0.400 )
11b3ae13	239	{ if( q<fCutKa->Eval(pm) )
e75c69cd	240	{return Pion();}
11b3ae13	241	else{ if( q<fCutPr->Eval(pm) )
e75c69cd	242	{return Kaon();}
e75c69cd	243	else{return Proton();}
11b3ae13	244	}
11b3ae13	245	}
5af4a2d0	246	if ( pm<=0.750 ){
	247	if ( q>fCutPr->Eval(pm) ){
	248	return Proton();
	249	}
	250	else {
	251	return Wpik(pm,q);
	252	}
	253	}
e75c69cd	254	if( pm<=1.10 ){ return Wpikp(pm,q); }
23efe5f1	255	return fPcode;
	256	}
	257	//-----------------------------------------------------------
	258	void AliITSPid::SetCut(Int_t n,Float_t pm,Float_t pilo,Float_t pihi,
	259	Float_t klo,Float_t khi,Float_t plo,Float_t phi)
	260	{
8057c644	261	// The cut-table initializer method.
e75c69cd	262	fCut[n][0]=pm;
	263	fCut[n][1]=pilo;
	264	fCut[n][2]=pihi;
	265	fCut[n][3]=klo;
	266	fCut[n][4]=khi;
	267	fCut[n][5]=plo;
	268	fCut[n][6]=phi;
23efe5f1	269	return ;
23efe5f1	270	}
38f7c306	271	//------------------------------------------------------------
e8d02863	272	void AliITSPid::SetVec(Int_t ntrack,const TVector& info) const
23efe5f1	273	{
8057c644	274	//Store track info in tracls table
e75c69cd	275	TClonesArray& arr=*fTrs;
23efe5f1	276	new( arr[ntrack] ) TVector(info);
	277	}
	278	//-----------------------------------------------------------
8057c644	279	TVector* AliITSPid::GetVec(Int_t ntrack) const
23efe5f1	280	{
8057c644	281	//Get given track from track table
e75c69cd	282	TClonesArray& arr=*fTrs;
23efe5f1	283	return (TVector*)arr[ntrack];
	284	}
	285	//-----------------------------------------------------------
	286	void AliITSPid::SetEdep(Int_t track,Float_t Edep)
	287	{
8057c644	288	//Set dEdx for given track
23efe5f1	289	TVector xx(0,11);
e75c69cd	290	if( ((TVector*)fTrs->At(track))->IsValid() )
e75c69cd	291	{TVector yy( ((TVector)fTrs->At(track)) );xx=yy; }
23efe5f1	292	Int_t j=(Int_t)xx(0); if(j>4)return;
23efe5f1	293	xx(++j)=Edep;xx(0)=j;
e75c69cd	294	TClonesArray &arr=*fTrs;
23efe5f1	295	new(arr[track])TVector(xx);
	296	}
	297	//-----------------------------------------------------------
	298	void AliITSPid::SetPmom(Int_t track,Float_t Pmom)
	299	{
8057c644	300	//Set momentum for given track
23efe5f1	301	TVector xx(0,11);
e75c69cd	302	if( ((TVector*)fTrs->At(track))->IsValid() )
e75c69cd	303	{TVector yy( ((TVector)fTrs->At(track)) );xx=yy; }
23efe5f1	304	xx(10)=Pmom;
e75c69cd	305	TClonesArray &arr=*fTrs;
23efe5f1	306	new(arr[track])TVector(xx);
	307	}
	308	//-----------------------------------------------------------
	309	void AliITSPid::SetPcod(Int_t track,Int_t partcode)
	310	{
8057c644	311	//Set particle code for given track
23efe5f1	312	TVector xx(0,11);
e75c69cd	313	if( ((TVector*)fTrs->At(track))->IsValid() )
e75c69cd	314	{TVector yy( ((TVector)fTrs->At(track)) );xx=yy; }
23efe5f1	315	if(xx(11)==0)
e75c69cd	316	{xx(11)=partcode; fMxtrs++;
e75c69cd	317	TClonesArray &arr=*fTrs;
23efe5f1	318	new(arr[track])TVector(xx);
	319	}
	320	}
	321	//-----------------------------------------------------------
	322	void AliITSPid::Print(Int_t track)
8057c644	323	{
	324	//Prints information for given track
	325	cout<<fMxtrs<<" tracks in AliITSPid obj."<<endl;
e75c69cd	326	if( ((TVector*)fTrs->At(track))->IsValid() )
e75c69cd	327	{TVector xx( ((TVector)fTrs->At(track)) );
23efe5f1	328	xx.Print();
	329	}
	330	else
	331	{cout<<"No data for track "<<track<<endl;return;}
	332	}
	333	//-----------------------------------------------------------
	334	void AliITSPid::Tab(void)
	335	{
8057c644	336	//Make PID for tracks stored in tracks table
e8f4b841	337	if(fTrs->GetEntries()==0){cout<<"No entries in TAB"<<endl;return;}
	338	cout<<"------------------------------------------------------------------------"<<endl;
	339	cout<<"Nq"<<" q1 "<<" q2 "<<" q3 "<<" q4 "<<" q5 "<<
	340	" Qtrm " <<" Wpi "<<" Wk "<<" Wp "<<"Pmom "<<endl;
	341	cout<<"------------------------------------------------------------------------"<<endl;
	342	for(Int_t i=0;i<fTrs->GetEntries();i++)
	343	{
	344	TVector xxx( ((TVector)fTrs->At(i)) );
	345	if( xxx.IsValid() && xxx(0)>0 )
23efe5f1	346	{
e8f4b841	347	TVector xx( ((TVector)fTrs->At(i)) );
	348	if(xx(0)>=2)
	349	{
	350	// 1)Calculate Qtrm
	351	xx(6)=(this->Qtrm(i));
23efe5f1	352
e8f4b841	353	}else{
	354	xx(6)=xx(1);
	355	}
	356	// 2)Calculate Wpi,Wk,Wp
	357	this->GetPcode(xx(6),xx(10)/1000.);
	358	xx(7)=GetWpi();
	359	xx(8)=GetWk();
	360	xx(9)=GetWp();
	361	// 3)Print table
	362	if(xx(0)>0){
	363	// cout<<xx(0)<<" ";
	364	for(Int_t j=1;j<11;j++){
	365	if(i<7){ cout.width(7);cout.precision(4);cout<<xx(j);}
	366	if(i>7){ cout.width(7);cout.precision(5);cout<<xx(j);}
	367	}
	368	cout<<endl;
	369	}
	370	// 4)Update data in TVector
	371	TClonesArray &arr=*fTrs;
	372	new(arr[i])TVector(xx);
23efe5f1	373	}
e8f4b841	374	else
	375	{/cout<<"No data for track "<<i<<endl;/}
	376	}// End loop for tracks
23efe5f1	377	}
	378	void AliITSPid::Reset(void)
	379	{
8057c644	380	//Reset tracks table
e75c69cd	381	for(Int_t i=0;i<fTrs->GetEntries();i++){
23efe5f1	382	TVector xx(0,11);
e75c69cd	383	TClonesArray &arr=*fTrs;
23efe5f1	384	new(arr[i])TVector(xx);
	385	}
	386	}
	387	//-----------------------------------------------------------
8221b41b	388	AliITSPid::AliITSPid(Int_t ntrack):
	389	fMxtrs(0),
	390	fTrs(0),
	391	fWpi(0),
	392	fWk(0),
	393	fWp(0),
	394	fRpik(0),
	395	fRppi(0),
	396	fRpka(0),
	397	fRp(0),
	398	fPcode(0),
	399	fSigmin(0.01),
	400	fSilent(0),
	401	fCutKa(0),
	402	fCutPr(0),
	403	fggpi(0),
	404	fggka(0),
	405	fggpr(0){
8057c644	406	//Constructor for AliITSPid class
8221b41b	407
e75c69cd	408	fTrs = new TClonesArray("TVector",ntrack);
e75c69cd	409	TClonesArray &arr=*fTrs;
23efe5f1	410	for(Int_t i=0;i<ntrack;i++)new(arr[i])TVector(0,11);
8221b41b	411	//
11b3ae13	412	fCutKa=new TF1("fcutka","pol4",0.05,0.4);
	413	Double_t ka[5]={25.616, -161.59, 408.97, -462.17, 192.86};
	414	fCutKa->SetParameters(ka);
	415	//
	416	fCutPr=new TF1("fcutpr","[0]/x/x+[1]",0.05,1.1);
	417	Double_t pr[2]={0.70675,0.4455};
	418	fCutPr->SetParameters(pr);
	419	//
	420	//---------- signal fit ----------
	421	{//Pions
	422	fGGpi[0]=new TF1("fp1pi","pol4",0.34,1.2);
e75c69cd	423	Double_t parpi0[10]={ -1.9096471071e+03, 4.5354331545e+04, -1.1860738840e+05,
11b3ae13	424	1.1405329025e+05, -3.8289694496e+04 };
e75c69cd	425	fGGpi[0]->SetParameters(parpi0);
11b3ae13	426	fGGpi[1]=new TF1("fp2pi","[0]/x/x+[1]",0.34,1.2);
e75c69cd	427	Double_t parpi1[10]={ 1.0791668283e-02, 9.7347716496e-01 };
e75c69cd	428	fGGpi[1]->SetParameters(parpi1);
11b3ae13	429	fGGpi[2]=new TF1("fp3pi","[0]/x/x+[1]",0.34,1.2);
e75c69cd	430	Double_t parpi2[10]={ 5.8191602279e-04, 9.7285601334e-02 };
e75c69cd	431	fGGpi[2]->SetParameters(parpi2);
11b3ae13	432	fGGpi[3]=new TF1("fp4pi","pol4",0.34,1.2);
e75c69cd	433	Double_t parpi3[10]={ 6.6267353195e+02, 7.1595101104e+02, -5.3095111914e+03,
11b3ae13	434	6.2900977606e+03, -2.2935862292e+03 };
e75c69cd	435	fGGpi[3]->SetParameters(parpi3);
11b3ae13	436	fGGpi[4]=new TF1("fp5pi","[0]/x/x+[1]",0.34,1.2);
e75c69cd	437	Double_t parpi4[10]={ 9.0419011783e-03, 1.1628922525e+00 };
e75c69cd	438	fGGpi[4]->SetParameters(parpi4);
11b3ae13	439	fGGpi[5]=new TF1("fp6pi","[0]/x/x+[1]",0.34,1.2);
e75c69cd	440	Double_t parpi5[10]={ 1.8324872519e-03, 2.1503968838e-01 };
e75c69cd	441	fGGpi[5]->SetParameters(parpi5);
11b3ae13	442	}//End Pions
	443	{//Kaons
	444	fGGka[0]=new TF1("fp1ka","pol4",0.24,1.2);
e75c69cd	445	Double_t parka0[20]={
11b3ae13	446	-1.1204243395e+02,4.6716191428e+01,2.2584059281e+03,
11b3ae13	447	-3.7123338009e+03,1.6003647641e+03 };
e75c69cd	448	fGGka[0]->SetParameters(parka0);
11b3ae13	449	fGGka[1]=new TF1("fp2ka","[0]/x/x+[1]",0.24,1.2);
e75c69cd	450	Double_t parka1[20]={
11b3ae13	451	2.5181172905e-01,8.7566001814e-01 };
e75c69cd	452	fGGka[1]->SetParameters(parka1);
11b3ae13	453	fGGka[2]=new TF1("fp3ka","pol6",0.24,1.2);
e75c69cd	454	Double_t parka2[20]={
11b3ae13	455	8.6236021573e+00,-7.0970427531e+01,2.4846827669e+02,
	456	-4.6094401290e+02,4.7546751408e+02,-2.5807112462e+02,
	457	5.7545491696e+01 };
e75c69cd	458	fGGka[2]->SetParameters(parka2);
11b3ae13	459	}//End Kaons
	460	{//Protons
	461	fGGpr[0]=new TF1("fp1pr","pol4",0.4,1.2);
e75c69cd	462	Double_t parpr0[10]={
11b3ae13	463	6.0150106543e+01,-8.8176206410e+02,3.1222644604e+03,
11b3ae13	464	-3.5269200901e+03,1.2859128345e+03 };
e75c69cd	465	fGGpr[0]->SetParameters(parpr0);
11b3ae13	466	fGGpr[1]=new TF1("fp2pr","[0]/x/x+[1]",0.4,1.2);
e75c69cd	467	Double_t parpr1[10]={
11b3ae13	468	9.4970837607e-01,7.3573504201e-01 };
e75c69cd	469	fGGpr[1]->SetParameters(parpr1);
11b3ae13	470	fGGpr[2]=new TF1("fp3pr","[0]/x/x+[1]",0.4,1.2);
e75c69cd	471	Double_t parpr2[10]={
11b3ae13	472	1.2498403757e-01,2.7845072306e-02 };
e75c69cd	473	fGGpr[2]->SetParameters(parpr2);
11b3ae13	474	}//End Protons
	475	//----------- end fit -----------
	476
	477	fggpr=new TF1("ggpr","gaus",0.4,1.2);
	478	fggpi=new TF1("ggpi","gaus+gaus(3)",0.4,1.2);
	479	fggka=new TF1("ggka","gaus",0.4,1.2);
	480
	481	//-------------------------------------------------
e75c69cd	482	const int kInf=10;
23efe5f1	483	// Ncut Pmom pilo pihi klo khi plo phi
	484	// cut[j] [0] [1] [2] [3] [4] [5] [6]
	485	//----------------------------------------------------------------
e75c69cd	486	SetCut( 1, 0.12 , 0. , 0. , kInf , kInf , kInf , kInf );
	487	SetCut( 2, 0.20 , 0. , 6.0 , 6.0 , kInf , kInf , kInf );
	488	SetCut( 3, 0.30 , 0. , 3.5 , 3.5 , 9.0 , 9.0 , kInf );
	489	SetCut( 4, 0.41 , 0. , 1.9 , 1.9 , 4.0 , 4.0 , kInf );
23efe5f1	490	//----------------------------------------------------------------
e75c69cd	491	SetCut( 5, 0.47 , 0.935, 0.139, 1.738 , 0.498 , 3.5 , kInf ); //410-470
e75c69cd	492	SetCut( 6, 0.53 , 0.914, 0.136, 1.493 , 0.436 , 3.0 , kInf ); //470-530
23efe5f1	493	//----------------------------------------------------------------
e75c69cd	494	SetCut( 7, 0.59 , 0.895, 0.131, 1.384 , 0.290 , 2.7 , kInf ); //530-590
	495	SetCut( 8, 0.65 , 0.887, 0.121, 1.167 , 0.287 , 2.5 , kInf ); //590-650
	496	SetCut( 9, 0.73 , 0.879, 0.120, 1.153 , 0.257 , 2.0 , kInf ); //650-730
23efe5f1	497	//----------------------------------------------------------------
38f7c306	498	SetCut( 10, 0.83 , 0.880, 0.126, 1.164 , 0.204 , 2.308 , 0.297 ); //730-830
	499	SetCut( 11, 0.93 , 0.918, 0.145, 1.164 , 0.204 , 2.00 , 0.168 ); //830-930
	500	SetCut( 12, 1.03 , 0.899, 0.128, 1.164 , 0.204 ,1.80 , 0.168);
	501	//------------------------ pi,K ---------------------
e75c69cd	502	fAprob[0][0]=1212; fAprob[1][0]=33.; // fAprob[0][i] - const for pions,cut[i+5]
e75c69cd	503	fAprob[0][1]=1022; fAprob[1][1]=46.2 ; // fAprob[1][i] - kaons
38f7c306	504	//---------------------------------------------------
e75c69cd	505	fAprob[0][2]= 889.7; fAprob[1][2]=66.58; fAprob[2][2]=14.53;
	506	fAprob[0][3]= 686.; fAprob[1][3]=88.8; fAprob[2][3]=19.27;
	507	fAprob[0][4]= 697.; fAprob[1][4]=125.6; fAprob[2][4]=28.67;
38f7c306	508	//------------------------ pi,K,p -------------------
e75c69cd	509	fAprob[0][5]= 633.7; fAprob[1][5]=100.1; fAprob[2][5]=37.99; // fAprob[2][i] - protons
	510	fAprob[0][6]= 469.5; fAprob[1][6]=20.74; fAprob[2][6]=25.43;
	511	fAprob[0][7]= 355.; fAprob[1][7]=
	512	355.*(20.74/469.5); fAprob[2][7]=34.08;
23efe5f1	513	}
11b3ae13	514	//End AliITSPid.cxx