+//------------------------------------------------------------//
+// 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 <iostream.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)
-//-----------------------------------------------------------
-Float_t AliITSPid::qcorr(Float_t xc)
-{
- assert(0);
- Float_t fcorr;
- fcorr=( 0.766 +0.9692*xc -1.267*xc*xc )*( 1.-TMath::Exp(-xc*64.75) );
- if(fcorr<=0.1)fcorr=0.1;
-return qtot/fcorr;
+
+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;
}
-//-----------------------------------------------------------
-#include "vector.h"
-#include <algorithm>
-Float_t AliITSPid::qtrm(Int_t track)
-{
+//
+ Float_t AliITSPid::Qtrm(Int_t track) {
+//
+// This calculates truncated mean signal for given track.
+//
TVector q(*( this->GetVec(track) ));
- Int_t nl=(Int_t)q(0);
- if(nl<1)return 0.;
- vector<float> vf(nl);
+ 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_t i=0;i<nl;i++){vf[i]=q(i+1);}
- sort(vf.begin(),vf.end());
+ 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));
}
-//-----------------------------------------------------------
-inline
-Int_t AliITSPid::wpik(Int_t nc,Float_t q)
-{
- return pion();
- Float_t qmpi,qmk,sigpi,sigk,dpi,dk,ppi,pk;
- qmpi =cut[nc][1];
- sigpi=cut[nc][2];
- qmk =cut[nc][3];
- sigk =cut[nc][4];
- Float_t dqpi=(q-qmpi)/sigpi;
- Float_t dqk =(q-qmk )/sigk;
- if( dqk<-1. )return pion();
- dpi =fabs(dqpi);
- dk =fabs(dqk);
- ppi =1.- TMath::Erf(dpi); // +0.5;
- pk =1.- TMath::Erf(dk); // +0.5;
- Float_t rpik=ppi/(pk+0.0000001);
- cout<<"q,dqpi,dqk, wpik: ppi,pk,rpik="
- <<q<<" "<<dqpi<<" "<<dqk<<" "<<ppi<<" "<<pk<<" "<<rpik<<endl;
- if( rpik>=1000. ){return pion();}
- if( rpik>0.001 ){fWk=1./(rpik+1.); fWpi=1.-fWk;return 0;}
- return pion();
+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;}
}
//-----------------------------------------------------------
-inline
-Int_t AliITSPid::wpikp(Int_t nc,Float_t q)
-{
- return pion();
+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);
- Float_t qmpi,qmk,qmp,sigpi,sigk,sigp,ppi,pka,ppr,rppi,rpka;
+ for(int i=0;i<3;i++){par[i]=fGGka[i]->Eval(pm);}
+ fggka->SetParameters(par);
- qmpi =cut[nc][1];
- sigpi=cut[nc][2];
- qmk =cut[nc][3];
- sigk =cut[nc][4];
- qmp =cut[nc][5];
- sigp =cut[nc][6];
+ for(int i=0;i<3;i++){par[i]=fGGpr[i]->Eval(pm);}
+ fggpr->SetParameters(par);
- ppi =TMath::Erf( TMath::Abs((q-qmpi)/sigpi) )+0.5;
- pka =TMath::Erf( TMath::Abs((q-qmk )/sigk) )+0.5;
- ppr =TMath::Erf( TMath::Abs((q-qmp )/sigp) )+0.5;
+ 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;
- rppi=ppr/ppi;
- rpka=ppr/pka;
- Float_t rp;
- if( rppi>rpka ) { rp=rpka; }
- else{ rp=rppi; }
- if( rp>=1000. ){ return proton(); }
- if( rp>0.001 ){ fWp=rp/(rp+1.);return proton(); }
-
- return 0;
+ 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)
{
- fWpi=fWk=fWp=-1.; fPcode=0;
-//1)
- if ( pm<=cut[1][0] )
- { return pion(); }
-//2)
- if ( pm<=cut[2][0] )
- { if( q<cut[2][2] ){ return pion(); } else { return kaon();} }
-//3)
- if ( pm<=cut[3][0] )
- if( q<cut[3][2])
- { return pion(); }
- else
- { if ( q<=cut[3][5] ) {return kaon();} else {return proton();}}
-//4)
- if ( pm<=cut[4][0] )
- if( q<cut[4][2] )
- { return pion(); }
- else
- { if( q<=cut[4][4] ) {return kaon();} else {return proton(); }}
-//5)
- if ( pm<=cut[5][0] )
- if ( q>cut[5][5] ) {return proton();} else {return wpik(5,q);};
-//6)
- if ( pm<=cut[6][0] )
- if ( q>cut[6][5] ) {return proton();} else {return wpik(6,q);};
-//7)
- if ( pm<=cut[7][0] )
- if ( q<=cut[7][5] ) {return fPcode;} else {return proton();};
-//8)
- if ( pm<=cut[8][0] )
- if ( q<=cut[8][5] ) {return fPcode;} else {return proton();};
-//9)
- if ( pm<=cut[9][0] )
- if ( q<=cut[9][5] ) {return fPcode;} else {return proton();};
-//10)
- if( pm<=cut[10][0] ){ return wpikp(10,q); }
-//11)
- if( pm<=cut[11][0] ){ return wpikp(11,q); }
-//12)
- if( pm<=cut[12][0] ){ return wpikp(12,q); }
+ //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)
{
- cut[n][0]=pm;
- cut[n][1]=pilo;
- cut[n][2]=pihi;
- cut[n][3]=klo;
- cut[n][4]=khi;
- cut[n][5]=plo;
- cut[n][6]=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,TVector info)
+//------------------------------------------------------------
+void AliITSPid::SetVec(Int_t ntrack,const TVector& info) const
{
-TClonesArray& arr=*trs;
+ //Store track info in tracls table
+TClonesArray& arr=*fTrs;
new( arr[ntrack] ) TVector(info);
}
//-----------------------------------------------------------
-TVector* AliITSPid::GetVec(Int_t ntrack)
+TVector* AliITSPid::GetVec(Int_t ntrack) const
{
-TClonesArray& arr=*trs;
+ //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*)trs->At(track))->IsValid() )
- {TVector yy( *((TVector*)trs->At(track)) );xx=yy; }
+ 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=*trs;
+ 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*)trs->At(track))->IsValid() )
- {TVector yy( *((TVector*)trs->At(track)) );xx=yy; }
+ if( ((TVector*)fTrs->At(track))->IsValid() )
+ {TVector yy( *((TVector*)fTrs->At(track)) );xx=yy; }
xx(10)=Pmom;
- TClonesArray &arr=*trs;
+ 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*)trs->At(track))->IsValid() )
- {TVector yy( *((TVector*)trs->At(track)) );xx=yy; }
+ if( ((TVector*)fTrs->At(track))->IsValid() )
+ {TVector yy( *((TVector*)fTrs->At(track)) );xx=yy; }
if(xx(11)==0)
- {xx(11)=partcode; mxtrs++;
- TClonesArray &arr=*trs;
+ {xx(11)=partcode; fMxtrs++;
+ TClonesArray &arr=*fTrs;
new(arr[track])TVector(xx);
}
}
//-----------------------------------------------------------
void AliITSPid::Print(Int_t track)
-{cout<<mxtrs<<" tracks in AliITSPid obj."<<endl;
- if( ((TVector*)trs->At(track))->IsValid() )
- {TVector xx( *((TVector*)trs->At(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
//-----------------------------------------------------------
void AliITSPid::Tab(void)
{
-if(trs->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<trs->GetEntries();i++)
-{
- TVector xx( *((TVector*)trs->At(i)) );
- if( xx.IsValid() && xx(0)>0 )
+ //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*)trs->At(i)) );
- if(xx(0)>=2)
- {
-// 1)Calculate Qtrm
- xx(6)=(this->qtrm(i));
+ 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++){
- cout.width(7);cout.precision(5);cout<<xx(j);
- }
- cout<<endl;
- }
-// 4)Update data in TVector
- TClonesArray &arr=*trs;
- new(arr[i])TVector(xx);
+ }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
+ else
+ {/*cout<<"No data for track "<<i<<endl;*/}
+ }// End loop for tracks
}
void AliITSPid::Reset(void)
{
- for(Int_t i=0;i<trs->GetEntries();i++){
+ //Reset tracks table
+ for(Int_t i=0;i<fTrs->GetEntries();i++){
TVector xx(0,11);
- TClonesArray &arr=*trs;
+ TClonesArray &arr=*fTrs;
new(arr[i])TVector(xx);
}
}
//-----------------------------------------------------------
-AliITSPid::AliITSPid(Int_t ntrack)
-{
- trs = new TClonesArray("TVector",ntrack);
- TClonesArray &arr=*trs;
+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);
- mxtrs=0;
-const int inf=10;
+ //
+ 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. , inf , inf , inf , inf );
- SetCut( 2, 0.20 , 0. , 6.0 , 6.0 , inf , inf , inf );
- SetCut( 3, 0.30 , 0. , 3.5 , 3.5 , 9.0 , 9.0 , inf );
- SetCut( 4, 0.41 , 0. , 1.9 , 1.9 , 4.0 , 4.0 , inf );
+ 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 , 1. , 0.12, 1.98 , 0.17 , 3.5 , inf );
- SetCut( 6, 0.53 , 1. , 0.12, 1.75 , 0.16 , 3.0 , inf );
+ 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. , 0. , 1.18 , 0.125 , 2.7 , inf );
- SetCut( 8, 0.65 , 0. , 0. , 1.18 , 0.125 , 2.5 , inf );
- SetCut( 9, 0.73 , 0. , 0. , 0. , 0.125 , 2.0 , inf );
+ 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. , 0. , 1.25 , 0.13 , 2.14 , 0.20 );
- SetCut( 11, 0.93 , 0. , 0. , 1.18 , 0.125 , 1.88 , 0.18 );
- SetCut( 12, 1.03 , 0. , 0. , 1.13 , 0.12 , 1.68 , 0.155);
+ 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