ITS/AliITSPidParItem.cxx

   1 /////////////////////////////////////////////////////////
   2 //Class for PID in the ITS                             //
   3 //                                                     //
   4 //                                                     //
   5 /////////////////////////////////////////////////////////
   6
   7 #include <Riostream.h>
   8 #include <TF1.h>
   9 #include "AliITSPidParItem.h"
  10
  11 ClassImp(AliITSPidParItem)
  12 //____________________________________________________________________
  13 AliITSPidParItem::AliITSPidParItem(){
  14   // default constructor
  15   fPCenter=0;
  16   fPWidth=0;
  17   for(Int_t i=0;i<39;i++){
  18   fBuff[i]=0;
  19   }
  20 }//____________________________________________________________________
  21 AliITSPidParItem::AliITSPidParItem(Float_t center,Float_t width,Double_t *buff){
  22   // standard constructor
  23   fPCenter=center;
  24   fPWidth=width;
  25   for (Int_t i=0;i<39;i++) fBuff[i]=buff[i];
  26
  27 }
  28
  29 //____________________________________________________________________
  30 void AliITSPidParItem::GetParameters(Double_t *buff) const{
  31   //get all the parameters
  32   for (Int_t i=0;i<39;i++) buff[i]=fBuff[i];
  33
  34 }
  35 //____________________________________________________________________
  36 void AliITSPidParItem::GetProtonPar(Double_t *buffp) const{
  37   //get the protons' parameters (Width Landau, Most Probable, Area, Width Gaussian, Chi2 fit, NDF fit, Integral fit)
  38   buffp[0]=fBuff[0];
  39   buffp[1]=fBuff[1];
  40   buffp[2]=fBuff[2];
  41   buffp[3]=fBuff[3];
  42   buffp[4]=fBuff[8];
  43   buffp[5]=fBuff[9];
  44   buffp[6]=fBuff[10];
  45 }
  46 //____________________________________________________________________
  47 void AliITSPidParItem::GetKaonPar(Double_t *buffk) const{
  48   //get the kaons' parameters (Width Landau, Most Probable, Area, Width Gaussian, Chi2 fit, NDF fit, Integral fit)
  49   buffk[0]=fBuff[13];
  50   buffk[1]=fBuff[14];
  51   buffk[2]=fBuff[15];
  52   buffk[3]=fBuff[16];
  53   buffk[4]=fBuff[21];
  54   buffk[5]=fBuff[22];
  55   buffk[6]=fBuff[23];
  56 }
  57 //____________________________________________________________________
  58 void AliITSPidParItem::GetPionPar(Double_t *buffpi) const{
  59   //get the pions' parameters (Width Landau, Most Probable, Area, Width Gaussian, Chi2 fit, NDF fit, Integral fit)
  60   buffpi[0]=fBuff[26];
  61   buffpi[1]=fBuff[27];
  62   buffpi[2]=fBuff[28];
  63   buffpi[3]=fBuff[29];
  64   buffpi[4]=fBuff[34];
  65   buffpi[5]=fBuff[35];
  66   buffpi[6]=fBuff[36];
  67 }
  68 //____________________________________________________________________
  69 void AliITSPidParItem::GetPar0(Double_t *buff0) const{
  70   //Width Landau for protons, kaons, pions.
  71   buff0[0]=fBuff[0];
  72   buff0[1]=fBuff[13];
  73   buff0[2]=fBuff[26];
  74 }
  75 //____________________________________________________________________
  76 void AliITSPidParItem::GetPar1(Double_t *buff1) const{
  77   //Most Probable for protons, kaons, pions.
  78   buff1[0]=fBuff[1];
  79   buff1[1]=fBuff[14];
  80   buff1[2]=fBuff[27];
  81 }
  82 //____________________________________________________________________
  83 void AliITSPidParItem::GetPar2(Double_t *buff2) const{
  84   //Area for protons, kaons, pions.
  85   buff2[0]=fBuff[2];
  86   buff2[1]=fBuff[15];
  87   buff2[2]=fBuff[28];
  88 }
  89 //____________________________________________________________________
  90 void AliITSPidParItem::GetPar3(Double_t *buff3) const{
  91   //Width Gaussian for protons, kaons, pions.
  92   buff3[0]=fBuff[3];
  93   buff3[1]=fBuff[16];
  94   buff3[2]=fBuff[29];
  95 }
  96 //____________________________________________________________________
  97 void AliITSPidParItem::GetChisquare(Double_t *buffchi) const{
  98   //Chi2 of the fit for protons, kaons, pions.
  99   buffchi[0]=fBuff[8];
 100   buffchi[1]=fBuff[21];
 101   buffchi[2]=fBuff[34];
 102 }
 103 //____________________________________________________________________
 104 void AliITSPidParItem::GetNDF(Double_t *buffndf) const{
 105   //NDF of the fit for protons, kaons, pions.
 106   buffndf[0]=fBuff[9];
 107   buffndf[1]=fBuff[22];
 108   buffndf[2]=fBuff[35];
 109 }
 110 //____________________________________________________________________
 111 void AliITSPidParItem::GetProParFun(Double_t *pfun) const{
 112   //some Protons parameters: Width Landau, Most Probable, Area, Width Gaussian, Integral fit
 113   pfun[0]=fBuff[0];
 114   pfun[1]=fBuff[1];
 115   pfun[2]=fBuff[2];
 116   pfun[3]=fBuff[3];
 117   pfun[4]=fBuff[10];
 118 }
 119 //____________________________________________________________________
 120 void AliITSPidParItem::GetKaoParFun(Double_t *kfun) const{
 121   //some Kaons parameters: Width Landau, Most Probable, Area, Width Gaussian, Integral fit
 122   kfun[0]=fBuff[13];
 123   kfun[1]=fBuff[14];
 124   kfun[2]=fBuff[15];
 125   kfun[3]=fBuff[16];
 126   kfun[4]=fBuff[23];
 127 }
 128 //____________________________________________________________________
 129 void AliITSPidParItem::GetPiParFun(Double_t *pifun) const{
 130   //some Pions parameters: Width Landau, Most Probable, Area, Width Gaussian, Integral fit
 131   pifun[0]=fBuff[26];
 132   pifun[1]=fBuff[27];
 133   pifun[2]=fBuff[28];
 134   pifun[3]=fBuff[29];
 135   pifun[4]=fBuff[36];
 136 }
 137 //____________________________________________________________________
 138 void AliITSPidParItem::GetRangeLim(Double_t *range) const{
 139   //Range limits for the response functions
 140   range[0]=fBuff[11];//proton low
 141   range[1]=fBuff[12];//proton high
 142   range[2]=fBuff[24];//kaon low
 143   range[3]=fBuff[25];//kaon high
 144   range[4]=fBuff[37];//pion low
 145   range[5]=fBuff[38];//pion high
 146 }
 147 //____________________________________________________________________
 148 void AliITSPidParItem::GetProtonParErr(Double_t *bufferp)const{
 149   //errors on the protons' parameters
 150   bufferp[0]=fBuff[4];
 151   bufferp[1]=fBuff[5];
 152   bufferp[2]=fBuff[6];
 153   bufferp[3]=fBuff[7];
 154 }
 155 //____________________________________________________________________
 156 void AliITSPidParItem::GetKaonParErr(Double_t *bufferk)const{
 157   //errors on the kaons' parameters
 158   bufferk[0]=fBuff[17];
 159   bufferk[1]=fBuff[18];
 160   bufferk[2]=fBuff[19];
 161   bufferk[3]=fBuff[20];
 162 }
 163 //____________________________________________________________________
 164 void AliITSPidParItem::GetPionParErr(Double_t *bufferpi)const{
 165   //errors on the pions' parameters
 166   bufferpi[0]=fBuff[30];
 167   bufferpi[1]=fBuff[31];
 168   bufferpi[2]=fBuff[32];
 169   bufferpi[3]=fBuff[33];
 170 }
 171 //____________________________________________________________________
 172 void AliITSPidParItem::PrintParameters() const {
 173   // Prints the data members of this class
 174 cout<<"==============================***************"<<endl;
 175 cout<<"Momentum (GeV/c) - center of the bin - "<<fPCenter<<endl;
 176 cout<<" Width of momentum bin (GeV/c) "<<fPWidth<<endl;
 177 for (Int_t i=0;i<39;i++) cout<<"Parameter"<<i<<" = "<<fBuff[i]<<endl;
 178
 179 }
 180
 181 //_______________________________________________________________________
 182 TF1* AliITSPidParItem::CookFunIts(TString namefun,Double_t *par,Double_t rangei,Double_t rangef,TString comment){
 183   //
 184   TF1 *fun;
 185     if (par[4]!=0) {
 186       fun=new TF1(comment,Langaufun2,rangei,rangef,5);
 187       fun->SetParameters(par);
 188
 189     }
 190     else {fun=new TF1(namefun,"0");}
 191     return fun;
 192 }
 193
 194 //_________________________________________________________________________
 195 Double_t AliITSPidParItem::Langaufun(Double_t *x, Double_t *par) {
 196
 197   //Fit parameters:
 198   //par[0]=Width (scale) parameter of Landau density
 199   //par[1]=Most Probable (MP, location) parameter of Landau density
 200   //par[2]=Total area (integral -inf to inf, normalization constant)
 201   //par[3]=Width (sigma) of convoluted Gaussian function
 202   //
 203   //In the Landau distribution (represented by the CERNLIB approximation),
 204   //the maximum is located at x=-0.22278298 with the location parameter=0.
 205   //This shift is corrected within this function, so that the actual
 206   //maximum is identical to the MP parameter.
 207
 208   // Numeric constants
 209   Double_t invsq2pi = 0.3989422804014;   // (2 pi)^(-1/2)
 210   Double_t mpshift  = -0.22278298;       // Landau maximum location
 211
 212   // Control constants
 213   Double_t np = 100.0;      // number of convolution steps
 214   Double_t sc =   5.0;      // convolution extends to +-sc Gaussian sigmas
 215
 216   // Variables
 217   Double_t xx;
 218   Double_t mpc;
 219   Double_t fland;
 220   Double_t sum = 0.0;
 221   Double_t xlow,xupp;
 222   Double_t step;
 223   Double_t i;
 224
 225
 226   // MP shift correction
 227   mpc = par[1] - mpshift * par[0];
 228
 229   // Range of convolution integral
 230   xlow = x[0] - sc * par[3];
 231   xupp = x[0] + sc * par[3];
 232
 233   step = (xupp-xlow) / np;
 234
 235   // Convolution integral of Landau and Gaussian by sum
 236   for(i=1.0; i<=np/2; i++) {
 237     xx = xlow + (i-.5) * step;
 238     fland = TMath::Landau(xx,mpc,par[0]) / par[0];
 239     sum += fland * TMath::Gaus(x[0],xx,par[3]);
 240
 241     xx = xupp - (i-.5) * step;
 242     fland = TMath::Landau(xx,mpc,par[0]) / par[0];
 243     sum += fland * TMath::Gaus(x[0],xx,par[3]);
 244   }
 245
 246   return (par[2] * step * sum * invsq2pi / par[3]);
 247 }
 248 //_______________________________________________________________________
 249 Double_t AliITSPidParItem::Langaufun2(Double_t *x, Double_t *par){
 250   //
 251   return 1/par[4]*Langaufun(x,par);
 252 }