TPC/AliTPCkineGrid.cxx

   1 /**************************************************************************
   2  * Copyright(c) 2001-2002, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id$ */
  17
  18 ////////////////////////////////////////////////////////////////////////
  19 // Class used by TPC tracking parameterization to handle to tracking
  20 // parameters (efficiencies, etc...) on a kinematic grid [pt,eta].
  21 // User has to provide the grid steps and the values of the parameter
  22 // in the points of the grid. The function GetValueAt(pt,eta) returns
  23 // the result of a linear interpolation at the point [pt,eta].
  24 //
  25 //  Origin: Andrea Dainese, Padova - e-mail: andrea.dainese@pd.infn.it
  26 ////////////////////////////////////////////////////////////////////////
  27
  28 //-- standard headers -----
  29 #include "Riostream.h"
  30 //--- Root headers --------
  31 #include <TMatrixD.h>
  32 #include <TArrayD.h>
  33 //-- AliRoot headers ------
  34 #include "AliTPCkineGrid.h"
  35 //-------------------------
  36
  37 ClassImp(AliTPCkineGrid)
  38
  39 //------------------------------------------------------------------------
  40   AliTPCkineGrid::AliTPCkineGrid()
  41     :TNamed(),
  42      fNpt(0),
  43      fNeta(0),
  44      fPt(0),
  45      fEta(0),
  46      fParams(0)
  47 {
  48 //------------------------------------------------------------------------
  49 // Default constructor
  50 //------------------------------------------------------------------------
  51
  52 }
  53 //------------------------------------------------------------------------
  54 AliTPCkineGrid::AliTPCkineGrid(Int_t npt,Int_t neta,
  55                                Double_t* pt,Double_t* eta)
  56   :TNamed(),
  57     fNpt(npt),
  58     fNeta(neta),
  59     fPt(0),
  60     fEta(0),
  61     fParams(0)
  62 {
  63 //------------------------------------------------------------------------
  64 // Standard constructor
  65 //------------------------------------------------------------------------
  66
  67
  68   fPt   = new TArrayD(fNpt);
  69   fEta  = new TArrayD(fNeta);
  70
  71   for(Int_t i=0; i<npt; i++)  (*fPt)[i]  = pt[i];
  72   for(Int_t i=0; i<neta; i++) (*fEta)[i] = eta[i];
  73
  74   fParams = new TMatrixD(fNpt,fNeta);
  75 }
  76 //-------------------------------------------------------------------------
  77 AliTPCkineGrid::AliTPCkineGrid(const AliTPCkineGrid& grid):TNamed(grid),
  78      fNpt(0),
  79      fNeta(0),
  80      fPt(0),
  81      fEta(0),
  82      fParams(0)
  83 {
  84 //-------------------------------------------------------------------------
  85 // Copy constructor
  86 //-------------------------------------------------------------------------
  87   fNpt = grid.fNpt;
  88   fNeta = grid.fNeta;
  89   fPt = new TArrayD(fNpt);
  90   for(Int_t i=0; i<fNpt; i++) (*fPt)[i] = grid.fPt->At(i);
  91   fEta = new TArrayD(fNeta);
  92   for(Int_t i=0; i<fNeta; i++) (*fEta)[i] = grid.fEta->At(i);
  93
  94   fParams = new TMatrixD(fNpt,fNeta);
  95   for(Int_t i=0; i<fNpt; i++) {
  96     for(Int_t j=0; j<fNeta; j++) (*fParams)(i,j)=(*grid.fParams)(i,j);
  97   }
  98 }
  99 //--------------------------------------------------------------------------
 100 AliTPCkineGrid::~AliTPCkineGrid() {
 101 //--------------------------------------------------------------------------
 102 // Destructor
 103 //--------------------------------------------------------------------------
 104   delete fPt;
 105   delete fEta;
 106   delete fParams;
 107 }
 108 //__________________________________________________________________________
 109 AliTPCkineGrid & AliTPCkineGrid::operator =(const AliTPCkineGrid & param)
 110 {
 111   //
 112   // assignment operator - dummy
 113   //
 114   fNpt=param.fNpt;
 115   return (*this);
 116 }
 117 //--------------------------------------------------------------------------
 118 void AliTPCkineGrid::GetArrayEta(Double_t* eta) const {
 119 //--------------------------------------------------------------------------
 120 // This functions returns an array with the eta points
 121 //--------------------------------------------------------------------------
 122   for(Int_t i=0;i<fNeta;i++) eta[i] = fEta->At(i);
 123   return;
 124 }
 125 //--------------------------------------------------------------------------
 126 void AliTPCkineGrid::GetArrayPt(Double_t* pt) const {
 127 //--------------------------------------------------------------------------
 128 // This functions returns an array with the pt points
 129 //--------------------------------------------------------------------------
 130   for(Int_t i=0;i<fNpt;i++) pt[i] = fPt->At(i);
 131   return;
 132 }
 133 //--------------------------------------------------------------------------
 134 Int_t AliTPCkineGrid::GetBin(Double_t pt,Double_t eta) const {
 135 //--------------------------------------------------------------------------
 136 // This functions tells in which bin of the grid a certain point falls
 137 //--------------------------------------------------------------------------
 138
 139   Int_t etaBin=0,ptBin=0,bin=0;
 140   eta = TMath::Abs(eta);
 141
 142   /*   this is how bins are numbered
 143
 144              ...  ... .
 145              ---+---+---
 146        ^      6 | 7 | 8
 147        |     ---+---+---
 148               3 | 4 | 5
 149        Pt    ---+---+---
 150               0 | 1 | 2
 151
 152                 eta ->
 153    */
 154
 155   if(eta < fEta->At(0)) {
 156     etaBin = 0;
 157   } else if(eta > fEta->At(fNeta-1)) {
 158     etaBin = fNeta;
 159   } else {
 160     for(Int_t i=0; i<fNeta; i++) {
 161       if(eta < fEta->At(i)) {
 162         etaBin = i;
 163         break;
 164       }
 165     }
 166   }
 167   if(pt < fPt->At(0)) {
 168     ptBin = 0;
 169   } else if(pt > fPt->At(fNpt-1)) {
 170     ptBin = fNpt;
 171   } else {
 172     for(Int_t i=0; i<fNpt; i++) {
 173       if(pt < fPt->At(i)) {
 174         ptBin = i;
 175         break;
 176       }
 177     }
 178   }
 179
 180   bin = ptBin*(fNeta+1) + etaBin;
 181
 182   return bin;
 183 }
 184 //--------------------------------------------------------------------------
 185 Double_t AliTPCkineGrid::GetParam(Int_t i) const {
 186 //--------------------------------------------------------------------------
 187 // This functions allows to get parameters using only one index
 188 //--------------------------------------------------------------------------
 189   Int_t ipt = (Int_t)i/fNeta;
 190   Int_t ieta = i-ipt*fNeta;
 191   return GetParam(ipt,ieta);
 192 }
 193 //--------------------------------------------------------------------------
 194 Double_t AliTPCkineGrid::GetValueAt(Double_t pt,Double_t eta) const {
 195 //--------------------------------------------------------------------------
 196 // This functions makes a linear interpolation at the point [pt,eta]
 197 //--------------------------------------------------------------------------
 198
 199   // determine the points to be used in the interpolation:
 200   //
 201   // eta
 202   eta = TMath::Abs(eta);
 203   Int_t etaLow=0,etaUp=0;
 204   if(eta < fEta->At(0)) {
 205     etaLow = 0;
 206     etaUp = 1;
 207   } else if(eta >= fEta->At(fNeta-1)) {
 208     etaLow = fNeta-2;
 209     etaUp = fNeta-1;
 210   } else {
 211     for(Int_t i=0; i<fNeta; i++) {
 212       if(eta < fEta->At(i)) {
 213         etaLow = i-1;
 214         etaUp = i;
 215         break;
 216       }
 217     }
 218   }
 219   //
 220   // pt
 221   Int_t ptLow=0,ptUp=0;
 222   if(pt < fPt->At(0)) {
 223     ptLow = 0;
 224     ptUp = 1;
 225   } else if(pt >= fPt->At(fNpt-1)) {
 226     ptLow = fNpt-2;
 227     ptUp = fNpt-1;
 228   } else {
 229     for(Int_t i=0; i<fNpt; i++) {
 230       if(pt < fPt->At(i)) {
 231         ptLow = i-1;
 232         ptUp = i;
 233         break;
 234       }
 235     }
 236   }
 237
 238   //cerr<<" Pt = ("<<ptLow<<","<<ptUp<<")   Eta = ("<<etaLow<<","<<etaUp<<")\n";
 239
 240   Double_t intValue=0,intValueEtaLow=0,intValueEtaUp=0;
 241   // interpolate, at etaLow, between ptLow and ptUp
 242   intValueEtaLow = (*fParams)(ptLow,etaLow)+
 243                   ((*fParams)(ptUp,etaLow)-(*fParams)(ptLow,etaLow))/
 244                   (fPt->At(ptUp)-fPt->At(ptLow))*(pt-fPt->At(ptLow));
 245   // interpolate, at etaUp, between ptLow and ptUp
 246   intValueEtaUp  = (*fParams)(ptLow,etaUp)+
 247                   ((*fParams)(ptUp,etaUp)-(*fParams)(ptLow,etaUp))/
 248                   (fPt->At(ptUp)-fPt->At(ptLow))*(pt-fPt->At(ptLow));
 249   // interpolate, at pt,  between etaLow and etaUp
 250   intValue       = intValueEtaLow+
 251                   (intValueEtaUp-intValueEtaLow)/
 252                   (fEta->At(etaUp)-fEta->At(etaLow))*(eta-fEta->At(etaLow));
 253
 254   if(intValue<0.) intValue=0.;
 255   return intValue;
 256 }
 257 //--------------------------------------------------------------------------
 258 void AliTPCkineGrid::SetParam(Int_t i,Double_t par) {
 259 //--------------------------------------------------------------------------
 260 // This functions allows to set parameters using only one index
 261 //--------------------------------------------------------------------------
 262   Int_t ipt = (Int_t)i/fNeta;
 263   Int_t ieta = i-ipt*fNeta;
 264   SetParam(ipt,ieta,par);
 265
 266   return;
 267 }
 268
 269
 270
 271
 272
 273
 274
 275
 276
 277
 278
 279
 280
 281