]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - ITS/AliITSPid.cxx
.so cleanup: removed from gSystem->Load()
[u/mrichter/AliRoot.git] / ITS / AliITSPid.cxx
index aa7c31d204ea9eeffa5498d5e20445fcdb5791e2..8343dbe73934f02650d375c56dab01cec36f0774 100644 (file)
@@ -1,16 +1,58 @@
+//------------------------------------------------------------//
+// 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 "TMath.h"
 #include <Riostream.h>
 #include <TClonesArray.h>
 #include <TVector.h>
-#include "AliTPCtrack.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.
 //
-Float_t AliITSPid::Qtrm(Int_t track)
-{
     TVector q(*( this->GetVec(track)  ));
     Int_t ml=(Int_t)q(0);
     if(ml<1)return 0.;
@@ -36,8 +78,10 @@ Float_t AliITSPid::Qtrm(Int_t track)
     return (q(6));
 }
 
-Float_t AliITSPid::Qtrm(Float_t qarr[6],Int_t narr)
-{
+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;};
@@ -60,8 +104,9 @@ Float_t AliITSPid::Qtrm(Float_t qarr[6],Int_t narr)
     return qm;
 }
 
-Int_t  AliITSPid::Wpik(Float_t pm,Float_t q)
-{
+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);
@@ -88,8 +133,10 @@ Int_t       AliITSPid::Wpik(Float_t pm,Float_t q)
   }else{return 0;}
 }
 //-----------------------------------------------------------
-Int_t  AliITSPid::Wpikp(Float_t pm,Float_t q)
-{
+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);
@@ -126,12 +173,14 @@ Int_t     AliITSPid::Wpikp(Float_t pm,Float_t q)
 //-----------------------------------------------------------
 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(AliTPCtrack *track)
+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();
@@ -139,7 +188,7 @@ Int_t   AliITSPid::GetPcode(AliTPCtrack *track)
       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->GetdEdx();
+      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;
@@ -147,6 +196,7 @@ Int_t   AliITSPid::GetPcode(AliTPCtrack *track)
 //------------------------------------------------------------
 Int_t   AliITSPid::GetPcode(AliITSIOTrack *track)
 {
+  //Returns particle code for given track(V1).
     Double_t px,py,pz;
     px=track->GetPx();
     py=track->GetPy();
@@ -163,10 +213,11 @@ 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();
+      //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]);
@@ -181,6 +232,7 @@ 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 )
@@ -191,9 +243,14 @@ Int_t      AliITSPid::GetPcode(Float_t q,Float_t pm)
            else{return Proton();}
            } 
        }
-    if ( pm<=0.750 )
-       if ( q>fCutPr->Eval(pm)  )
-           {return Proton();} else {return Wpik(pm,q);};
+    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;    
 }
@@ -201,6 +258,7 @@ Int_t       AliITSPid::GetPcode(Float_t q,Float_t pm)
 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;
@@ -211,20 +269,23 @@ void      AliITSPid::SetCut(Int_t n,Float_t pm,Float_t pilo,Float_t pihi,
     return ;    
 }
 //------------------------------------------------------------
-void AliITSPid::SetVec(Int_t ntrack,TVector info)
+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)
+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; }
@@ -236,6 +297,7 @@ void AliITSPid::SetEdep(Int_t track,Float_t Edep)
 //-----------------------------------------------------------
 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; }
@@ -246,6 +308,7 @@ void AliITSPid::SetPmom(Int_t track,Float_t Pmom)
 //-----------------------------------------------------------
 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; }
@@ -257,7 +320,9 @@ void AliITSPid::SetPcod(Int_t track,Int_t partcode)
 }
 //-----------------------------------------------------------
 void AliITSPid::Print(Int_t track)
-{cout<<fMxtrs<<" tracks in AliITSPid obj."<<endl;
+{
+  //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();
@@ -268,49 +333,51 @@ void AliITSPid::Print(Int_t track)
 //-----------------------------------------------------------
 void AliITSPid::Tab(void)
 {
-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 )
+  //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));
+         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{
+             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)
 {
+  //Reset tracks table
   for(Int_t i=0;i<fTrs->GetEntries();i++){
     TVector xx(0,11);
     TClonesArray &arr=*fTrs;
@@ -318,14 +385,30 @@ void AliITSPid::Reset(void)
   }
 }
 //-----------------------------------------------------------
-AliITSPid::AliITSPid(Int_t ntrack)
-{
-  fSigmin=0.01;
+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);
-    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);