]>
Commit | Line | Data |
---|---|---|
ddcce466 | 1 | ///************************************************************************** |
2 | // * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * | |
3 | // * * | |
4 | // * Author: Alberto Pulvirenti. * | |
5 | // * * | |
6 | // * Permission to use, copy, modify and distribute this software and its * | |
7 | // * documentation strictly for non-commercial purposes is hereby granted * | |
8 | // * without fee, provided that the above copyright notice appears in all * | |
9 | // * copies and that both the copyright notice and this permission notice * | |
10 | // * appear in the supporting documentation. The authors make no claims * | |
11 | // * about the suitability of this software for any purpose. * | |
12 | // * It is provided "as is" without express or implied warranty. * | |
13 | // * * | |
14 | // * AN ITS STAND-ALONE "NEURAL" TRACK FINDER * | |
15 | // * ---------------------------------------- * | |
16 | // * This class implements the Denby-Peterson algorithm for track finding * | |
17 | // * in the ITS stand-alone, by means of a neural network simulation. * | |
18 | // * Many parameters have to be set for the neural network to operate * | |
19 | // * correctly and with a good efficiency. * | |
20 | // * The neural tracker must be feeded with a TTree filled with objects * | |
21 | // * of the class "AliITSNeuralPoint", into a single branch called * | |
22 | // * "Points". * | |
23 | // **************************************************************************/ | |
24 | ||
25 | //#include <fstream> | |
b9d722bc | 26 | #include <Riostream.h> |
27 | #include <stdlib.h> | |
28 | ||
ddcce466 | 29 | //#include <TROOT.h> |
b9d722bc | 30 | #include <TFile.h> |
31 | #include <TTree.h> | |
32 | #include <TMath.h> | |
33 | #include <TLine.h> | |
ddcce466 | 34 | //#include <TMarker.h> |
b9d722bc | 35 | #include <TRandom.h> |
36 | #include <TString.h> | |
37 | #include <TCanvas.h> | |
38 | #include <TVector3.h> | |
ddcce466 | 39 | //#include <TParticle.h> |
b9d722bc | 40 | #include <TObjArray.h> |
ddcce466 | 41 | //#include <TList.h> |
b9d722bc | 42 | |
b9d722bc | 43 | #include "AliITSNeuralTracker.h" |
44 | ||
b9f05b32 | 45 | //using namespace std; |
b9d722bc | 46 | ClassImp(AliITSNeuralTracker) |
47 | ||
48 | //-------------------------------------------------------------------------------------------- | |
49 | ||
50 | AliITSNeuralTracker::AliITSNeuralTracker() | |
51 | { | |
52 | // CONSTRUCTOR | |
53 | // | |
54 | // Initializes some data-members: | |
55 | // - all pointers to NULL | |
56 | // - theta cut to 180 deg. (= no theta cut) | |
57 | // - initial choice of only 1 azimuthal sector | |
58 | // - initial values for the neural network parameters. | |
59 | // | |
60 | // With these settings the neural tracker can't work | |
61 | // because it has not any curvature cut set. | |
62 | ||
63 | fCurvNum = 0; | |
64 | fCurvCut = 0; | |
65 | ||
66 | fSectorNum = 1; | |
67 | fSectorWidth = TMath::Pi() * 2.0; | |
68 | fPolarInterval = 45.0; | |
69 | ||
70 | fStructureOK = kFALSE; | |
71 | ||
72 | fVX = fVY = fVZ = 0.0; | |
73 | ||
74 | Int_t ilayer, itheta; | |
75 | for (ilayer = 0; ilayer < 6; ilayer++) { | |
76 | for (itheta = 0; itheta < 180; itheta++) fPoints[ilayer][itheta] = 0; | |
77 | fThetaCut2DMin[ilayer] = 0.0; | |
78 | fThetaCut2DMax[ilayer] = TMath::Pi(); | |
79 | fThetaCut3DMin[ilayer] = 0.0; | |
80 | fThetaCut3DMax[ilayer] = TMath::Pi(); | |
81 | fHelixMatchCutMin[ilayer] = 1.0; | |
82 | fHelixMatchCutMax[ilayer] = 1.0; | |
83 | } | |
84 | ||
85 | fEdge1 = 0.3; | |
86 | fEdge2= 0.7; | |
87 | ||
88 | fTemperature = 1.0; | |
89 | fStabThreshold = 0.001; | |
90 | fGain2CostRatio = 1.0; | |
91 | fAlignExponent = 1.0; | |
92 | fActMinimum = 0.5; | |
93 | ||
94 | fNeurons = 0; | |
95 | ||
96 | fChains = new TTree("TreeC", "Sextines of points"); | |
97 | fChains->Branch("l0", &fPoint[0], "l0/I"); | |
98 | fChains->Branch("l1", &fPoint[1], "l1/I"); | |
99 | fChains->Branch("l2", &fPoint[2], "l2/I"); | |
100 | fChains->Branch("l3", &fPoint[3], "l3/I"); | |
101 | fChains->Branch("l4", &fPoint[4], "l4/I"); | |
102 | fChains->Branch("l5", &fPoint[5], "l5/I"); | |
103 | } | |
104 | // | |
105 | //-------------------------------------------------------------------------------------------- | |
106 | // | |
107 | AliITSNeuralTracker::~AliITSNeuralTracker() | |
108 | { | |
109 | // DESTRUCTOR | |
110 | // | |
111 | // It Destroys all the dynamic arrays and | |
112 | // clears the TCollections and the points tree | |
113 | ||
b9d722bc | 114 | delete [] fCurvCut; |
115 | ||
116 | Int_t ilayer, itheta; | |
117 | if (fStructureOK) { | |
b9d722bc | 118 | for (ilayer = 0; ilayer < 6; ilayer++) { |
119 | for (itheta = 0; itheta < 180; itheta++) { | |
120 | fPoints[ilayer][itheta]->SetOwner(); | |
121 | delete fPoints[ilayer][itheta]; | |
122 | } | |
123 | } | |
b9d722bc | 124 | fNeurons->SetOwner(); |
125 | delete fNeurons; | |
126 | } | |
127 | ||
b9d722bc | 128 | } |
129 | // | |
130 | //-------------------------------------------------------------------------------------------- | |
131 | // | |
853a0f19 | 132 | void AliITSNeuralTracker::Display(TCanvas*& canv) const |
b9d722bc | 133 | { |
853a0f19 | 134 | // Displays the neural segments |
135 | ||
b9d722bc | 136 | Double_t x1, y1, x2, y2; |
137 | canv->Clear(); | |
138 | TObjArrayIter iter(fNeurons); | |
139 | for (;;) { | |
140 | AliITSneuron *unit = (AliITSneuron*)iter.Next(); | |
141 | if (!unit) break; | |
ddcce466 | 142 | if (unit->Activation() < fActMinimum) continue; |
143 | x1 = unit->Inner()->X(); | |
144 | x2 = unit->Outer()->X(); | |
145 | y1 = unit->Inner()->Y(); | |
146 | y2 = unit->Outer()->Y(); | |
b9d722bc | 147 | TLine *line = new TLine(x1, y1, x2, y2); |
148 | canv->cd(); | |
149 | line->Draw(); | |
150 | } | |
151 | canv->Update(); | |
152 | } | |
153 | // | |
154 | //-------------------------------------------------------------------------------------------- | |
155 | // | |
156 | void AliITSNeuralTracker::SetThetaCuts2D(Double_t *min, Double_t *max) | |
157 | { | |
853a0f19 | 158 | // Cut setter |
159 | ||
b9d722bc | 160 | Int_t i; |
161 | Double_t temp; | |
162 | for (i = 0; i < 5; i++) { | |
163 | if (min[i] > max[i]) { | |
164 | temp = min[i]; | |
165 | min[i] = max[i]; | |
166 | max[i] = temp; | |
167 | } | |
168 | fThetaCut2DMin[i] = min[i]; | |
169 | fThetaCut2DMax[i] = max[i]; | |
170 | } | |
171 | for (i = 0; i < 5; i++) { | |
172 | cout << "Theta 2D cut for layer " << i << " = " << fThetaCut2DMin[i] << " to " << fThetaCut2DMax[i] << endl; | |
173 | } | |
174 | } | |
175 | // | |
176 | //-------------------------------------------------------------------------------------------- | |
177 | // | |
178 | void AliITSNeuralTracker::SetThetaCuts3D(Double_t *min, Double_t *max) | |
179 | { | |
853a0f19 | 180 | // Cut setter |
181 | ||
b9d722bc | 182 | Int_t i; |
183 | Double_t temp; | |
184 | for (i = 0; i < 5; i++) { | |
185 | if (min[i] > max[i]) { | |
186 | temp = min[i]; | |
187 | min[i] = max[i]; | |
188 | max[i] = temp; | |
189 | } | |
190 | fThetaCut3DMin[i] = min[i]; | |
191 | fThetaCut3DMax[i] = max[i]; | |
192 | } | |
193 | for (i = 0; i < 5; i++) { | |
194 | cout << "Theta 3D cut for layer " << i << " = " << fThetaCut3DMin[i] << " to " << fThetaCut3DMax[i] << endl; | |
195 | } | |
196 | } | |
197 | // | |
198 | //-------------------------------------------------------------------------------------------- | |
199 | // | |
200 | void AliITSNeuralTracker::SetHelixMatchCuts(Double_t *min, Double_t *max) | |
201 | { | |
853a0f19 | 202 | // Cut setter |
203 | ||
b9d722bc | 204 | Int_t i; |
205 | Double_t temp; | |
206 | for (i = 0; i < 5; i++) { | |
207 | if (min[i] > max[i]) { | |
208 | temp = min[i]; | |
209 | min[i] = max[i]; | |
210 | max[i] = temp; | |
211 | } | |
212 | fHelixMatchCutMin[i] = min[i]; | |
213 | fHelixMatchCutMax[i] = max[i]; | |
214 | } | |
215 | for (i = 0; i < 5; i++) { | |
216 | cout << "Helix-match cut for layer " << i << " = " << fHelixMatchCutMin[i] << " to " << fHelixMatchCutMax[i] << endl; | |
217 | } | |
218 | } | |
219 | // | |
220 | //-------------------------------------------------------------------------------------------- | |
221 | // | |
222 | void AliITSNeuralTracker::SetCurvatureCuts(Int_t ncuts, Double_t *curv) | |
223 | { | |
224 | // CURVATURE CUTS SETTER | |
225 | // | |
226 | // Requires an array of double values and its dimension | |
227 | // After sorting it in increasing order, the array of curvature cuts | |
228 | // is dinamically allocated, and filled with the sorted cuts array. | |
229 | // A message is shown which lists all the curvature cuts. | |
230 | ||
231 | Int_t i, *ind = new Int_t[ncuts]; | |
232 | TMath::Sort(ncuts, curv, ind, kFALSE); | |
233 | fCurvCut = new Double_t[ncuts]; | |
234 | cout << "\n" << ncuts << " curvature cuts defined" << endl << "-----" << endl; | |
235 | for (i = 0; i < ncuts; i++) { | |
236 | fCurvCut[i] = curv[ind[i]]; | |
237 | cout << "Cut #" << i + 1 << " : " << fCurvCut[i] << endl; | |
238 | } | |
239 | cout << "-----" << endl; | |
240 | fCurvNum = ncuts; | |
241 | } | |
242 | // | |
243 | //-------------------------------------------------------------------------------------------- | |
244 | // | |
245 | void AliITSNeuralTracker::CreateArrayStructure(Int_t nsectors) | |
246 | { | |
247 | // ARRAY CREATOR | |
248 | // | |
249 | // Organizes the array structure to store all points in. | |
250 | // | |
251 | // The array is organized into a "multi-level" TCollection: | |
252 | // - 6 fPoints[] TObjArray containing a TObjArray for each layer | |
253 | // - each TObject contains a TObjArray for each sector. | |
254 | ||
255 | // sets the number of sectors and their width. | |
256 | fSectorNum = nsectors; | |
257 | fSectorWidth = TMath::Pi() * 2.0 / (Double_t)fSectorNum; | |
258 | ||
259 | // creates the TCollection structure | |
260 | Int_t ilayer, isector, itheta; | |
261 | TObjArray *sector = 0; | |
262 | for (ilayer = 0; ilayer < 6; ilayer++) { | |
263 | for (itheta = 0; itheta < 180; itheta++) { | |
264 | if (fPoints[ilayer][itheta]) { | |
265 | fPoints[ilayer][itheta]->SetOwner(); | |
266 | delete fPoints[ilayer][itheta]; | |
267 | } | |
268 | fPoints[ilayer][itheta] = new TObjArray(nsectors); | |
269 | for (isector = 0; isector < nsectors; isector++) { | |
270 | sector = new TObjArray; | |
271 | sector->SetOwner(); | |
272 | fPoints[ilayer][itheta]->AddAt(sector, isector); | |
273 | } | |
274 | } | |
275 | } | |
276 | ||
277 | // Sets a checking flag to TRUE. | |
278 | // This flag is checked before filling up the arrays with the points. | |
279 | fStructureOK = kTRUE; | |
280 | } | |
281 | // | |
282 | //-------------------------------------------------------------------------------------------- | |
283 | // | |
853a0f19 | 284 | Int_t AliITSNeuralTracker::ArrangePoints(TTree* ptstree) |
b9d722bc | 285 | { |
286 | // POINTS STORAGE INTO ARRAY | |
287 | // | |
288 | // Reads points from the tree and creates AliITSNode objects for each one, | |
289 | // storing them into the array structure defined above. | |
290 | // Returns the number of points collected (if successful) or 0 (otherwise) | |
291 | ||
292 | // check: if the points tree is NULL or empty, there is nothing to do... | |
853a0f19 | 293 | if ( !ptstree || (ptstree && !(Int_t)ptstree->GetEntries()) ) { |
b9d722bc | 294 | Error("ArrangePoints", "Points tree is NULL or empty: no points to arrange"); |
295 | return 0; | |
296 | } | |
297 | ||
298 | if (!fStructureOK) { | |
299 | Error("ArrangePoints", "Structure NOT defined. Call CreateArrayStructure() first"); | |
300 | return 0; | |
301 | } | |
302 | ||
303 | Int_t isector, itheta, ientry, ilayer, nentries, pos; | |
304 | TObjArray *sector = 0; | |
305 | AliITSNode *created = 0; | |
306 | AliITSNeuralPoint *cursor = 0; | |
307 | ||
853a0f19 | 308 | ptstree->SetBranchAddress("pos", &pos); |
309 | ptstree->SetBranchAddress("Points", &cursor); | |
310 | nentries = (Int_t)ptstree->GetEntries(); | |
b9d722bc | 311 | |
312 | for (ientry = 0; ientry < nentries; ientry++) { | |
853a0f19 | 313 | ptstree->GetEntry(ientry); |
b9d722bc | 314 | // creates the object |
315 | created = new AliITSNode(cursor, kTRUE); | |
316 | created->SetUser(-1); | |
ddcce466 | 317 | created->PosInTree() = pos; |
b9d722bc | 318 | // finds the sector in phi |
319 | isector = created->GetSector(fSectorWidth); | |
320 | itheta = created->GetThetaCell(); | |
321 | ilayer = created->GetLayer(); | |
322 | if (ilayer < 0 || ilayer > 5) { | |
323 | Error("ArrangePoints", "Layer value %d not allowed. Aborting...", ilayer); | |
324 | return 0; | |
325 | } | |
326 | // selects the right TObjArray to store object into | |
327 | sector = (TObjArray*)fPoints[ilayer][itheta]->At(isector); | |
328 | sector->AddLast(created); | |
329 | } | |
330 | ||
331 | // returns the number of points processed | |
332 | return ientry; | |
333 | } | |
334 | // | |
335 | //-------------------------------------------------------------------------------------------- | |
336 | // | |
337 | void AliITSNeuralTracker::PrintPoints() | |
338 | { | |
339 | // creates the TCollection structure | |
340 | TObjArray *sector = 0; | |
341 | Int_t ilayer, isector, itheta; | |
342 | fstream file("test.log", ios::out); | |
343 | for (ilayer = 0; ilayer < 6; ilayer++) { | |
344 | for (isector = 0; isector < fSectorNum; isector++) { | |
345 | for (itheta = 0; itheta < 180; itheta++) { | |
346 | sector = (TObjArray*)fPoints[ilayer][itheta]->At(isector); | |
347 | file << ilayer << " " << isector << " " << itheta; | |
348 | file << " " << sector->GetSize() << " points" << endl; | |
349 | } | |
350 | } | |
351 | } | |
352 | file.close(); | |
353 | } | |
354 | // | |
355 | //-------------------------------------------------------------------------------------------- | |
356 | // | |
357 | Bool_t AliITSNeuralTracker::PassCurvCut | |
853a0f19 | 358 | (AliITSNode *p1, AliITSNode *p2, Int_t curvindex, Double_t vx, Double_t vy, Double_t vz) |
b9d722bc | 359 | { |
360 | // CURVATURE CUT EVALUATOR | |
361 | // | |
362 | // Checks the passsed point pair w.r. to the current curvature cut | |
363 | // Returns the result of the check. | |
364 | ||
853a0f19 | 365 | if (curvindex < 0 || curvindex >= fCurvNum) { |
366 | Error("PassCurvCut", "Curv index %d out of range", curvindex); | |
b9d722bc | 367 | return kFALSE; |
368 | } | |
369 | ||
370 | // Find the reference layer | |
371 | Int_t lay1 = p1->GetLayer(); | |
372 | Int_t lay2 = p2->GetLayer(); | |
853a0f19 | 373 | Int_t reflayer = (lay1 < lay2) ? lay1 : lay2; |
b9d722bc | 374 | |
375 | Double_t x1 = p1->X() - vx; | |
376 | Double_t x2 = p2->X() - vx; | |
377 | Double_t y1 = p1->Y() - vy; | |
378 | Double_t y2 = p2->Y() - vy; | |
379 | Double_t z1 = p1->Z() - vz; | |
380 | Double_t z2 = p2->Z() - vz; | |
381 | Double_t r1 = sqrt(x1*x1 + y1*y1); | |
382 | Double_t r2 = sqrt(x2*x2 + y2*y2); | |
383 | ||
384 | // calculation of curvature | |
385 | Double_t dx = p1->X() - p2->X(), dy = p1->Y() - p2->Y(); | |
386 | Double_t num = 2 * (x1*y2 - x2*y1); | |
387 | Double_t den = r1*r2*sqrt(dx*dx + dy*dy); | |
388 | Double_t curv = 0.; | |
389 | /* FOR OLD VERSION | |
390 | if (den != 0.) { | |
391 | curv = fabs(num / den); | |
853a0f19 | 392 | if (curv > fCurvCut[curvindex]) return kFALSE; |
b9d722bc | 393 | return kTRUE; |
394 | } | |
395 | else | |
396 | return kFALSE; | |
397 | */ | |
398 | // NEW VERSION | |
399 | if (den != 0.) { | |
b9f05b32 | 400 | curv = TMath::Abs(num / den); |
853a0f19 | 401 | if (curv > fCurvCut[curvindex]) return kFALSE; |
b9d722bc | 402 | } |
403 | else | |
404 | return kFALSE; | |
405 | // calculation of helix matching | |
406 | Double_t arc1 = 2.0 * r1 * curv; | |
407 | Double_t arc2 = 2.0 * r2 * curv; | |
853a0f19 | 408 | Double_t helmatch = 0.0; |
b9d722bc | 409 | if (arc1 > -1.0 && arc1 < 1.0) arc1 = asin(arc1); |
410 | else arc1 = ((arc1 > 0.0) ? 0.5 : 1.5) * TMath::Pi(); | |
411 | if (arc2 > -1.0 && arc2 < 1.0) arc2 = asin(arc2); | |
412 | else arc2 = ((arc2 > 0.0) ? 0.5 : 1.5) * TMath::Pi(); | |
413 | arc1 /= 2.0 * curv; | |
414 | arc2 /= 2.0 * curv; | |
415 | if (arc1 == 0.0 || arc2 == 0.0) return kFALSE; | |
b9f05b32 | 416 | helmatch = TMath::Abs(z1 / arc1 - z2 / arc2); |
853a0f19 | 417 | return (helmatch >= fHelixMatchCutMin[reflayer] && helmatch <= fHelixMatchCutMax[reflayer]); |
b9d722bc | 418 | // END NEW VERSION |
419 | } | |
420 | // | |
421 | //-------------------------------------------------------------------------------------------- | |
422 | // | |
423 | Int_t AliITSNeuralTracker::PassAllCuts | |
853a0f19 | 424 | (AliITSNode *p1, AliITSNode *p2, Int_t curvindex, Double_t vx, Double_t vy, Double_t vz) |
b9d722bc | 425 | { |
426 | // GLOBAL CUT EVALUATOR | |
427 | // | |
428 | // Checks all cuts for the passed point pair. | |
429 | // Return values: | |
430 | // | |
431 | // 0 - All cuts passed | |
432 | // 1 - theta 2D cut not passed | |
433 | // 2 - theta 3D cut not passed | |
434 | // 3 - curvature calculated but cut not passed | |
435 | // 4 - curvature not calculated (division by zero) | |
436 | // 5 - helix cut not passed | |
437 | // 6 - curvature inxed out of range | |
438 | ||
853a0f19 | 439 | if (curvindex < 0 || curvindex >= fCurvNum) return 6; |
b9d722bc | 440 | |
441 | // Find the reference layer | |
442 | Int_t lay1 = p1->GetLayer(); | |
443 | Int_t lay2 = p2->GetLayer(); | |
853a0f19 | 444 | Int_t reflayer = (lay1 < lay2) ? lay1 : lay2; |
b9d722bc | 445 | |
446 | // Swap points in order that r1 < r2 | |
447 | AliITSNode *temp = 0; | |
448 | if (p2->GetLayer() < p1->GetLayer()) { | |
449 | temp = p2; | |
450 | p2 = p1; | |
451 | p1 = temp; | |
452 | } | |
453 | ||
454 | // shift XY coords to the reference to the vertex position, | |
455 | // for easier calculus of quantities. | |
456 | Double_t x1 = p1->X() - vx; | |
457 | Double_t x2 = p2->X() - vx; | |
458 | Double_t y1 = p1->Y() - vy; | |
459 | Double_t y2 = p2->Y() - vy; | |
460 | Double_t z1 = p1->Z() - vz; | |
461 | Double_t z2 = p2->Z() - vz; | |
462 | Double_t r1 = sqrt(x1*x1 + y1*y1); | |
463 | Double_t r2 = sqrt(x2*x2 + y2*y2); | |
464 | ||
465 | // Check for theta cut | |
466 | Double_t dtheta, dtheta3; | |
467 | TVector3 v01(z1, r1, 0.0); | |
468 | TVector3 v12(z2 - z1, r2 - r1, 0.0); | |
469 | dtheta = v01.Angle(v12) * 180.0 / TMath::Pi(); | |
853a0f19 | 470 | TVector3 vv01(x1, y1, z1); |
471 | TVector3 vv12(x2 - x1, y2 - y1, z2 - z1); | |
472 | dtheta3 = vv01.Angle(vv12) * 180.0 / TMath::Pi(); | |
473 | if (dtheta < fThetaCut2DMin[reflayer] || dtheta > fThetaCut2DMax[reflayer]) return 1; | |
474 | if (dtheta3 < fThetaCut3DMin[reflayer] || dtheta3 > fThetaCut3DMax[reflayer]) return 2; | |
b9d722bc | 475 | |
476 | // calculation of curvature | |
477 | Double_t dx = x1 - x2, dy = y1 - y2; | |
478 | Double_t num = 2 * (x1*y2 - x2*y1); | |
479 | Double_t den = r1*r2*sqrt(dx*dx + dy*dy); | |
480 | Double_t curv = 0.; | |
481 | if (den != 0.) { | |
b9f05b32 | 482 | curv = TMath::Abs(num / den); |
853a0f19 | 483 | if (curv > fCurvCut[curvindex]) return 3; |
b9d722bc | 484 | } |
485 | else | |
486 | return 4; | |
487 | ||
488 | // calculation of helix matching | |
489 | Double_t arc1 = 2.0 * r1 * curv; | |
490 | Double_t arc2 = 2.0 * r2 * curv; | |
853a0f19 | 491 | Double_t helmatch = 0.0; |
b9d722bc | 492 | if (arc1 > -1.0 && arc1 < 1.0) arc1 = asin(arc1); |
493 | else arc1 = ((arc1 > 0.0) ? 0.5 : 1.5) * TMath::Pi(); | |
494 | if (arc2 > -1.0 && arc2 < 1.0) arc2 = asin(arc2); | |
495 | else arc2 = ((arc2 > 0.0) ? 0.5 : 1.5) * TMath::Pi(); | |
496 | arc1 /= 2.0 * curv; | |
497 | arc2 /= 2.0 * curv; | |
498 | if (arc1 == 0.0 || arc2 == 0.0) return kFALSE; | |
b9f05b32 | 499 | helmatch = TMath::Abs(z1 / arc1 - z2 / arc2); |
853a0f19 | 500 | if (helmatch < fHelixMatchCutMin[reflayer] || helmatch > fHelixMatchCutMax[reflayer]) return 5; |
b9d722bc | 501 | |
502 | return 0; | |
503 | } | |
504 | // | |
505 | //-------------------------------------------------------------------------------------------- | |
506 | // | |
507 | void AliITSNeuralTracker::StoreAbsoluteMatches() | |
508 | { | |
509 | // Stores in the 'fMatches' array of each node all the points in the | |
510 | // adjacent layers which allow to create neurons accordin to the | |
511 | // helix and theta cut, and also to the largest curvature cut | |
512 | ||
513 | Int_t ilayer, isector, itheta1, itheta2, check; | |
514 | TObjArray *list1 = 0, *list2 = 0; | |
515 | AliITSNode *node1 = 0, *node2 = 0; | |
853a0f19 | 516 | Double_t thetamin, thetamax; |
b9d722bc | 517 | Int_t imin, imax; |
518 | ||
519 | for (isector = 0; isector < fSectorNum; isector++) { | |
520 | // sector is chosen once for both lists | |
521 | for (ilayer = 0; ilayer < 5; ilayer++) { | |
522 | for (itheta1 = 0; itheta1 < 180; itheta1++) { | |
523 | list1 = (TObjArray*)fPoints[ilayer][itheta1]->At(isector); | |
524 | TObjArrayIter iter1(list1); | |
525 | while ( (node1 = (AliITSNode*)iter1.Next()) ) { | |
526 | if (node1->GetUser() >= 0) continue; | |
ddcce466 | 527 | node1->Matches()->Clear(); |
853a0f19 | 528 | thetamin = node1->ThetaDeg() - fPolarInterval; |
529 | thetamax = node1->ThetaDeg() + fPolarInterval; | |
530 | imin = (Int_t)thetamin; | |
531 | imax = (Int_t)thetamax; | |
b9d722bc | 532 | if (imin < 0) imin = 0; |
533 | if (imax > 179) imax = 179; | |
534 | for (itheta2 = imin; itheta2 <= imax; itheta2++) { | |
535 | list2 = (TObjArray*)fPoints[ilayer + 1][itheta2]->At(isector); | |
536 | TObjArrayIter iter2(list2); | |
537 | while ( (node2 = (AliITSNode*)iter2.Next()) ) { | |
538 | check = PassAllCuts(node1, node2, fCurvNum - 1, fVX, fVY, fVZ); | |
539 | switch (check) { | |
540 | case 0: | |
ddcce466 | 541 | node1->Matches()->AddLast(node2); |
b9d722bc | 542 | break; |
543 | case 1: | |
544 | //Info("StoreAbsoluteMatches", "Layer %d: THETA 2D cut not passed", ilayer); | |
545 | break; | |
546 | case 2: | |
547 | //Info("StoreAbsoluteMatches", "Layer %d: THETA 3D cut not passed", ilayer); | |
548 | break; | |
549 | case 3: | |
550 | //Info("StoreAbsoluteMatches", "Layer %d: CURVATURE cut not passed", ilayer); | |
551 | break; | |
552 | case 4: | |
553 | //Info("StoreAbsoluteMatches", "Layer %d: Division by ZERO in curvature evaluation", ilayer); | |
554 | break; | |
555 | case 5: | |
556 | //Info("StoreAbsoluteMatches", "Layer %d: HELIX-MATCH cut not passed", ilayer); | |
557 | break; | |
558 | case 6: | |
559 | //Error("PassAllCuts", "Layer %d: Curv index out of range", ilayer); | |
560 | break; | |
561 | default: | |
562 | Warning("StoreAbsoluteMatches", "Layer %d: %d: unrecognized return value", ilayer, check); | |
563 | } | |
564 | } // while (node2...) | |
565 | } // for (itheta2...) | |
566 | } // while (node1...) | |
567 | } // for (itheta...) | |
568 | } // for (ilayer...) | |
569 | } // for (isector...) | |
570 | } | |
571 | // | |
572 | //-------------------------------------------------------------------------------------------- | |
573 | // | |
574 | void AliITSNeuralTracker::PrintMatches(Bool_t stop) | |
575 | { | |
853a0f19 | 576 | // Prints the matches for each point |
577 | ||
b9d722bc | 578 | TFile *ft = new TFile("its_findables_v1.root"); |
579 | TTree *tt = (TTree*)ft->Get("Tracks"); | |
580 | Int_t it, nP, nU, lab, nF = (Int_t)tt->GetEntries(); | |
581 | tt->SetBranchAddress("nhitsP", &nP); | |
582 | tt->SetBranchAddress("nhitsU", &nU); | |
583 | tt->SetBranchAddress("label", &lab); | |
584 | TString strP("|"), strU("|"); | |
585 | for (it = 0; it < nF; it++) { | |
586 | tt->GetEntry(it); | |
587 | if (nP >= 5) strP.Append(Form("%d|", lab)); | |
588 | if (nU >= 5) strU.Append(Form("%d|", lab)); | |
589 | } | |
590 | ||
591 | TObjArray *sector = 0; | |
592 | Int_t ilayer, isector, itheta, id[3]; | |
593 | AliITSNode *node1 = 0, *node2 = 0; | |
594 | ||
595 | for (ilayer = 0; ilayer < 6; ilayer++) { | |
596 | for (isector = 0; isector < fSectorNum; isector++) { | |
597 | for (itheta = 0; itheta < 180; itheta++) { | |
598 | sector = (TObjArray*)fPoints[ilayer][itheta]->At(isector); | |
599 | TObjArrayIter points(sector); | |
600 | while ( (node1 = (AliITSNode*)points.Next()) ) { | |
601 | for (it = 0; it < 3; it++) id[it] = node1->GetLabel(it); | |
ddcce466 | 602 | nF = (Int_t)node1->Matches()->GetSize(); |
b9d722bc | 603 | cout << "Node layer: " << node1->GetLayer() << ", labels: "; |
604 | cout << id[0] << " " << id[1] << " " << id[2] << " --> "; | |
605 | if (!nF) { | |
606 | cout << "NO MatchES!!!" << endl; | |
607 | continue; | |
608 | } | |
609 | cout << nF << " Matches" << endl; | |
610 | for (it = 0; it < 3; it++) { | |
611 | if (strP.Contains(Form("|%d|", id[it]))) | |
612 | cout << "Belongs to findable (originary) track " << id[it] << endl; | |
613 | if (strU.Contains(Form("|%d|", id[it]))) | |
614 | cout << "Belongs to findable (post-Kalman) track " << id[it] << endl; | |
615 | } | |
ddcce466 | 616 | TObjArrayIter matches(node1->Matches()); |
853a0f19 | 617 | while ( (node2 = (AliITSNode*)matches.Next()) ) { |
b9d722bc | 618 | cout << "Match with " << node2; |
619 | Int_t *sh = node1->SharedID(node2); | |
620 | for (Int_t k = 0; k < 3; k++) | |
621 | if (sh[k] > -1) cout << " " << sh[k]; | |
622 | cout << endl; | |
623 | } | |
624 | if (stop) { | |
625 | cout << "Press a key" << endl; | |
626 | cin.get(); | |
627 | } | |
628 | } | |
629 | } | |
630 | } | |
631 | } | |
632 | } | |
633 | // | |
634 | //-------------------------------------------------------------------------------------------- | |
635 | // | |
636 | void AliITSNeuralTracker::ResetNodes(Int_t isector) | |
637 | { | |
853a0f19 | 638 | // Clears some utility arrays for each node |
639 | ||
b9d722bc | 640 | Int_t ilayer, itheta; |
641 | TObjArray *sector = 0; | |
642 | AliITSNode *node = 0; | |
643 | for (ilayer = 0; ilayer < 6; ilayer++) { | |
644 | for (itheta = 0; itheta < 180; itheta++) { | |
645 | sector = (TObjArray*)fPoints[ilayer][itheta]->At(isector); | |
646 | TObjArrayIter iter(sector); | |
647 | for (;;) { | |
648 | node = (AliITSNode*)iter.Next(); | |
649 | if (!node) break; | |
ddcce466 | 650 | node->InnerOf()->Clear(); |
651 | node->OuterOf()->Clear(); | |
652 | delete node->InnerOf(); | |
653 | delete node->OuterOf(); | |
654 | node->InnerOf() = new TObjArray; | |
655 | node->OuterOf() = new TObjArray; | |
b9d722bc | 656 | } |
657 | } | |
658 | } | |
659 | } | |
660 | // | |
661 | //-------------------------------------------------------------------------------------------- | |
662 | // | |
663 | void AliITSNeuralTracker::NeuralTracking(const char* filesave, TCanvas*& display) | |
853a0f19 | 664 | { |
b9d722bc | 665 | // This is the public method that operates the tracking. |
666 | // It works sector by sector, and at the end saves the found tracks. | |
667 | // Other methods are privare because they have no meaning id used alone, | |
668 | // and sometimes they could get segmentation faults due to uninitialized | |
669 | // datamembert they require to work on. | |
670 | // The argument is a file where the final results have to be stored. | |
853a0f19 | 671 | |
b9d722bc | 672 | Bool_t isStable = kFALSE; |
673 | Int_t i, nUnits = 0, nLinks = 0, nTracks = 0, sectTracks = 0, totTracks = 0; | |
674 | ||
675 | // tracking through sectors | |
676 | cout << endl; | |
677 | Int_t sector, curv; | |
678 | for (sector = 0; sector < fSectorNum; sector++) { | |
679 | cout << "\rSector " << sector << ": " << endl; | |
680 | sectTracks = 0; | |
681 | for(curv = 0; curv < fCurvNum; curv++) { | |
682 | cout << "- curvature step " << curv + 1; | |
683 | cout << " (" << fCurvCut[curv] << "): "; | |
684 | // units creation | |
685 | nUnits = CreateNeurons(sector, curv); | |
686 | if (!nUnits) { | |
687 | cout << "no neurons --> skipped" << endl; | |
688 | continue; | |
689 | } | |
690 | cout << nUnits << " units, " << flush; | |
691 | // units connection | |
692 | nLinks = LinkNeurons(); | |
693 | if (!nLinks) { | |
694 | cout << "no connections --> skipped" << endl; | |
695 | continue; | |
696 | } | |
697 | cout << nLinks << " connections, " << flush; | |
698 | // neural updating | |
699 | for (i = 0;; i++) { | |
700 | isStable = Update(); | |
701 | if (display) Display(display); | |
702 | TObjArrayIter iter(fNeurons); | |
703 | for (;;) { | |
704 | AliITSneuron *n = (AliITSneuron*)iter.Next(); | |
705 | if (!n) break; | |
706 | } | |
707 | if (isStable) break; | |
708 | } | |
709 | cout << i << " cycles --> " << flush; | |
710 | // tracks saving | |
711 | CleanNetwork(); | |
712 | nTracks = Save(sector); | |
713 | cout << nTracks << " tracks" << endl; | |
714 | sectTracks += nTracks; | |
715 | } | |
716 | totTracks += sectTracks; | |
717 | //cout << sectTracks << " tracks found (total = " << totTracks << ") " << flush; | |
718 | } | |
719 | ||
720 | cout << endl << totTracks << " tracks found!!!" << endl; | |
721 | cout << endl << "Saving results in file " << filesave << "..." << flush; | |
722 | TFile *f = new TFile(filesave, "recreate"); | |
723 | fChains->Write("TreeC"); | |
724 | f->Close(); | |
725 | } | |
726 | // | |
727 | //-------------------------------------------------------------------------------------------- | |
728 | // | |
853a0f19 | 729 | Int_t AliITSNeuralTracker::CreateNeurons(Int_t sectoridx, Int_t curvidx) |
730 | { | |
b9d722bc | 731 | // Fills the neuron arrays, following the cut criteria for the selected step |
732 | // (secnum = sector to analyze, curvnum = curvature cut step to use) | |
733 | // It also sets the initial random activation. | |
734 | // In order to avoid a large number of 'new' operations, all existing neurons | |
735 | // are reset and initialized with the new values, and are created new unit only if | |
736 | // it turns out to be necessary | |
737 | // the 'flag' argument is used to decide if the lower edge in the intevral | |
738 | // of the accepted curvatures is given by zero (kFALSE) or by the preceding used cut (kTRUE) | |
739 | // (if this is the first step, anyway, the minimum is always zero) | |
853a0f19 | 740 | |
741 | ResetNodes(sectoridx); | |
b9d722bc | 742 | |
743 | if (fNeurons) delete fNeurons; | |
744 | fNeurons = new TObjArray; | |
745 | ||
746 | AliITSneuron *unit = 0; | |
747 | Int_t itheta, neurons = 0; | |
853a0f19 | 748 | TObjArray *lstsector = 0; |
b9d722bc | 749 | |
750 | // NEW VERSION | |
751 | Double_t vx[6], vy[6], vz[6]; | |
752 | AliITSNode *p[6] = {0, 0, 0, 0, 0, 0}; | |
753 | for (itheta = 0; itheta < 180; itheta++) { | |
853a0f19 | 754 | lstsector = (TObjArray*)fPoints[0][itheta]->At(sectoridx); |
755 | TObjArrayIter lay0(lstsector); | |
b9d722bc | 756 | while ( (p[0] = (AliITSNode*)lay0.Next()) ) { |
757 | if (p[0]->GetUser() >= 0) continue; | |
758 | vx[0] = fVX; | |
759 | vy[0] = fVY; | |
760 | vz[0] = fVZ; | |
ddcce466 | 761 | TObjArrayIter lay1(p[0]->Matches()); |
b9d722bc | 762 | while ( (p[1] = (AliITSNode*)lay1.Next()) ) { |
763 | if (p[1]->GetUser() >= 0) continue; | |
853a0f19 | 764 | if (!PassCurvCut(p[0], p[1], curvidx, fVX, fVY, fVZ)) continue; |
b9d722bc | 765 | unit = new AliITSneuron; |
ddcce466 | 766 | unit->Inner() = p[0]; |
767 | unit->Outer() = p[1]; | |
768 | unit->Activation() = gRandom->Rndm() * (fEdge1 - fEdge2) + fEdge2; | |
769 | unit->Gain() = new TObjArray; | |
b9d722bc | 770 | fNeurons->AddLast(unit); |
ddcce466 | 771 | p[0]->InnerOf()->AddLast(unit); |
772 | p[1]->OuterOf()->AddLast(unit); | |
b9d722bc | 773 | neurons++; |
774 | vx[1] = p[0]->X(); | |
775 | vy[1] = p[0]->Y(); | |
776 | vz[1] = p[0]->Z(); | |
ddcce466 | 777 | TObjArrayIter lay2(p[1]->Matches()); |
b9d722bc | 778 | while ( (p[2] = (AliITSNode*)lay2.Next()) ) { |
779 | if (p[2]->GetUser() >= 0) continue; | |
853a0f19 | 780 | if (!PassCurvCut(p[1], p[2], curvidx, vx[1], vy[1], vz[1])) continue; |
b9d722bc | 781 | unit = new AliITSneuron; |
ddcce466 | 782 | unit->Inner() = p[1]; |
783 | unit->Outer() = p[2]; | |
784 | unit->Activation() = gRandom->Rndm() * (fEdge1 - fEdge2) + fEdge2; | |
785 | unit->Gain() = new TObjArray; | |
b9d722bc | 786 | fNeurons->AddLast(unit); |
ddcce466 | 787 | p[1]->InnerOf()->AddLast(unit); |
788 | p[2]->OuterOf()->AddLast(unit); | |
b9d722bc | 789 | neurons++; |
790 | vx[2] = p[1]->X(); | |
791 | vy[2] = p[1]->Y(); | |
792 | vz[2] = p[1]->Z(); | |
ddcce466 | 793 | TObjArrayIter lay3(p[2]->Matches()); |
b9d722bc | 794 | while ( (p[3] = (AliITSNode*)lay3.Next()) ) { |
795 | if (p[3]->GetUser() >= 0) continue; | |
853a0f19 | 796 | if (!PassCurvCut(p[2], p[3], curvidx, vx[2], vy[2], vz[2])) continue; |
b9d722bc | 797 | unit = new AliITSneuron; |
ddcce466 | 798 | unit->Inner() = p[2]; |
799 | unit->Outer() = p[3]; | |
800 | unit->Activation() = gRandom->Rndm() * (fEdge1 - fEdge2) + fEdge2; | |
801 | unit->Gain() = new TObjArray; | |
b9d722bc | 802 | fNeurons->AddLast(unit); |
ddcce466 | 803 | p[2]->InnerOf()->AddLast(unit); |
804 | p[3]->OuterOf()->AddLast(unit); | |
b9d722bc | 805 | neurons++; |
806 | vx[3] = p[2]->X(); | |
807 | vy[3] = p[2]->Y(); | |
808 | vz[3] = p[2]->Z(); | |
ddcce466 | 809 | TObjArrayIter lay4(p[3]->Matches()); |
b9d722bc | 810 | while ( (p[4] = (AliITSNode*)lay4.Next()) ) { |
811 | if (p[4]->GetUser() >= 0) continue; | |
853a0f19 | 812 | if (!PassCurvCut(p[3], p[4], curvidx, vx[3], vy[3], vz[3])) continue; |
b9d722bc | 813 | unit = new AliITSneuron; |
ddcce466 | 814 | unit->Inner() = p[3]; |
815 | unit->Outer() = p[4]; | |
816 | unit->Activation() = gRandom->Rndm() * (fEdge1 - fEdge2) + fEdge2; | |
817 | unit->Gain() = new TObjArray; | |
b9d722bc | 818 | fNeurons->AddLast(unit); |
ddcce466 | 819 | p[3]->InnerOf()->AddLast(unit); |
820 | p[4]->OuterOf()->AddLast(unit); | |
b9d722bc | 821 | neurons++; |
822 | vx[4] = p[3]->X(); | |
823 | vy[4] = p[3]->Y(); | |
824 | vz[4] = p[3]->Z(); | |
ddcce466 | 825 | TObjArrayIter lay5(p[4]->Matches()); |
b9d722bc | 826 | while ( (p[5] = (AliITSNode*)lay5.Next()) ) { |
827 | if (p[5]->GetUser() >= 0) continue; | |
853a0f19 | 828 | if (!PassCurvCut(p[4], p[5], curvidx, vx[4], vy[4], vz[4])) continue; |
b9d722bc | 829 | unit = new AliITSneuron; |
ddcce466 | 830 | unit->Inner() = p[4]; |
831 | unit->Outer() = p[5]; | |
832 | unit->Activation() = gRandom->Rndm() * (fEdge1 - fEdge2) + fEdge2; | |
833 | unit->Gain() = new TObjArray; | |
b9d722bc | 834 | fNeurons->AddLast(unit); |
ddcce466 | 835 | p[4]->InnerOf()->AddLast(unit); |
836 | p[5]->OuterOf()->AddLast(unit); | |
b9d722bc | 837 | neurons++; |
838 | } // while (p[5]) | |
839 | } // while (p[4]) | |
840 | } // while (p[3]) | |
841 | } // while (p[2]) | |
842 | } // while (p[1]) | |
843 | } // while (p[0]) | |
844 | } // for (itheta...) | |
845 | // END OF NEW VERSION | |
846 | ||
847 | /* OLD VERSION | |
848 | for (ilayer = 0; ilayer < 6; ilayer++) { | |
849 | for (itheta = 0; itheta < 180; itheta++) { | |
853a0f19 | 850 | lstsector = (TObjArray*)fPoints[ilayer][itheta]->At(sectoridx); |
851 | TObjArrayIter inners(lstsector); | |
b9d722bc | 852 | while ( (inner = (AliITSNode*)inners.Next()) ) { |
853 | if (inner->GetUser() >= 0) continue; | |
ddcce466 | 854 | TObjArrayIter outers(inner->Matches()); |
b9d722bc | 855 | while ( (outer = (AliITSNode*)outers.Next()) ) { |
856 | if (outer->GetUser() >= 0) continue; | |
853a0f19 | 857 | if (!PassCurvCut(inner, outer, curvidx, fVX, fVY, fVZ)) continue; |
b9d722bc | 858 | unit = new AliITSneuron; |
ddcce466 | 859 | unit->Inner() = inner; |
860 | unit->Outer() = outer; | |
861 | unit->Activation() = gRandom->Rndm() * (fEdge1 - fEdge2) + fEdge2; | |
862 | unit->Gain() = new TObjArray; | |
b9d722bc | 863 | fNeurons->AddLast(unit); |
ddcce466 | 864 | inner->InnerOf()->AddLast(unit); |
865 | outer->OuterOf()->AddLast(unit); | |
b9d722bc | 866 | neurons++; |
867 | } // for (;;) | |
868 | } // for (;;) | |
869 | } // for (itheta...) | |
870 | } // for (ilayer...) | |
871 | */ | |
872 | ||
873 | fNeurons->SetOwner(); | |
874 | return neurons; | |
875 | } | |
876 | // | |
877 | // | |
878 | // | |
853a0f19 | 879 | Int_t AliITSNeuralTracker::LinkNeurons() const |
880 | { | |
b9d722bc | 881 | // Creates the neural synapses among all neurons |
882 | // which share a point. | |
853a0f19 | 883 | |
b9d722bc | 884 | Int_t total = 0; |
853a0f19 | 885 | TObjArrayIter iter(fNeurons), *testiter; |
b9d722bc | 886 | AliITSneuron *neuron = 0, *test = 0; |
887 | for (;;) { | |
888 | neuron = (AliITSneuron*)iter.Next(); | |
889 | if (!neuron) break; | |
890 | // gain contributors | |
ddcce466 | 891 | testiter = (TObjArrayIter*)neuron->Inner()->OuterOf()->MakeIterator(); |
b9d722bc | 892 | for (;;) { |
853a0f19 | 893 | test = (AliITSneuron*)testiter->Next(); |
b9d722bc | 894 | if (!test) break; |
895 | neuron->Add2Gain(test, fGain2CostRatio, fAlignExponent); | |
896 | total++; | |
897 | } | |
853a0f19 | 898 | delete testiter; |
ddcce466 | 899 | testiter = (TObjArrayIter*)neuron->Outer()->InnerOf()->MakeIterator(); |
b9d722bc | 900 | for (;;) { |
853a0f19 | 901 | test = (AliITSneuron*)testiter->Next(); |
b9d722bc | 902 | if (!test) break; |
903 | neuron->Add2Gain(test, fGain2CostRatio, fAlignExponent); | |
904 | total++; | |
905 | } | |
853a0f19 | 906 | delete testiter; |
b9d722bc | 907 | } |
908 | return total; | |
909 | } | |
910 | // | |
911 | // | |
912 | // | |
913 | Bool_t AliITSNeuralTracker::Update() | |
853a0f19 | 914 | { |
b9d722bc | 915 | // A complete updating cycle with the asynchronous method |
916 | // when the neural network is stable, kTRUE is returned. | |
853a0f19 | 917 | |
b9d722bc | 918 | Double_t actVar = 0.0, totDiff = 0.0; |
919 | TObjArrayIter iter(fNeurons); | |
920 | AliITSneuron *unit; | |
921 | for (;;) { | |
922 | unit = (AliITSneuron*)iter.Next(); | |
923 | if (!unit) break; | |
924 | actVar = Activate(unit); | |
925 | // calculation the relative activation variation | |
926 | totDiff += actVar; | |
927 | } | |
928 | totDiff /= fNeurons->GetSize(); | |
929 | return (totDiff < fStabThreshold); | |
930 | } | |
931 | // | |
932 | // | |
933 | // | |
934 | void AliITSNeuralTracker::CleanNetwork() | |
853a0f19 | 935 | { |
b9d722bc | 936 | // Removes all deactivated neurons, and resolves the competitions |
937 | // to the advantage of the most active unit | |
853a0f19 | 938 | |
b9d722bc | 939 | AliITSneuron *unit = 0; |
940 | TObjArrayIter neurons(fNeurons); | |
941 | while ( (unit = (AliITSneuron*)neurons.Next()) ) { | |
ddcce466 | 942 | if (unit->Activation() < fActMinimum) { |
943 | unit->Inner()->InnerOf()->Remove(unit); | |
944 | unit->Outer()->OuterOf()->Remove(unit); | |
b9d722bc | 945 | delete fNeurons->Remove(unit); |
946 | } | |
947 | } | |
948 | return; | |
949 | Bool_t removed; | |
950 | Int_t nIn, nOut; | |
951 | AliITSneuron *enemy = 0; | |
952 | neurons.Reset(); | |
953 | while ( (unit = (AliITSneuron*)neurons.Next()) ) { | |
ddcce466 | 954 | nIn = (Int_t)unit->Inner()->InnerOf()->GetSize(); |
955 | nOut = (Int_t)unit->Outer()->OuterOf()->GetSize(); | |
b9d722bc | 956 | if (nIn < 2 && nOut < 2) continue; |
957 | removed = kFALSE; | |
958 | if (nIn > 1) { | |
ddcce466 | 959 | TObjArrayIter competing(unit->Inner()->InnerOf()); |
b9d722bc | 960 | while ( (enemy = (AliITSneuron*)competing.Next()) ) { |
ddcce466 | 961 | if (unit->Activation() > enemy->Activation()) { |
962 | enemy->Inner()->InnerOf()->Remove(enemy); | |
963 | enemy->Outer()->OuterOf()->Remove(enemy); | |
b9d722bc | 964 | delete fNeurons->Remove(enemy); |
965 | } | |
966 | else { | |
ddcce466 | 967 | unit->Inner()->InnerOf()->Remove(unit); |
968 | unit->Outer()->OuterOf()->Remove(unit); | |
b9d722bc | 969 | delete fNeurons->Remove(unit); |
970 | removed = kTRUE; | |
971 | break; | |
972 | } | |
973 | } | |
974 | if (removed) continue; | |
975 | } | |
976 | if (nOut > 1) { | |
ddcce466 | 977 | TObjArrayIter competing(unit->Outer()->OuterOf()); |
b9d722bc | 978 | while ( (enemy = (AliITSneuron*)competing.Next()) ) { |
ddcce466 | 979 | if (unit->Activation() > enemy->Activation()) { |
980 | enemy->Inner()->InnerOf()->Remove(enemy); | |
981 | enemy->Outer()->OuterOf()->Remove(enemy); | |
b9d722bc | 982 | delete fNeurons->Remove(enemy); |
983 | } | |
984 | else { | |
ddcce466 | 985 | unit->Inner()->InnerOf()->Remove(unit); |
986 | unit->Outer()->OuterOf()->Remove(unit); | |
b9d722bc | 987 | delete fNeurons->Remove(unit); |
988 | removed = kTRUE; | |
989 | break; | |
990 | } | |
991 | } | |
992 | } | |
993 | } | |
994 | } | |
995 | // | |
996 | // | |
997 | // | |
853a0f19 | 998 | Bool_t AliITSNeuralTracker::CheckOccupation() const |
b9d722bc | 999 | { |
853a0f19 | 1000 | // checks if a track recognized has a point for each layer |
b9d722bc | 1001 | Int_t i; |
1002 | for (i = 0; i < 6; i++) if (fPoint[i] < 0) return kFALSE; | |
1003 | return kTRUE; | |
1004 | } | |
1005 | // | |
1006 | // | |
1007 | // | |
853a0f19 | 1008 | Int_t AliITSNeuralTracker::Save(Int_t sectorid) |
1009 | { | |
b9d722bc | 1010 | // Creates chains of active neurons, in order to |
1011 | // find the tracks obtained as the neural network output. | |
853a0f19 | 1012 | |
b9d722bc | 1013 | // every chain is started from the neurons in the first 2 layers |
1014 | Int_t i, check, stored = 0; | |
853a0f19 | 1015 | Int_t a = sectorid; a = 0; |
1016 | Double_t testact = 0; | |
b9d722bc | 1017 | AliITSneuron *unit = 0, *cursor = 0, *fwd = 0; |
1018 | AliITSNode *node = 0; | |
853a0f19 | 1019 | TObjArrayIter iter(fNeurons), *fwditer; |
b9d722bc | 1020 | TObjArray *removedUnits = new TObjArray(0); |
1021 | TObjArray *removedPoints = new TObjArray(6); | |
1022 | removedUnits->SetOwner(kFALSE); | |
1023 | removedPoints->SetOwner(kFALSE); | |
1024 | for (;;) { | |
1025 | for (i = 0; i < 6; i++) fPoint[i] = -1; | |
1026 | unit = (AliITSneuron*)iter.Next(); | |
1027 | if (!unit) break; | |
ddcce466 | 1028 | if (unit->Inner()->GetLayer() > 0) continue; |
1029 | fPoint[unit->Inner()->GetLayer()] = unit->Inner()->PosInTree(); | |
1030 | fPoint[unit->Outer()->GetLayer()] = unit->Outer()->PosInTree(); | |
1031 | node = unit->Outer(); | |
b9d722bc | 1032 | removedUnits->AddLast(unit); |
ddcce466 | 1033 | removedPoints->AddAt(unit->Inner(), unit->Inner()->GetLayer()); |
1034 | removedPoints->AddAt(unit->Outer(), unit->Outer()->GetLayer()); | |
b9d722bc | 1035 | while (node) { |
853a0f19 | 1036 | testact = fActMinimum; |
ddcce466 | 1037 | fwditer = (TObjArrayIter*)node->InnerOf()->MakeIterator(); |
b9d722bc | 1038 | fwd = 0; |
1039 | for (;;) { | |
853a0f19 | 1040 | cursor = (AliITSneuron*)fwditer->Next(); |
b9d722bc | 1041 | if (!cursor) break; |
ddcce466 | 1042 | if (cursor->Used()) continue; |
1043 | if (cursor->Activation() >= testact) { | |
1044 | testact = cursor->Activation(); | |
b9d722bc | 1045 | fwd = cursor; |
1046 | } | |
1047 | } | |
1048 | if (!fwd) break; | |
1049 | removedUnits->AddLast(fwd); | |
ddcce466 | 1050 | node = fwd->Outer(); |
1051 | fPoint[node->GetLayer()] = node->PosInTree(); | |
b9d722bc | 1052 | removedPoints->AddAt(node, node->GetLayer()); |
1053 | } | |
1054 | check = 0; | |
1055 | for (i = 0; i < 6; i++) if (fPoint[i] >= 0) check++; | |
1056 | if (check >= 6) { | |
1057 | stored++; | |
1058 | fChains->Fill(); | |
1059 | for (i = 0; i < 6; i++) { | |
1060 | node = (AliITSNode*)removedPoints->At(i); | |
1061 | node->SetUser(1); | |
1062 | } | |
853a0f19 | 1063 | fwditer = (TObjArrayIter*)removedUnits->MakeIterator(); |
b9d722bc | 1064 | for (;;) { |
853a0f19 | 1065 | cursor = (AliITSneuron*)fwditer->Next(); |
b9d722bc | 1066 | if(!cursor) break; |
ddcce466 | 1067 | cursor->Used() = 1; |
b9d722bc | 1068 | } |
1069 | } | |
1070 | } | |
1071 | ||
1072 | return stored; | |
1073 | } | |
1074 | /* | |
853a0f19 | 1075 | Int_t AliITSNeuralTracker::Save(Int_t sectoridx) |
b9d722bc | 1076 | // Creates chains of active neurons, in order to |
1077 | // find the tracks obtained as the neural network output. | |
1078 | { | |
1079 | // Reads the final map of neurons and removes all | |
1080 | // units with too small activation | |
1081 | cout << "saving..." << flush; | |
1082 | ||
1083 | // every chain is started from the neurons in the first 2 layers | |
1084 | // 00111111 00011111 00101111 00110111 | |
853a0f19 | 1085 | Int_t allowmask[] = { 0x3F, 0x1F, 0x2F, 0x37, |
b9d722bc | 1086 | 0x3B, 0x3D, 0x3E }; |
1087 | // 00111011 00111101 00111110 | |
1088 | ||
853a0f19 | 1089 | Double_t test_fwd = 0., testback = 0.; |
b9d722bc | 1090 | Int_t ilayer, itheta; |
1091 | AliITSNode *node = 0; | |
1092 | AliITSneuron *fwd = 0, *back = 0; | |
853a0f19 | 1093 | TList *listsector = 0; |
b9d722bc | 1094 | |
1095 | cout << "A -" << fActMinimum << "-" << flush; | |
1096 | for (ilayer = 0; ilayer < 6; ilayer++) { | |
1097 | for (itheta = 0; itheta < 180; itheta++) { | |
853a0f19 | 1098 | listsector = (TList*)fPoints[ilayer][itheta]->At(sectoridx); |
1099 | TListIter iter(listsector); | |
b9d722bc | 1100 | while ( (node = (AliITSNode*)iter.Next()) ) { |
ddcce466 | 1101 | TListIter fwditer(node->InnerOf()); |
1102 | TListIter backiter(node->OuterOf()); | |
853a0f19 | 1103 | testfwd = testback = fActMinimum; |
1104 | while ( (fwd = (AliITSneuron*)fwditer.Next()) ) { | |
ddcce466 | 1105 | if (fwd->Activation() > testfwd) { |
1106 | testfwd = fwd->Activation(); | |
1107 | node->fNext = fwd->Outer(); | |
b9d722bc | 1108 | } |
1109 | } | |
853a0f19 | 1110 | while ( (back = (AliITSneuron*)backiter.Next()) ) { |
ddcce466 | 1111 | if (back->Activation() > testback) { |
1112 | testback = back->Activation(); | |
1113 | node->fPrev = back->Inner(); | |
b9d722bc | 1114 | } |
1115 | } | |
1116 | } | |
1117 | } | |
1118 | } | |
1119 | ||
1120 | cout << "B" << flush; | |
1121 | for (ilayer = 0; ilayer < 5; ilayer++) { | |
1122 | for (itheta = 0; itheta < 180; itheta++) { | |
853a0f19 | 1123 | listsector = (TList*)fPoints[ilayer][itheta]->At(sectoridx); |
1124 | TListIter iter(listsector); | |
b9d722bc | 1125 | while ( (node = (AliITSNode*)iter.Next()) ) { |
1126 | if (node->fNext) { | |
1127 | if (node->fNext->fPrev != node) node->fNext = 0; | |
1128 | } | |
1129 | } | |
1130 | } | |
1131 | } | |
1132 | ||
1133 | cout << "C" << flush; | |
1134 | Bool_t ok; | |
1135 | Int_t stored = 0, l1, l2, i, mask, im; | |
1136 | AliITSNode *temp = 0; | |
1137 | TList *list[2]; | |
1138 | for (itheta = 0; itheta < 180; itheta++) { | |
853a0f19 | 1139 | list[0] = (TList*)fPoints[0][itheta]->At(sectoridx); |
1140 | list[1] = (TList*)fPoints[1][itheta]->At(sectoridx); | |
b9d722bc | 1141 | for (i = 0; i < 2; i++) { |
1142 | TListIter iter(list[i]); | |
1143 | while ( (node = (AliITSNode*)iter.Next()) ) { | |
1144 | //cout << endl << node << flush; | |
1145 | AliITSneuralChain *chain = new AliITSneuralChain; | |
1146 | for (temp = node; temp; temp = temp->fNext) { | |
1147 | chain->Insert((AliITSNeuralPoint*)temp); | |
1148 | } | |
1149 | ok = kFALSE; | |
1150 | mask = chain->OccupationMask(); | |
853a0f19 | 1151 | for (im = 0; im < 7; im++) ok = ok || (mask == allowmask[im]); |
b9d722bc | 1152 | if (!ok) { |
1153 | delete chain; | |
1154 | continue; | |
1155 | } | |
1156 | //cout << " lab " << flush; | |
1157 | chain->AssignLabel(); | |
1158 | //cout << " add " << flush; | |
1159 | fTracks->AddLast(chain); | |
1160 | stored++; | |
1161 | //cout << " " << stored << flush; | |
1162 | } | |
1163 | } | |
1164 | } | |
1165 | ||
1166 | cout << "D" << flush; | |
1167 | return stored; | |
1168 | } | |
1169 | */ | |
1170 | // | |
1171 | // | |
1172 | // | |
1173 | Double_t AliITSNeuralTracker::Activate(AliITSneuron* &unit) | |
853a0f19 | 1174 | { |
b9d722bc | 1175 | // Computes the new neural activation and |
1176 | // returns the variation w.r.t the previous one | |
853a0f19 | 1177 | |
b9d722bc | 1178 | if (!unit) return 0.0; |
853a0f19 | 1179 | Double_t sumgain = 0.0, sumcost = 0.0, input, actOld, actNew; |
b9d722bc | 1180 | |
1181 | // sum gain contributions | |
ddcce466 | 1182 | TObjArrayIter *iter = (TObjArrayIter*)unit->Gain()->MakeIterator(); |
b9d722bc | 1183 | AliITSlink *link; |
1184 | for(;;) { | |
1185 | link = (AliITSlink*)iter->Next(); | |
1186 | if (!link) break; | |
ddcce466 | 1187 | sumgain += link->Linked()->Activation() * link->Weight(); |
b9d722bc | 1188 | } |
1189 | ||
1190 | // sum cost contributions | |
ddcce466 | 1191 | TObjArrayIter *testiter = (TObjArrayIter*)unit->Inner()->InnerOf()->MakeIterator(); |
b9d722bc | 1192 | AliITSneuron *linked = 0; |
1193 | for (;;) { | |
853a0f19 | 1194 | linked = (AliITSneuron*)testiter->Next(); |
b9d722bc | 1195 | if (!linked) break; |
1196 | if (linked == unit) continue; | |
ddcce466 | 1197 | sumcost += linked->Activation(); |
b9d722bc | 1198 | } |
853a0f19 | 1199 | delete testiter; |
ddcce466 | 1200 | testiter = (TObjArrayIter*)unit->Outer()->OuterOf()->MakeIterator(); |
b9d722bc | 1201 | for (;;) { |
853a0f19 | 1202 | linked = (AliITSneuron*)testiter->Next(); |
b9d722bc | 1203 | if (!linked) break; |
1204 | if (linked == unit) continue; | |
ddcce466 | 1205 | sumcost += linked->Activation(); |
b9d722bc | 1206 | } |
1207 | ||
853a0f19 | 1208 | //cout << "gain = " << sumgain << ", cost = " << sumcost << endl; |
b9d722bc | 1209 | |
853a0f19 | 1210 | input = (sumgain - sumcost) / fTemperature; |
ddcce466 | 1211 | actOld = unit->Activation(); |
b9d722bc | 1212 | #ifdef NEURAL_LINEAR |
1213 | if (input <= -2.0 * fTemperature) | |
1214 | actNew = 0.0; | |
1215 | else if (input >= 2.0 * fTemperature) | |
1216 | actNew = 1.0; | |
1217 | else | |
1218 | actNew = input / (4.0 * fTemperature) + 0.5; | |
1219 | #else | |
1220 | actNew = 1.0 / (1.0 + TMath::Exp(-input)); | |
1221 | #endif | |
ddcce466 | 1222 | unit->Activation() = actNew; |
b9d722bc | 1223 | return TMath::Abs(actNew - actOld); |
1224 | } | |
1225 | // | |
1226 | // | |
1227 | // | |
853a0f19 | 1228 | void AliITSNeuralTracker::AliITSneuron::Add2Gain(AliITSneuron *n, Double_t multconst, Double_t exponent) |
b9d722bc | 1229 | { |
853a0f19 | 1230 | // Adds a neuron to the collection of sequenced ones of another neuron |
1231 | ||
b9d722bc | 1232 | AliITSlink *link = new AliITSlink; |
ddcce466 | 1233 | link->Linked() = n; |
b9d722bc | 1234 | Double_t weight = Weight(n); |
ddcce466 | 1235 | link->Weight() = multconst * TMath::Power(weight, exponent); |
b9d722bc | 1236 | fGain->AddLast(link); |
1237 | } | |
1238 | // | |
1239 | // | |
1240 | // | |
1241 | Double_t AliITSNeuralTracker::AliITSneuron::Weight(AliITSneuron *n) | |
1242 | { | |
853a0f19 | 1243 | // computes synaptic weight |
ddcce466 | 1244 | TVector3 me(Outer()->X() - Inner()->X(), Outer()->Y() - Inner()->Y(), Outer()->Z() - Inner()->Z()); |
1245 | TVector3 it(n->Outer()->X() - n->Inner()->X(), n->Outer()->Y() - n->Inner()->Y(), n->Outer()->Z() - n->Inner()->Z()); | |
b9d722bc | 1246 | |
1247 | Double_t angle = me.Angle(it); | |
1248 | Double_t weight = 1.0 - sin(angle); | |
1249 | return weight; | |
1250 | } | |
ddcce466 | 1251 | // |
1252 | // | |
1253 | // | |
1254 | // | |
1255 | AliITSNeuralTracker::AliITSNeuralTracker(const AliITSNeuralTracker &t) | |
1256 | : TObject((TObject&)t) | |
1257 | { | |
1258 | // NO copy constructor | |
1259 | Fatal("AliITSNeuralTracker", "No copy constructor allowed!"); | |
1260 | exit(0); | |
1261 | } | |
1262 | AliITSNeuralTracker& AliITSNeuralTracker::operator=(const AliITSNeuralTracker&) | |
1263 | { | |
1264 | // NO assignment operator | |
1265 | Fatal("AliITSNeuralTracker", "No assignment allowed!"); | |
1266 | return *this; | |
1267 | } | |
1268 | // | |
1269 | AliITSNeuralTracker::AliITSNode::AliITSNode(const AliITSNode &t) | |
1270 | : AliITSNeuralPoint((AliITSNeuralPoint&)t) | |
1271 | { | |
1272 | // NO copy constructor | |
1273 | Fatal("AliITSNode", "No copy constructor allowed!"); | |
1274 | exit(0); | |
1275 | } | |
1276 | AliITSNeuralTracker::AliITSNode& AliITSNeuralTracker::AliITSNode::operator=(const AliITSNode&) | |
1277 | { | |
1278 | // NO assignment operator | |
1279 | Fatal("AliITSNode", "No assignment allowed!"); | |
1280 | return *this; | |
1281 | } | |
1282 | // | |
1283 | AliITSNeuralTracker::AliITSneuron::AliITSneuron(const AliITSneuron &t) | |
1284 | : TObject((TObject&)t) | |
1285 | { | |
1286 | // NO copy constructor | |
1287 | Fatal("AliITSneuron", "No copy constructor allowed!"); | |
1288 | exit(0); | |
1289 | } | |
1290 | AliITSNeuralTracker::AliITSneuron& AliITSNeuralTracker::AliITSneuron::operator=(const AliITSneuron&) | |
1291 | { | |
1292 | // NO assignment operator | |
1293 | Fatal("AliITSneuron", "No assignment allowed!"); | |
1294 | return *this; | |
1295 | } | |
1296 | // | |
1297 | AliITSNeuralTracker::AliITSlink::AliITSlink(const AliITSlink &t) | |
1298 | : TObject((TObject&)t) | |
1299 | { | |
1300 | // NO copy constructor | |
1301 | Fatal("AliITSlink", "No copy constructor allowed!"); | |
1302 | exit(0); | |
1303 | } | |
1304 | AliITSNeuralTracker::AliITSlink& AliITSNeuralTracker::AliITSlink::operator=(const AliITSlink&) | |
1305 | { | |
1306 | // NO assignment operator | |
1307 | Fatal("AliITSlink", "No assignment allowed!"); | |
1308 | return *this; | |
1309 | } | |
1310 | ||
1311 | ||
1312 |