CbmRoot
CbmTrdSetTracksPidANN.cxx
Go to the documentation of this file.
1 // -------------------------------------------------------------------------
2 // ----- CbmTrdSetTracksPidANN source file -----
3 // ----- Created 06/03/07 by Simeon Lebedev -----
4 // -------------------------------------------------------------------------
6 
7 #include "CbmTrdHit.h"
8 #include "CbmTrdTrack.h"
9 
10 #include "FairRootManager.h"
11 
12 #include "TClonesArray.h"
13 #include "TMath.h"
14 #include "TMultiLayerPerceptron.h"
15 #include "TString.h"
16 #include "TSystem.h"
17 #include "TTree.h"
18 
19 #include <fstream>
20 #include <iostream>
21 #include <sstream>
22 #include <vector>
23 using namespace std;
24 
26  : CbmTrdSetTracksPidANN("TrdSetTracksPidANN", "") {}
27 
28 CbmTrdSetTracksPidANN::CbmTrdSetTracksPidANN(const char* name, const char*)
29  : FairTask(name)
30  , fTrackArray(NULL)
31  , fTrdHitArray(NULL)
32  , fNofTracks(0)
33  , fANNPar1(-1.)
34  , fANNPar2(-1.)
35  , fNN()
36  , fTRDGeometryType("h++") {}
37 
38 
40 
42 
44  string fileName = string(gSystem->Getenv("VMCWORKDIR"));
45  vector<string> weightsFilesANN;
46 
47  if (fTRDGeometryType != "h++") {
48  cout << "-E- CbmTrdSetTracksPidANN::Init: " << fTRDGeometryType
49  << " is wrong geometry type." << endl;
50  return kFALSE;
51  }
52  fileName += "/parameters/trd/elid/ann/" + string(fTRDGeometryType) + "/";
53 
54  for (int i = 0; i < 10; i++) {
55  stringstream ss;
56  ss << fileName << "ann_weights_" << (i + 1) << ".txt";
57  weightsFilesANN.push_back(ss.str());
58  }
59 
60  if (fTRDGeometryType == "h++") {
61  fANNPar1 = 1.06;
62  fANNPar2 = 0.57;
63  }
64 
65  for (UInt_t i = 0; i < weightsFilesANN.size(); i++) {
66  ifstream myfile(weightsFilesANN[i].c_str());
67  if (!myfile.is_open()) {
68  cout << "-E- CbmTrdSetTracksPidANN::Init: "
69  << "Could not open input file:" << weightsFilesANN[i] << endl;
70  weightsFilesANN[i] = "";
71  //return kFALSE;
72  }
73  myfile.close();
74  }
75 
76  Float_t inVector[10];
77  Int_t xOut; //output value
78 
79  //init TTree as input data to neural net
80  TTree* simu = new TTree("MonteCarlo", "MontecarloData");
81  simu->Branch("x1", &inVector[0], "x1/F");
82  simu->Branch("x2", &inVector[1], "x2/F");
83  simu->Branch("x3", &inVector[2], "x3/F");
84  simu->Branch("x4", &inVector[3], "x4/F");
85  simu->Branch("x5", &inVector[4], "x5/F");
86  simu->Branch("x6", &inVector[5], "x6/F");
87  simu->Branch("x7", &inVector[6], "x7/F");
88  simu->Branch("x8", &inVector[7], "x8/F");
89  simu->Branch("x9", &inVector[8], "x9/F");
90  simu->Branch("x10", &inVector[9], "x10/F");
91  simu->Branch("xOut", &xOut, "xOut/I");
92 
93  fNN.clear();
94 
95  for (int iH = 0; iH < 10; iH++) {
96  if (weightsFilesANN[iH] == "") {
97  fNN.push_back(NULL);
98  continue;
99  }
100  stringstream ss;
101  for (int iL = 0; iL <= iH; iL++) {
102  if (iL != iH) ss << "x" << (iL + 1) << ",";
103  int nofHidden = 2 * (iH + 1);
104  if (iL == iH) ss << "x" << (iL + 1) << ":" << nofHidden << ":xOut";
105  }
106  TMultiLayerPerceptron* ann =
107  new TMultiLayerPerceptron(ss.str().c_str(), simu);
108  ann->LoadWeights(weightsFilesANN[iH].c_str());
109  fNN.push_back(ann);
110  }
111 
112  return kTRUE;
113 }
114 
116  if (!ReadData()) { Fatal("CbmTrdSetTracksPidANN::Init", "ReadData"); }
117 
118  FairRootManager* ioman = FairRootManager::Instance();
119  if (NULL == ioman) {
120  Fatal("CbmTrdSetTracksPidANN::Init", "RootManager not instantised");
121  }
122 
123  fTrackArray = (TClonesArray*) ioman->GetObject("TrdTrack");
124  if (NULL == fTrackArray) {
125  Fatal("CbmTrdSetTracksPidANN::Init", "No TrdTrack array");
126  }
127 
128  fTrdHitArray = (TClonesArray*) ioman->GetObject("TrdHit");
129  if (NULL == fTrdHitArray) {
130  Fatal("CbmTrdSetTracksPidANN::Init", "No TrdHit array");
131  }
132 
133  return kSUCCESS;
134 }
135 
137  Int_t nTracks = fTrackArray->GetEntriesFast();
138  std::vector<Double_t> eLossVector;
139 
140  for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
141  eLossVector.clear();
142  CbmTrdTrack* pTrack = (CbmTrdTrack*) fTrackArray->At(iTrack);
143  Int_t nofHits = pTrack->GetNofHits();
144  for (Int_t iTRD = 0; iTRD < nofHits; iTRD++) {
145  Int_t index = pTrack->GetHitIndex(iTRD);
146  CbmTrdHit* trdHit = (CbmTrdHit*) fTrdHitArray->At(index);
147  eLossVector.push_back(trdHit->GetELoss());
148  }
149 
150  //------------------transform Data BEGIN--------------------------
151  for (UInt_t j = 0; j < eLossVector.size(); j++) {
152  eLossVector[j] = eLossVector[j] * 1.e6;
153  eLossVector[j] = (eLossVector[j] - fANNPar1) / fANNPar2 - 0.225;
154  }
155  sort(eLossVector.begin(), eLossVector.end());
156  for (UInt_t j = 0; j < eLossVector.size(); j++) {
157  eLossVector[j] = TMath::LandauI(eLossVector[j]);
158  }
159  //------------------transform Data END----------------------------
160 
161  Int_t iANN = nofHits - 1;
162  Double_t nnEval;
163  if (iANN < 0 || iANN >= 12 || fNN[iANN] == NULL) {
164  nnEval = -2.;
165  } else {
166  nnEval = fNN[iANN]->Evaluate(0, &eLossVector[0]);
167  if (TMath::IsNaN(nnEval) == 1) {
168  cout << " -W- CbmTrdSetTracksPidANN: nnEval nan " << endl;
169  nnEval = -2;
170  }
171  }
172  pTrack->SetPidANN(nnEval);
173  }
174 }
175 
177 
CbmTrdSetTracksPidANN::CbmTrdSetTracksPidANN
CbmTrdSetTracksPidANN()
Definition: CbmTrdSetTracksPidANN.cxx:25
CbmTrack::GetNofHits
virtual Int_t GetNofHits() const
Definition: CbmTrack.h:53
CbmTrdHit::GetELoss
Double_t GetELoss() const
Definition: CbmTrdHit.h:79
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmTrdHit
data class for a reconstructed Energy-4D measurement in the TRD
Definition: CbmTrdHit.h:35
CbmTrdSetTracksPidANN::~CbmTrdSetTracksPidANN
virtual ~CbmTrdSetTracksPidANN()
Definition: CbmTrdSetTracksPidANN.cxx:39
CbmTrdSetTracksPidANN::fANNPar1
Double_t fANNPar1
Definition: CbmTrdSetTracksPidANN.h:95
CbmTrdSetTracksPidANN::fNN
std::vector< TMultiLayerPerceptron * > fNN
Definition: CbmTrdSetTracksPidANN.h:98
CbmTrdSetTracksPidANN.h
CbmTrdTrack::SetPidANN
void SetPidANN(Double_t pid)
Definition: CbmTrdTrack.h:42
CbmTrack::GetHitIndex
Int_t GetHitIndex(Int_t iHit) const
Definition: CbmTrack.h:54
CbmTrdSetTracksPidANN::fTrackArray
TClonesArray * fTrackArray
Definition: CbmTrdSetTracksPidANN.h:91
CbmTrdHit.h
Class for hits in TRD detector.
ClassImp
ClassImp(CbmConverterManager) InitStatus CbmConverterManager
Definition: CbmConverterManager.cxx:12
CbmTrdSetTracksPidANN::ReadData
Bool_t ReadData()
Definition: CbmTrdSetTracksPidANN.cxx:43
CbmTrdTrack
Definition: CbmTrdTrack.h:22
CbmTrdSetTracksPidANN::SetParContainers
virtual void SetParContainers()
Definition: CbmTrdSetTracksPidANN.cxx:41
CbmTrdSetTracksPidANN::Init
virtual InitStatus Init()
Definition: CbmTrdSetTracksPidANN.cxx:115
CbmTrdSetTracksPidANN::fTRDGeometryType
TString fTRDGeometryType
Definition: CbmTrdSetTracksPidANN.h:100
CbmTrdSetTracksPidANN::fTrdHitArray
TClonesArray * fTrdHitArray
Definition: CbmTrdSetTracksPidANN.h:92
CbmTrdTrack.h
CbmTrdSetTracksPidANN
Definition: CbmTrdSetTracksPidANN.h:32
CbmTrdSetTracksPidANN::fANNPar2
Double_t fANNPar2
Definition: CbmTrdSetTracksPidANN.h:96
CbmTrdSetTracksPidANN::Finish
virtual void Finish()
Definition: CbmTrdSetTracksPidANN.cxx:176
CbmTrdSetTracksPidANN::Exec
virtual void Exec(Option_t *opt)
Definition: CbmTrdSetTracksPidANN.cxx:136