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