[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){
    // Copy constructor. This is a function which is not allowed to be
    // done: protected It exits with an error.
    // Inputs:
    //      AliITSPid &source  An AliITSPid class.
    // Outputs:
    //      none.
    // Return:
    //      none.

  Error("AliITSPid","You are not allowed to make a copy of the AliITSPid\n");
}
//______________________________________________________________________
AliITSPid& AliITSPid::operator=(const AliITSPid& /* source */){
    // Assignment operator. This is a function which is not allowed to be
    // done to the ITS. It exits with an error.
    // Inputs:
    //      AliITSPid &source  An AliITSPid class.
    // Outputs:
    //      none.
    // Return:
    //      none.

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