TRD/Cal/AliTRDCalPIDLQ.cxx

   1 /**************************************************************************
   2 * Copyright(c) 1998-1999, 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 //                                                                      //
  20 // Container for the distributions of dE/dx and the time bin of the     //
  21 // max. cluster for electrons and pions                                 //
  22 //                                                                      //
  23 // Authors:                                                             //
  24 //   Prashant Shukla <shukla@pi0.physi.uni-heidelberg.de>               //
  25 //   Alex Bercuci <a.bercuci@gsi.de>                                    //
  26 //                                                                      //
  27 //////////////////////////////////////////////////////////////////////////
  28
  29 #include <TFile.h>
  30 #include <TROOT.h>
  31 #include <TSystem.h>
  32
  33 #include "AliLog.h"
  34 #include "AliPID.h"
  35
  36 #include "TKDPDF.h"
  37
  38 #include "AliTRDCalPIDLQ.h"
  39 #include "AliTRDcalibDB.h"
  40
  41 ClassImp(AliTRDCalPIDLQ)
  42
  43 Float_t AliTRDCalPIDLQ::fgTrackSegLength[kNLength] = { 3.7, 3.9, 4.2, 5.0 };
  44
  45 //_________________________________________________________________________
  46 AliTRDCalPIDLQ::AliTRDCalPIDLQ()
  47   :AliTRDCalPID("pid", "LQ PID references for TRD")
  48 {
  49   //
  50   //  The Default constructor
  51   //
  52
  53   //Init();
  54
  55 }
  56
  57 //_________________________________________________________________________
  58 AliTRDCalPIDLQ::AliTRDCalPIDLQ(const Text_t *name, const Text_t *title)
  59   :AliTRDCalPID(name,title)
  60 {
  61   //
  62   //  The main constructor
  63   //
  64
  65   Init();
  66
  67 }
  68
  69 //_________________________________________________________________________
  70 Bool_t AliTRDCalPIDLQ::LoadReferences(Char_t *refFile)
  71 {
  72   //
  73   // Read the TRD dEdx histograms.
  74   //
  75
  76   if(gSystem->AccessPathName(refFile, kReadPermission)){
  77     AliError(Form("File %s.root not readable", refFile));
  78     return kFALSE;
  79   }
  80   if(!TFile::Open(refFile)){
  81     AliError(Form("File %s corrupted", refFile));
  82     return kFALSE;
  83   }
  84   TObjArray *pdf(NULL);
  85   if (!( pdf = dynamic_cast<TObjArray*>(gFile->Get("PDF_LQ")))) {
  86     AliError("PID data not available");
  87     return kFALSE;
  88   }
  89   fModel=(TObjArray*)pdf->Clone("PDF");
  90   gFile->Close();
  91   return kTRUE;
  92 }
  93
  94 //_________________________________________________________________________
  95 TObject* AliTRDCalPIDLQ::GetModel(Int_t ip, Int_t is, Int_t slices) const
  96 {
  97   //
  98   // Returns one selected dEdx histogram
  99   //
 100
 101   if (is < 0 || is >= AliPID::kSPECIES) return NULL;
 102   if(ip<0 || ip>= kNMom ) return NULL;
 103
 104   AliDebug(2, Form("Retrive dEdx distribution for %s @ p=%5.2f [GeV/c].", AliPID::ParticleShortName(is), fgTrackMomentum[ip]));
 105   return fModel->At(GetModelID(ip, is, slices));
 106 }
 107
 108 //_________________________________________________________________________
 109 Int_t AliTRDCalPIDLQ::GetNrefs()
 110 {
 111 // Returns the number of PDF distribution
 112   return AliTRDCalPID::kNMom*AliPID::kSPECIES*2;
 113 }
 114
 115 //_________________________________________________________________________
 116 Double_t AliTRDCalPIDLQ::GetProbability(Int_t spec, Float_t mom
 117                                       , const Float_t * const dedx
 118                                       , Float_t length
 119                                       , Int_t slices) const
 120 {
 121 //
 122 // Core function of AliTRDCalPID class for calculating the
 123 // likelihood for species "spec" (see AliTRDtrack::kNspecie) of a
 124 // given momentum "mom" and a given dE/dx (slice "dedx") yield per
 125 // layer
 126 //
 127
 128   if (spec < 0 || spec >= AliPID::kSPECIES) return 0.;
 129
 130   Bool_t k2D(slices==2);
 131   Double_t x[]={0., 0.};
 132   if(!CookdEdx(dedx, x, k2D)) return 0.;
 133
 134   // find the interval in momentum and track segment length which applies for this data
 135 /*  Int_t ilength = 1;
 136   while(ilength<kNLength-1 && length>fgTrackSegLength[ilength]){
 137     ilength++;
 138   }*/
 139   Int_t imom = 1;
 140   while(imom<kNMom-1 && mom>fgTrackMomentum[imom]) imom++;
 141
 142   Double_t p[2], e[2], r;
 143   TKDPDF *pdf(NULL);
 144
 145   AliDebug(2, Form("Looking %cD for %s@p=%6.4f[GeV/c] log(dEdx)={%7.2f %7.2f}[a.u.] l=%4.2f[cm].", k2D?'2':'1', AliPID::ParticleShortName(spec), mom, x[0], x[1], length));
 146   if(!(pdf = dynamic_cast<TKDPDF*>(fModel->At(GetModelID(imom-1, spec, slices))))) {
 147     AliError(Form("%cD Ref data @ idx[%d] not found in DB.", k2D?'2':'1', GetModelID(imom-1, spec, slices)));
 148     return 0.;
 149   }
 150   if(!pdf->GetSize()) pdf->Bootstrap();
 151   pdf->Eval(x, r, e[0], kFALSE);
 152   p[0]=TMath::Abs(r); // conversion from interpolation to PDF
 153   AliDebug(2, Form("LQ=%6.3f+-%5.2f%% @ %4.1f[GeV/c]", p[0], 1.E2*e[0]/p[0], fgTrackMomentum[imom-1]));
 154
 155   if(!(pdf = dynamic_cast<TKDPDF*>(fModel->At(GetModelID(imom, spec, slices))))){
 156     AliError(Form("%cD Ref data @ idx[%d] not found in DB.", k2D?'2':'1', GetModelID(imom, spec, slices)));
 157     return p[0];
 158   }
 159   if(!pdf->GetSize()) pdf->Bootstrap();
 160   pdf->Eval(x, r, e[1], kFALSE);
 161   p[1]=TMath::Abs(r); // conversion from interpolation to PDF
 162   AliDebug(2, Form("LQ=%6.3f+-%5.2f%% @ %4.1f[GeV/c]", p[1], 1.E2*e[1]/p[1], fgTrackMomentum[imom]));
 163
 164   // return interpolation over momentum binning
 165   if(mom < fgTrackMomentum[0]) return p[0];
 166   else if(mom > fgTrackMomentum[kNMom-1]) return p[1];
 167   else{
 168     Double_t lmom[2] = {fgTrackMomentum[imom-1],  fgTrackMomentum[imom]};
 169     return p[0] + (p[1] - p[0])*(mom - lmom[0])/(lmom[1] - lmom[0]);
 170   }
 171 }
 172
 173 //_________________________________________________________________________
 174 void AliTRDCalPIDLQ::Init()
 175 {
 176 //
 177 // PID LQ list initialization
 178 //
 179   fModel = new TObjArray(AliPID::kSPECIES  * kNMom * 2);
 180   fModel->SetOwner();
 181 }