Line data Source code
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 16 : ClassImp(AliTRDCalPIDLQ)
42 :
43 : Float_t AliTRDCalPIDLQ::fgTrackSegLength[kNLength] = { 3.7, 3.9, 4.2, 5.0 };
44 :
45 : //_________________________________________________________________________
46 : AliTRDCalPIDLQ::AliTRDCalPIDLQ()
47 2 : :AliTRDCalPID("pid", "LQ PID references for TRD")
48 10 : {
49 : //
50 : // The Default constructor
51 : //
52 :
53 : //Init();
54 :
55 4 : }
56 :
57 : //_________________________________________________________________________
58 : AliTRDCalPIDLQ::AliTRDCalPIDLQ(const Text_t *name, const Text_t *title)
59 0 : :AliTRDCalPID(name,title)
60 0 : {
61 : //
62 : // The main constructor
63 : //
64 :
65 0 : Init();
66 :
67 0 : }
68 :
69 : //_________________________________________________________________________
70 : Bool_t AliTRDCalPIDLQ::LoadReferences(Char_t *refFile)
71 : {
72 : //
73 : // Read the TRD dEdx histograms.
74 : //
75 :
76 0 : if(gSystem->AccessPathName(refFile, kReadPermission)){
77 0 : AliError(Form("File %s.root not readable", refFile));
78 0 : return kFALSE;
79 : }
80 0 : if(!TFile::Open(refFile)){
81 0 : AliError(Form("File %s corrupted", refFile));
82 0 : return kFALSE;
83 : }
84 : TObjArray *pdf(NULL);
85 0 : if (!( pdf = dynamic_cast<TObjArray*>(gFile->Get("PDF_LQ")))) {
86 0 : AliError("PID data not available");
87 0 : return kFALSE;
88 : }
89 0 : fModel=(TObjArray*)pdf->Clone("PDF");
90 0 : gFile->Close();
91 0 : return kTRUE;
92 0 : }
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 0 : if (is < 0 || is >= AliPID::kSPECIES) return NULL;
102 0 : if(ip<0 || ip>= kNMom ) return NULL;
103 :
104 0 : AliDebug(2, Form("Retrive dEdx distribution for %s @ p=%5.2f [GeV/c].", AliPID::ParticleShortName(is), fgTrackMomentum[ip]));
105 0 : return fModel->At(GetModelID(ip, is, slices));
106 0 : }
107 :
108 : //_________________________________________________________________________
109 : Int_t AliTRDCalPIDLQ::GetNrefs()
110 : {
111 : // Returns the number of PDF distribution
112 0 : 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 0 : if (spec < 0 || spec >= AliPID::kSPECIES) return 0.;
129 :
130 0 : Bool_t k2D(slices==2);
131 0 : Double_t x[]={0., 0.};
132 0 : 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 0 : while(imom<kNMom-1 && mom>fgTrackMomentum[imom]) imom++;
141 :
142 0 : Double_t p[2], e[2], r;
143 : TKDPDF *pdf(NULL);
144 :
145 0 : 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 0 : if(!(pdf = dynamic_cast<TKDPDF*>(fModel->At(GetModelID(imom-1, spec, slices))))) {
147 0 : AliError(Form("%cD Ref data @ idx[%d] not found in DB.", k2D?'2':'1', GetModelID(imom-1, spec, slices)));
148 0 : return 0.;
149 : }
150 0 : if(!pdf->GetSize()) pdf->Bootstrap();
151 0 : pdf->Eval(x, r, e[0], kFALSE);
152 0 : p[0]=TMath::Abs(r); // conversion from interpolation to PDF
153 0 : 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 0 : if(!(pdf = dynamic_cast<TKDPDF*>(fModel->At(GetModelID(imom, spec, slices))))){
156 0 : AliError(Form("%cD Ref data @ idx[%d] not found in DB.", k2D?'2':'1', GetModelID(imom, spec, slices)));
157 0 : return p[0];
158 : }
159 0 : if(!pdf->GetSize()) pdf->Bootstrap();
160 0 : pdf->Eval(x, r, e[1], kFALSE);
161 0 : p[1]=TMath::Abs(r); // conversion from interpolation to PDF
162 0 : 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 0 : if(mom < fgTrackMomentum[0]) return p[0];
166 0 : else if(mom > fgTrackMomentum[kNMom-1]) return p[1];
167 : else{
168 0 : Double_t lmom[2] = {fgTrackMomentum[imom-1], fgTrackMomentum[imom]};
169 0 : return p[0] + (p[1] - p[0])*(mom - lmom[0])/(lmom[1] - lmom[0]);
170 : }
171 0 : }
172 :
173 : //_________________________________________________________________________
174 : void AliTRDCalPIDLQ::Init()
175 : {
176 : //
177 : // PID LQ list initialization
178 : //
179 0 : fModel = new TObjArray(AliPID::kSPECIES * kNMom * 2);
180 0 : fModel->SetOwner();
181 0 : }
|
