1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 ////////////////////////////////////////////////////////////////////////////
20 // Container of the distributions for the neural network //
23 // Alex Wilk <wilka@uni-muenster.de> //
25 ////////////////////////////////////////////////////////////////////////////
29 #include <TMultiLayerPerceptron.h>
34 #include "AliTRDgeometry.h"
35 #include "AliTRDCalPIDNN.h"
36 #include "AliTRDcalibDB.h"
38 ClassImp(AliTRDCalPIDNN)
40 //_________________________________________________________________________
41 AliTRDCalPIDNN::AliTRDCalPIDNN()
42 :AliTRDCalPID("pid", "NN PID references for TRD")
45 // The Default constructor
52 //_________________________________________________________________________
53 AliTRDCalPIDNN::AliTRDCalPIDNN(const Text_t *name, const Text_t *title)
54 :AliTRDCalPID(name,title)
57 // The main constructor
64 // //_____________________________________________________________________________
65 // AliTRDCalPIDNN::AliTRDCalPIDNN(const AliTRDCalPIDNN &c)
67 // ,fMeanChargeRatio(c.fMeanChargeRatio)
71 // // Copy constructor
74 // if (this != &c) ((AliTRDCalPIDNN &) c).Copy(*this);
78 //_________________________________________________________________________
79 AliTRDCalPIDNN::~AliTRDCalPIDNN()
85 //if (fModel) delete fModel;
89 //_________________________________________________________________________
90 Bool_t AliTRDCalPIDNN::LoadReferences(Char_t *refFile)
93 // Read the TRD Neural Networks
97 TFile *nnFile = TFile::Open(refFile,"READ");
98 if (!nnFile || !nnFile->IsOpen()) {
99 AliError(Form("Opening TRD histgram file %s failed",refFile));
105 for (Int_t imom = 0; imom < kNMom; imom++) {
106 for (Int_t iplane = 0; iplane < AliTRDgeometry::kNlayer; iplane++) {
107 TMultiLayerPerceptron *nn = (TMultiLayerPerceptron *)
108 nnFile->Get(Form("NN_Mom%d_Plane%d",imom,iplane));
109 fModel->AddAt(nn,GetModelID(imom,0,iplane));
120 //_________________________________________________________________________
121 TObject *AliTRDCalPIDNN::GetModel(Int_t ip, Int_t, Int_t iplane) const
124 // Returns one selected TMultiLayerPerceptron. iType not used.
127 if (ip<0 || ip>= kNMom) return 0x0;
129 AliInfo(Form("Retrive MultiLayerPerceptron for %5.2f GeV/c for plane %d"
133 return fModel->At(GetModelID(ip, 0, iplane));
137 //_________________________________________________________________________
138 Double_t AliTRDCalPIDNN::GetProbability(Int_t spec, Float_t mom
139 , const Float_t * const dedx
140 , Float_t, Int_t iplane) const
143 // Core function of AliTRDCalPID class for calculating the
144 // likelihood for species "spec" (see AliTRDtrack::kNspecie) of a
145 // given momentum "mom", a given dE/dx (slice "dedx") yield per
146 // layer in a given layer (iplane)
149 if (spec < 0 || spec >= AliPID::kSPECIES) return 0.;
151 // find the interval in momentum and track segment length which applies for this data
154 while (imom<AliTRDCalPID::kNMom-1 && mom>fgTrackMomentum[imom]) imom++;
156 Double_t mom1 = fgTrackMomentum[imom-1], mom2 = fgTrackMomentum[imom];
158 TMultiLayerPerceptron *nn = 0x0;
159 if(!(nn = (TMultiLayerPerceptron *) fModel->At(GetModelID(imom-1, spec, iplane/*, ilength*/)))){
160 AliInfo(Form("Looking for mom(%f) plane(%d)", mom-1, iplane));
161 AliError(Form("NN id %d not found in DB.", GetModelID(imom-1, spec, iplane)));
165 Double_t ddedx[AliTRDCalPID::kNSlicesNN];
167 for (int inode=0; inode<AliTRDCalPID::kNSlicesNN; inode++) {
168 ddedx[inode] = (((Double_t) dedx[inode]/kMLPscale)*3) // Bug fix! Needs new reference data or different calculation of dedx!!!!
169 / (AliTRDcalibDB::Instance()->GetNumberOfTimeBins()/AliTRDCalPID::kNSlicesNN);
172 lNN1 = nn->Evaluate(spec, ddedx);
174 if(!(nn = (TMultiLayerPerceptron*)fModel->At(GetModelID(imom, spec, iplane/*, ilength*/)))){
175 AliInfo(Form("Looking for mom(%f) plane(%d)", mom, iplane));
176 AliError(Form("NN id %d not found in DB.", GetModelID(imom, spec, iplane)));
179 lNN2 = nn->Evaluate(spec, ddedx);
181 // return interpolation over momentum binning
182 if (mom < fgTrackMomentum[0]) {
185 else if (mom > fgTrackMomentum[AliTRDCalPID::kNMom-1]) {
189 return lNN1 + (lNN2 - lNN1)*(mom - mom1)/(mom2 - mom1);
194 //_________________________________________________________________________
195 void AliTRDCalPIDNN::Init()
201 fModel = new TObjArray(AliTRDgeometry::kNlayer * AliTRDCalPID::kNMom);
206 //_________________________________________________________________________
207 Int_t AliTRDCalPIDNN::GetModelID(Int_t mom, Int_t /*ii*/, Int_t plane) const
210 // returns the ID of the NN distribution (66 MLPs, ID from 56 to 121)
212 return /*AliPID::kSPECIES * AliTRDCalPID::kNMom + */
213 plane * AliTRDCalPID::kNMom + mom;