]>
Commit | Line | Data |
---|---|---|
cec46807 | 1 | #if !defined(__CINT__) || defined(__MAKECINT__) |
2 | #include "Riostream.h" | |
8b462fd8 | 3 | #include "TFile.h" |
4 | #include "TTree.h" | |
cec46807 | 5 | #include "TStopwatch.h" |
6 | ||
7 | #include "AliRun.h" | |
8 | #include "AliMagF.h" | |
9 | #include "AliRunLoader.h" | |
10 | #include "AliTPCLoader.h" | |
11 | #include "AliITSLoader.h" | |
12 | #include "AliITS.h" | |
13 | #include "AliITSgeom.h" | |
8b462fd8 | 14 | #include "AliITStrackerANN.h" |
cec46807 | 15 | #include "AliESD.h" |
16 | #endif | |
17 | ||
18 | extern AliRun *gAlice; | |
19 | ||
20 | Int_t AliITSFindTracksANN | |
21 | (Int_t nsectors = 10, // number of azimutal sectors | |
22 | Int_t nev = 5) // number of events to process | |
23 | { | |
24 | ||
25 | // ================================== | |
26 | // ==== EVENT READING =============== | |
27 | // ================================== | |
28 | ||
29 | // Remove an already existing Run Loader | |
30 | if (gAlice) { | |
31 | delete gAlice->GetRunLoader(); | |
32 | delete gAlice; | |
33 | gAlice=0; | |
34 | } | |
35 | ||
36 | // Instance the new Run Loader | |
37 | AliRunLoader* rl = AliRunLoader::Open("galice.root"); | |
38 | if (rl == 0x0) { | |
39 | cerr<<"AliITSFindTracks.C : Can not open session RL=NULL"<< endl; | |
40 | return 3; | |
41 | } | |
42 | ||
43 | // Instance the ITS Loader | |
44 | AliITSLoader* itsl = (AliITSLoader*)rl->GetLoader("ITSLoader"); | |
45 | if (itsl == 0x0) { | |
46 | cerr<<"AliITSFindTracksV2.C : Can not get ITS loader"<<endl; | |
47 | return 4; | |
48 | } | |
49 | ||
50 | // Load the gAlice | |
51 | if (rl->LoadgAlice()) { | |
52 | cerr<<"AliITSFindTracksV2.C : LoadgAlice returned error"<<endl; | |
53 | delete rl; | |
54 | return 3; | |
55 | } | |
56 | ||
57 | // Set NECESSARY conversion constant for magnetic field | |
58 | AliKalmanTrack::SetConvConst | |
59 | ( | |
60 | 1000/0.299792458/rl->GetAliRun()->Field()->SolenoidField() | |
61 | ); | |
62 | ||
63 | // Get the ITS data & geometry & recpoints | |
64 | AliITS *dITS = (AliITS*)rl->GetAliRun()->GetDetector("ITS"); | |
65 | if (!dITS) { | |
66 | cerr<<"AliITSFindClusters.C : Can not find the ITS detector !"<<endl; | |
67 | return 6; | |
68 | } | |
69 | AliITSgeom *geom = dITS->GetITSgeom(); | |
70 | itsl->LoadRecPoints("read"); | |
71 | ||
72 | // ================================== | |
73 | // ==== CURVATURE CUT DEFINITION ==== | |
74 | // ================================== | |
75 | ||
76 | // These values define the curvature cuts for all steps | |
77 | // within a sector. | |
78 | // For a greater clarity, the cuts are indicated in units | |
79 | // of transverse momentum (GeV/c) but these value have no | |
80 | // exact physical meaning, but are useful to understand | |
81 | // well what means a choice in the value of a certain | |
82 | // curvature constraint | |
83 | // NOTE: be careful to make sure that the 'ncuts' variable | |
84 | // have the same value of the dimension of the allocated arrays | |
85 | ||
86 | Int_t ncuts; | |
87 | Double_t *p, *cut; | |
88 | ||
89 | ncuts = 6; | |
90 | p = new Double_t[6]; | |
91 | cut = new Double_t[6]; | |
92 | cut[ 0] = 0.0006; | |
93 | cut[ 1] = 0.0010; | |
94 | cut[ 2] = 0.0013; | |
95 | cut[ 3] = 0.0018; | |
96 | cut[ 4] = 0.0020; | |
97 | cut[ 5] = 0.0030; | |
98 | ||
99 | ||
100 | // ========================== | |
101 | // ==== OTHER PARAMETERS ==== | |
102 | // ========================== | |
103 | ||
104 | ||
105 | Double_t helix[5] = { 0.03, 0.02, 0.01, 0.03, 0.2 }; | |
106 | Double_t theta2D[5] = { 1.5, 1.0, 1.0, 3.0, 10.0 }; | |
107 | Double_t theta3D[5] = { 200., 200., 200., 200.0, 200.0 }; | |
108 | ||
109 | ||
110 | /* | |
111 | Double_t helix[5] = { 0.3, 0.3, 0.3, 0.5, 0.5 }; | |
112 | Double_t theta2D[5] = { 5.0, 5.0, 5.0, 5.0, 5.0 }; | |
113 | Double_t theta3D[5] = { 5.0, 5.0, 5.0, 5.0, 5.0 }; | |
114 | */ | |
115 | ||
116 | Double_t temp = 1.0; // temperature parameter | |
117 | Double_t var = 0.00001; // stabilization threshold | |
118 | ||
119 | Double_t exp = 20.0; // straight-line excitator | |
120 | Double_t gtoc = 6.0; // gain/cost contribution ratio | |
121 | ||
122 | Double_t min = 0.4; // minimum in random activation initialization | |
123 | Double_t max = 0.6; // maximum in random activation initialization | |
124 | Double_t actmin = 0.55; // activation threshold for binary map conversion | |
125 | ||
126 | // Instance the Tracker | |
127 | AliITStrackerANN tracker(geom, 2); | |
128 | ||
129 | // Set cuts | |
130 | tracker.SetVertex(0.0, 0.0, 0.0); | |
131 | ||
132 | tracker.SetCuts(ncuts, cut, theta2D, theta3D, helix); | |
133 | tracker.SetTemperature(temp); | |
134 | tracker.SetVariationLimit(var); | |
135 | tracker.SetGainToCostRatio(gtoc); | |
136 | tracker.SetWeightExponent(exp); | |
137 | tracker.SetInitInterval(min, max); | |
138 | tracker.SetActThreshold(actmin); | |
139 | ||
140 | tracker.SetPolarInterval(45.0); | |
141 | ||
142 | // Set number of events | |
143 | if (nev > rl->GetNumberOfEvents()) nev = rl->GetNumberOfEvents(); | |
144 | Int_t rc = 0; | |
145 | ||
146 | // Get ESD files | |
147 | TFile *itsf=TFile::Open("AliESDits.root","RECREATE"); | |
148 | if ((!itsf)||(!itsf->IsOpen())) { | |
149 | cerr << "Can't AliESDits.root !\n"; | |
150 | return 1; | |
151 | } | |
152 | TFile *tpcf=TFile::Open("AliESDtpc.root"); | |
153 | if ((!tpcf)||(!tpcf->IsOpen())) { | |
154 | cerr<<"Can't AliESDtpc.root !\n"; | |
155 | return 1; | |
156 | } | |
8b462fd8 | 157 | AliESD* event = new AliESD; |
158 | TTree* esdTree = (TTree*) tpcf->Get("esdTree"); | |
159 | if (!esdTree) { | |
160 | cerr<<"no ESD tree found !\n"; | |
161 | return 1; | |
162 | } | |
163 | esdTree->SetBranchAddress("ESD", &event); | |
cec46807 | 164 | |
165 | // Loop on events | |
cec46807 | 166 | TStopwatch timer; |
167 | for (Int_t i = 0; i < nev; i++) { | |
cec46807 | 168 | cerr << "Processing event number: " << i << endl; |
8b462fd8 | 169 | esdTree->GetEvent(i); |
cec46807 | 170 | |
171 | rl->GetEvent(i); | |
172 | ||
173 | // Get clusters from file | |
174 | TTree *cTree=itsl->TreeR(); | |
175 | if (!cTree) { | |
176 | cerr<<"AliITSFindTracksANN.C : NO clusters tree!" << endl; | |
177 | return 4; | |
178 | } | |
179 | ||
180 | // Load clusters in tracker | |
181 | tracker.LoadClusters(cTree); | |
182 | ||
183 | // Create array structure and arrange points in it | |
184 | tracker.CreateArrayStructure(nsectors); | |
185 | tracker.ArrangePoints("its.recpoints.txt"); | |
186 | tracker.StoreOverallMatches(); | |
187 | //tracker.PrintMatches(kFALSE); | |
188 | ||
189 | // *************************** | |
190 | // NEURAL TRACKING OPERATION | |
191 | // *************************** | |
192 | ||
193 | Bool_t isStable = kFALSE; | |
194 | Int_t i, nUnits = 0, nLinks = 0, nTracks = 0; | |
195 | Int_t sectTracks = 0, totTracks = 0; | |
196 | ||
197 | // tracking through sectors | |
198 | cout << endl; | |
199 | Int_t sector, curv; | |
200 | // nUnits = tracker.CreateNetwork(0, ncuts - 1); | |
201 | // cout << endl << nUnits << " NEURONS CREATED" << endl; | |
202 | // return; | |
203 | ||
204 | for (sector = 0; sector < nsectors; sector++) { | |
205 | sectTracks = 0; | |
206 | for(curv = 0; curv < ncuts; curv++) { | |
207 | // units creation | |
208 | nUnits = tracker.CreateNetwork(sector, curv); | |
209 | if (!nUnits) { | |
210 | cout << "no neurons --> skipped" << endl; | |
211 | continue; | |
212 | } | |
213 | // units connection | |
214 | nLinks = tracker.LinkNeurons(); | |
215 | if (!nLinks) { | |
216 | cout << "no connections --> skipped" << endl; | |
217 | continue; | |
218 | } | |
219 | // neural updating | |
220 | for (i = 0;; i++) { | |
221 | isStable = tracker.Update(); | |
222 | if (isStable) break; | |
223 | } | |
224 | // tracks saving | |
225 | tracker.CleanNetwork(); | |
226 | nTracks = tracker.StoreTracks(); | |
227 | cout << nUnits << " units, "; | |
228 | cout << nLinks << " connections, "; | |
229 | cout << i << " cycles --> "; | |
230 | cout << nTracks << " tracks" << endl; | |
231 | sectTracks += nTracks; | |
232 | } | |
233 | cout << "\n >>> Total tracks found in sector: " << sectTracks << endl; | |
234 | totTracks += sectTracks; | |
235 | } | |
236 | tracker.ExportTracks("ITS.Neural.root"); | |
237 | cout << "\n\n--- Totally found " << totTracks << " tracks ---\n\n" << endl; | |
238 | ||
239 | // *************************** | |
240 | ||
241 | // End of operations: unload clusters | |
242 | tracker.UnloadClusters(); | |
243 | ||
cec46807 | 244 | } |
245 | timer.Stop(); | |
246 | timer.Print(); | |
247 | ||
248 | // Close files & delete objects | |
8b462fd8 | 249 | delete event; |
cec46807 | 250 | tpcf->Close(); |
251 | itsf->Close(); | |
252 | delete rl; | |
253 | return rc; | |
254 | } | |
255 | ||
256 |