CbmRoot
L1CATrackFinder.cxx
Go to the documentation of this file.
1 /*
2  *=====================================================
3  *
4  * CBM Level 1 4D Reconstruction
5  *
6  * Authors: V.Akishina, I.Kisel, S.Gorbunov, I.Kulakov, M.Zyzak
7  * Documentation: V.Akishina
8  *
9  * e-mail : v.akishina@gsi.de
10  *
11  *=====================================================
12  *
13  * Finds tracks using the Cellular Automaton algorithm
14  *
15  */
16 
17 #include "L1AddMaterial.h"
18 #include "L1Algo.h"
19 #include "L1Branch.h"
20 #include "L1Extrapolation.h"
21 #include "L1Filtration.h"
22 #include "L1HitPoint.h"
23 #include "L1Track.h"
24 #include "L1TrackPar.h"
25 //#include "TDHelper.h"
26 #include "L1Grid.h"
27 #include "L1HitArea.h"
28 #include "L1Portion.h"
29 
30 #ifdef _OPENMP
31 #include "omp.h"
32 #include "pthread.h"
33 #endif
34 
35 #include "L1HitsSortHelper.h"
36 
37 
38 #include "L1Timer.h"
39 #include "TRandom.h"
40 #ifdef DRAW
41 #include "utils/L1AlgoDraw.h"
42 #endif
43 #ifdef PULLS
44 #include "utils/L1AlgoPulls.h"
45 #endif
46 #ifdef TRIP_PERFORMANCE
47 #include "utils/L1AlgoEfficiencyPerformance.h"
48 #endif
49 #ifdef DOUB_PERFORMANCE
50 #include "utils/L1AlgoEfficiencyPerformance.h"
51 #endif
52 
53 #ifdef TBB
54 #include "L1AlgoTBB.h"
55 #endif
56 
57 #include <algorithm>
58 #include <fstream>
59 #include <iostream>
60 #include <list>
61 #include <map>
62 #include <stdio.h>
63 
64 
65 // using namespace std;
66 using std::cout;
67 using std::endl;
68 using std::vector;
69 
70 
73 inline void L1Algo::f10( // input
74  Tindex start_lh,
75  Tindex n1_l,
76  L1HitPoint* StsHits_l,
77  // output
78  fvec* u_front_l,
79  fvec* u_back_l,
80  fvec* zPos_l,
81  THitI* hitsl,
82  fvec* HitTime_l,
83  fvec* HitTimeEr,
84  fvec* Event_l,
85  fvec* d_x,
86  fvec* d_y,
87  fvec* d_xy,
88  fvec* d_u,
89  fvec* d_v) {
90  const Tindex& end_lh = start_lh + n1_l;
91 
92 
93  for (Tindex ilh = start_lh, i1 = 0; ilh < end_lh; ++ilh, ++i1) {
94  Tindex i1_V = i1 / fvecLen;
95  Tindex i1_4 = i1 % fvecLen;
96  L1HitPoint& hitl = StsHits_l[ilh];
97 
98 
99 #ifdef USE_EVENT_NUMBER
100  Event_l[i1_V][i1_4] = hitl.track;
101 #endif
102 
103  HitTime_l[i1_V][i1_4] = hitl.time;
104  HitTimeEr[i1_V][i1_4] = hitl.timeEr;
105 
106  hitsl[i1] = ilh;
107  u_front_l[i1_V][i1_4] = hitl.U();
108  u_back_l[i1_V][i1_4] = hitl.V();
109 
110  if (fUseHitErrors) {
111  // d_x[i1_V][i1_4] = hitl.dX();
112  // d_y[i1_V][i1_4] = hitl.dY();
113  // d_xy[i1_V][i1_4] = hitl.dXY();
114  d_u[i1_V][i1_4] = hitl.dU();
115  d_v[i1_V][i1_4] = hitl.dV();
116  }
117 
118  zPos_l[i1_V][i1_4] = hitl.Z();
119  }
120 }
121 
122 
124 inline void L1Algo::f11(
125  int istal,
126  int istam,
127  Tindex n1_V,
128 
129  fvec* u_front_l,
130  fvec* u_back_l,
131  fvec* zPos_l,
132  fvec* HitTime_l,
133  fvec* HitTimeEr,
134  // output
135  // L1TrackPar *T_1,
136  L1TrackPar* T_1,
137  L1FieldRegion* fld_1,
138  fvec* d_x,
139  fvec* d_y,
140  fvec* d_xy,
141  fvec* d_u,
142  fvec* d_v) {
143  L1Station& stal = vStations[istal];
144  L1Station& stam = vStations[istam];
145  fvec zstal = stal.z;
146  fvec zstam = stam.z;
147 
148  int istal_global = 5;
149  int istam_global = 6;
150  L1Station& stal_global = vStations[istal_global];
151  L1Station& stam_global = vStations[istam_global];
152  fvec zstal_global = stal_global.z;
153  fvec zstam_global = stam_global.z;
154 
155  L1FieldRegion
156  fld0; // by left hit, target and "center" (station in-between). Used here for extrapolation to target and to middle hit
157  L1FieldValue l_B_global, centB, centB_global,
158  l_B _fvecalignment; // field for singlet creation
159  L1FieldValue m_B,
160  m_B_global _fvecalignment; // field for the next extrapolations
161 
162  for (int i1_V = 0; i1_V < n1_V; i1_V++) {
163  L1TrackPar& T = T_1[i1_V];
164  L1FieldRegion& fld1 = fld_1
165  [i1_V]; // by middle hit, left hit and "center" . Will be used for extrapolation to right hit
166 
167  // get the magnetic field along the track.
168  // (suppose track is straight line with origin in the target)
169  fvec& u = u_front_l[i1_V];
170  fvec& v = u_back_l[i1_V];
171  fvec xl, yl; // left(1-st) hit coor
172  fvec zl = zPos_l[i1_V];
173  fvec& time = HitTime_l[i1_V];
174  fvec& timeEr = HitTimeEr[i1_V];
175  const fvec& dzli = 1. / (zl - targZ);
176 
177  fvec dx1, dy1, dxy1 = 0;
178 
179  dUdV_to_dX(d_u[i1_V], d_v[i1_V], dx1, stal);
180  dUdV_to_dY(d_u[i1_V], d_v[i1_V], dy1, stal);
181  dUdV_to_dXdY(d_u[i1_V], d_v[i1_V], dxy1, stal);
182 
183  StripsToCoor(u, v, xl, yl, stal);
184 
185  const fvec& tx = (xl - targX) * dzli;
186  const fvec& ty = (yl - targY) * dzli;
187 
188  // estimate field for singlet creation
189  int istac = istal / 2; // "center" station
190  // if (istal >= NMvdStations) // to make it in the same way for both with and w\o mvd
191  // istac = (istal-NMvdStations)/2+NMvdStations;
192  L1Station& stac = vStations[istac];
193  fvec zstac = stac.z;
194 
195  int istac_global = istal_global / 2; // "center" station
196 
197  L1Station& stac_global = vStations[istac_global];
198  fvec zstac_global = stac.z;
199 
200  stac.fieldSlice.GetFieldValue(
201  targX + tx * (zstac - targZ), targY + ty * (zstac - targZ), centB);
202  stal.fieldSlice.GetFieldValue(
203  targX + tx * (zstal - targZ), targY + ty * (zstal - targZ), l_B);
204 
205  stam_global.fieldSlice.GetFieldValue(targX + tx * (zstam_global - targZ),
206  targY + ty * (zstam_global - targZ),
207  m_B_global);
208  stal_global.fieldSlice.GetFieldValue(targX + tx * (zstal_global - targZ),
209  targY + ty * (zstal_global - targZ),
210  l_B_global);
211  stac_global.fieldSlice.GetFieldValue(targX + tx * (zstac_global - targZ),
212  targY + ty * (zstac_global - targZ),
213  centB_global);
214 
215  if (istac != istal)
216  fld0.Set(l_B, stal.z, centB, stac.z, targB, targZ);
217  else
218  fld0.Set(l_B, zstal, targB, targZ);
219  // estimate field for the next extrapolations
220  stam.fieldSlice.GetFieldValue(
221  targX + tx * (zstam - targZ), targY + ty * (zstam - targZ), m_B);
222 #define USE_3HITS
223 #ifndef USE_3HITS
224  if (istac != istal)
225  fld1.Set(m_B, zstam, l_B, zstal, centB, zstac);
226  else
227  fld1.Set(m_B, zstam, l_B, zstal, targB, targZ);
228 #else // if USE_3HITS // the best now
229  L1FieldValue r_B _fvecalignment;
230  L1Station& star = vStations[istam + 1];
231  fvec zstar = star.z;
232  star.fieldSlice.GetFieldValue(
233  targX + tx * (zstar - targZ), targY + ty * (zstar - targZ), r_B);
234  fld1.Set(r_B, zstar, m_B, zstam, l_B, zstal);
235  if ((istam + 1) >= NFStations)
236  fld1.Set(m_B_global,
237  zstam_global,
238  l_B_global,
239  zstal_global,
240  centB_global,
241  zstac_global);
242 #endif // USE_3HITS
243 
244  T.chi2 = 0.;
245  T.NDF = 2.;
246  if ((isec == kAllSecIter) || (isec == kAllSecEIter)
247  || (isec == kAllSecJumpIter))
248  T.NDF = 0;
249  T.tx = tx;
250  T.ty = ty;
251 
252  T.qp = 0.;
253  T.C20 = T.C21 = 0;
254  T.C30 = T.C31 = T.C32 = 0;
255  T.C40 = T.C41 = T.C42 = T.C43 = 0;
256  T.C50 = T.C51 = T.C52 = T.C53 = T.C54 = 0;
257  T.C22 = T.C33 = MaxSlope * MaxSlope / 9.;
258  T.C44 = MaxInvMom / 3. * MaxInvMom / 3.;
259  T.C55 = timeEr * timeEr;
260 
261  T.t = time;
262 
263  // #define BEGIN_FROM_TARGET
264 #ifndef BEGIN_FROM_TARGET // the best now
265 
266  T.x = xl;
267  T.y = yl;
268  T.z = zl;
269  T.C00 = stal.XYInfo.C00;
270  T.C10 = stal.XYInfo.C10;
271  T.C11 = stal.XYInfo.C11;
272 
273  if (fUseHitErrors) {
274  T.C00 = dx1 * dx1;
275  T.C10 = dxy1;
276  T.C11 = dy1 * dy1;
277  }
278 
279 
280  if (fGlobal || fmCBMmode)
281  // add the target
282  {
283 
284  if (istal < NFStations) {
285  fvec eX, eY, J04, J14;
286  fvec dz;
287  dz = targZ - zl;
288  L1ExtrapolateJXY0(T.tx, T.ty, dz, fld0, eX, eY, J04, J14);
289  fvec J[6];
290  J[0] = dz;
291  J[1] = 0;
292  J[2] = J04;
293  J[3] = 0;
294  J[4] = dz;
295  J[5] = J14;
296  L1FilterVtx(T, targX, targY, TargetXYInfo, eX, eY, J);
297  } else {
298  L1ExtrapolateLine(T, targZ);
299  L1FilterXY(T, targX, targY, TargetXYInfo);
300  }
301  }
302 
303  else
304 
305  //add the target
306  {
307  fvec eX, eY, J04, J14;
308  fvec dz;
309  dz = targZ - zl;
310  L1ExtrapolateJXY0(T.tx, T.ty, dz, fld0, eX, eY, J04, J14);
311  fvec J[6];
312  J[0] = dz;
313  J[1] = 0;
314  J[2] = J04;
315  J[3] = 0;
316  J[4] = dz;
317  J[5] = J14;
318  L1FilterVtx(T, targX, targY, TargetXYInfo, eX, eY, J);
319  }
320 
321 
322 #else // TODO: doesn't work. Why?
323 
324  T.x = targX;
325  T.y = targY;
326 
327  T.z = targZ;
328  T.C00 = TargetXYInfo.C00;
329  T.C10 = TargetXYInfo.C10;
330  T.C11 = TargetXYInfo.C11;
331 
332  if (fGlobal || fmCBMmode) // extrapolate to left hit
333  {
334  if (istal <= NFStations)
335  L1Extrapolate0(T, zl, fld0);
336  else
337  L1ExtrapolateLine(T, zl);
338  } else
339  L1Extrapolate0(T, zl, fld0);
340 
341 
342  for (int ista = 0; ista <= istal - 1; ista++) {
343  if ((isec != kAllPrimEIter) && (isec != kAllSecEIter)) {
344 #ifdef USE_RL_TABLE
345  L1AddMaterial(T, fRadThick[ista].GetRadThick(T.x, T.y), MaxInvMom);
346 #else
347  L1AddMaterial(T, vStations[ista].materialInfo, MaxInvMom);
348 #endif
349  if (ista == NMvdStations - 1) L1AddPipeMaterial(T, MaxInvMom);
350  } else {
351 #ifdef USE_RL_TABLE
352  L1AddMaterial(T,
353  fRadThick[ista].GetRadThick(T.x, T.y),
354  MaxInvMom,
355  1,
356  0.000511f * 0.000511f);
357 
358 #else
359  L1AddMaterial(
360  T, vStations[ista].materialInfo, MaxInvMom, 1, 0.000511f * 0.000511f);
361 #endif
362  if (ista == NMvdStations - 1)
363  L1AddPipeMaterial(T, MaxInvMom, 1, 0.000511f * 0.000511f);
364  }
365  }
366 
367  // add left hit
368  L1UMeasurementInfo info = stal.frontInfo;
369 
370  if (fUseHitErrors) info.sigma2 = du1 * du1;
371 
372  if (fGlobal || fmCBMmode) {
373  if (istal < NFStations)
374  L1Filter(T, info, u);
375  else
376  L1FilterNoField(T, info, u);
377  } else
378  L1Filter(T, info, u);
379 
380  info = stal.backInfo;
381 
382  if (fUseHitErrors) info.sigma2 = dv1 * dv1;
383 
384  if (fGlobal || fmCBMmode) {
385  if (istal < NFStations)
386  L1Filter(T, info, v);
387  else
388  L1FilterNoField(T, info, v);
389  } else
390  L1Filter(T, info, v);
391 
392 #endif
393 
394  FilterTime(T, time, timeEr);
395 
396  if ((isec != kAllPrimEIter) && (isec != kAllSecEIter)) {
397 #ifdef USE_RL_TABLE
398  L1AddMaterial(T, fRadThick[istal].GetRadThick(T.x, T.y), MaxInvMom);
399 #else
400  L1AddMaterial(T, stal.materialInfo, MaxInvMom);
401 #endif
402  if ((istam >= NMvdStations) && (istal <= NMvdStations - 1))
403  L1AddPipeMaterial(T, MaxInvMom);
404  } else {
405 #ifdef USE_RL_TABLE
406  L1AddMaterial(T,
407  fRadThick[istal].GetRadThick(T.x, T.y),
408  MaxInvMom,
409  1,
410  0.000511f * 0.000511f);
411 #else
412  L1AddMaterial(T, stal.materialInfo, MaxInvMom, 1, 0.000511f * 0.000511f);
413 #endif
414  if ((istam >= NMvdStations) && (istal <= NMvdStations - 1))
415  L1AddPipeMaterial(T, MaxInvMom, 1, 0.000511f * 0.000511f);
416  }
417  fvec dz = zstam - T.z;
418  L1ExtrapolateTime(T, dz);
419 
420  if (fGlobal || fmCBMmode)
421  // extrapolate to left hit
422  {
423  if (istam < NFStations)
424  L1Extrapolate0(T, zstam, fld0);
425  else
426  L1ExtrapolateLine(T, zstam); // TODO: fld1 doesn't work!
427  } else
428  L1Extrapolate0(T, zstam, fld0); // TODO: fld1 doesn't work!
429 
430  } // i1_V
431 }
432 
433 
435 inline void L1Algo::f20( // input
436  Tindex n1,
437  L1Station& stam,
438 
439  L1HitPoint* vStsHits_m,
440  L1TrackPar* T_1,
441  THitI* hitsl_1,
442 
443  // output
444  Tindex& n2,
445  vector<THitI>& i1_2,
446 #ifdef DOUB_PERFORMANCE
447  vector<THitI>& hitsl_2,
448 #endif // DOUB_PERFORMANCE
449  vector<THitI>& hitsm_2,
450  fvec* Event,
451  vector<bool>& lmDuplets) {
452  n2 = 0; // number of doublet
453  for (Tindex i1 = 0; i1 < n1; ++i1) // for each singlet
454  {
455  const Tindex& i1_V = i1 / fvecLen;
456  const Tindex& i1_4 = i1 % fvecLen;
457  L1TrackPar& T1 = T_1[i1_V];
458 
459  const int n2Saved = n2;
460 
461  const fvec& Pick_m22 =
462  (DOUBLET_CHI2_CUT
463  - T1.chi2); // if make it bigger the found hits will be rejected later because of the chi2 cut.
464  // Pick_m22 is not used, search for mean squared, 2nd version
465 
466  // -- collect possible doublets --
467  const fscal& iz = 1 / T1.z[i1_4];
468  const float& timeError = T1.C55[i1_4];
469  const float& time = T1.t[i1_4];
470 
471  L1HitAreaTime areaTime(
472  vGridTime[&stam - vStations],
473  T1.x[i1_4] * iz,
474  T1.y[i1_4] * iz,
475  (sqrt(Pick_m22 * (T1.C00 + stam.XYInfo.C00)) + MaxDZ * fabs(T1.tx))[i1_4]
476  * iz,
477  (sqrt(Pick_m22 * (T1.C11 + stam.XYInfo.C11)) + MaxDZ * fabs(T1.ty))[i1_4]
478  * iz,
479  time,
480  sqrt(timeError) * 5);
481 
482 
483  THitI imh = 0;
484 
485 
486  while (areaTime.GetNext(imh))
487 
488  {
489  const L1HitPoint& hitm = vStsHits_m[imh];
490 
491 
492  // check y-boundaries
493  if (fabs(time - hitm.time)
494  > sqrt(timeError + hitm.timeEr * hitm.timeEr) * 5)
495  continue;
496  if (fabs(time - hitm.time) > 40) continue;
497 
498 #ifdef USE_EVENT_NUMBER
499  if ((Event[i1_V][i1_4] != hitm.n)) continue;
500 #endif
501  // - check whether hit belong to the window ( track position +\- errors ) -
502  const fscal& zm = hitm.Z();
503  L1TrackPar T1_new = T1;
504  fvec dz = fvec(zm) - T1.z;
505 
506  L1ExtrapolateTime(T1_new, dz);
507 
508  // if (fabs(T1_new.t[i1_4]-hitm.time)>sqrt(T1_new.C55[i1_4]+hitm.timeEr*hitm.timeEr)*4) continue;
509  // if (fabs(T1_new.t[i1_4]-hitm.time)>sqrt(2.9*2.9)*5) continue;
510  // const fscal &dt_est2 = Pick_m22[i1_4] * fabs(T1_new.C55[i1_4] + hitm.timeEr*hitm.timeEr);
511  // const fscal &dt = hitm.time - T1_new.t[i1_4];
512  // if ( dt*dt > dt_est2 && dt < 0 ) continue;
513 
514 
515  fvec y, C11;
516  L1ExtrapolateYC11Line(T1, zm, y, C11);
517 
518  fscal dy_est2 = Pick_m22[i1_4] * fabs(C11[i1_4] + stam.XYInfo.C11[i1_4]);
519 
520  if (fUseHitErrors) {
521  fvec dym = 0;
522  dUdV_to_dY(hitm.dU(), hitm.dV(), dym, stam);
523  dy_est2 = Pick_m22[i1_4] * fabs(C11[i1_4] + dym[0] * dym[0]);
524  }
525 
526  fvec xm, ym = 0;
527 
528  StripsToCoor(hitm.U(), hitm.V(), xm, ym, stam);
529 
530  const fscal& dY = ym[i1_4] - y[i1_4];
531 
532  if (dY * dY > dy_est2 && dY < 0) continue;
533 
534  // check x-boundaries
535  fvec x, C00;
536  L1ExtrapolateXC00Line(T1, zm, x, C00);
537 
538  fscal dx_est2 = Pick_m22[i1_4] * fabs(C00[i1_4] + stam.XYInfo.C00[i1_4]);
539 
540  if (fUseHitErrors) {
541  fvec dxm = 0;
542  dUdV_to_dX(hitm.dU(), hitm.dV(), dxm, stam);
543  dx_est2 = Pick_m22[i1_4] * fabs(C00[i1_4] + dxm[0] * dxm[0]);
544  }
545 
546  const fscal& dX = xm[i1_4] - x[i1_4];
547  if (dX * dX > dx_est2) continue;
548 
549  // check chi2
550  fvec C10;
551  L1ExtrapolateC10Line(T1, zm, C10);
552  fvec chi2 = T1.chi2;
553 
554  L1UMeasurementInfo info = stam.frontInfo;
555 
556  if (fUseHitErrors) info.sigma2 = hitm.dU() * hitm.dU();
557 
558  L1FilterChi2XYC00C10C11(info, x, y, C00, C10, C11, chi2, hitm.U());
559 #ifdef DO_NOT_SELECT_TRIPLETS
560  if (isec != TRACKS_FROM_TRIPLETS_ITERATION)
561 #endif
562  if (chi2[i1_4] > DOUBLET_CHI2_CUT) continue;
563  // T1.t[i1_4] = hitm.time;
564 
565 #ifdef USE_EVENT_NUMBER
566  T1.n[i1_4] = hitm.n;
567 #endif
568 
569  info = stam.backInfo;
570 
571  if (fUseHitErrors) info.sigma2 = hitm.dV() * hitm.dV();
572 
573  L1FilterChi2(info, x, y, C00, C10, C11, chi2, hitm.V());
574 
575  FilterTime(T1_new, hitm.time, hitm.timeEr);
576 
577 #ifdef DO_NOT_SELECT_TRIPLETS
578  if (isec != TRACKS_FROM_TRIPLETS_ITERATION)
579 #endif
580  if (chi2[i1_4] > DOUBLET_CHI2_CUT) continue;
581 
582  i1_2.push_back(i1);
583 #ifdef DOUB_PERFORMANCE
584  hitsl_2.push_back(hitsl_1[i1]);
585 #endif // DOUB_PERFORMANCE
586  hitsm_2.push_back(imh);
587 
588  TripForHit[0][hitsl_1[i1]
589  + StsHitsUnusedStartIndex[&stam - vStations - 1]] = 0;
590  TripForHit[1][hitsl_1[i1]
591  + StsHitsUnusedStartIndex[&stam - vStations - 1]] = 0;
592 
593  TripForHit[0][hitsl_1[i1]
594  + StsHitsUnusedStartIndex[&stam - vStations - 2]] = 0;
595  TripForHit[1][hitsl_1[i1]
596  + StsHitsUnusedStartIndex[&stam - vStations - 2]] = 0;
597 
598  if (n2 > 8000) return;
599 
600  n2++;
601  }
602  lmDuplets[hitsl_1[i1]] = (n2Saved < n2);
603  } // for i1
604 }
605 
606 
609 inline void L1Algo::f30( // input
610  L1HitPoint* vStsHits_r,
611  L1Station& stam,
612  L1Station& star,
613  int istam,
614  int istar,
615  L1HitPoint* vStsHits_m,
616  L1TrackPar* T_1,
617  L1FieldRegion* fld_1,
618  THitI* hitsl_1,
619  Tindex n2,
620  vector<THitI>& hitsm_2,
621  vector<THitI>& i1_2,
622  const vector<bool>& /*mrDuplets*/,
623  // output
624  Tindex& n3,
625  nsL1::vector<L1TrackPar>::TSimd& T_3,
626  vector<THitI>& hitsl_3,
627  vector<THitI>& hitsm_3,
628  vector<THitI>& hitsr_3,
629  nsL1::vector<fvec>::TSimd& u_front_3,
630  nsL1::vector<fvec>::TSimd& u_back_3,
631  nsL1::vector<fvec>::TSimd& z_Pos_3,
632  // nsL1::vector<fvec>::TSimd& dx_,
633  // nsL1::vector<fvec>::TSimd& dy_,
634  nsL1::vector<fvec>::TSimd& dv_,
635  nsL1::vector<fvec>::TSimd& du_,
636  nsL1::vector<fvec>::TSimd& timeR,
637  nsL1::vector<fvec>::TSimd& timeER) {
638  THitI hitsl_2[fvecLen];
639  THitI hitsm_2_tmp[fvecLen];
640  fvec fvec_0;
641  L1TrackPar L1TrackPar_0;
642 
643  Tindex n3_V = 0, n3_4 = 0;
644 
645  T_3.push_back(L1TrackPar_0);
646  u_front_3.push_back(fvec_0);
647  u_back_3.push_back(fvec_0);
648  z_Pos_3.push_back(fvec_0);
649  // dx_.push_back(fvec_0);
650  // dy_.push_back(fvec_0);
651  du_.push_back(fvec_0);
652  dv_.push_back(fvec_0);
653  timeR.push_back(fvec_0);
654  timeER.push_back(fvec_0);
655 
656 
657  L1TrackPar T2;
658  L1FieldRegion f2;
659  // pack the data
660  fvec u_front_2;
661  fvec u_back_2;
662  fvec dx2;
663  fvec dy2;
664  fvec du2;
665  fvec dv2;
666  fvec zPos_2;
667  fvec timeM;
668  fvec timeMEr;
669  fvec num;
670 
671 
672  // ---- Add the middle hits to parameters estimation. Propagate to right station. ----
673  if (istar < NStations) {
674  for (Tindex i2 = 0; i2 < n2;) {
675  Tindex n2_4 = 0;
676  for (; n2_4 < fvecLen && i2 < n2; n2_4++, i2++) {
677  // if (!mrDuplets[hitsm_2[i2]]) {
678  // n2_4--;
679  // continue;
680  // }
681  const Tindex& i1 = i1_2[i2];
682  const Tindex i1_V = i1 / fvecLen;
683  const Tindex i1_4 = i1 % fvecLen;
684 
685 
686  const L1TrackPar& T1 = T_1[i1_V];
687  const L1FieldRegion& f1 = fld_1[i1_V];
688  T2.SetOneEntry(n2_4, T1, i1_4);
689  f2.SetOneEntry(n2_4, f1, i1_4);
690 
691  const Tindex& imh = hitsm_2[i2];
692  const L1HitPoint& hitm = vStsHits_m[imh];
693  u_front_2[n2_4] = hitm.U();
694  u_back_2[n2_4] = hitm.V();
695  zPos_2[n2_4] = hitm.Z();
696  timeM[n2_4] = hitm.time;
697  timeMEr[n2_4] = hitm.timeEr;
698  // num[n2_4] = hitm.track;
699  du2[n2_4] = hitm.dU();
700  dv2[n2_4] = hitm.dV();
701 
702  hitsl_2[n2_4] = hitsl_1[i1];
703  hitsm_2_tmp[n2_4] = hitsm_2[i2];
704  } // n2_4
705 
706  dUdV_to_dX(du2, dv2, dx2, stam);
707  dUdV_to_dY(du2, dv2, dy2, stam);
708 
709  fvec dz = zPos_2 - T2.z;
710 
711  L1ExtrapolateTime(T2, dz);
712  // add middle hit
713  L1ExtrapolateLine(T2, zPos_2);
714 
715  L1UMeasurementInfo info = stam.frontInfo;
716 
717  if (fUseHitErrors) info.sigma2 = du2 * du2;
718 
719  if (fGlobal || fmCBMmode) {
720  if (istam < NFStations)
721  L1Filter(T2, info, u_front_2);
722  else
723  L1FilterNoField(T2, info, u_front_2);
724  } else
725  L1Filter(T2, info, u_front_2);
726 
727  info = stam.backInfo;
728  if (fUseHitErrors) info.sigma2 = dv2 * dv2;
729 
730  if (fGlobal || fmCBMmode) {
731  if (istam < NFStations)
732  L1Filter(T2, info, u_back_2);
733  else
734  L1FilterNoField(T2, info, u_back_2);
735  } else
736  L1Filter(T2, info, u_back_2);
737 
738  FilterTime(T2, timeM, timeMEr);
739 
740  if ((isec != kAllPrimEIter) && (isec != kAllSecEIter)) {
741 #ifdef USE_RL_TABLE
742  L1AddMaterial(T2, fRadThick[istam].GetRadThick(T2.x, T2.y), T2.qp);
743 #else
744  L1AddMaterial(T2, stam.materialInfo, T2.qp);
745 #endif
746  if ((istar >= NMvdStations) && (istam <= NMvdStations - 1))
747  L1AddPipeMaterial(T2, T2.qp);
748  } else {
749 #ifdef USE_RL_TABLE
750  L1AddMaterial(T2,
751  fRadThick[istam].GetRadThick(T2.x, T2.y),
752  T2.qp,
753  1,
754  0.000511f * 0.000511f);
755 #else
756  L1AddMaterial(T2, stam.materialInfo, T2.qp, 1, 0.000511f * 0.000511f);
757 #endif
758  if ((istar >= NMvdStations) && (istam <= NMvdStations - 1))
759  L1AddPipeMaterial(T2, T2.qp, 1, 0.000511f * 0.000511f);
760  }
761 
762  fvec dz2 = star.z - T2.z;
763  L1ExtrapolateTime(T2, dz2);
764  // extrapolate to the right hit station
765 
766  if (fGlobal || fmCBMmode)
767  // extrapolate to the right hit station
768  {
769  if (istar <= NFStations)
770  L1Extrapolate(T2, star.z, T2.qp, f2);
771  else
772  L1ExtrapolateLine(T2, star.z);
773  } else
774  L1Extrapolate(T2, star.z, T2.qp, f2);
775 
776  // ---- Find the triplets(right hit). Reformat data in the portion of triplets. ----
777  for (Tindex i2_4 = 0; i2_4 < n2_4; ++i2_4) {
778  // if (T2.C00[i2_4] < 0 || T2.C11[i2_4] < 0 || T2.C22[i2_4] < 0
779  // || T2.C33[i2_4] < 0 || T2.C44[i2_4] < 0 || T2.C55[i2_4] < 0)
780  // continue;
781  if (T2.C00[i2_4] < 0 || T2.C11[i2_4] < 0 || T2.C22[i2_4] < 0
782  || T2.C33[i2_4] < 0)
783  continue;
784 
785  const fvec& Pick_r22 = (TRIPLET_CHI2_CUT - T2.chi2);
786  const float& timeError = T2.C55[i2_4];
787  const float& time = T2.t[i2_4];
788 
789  // find first possible hit
790 
791 #ifdef DO_NOT_SELECT_TRIPLETS
792  if (isec == TRACKS_FROM_TRIPLETS_ITERATION) Pick_r22 = Pick_r2 + 1;
793 #endif
794  const fscal& iz = 1 / T2.z[i2_4];
795  L1HitAreaTime area(vGridTime[&star - vStations],
796  T2.x[i2_4] * iz,
797  T2.y[i2_4] * iz,
798  (sqrt(Pick_r22 * (T2.C00 + stam.XYInfo.C00))
799  + MaxDZ * fabs(T2.tx))[i2_4]
800  * iz,
801  (sqrt(Pick_r22 * (T2.C11 + stam.XYInfo.C11))
802  + MaxDZ * fabs(T2.ty))[i2_4]
803  * iz,
804  time,
805  sqrt(timeError) * 5);
806 
807  THitI irh = 0;
808 
809  while (area.GetNext(irh)) {
810  const L1HitPoint& hitr = vStsHits_r[irh];
811 
812 #ifdef USE_EVENT_NUMBER
813  if ((T2.n[i2_4] != hitr.n)) continue;
814 #endif
815  const fscal& zr = hitr.Z();
816  // const fscal& yr = hitr.Y();
817 
818  fvec xr, yr = 0;
819 
820  StripsToCoor(hitr.U(), hitr.V(), xr, yr, star);
821 
822  fvec dz3 = zr - T2.z;
823  L1ExtrapolateTime(T2, dz3);
824 
825 
826  if (fabs(T2.t[i2_4] - hitr.time)
827  > sqrt(T2.C55[i2_4] + hitr.timeEr) * 5)
828  continue;
829  if (fabs(T2.t[i2_4] - hitr.time) > 40) continue;
830  //
831  // if (fabs(T2.t[i2_4]-hitr.time)>sqrt(2.9*2.9)*5) continue;
832 
833  // - check whether hit belong to the window ( track position +\- errors ) -
834  // check lower boundary
835  fvec y, C11;
836  L1ExtrapolateYC11Line(T2, zr, y, C11);
837 
838  fscal dy_est2 =
839  (Pick_r22[i2_4]
840  * (fabs(
841  C11[i2_4]
842  + star.XYInfo.C11
843  [i2_4]))); // TODO for FastPrim dx < dy - other sort is optimal. But not for doublets
844 
845  if (fUseHitErrors) {
846  fvec dyr = 0;
847  dUdV_to_dY(hitr.dU(), hitr.dV(), dyr, star);
848  dy_est2 = (Pick_r22[i2_4] * (fabs(C11[i2_4] + dyr[0] * dyr[0])));
849  }
850 
851  const fscal& dY = yr[i2_4] - y[i2_4];
852  const fscal& dY2 = dY * dY;
853  if (dY2 > dy_est2 && dY2 < 0)
854  continue; // if (yr < y_minus_new[i2_4]) continue;
855 
856  // check upper boundary
857  if (dY2 > dy_est2)
858  continue; // if (yr > y_plus_new [i2_4] ) continue;
859  // check x-boundaries
860  fvec x, C00;
861  L1ExtrapolateXC00Line(T2, zr, x, C00);
862 
863  fscal dx_est2 =
864  (Pick_r22[i2_4] * (fabs(C00[i2_4] + star.XYInfo.C00[i2_4])));
865 
866  if (fUseHitErrors) {
867  fvec dxr = 0;
868  dUdV_to_dX(hitr.dU(), hitr.dV(), dxr, star);
869  dx_est2 = (Pick_r22[i2_4] * (fabs(C00[i2_4] + dxr[0] * dxr[0])));
870  }
871 
872  const fscal& dX = xr[i2_4] - x[i2_4];
873  if (dX * dX > dx_est2) continue;
874  // check chi2 // not effective
875  fvec C10;
876  L1ExtrapolateC10Line(T2, zr, C10);
877  fvec chi2 = T2.chi2;
878 
879  info = star.frontInfo;
880 
881  if (fUseHitErrors) info.sigma2 = hitr.dU() * hitr.dU();
882 
883  L1FilterChi2XYC00C10C11(info, x, y, C00, C10, C11, chi2, hitr.U());
884  info = star.backInfo;
885 
886  if (fUseHitErrors) info.sigma2 = hitr.dV() * hitr.dV();
887 
888  L1FilterChi2(info, x, y, C00, C10, C11, chi2, hitr.V());
889 
890  L1TrackPar T = T2;
891 
892  FilterTime(T, hitr.time, hitr.timeEr);
893 #ifdef DO_NOT_SELECT_TRIPLETS
894  if (isec != TRACKS_FROM_TRIPLETS_ITERATION)
895 #endif
896 
897  if (fGlobal || fmCBMmode)
898  if (chi2[i2_4] > TRIPLET_CHI2_CUT || C00[i2_4] < 0
899  || C11[i2_4] < 0)
900  continue;
901  else if (chi2[i2_4] > TRIPLET_CHI2_CUT || C00[i2_4] < 0
902  || C11[i2_4] < 0 || T.C55[i2_4] < 0)
903  continue; // chi2_triplet < CHI2_CUT
904 
905 
906  // pack triplet
907  L1TrackPar& T3 = T_3[n3_V];
908 
909  hitsl_3.push_back(hitsl_2[i2_4]);
910  hitsm_3.push_back(hitsm_2_tmp[i2_4]);
911  hitsr_3.push_back(irh);
912 
913 
914  T3.SetOneEntry(n3_4, T2, i2_4);
915  u_front_3[n3_V][n3_4] = hitr.U();
916  u_back_3[n3_V][n3_4] = hitr.V();
917  //dx_[n3_V][n3_4] = hitr.dX();
918  // dy_[n3_V][n3_4] = hitr.dY();
919  du_[n3_V][n3_4] = hitr.dU();
920  dv_[n3_V][n3_4] = hitr.dV();
921  z_Pos_3[n3_V][n3_4] = zr;
922  timeR[n3_V][n3_4] = hitr.time;
923  timeER[n3_V][n3_4] = hitr.timeEr;
924 
925  n3++;
926  n3_V = n3 / fvecLen;
927  n3_4 = n3 % fvecLen;
928 
929  if (!n3_4) {
930  T_3.push_back(L1TrackPar_0);
931  u_front_3.push_back(fvec_0);
932  u_back_3.push_back(fvec_0);
933  z_Pos_3.push_back(fvec_0);
934  // dx_.push_back(fvec_0);
935  // dy_.push_back(fvec_0);
936  du_.push_back(fvec_0);
937  dv_.push_back(fvec_0);
938  timeR.push_back(fvec_0);
939  timeER.push_back(fvec_0);
940  }
941 
942  if (n3 > 4000) return;
943  }
944  } // i2_4
945  } // i2_V
946  } // if istar
947 }
948 
950 inline void L1Algo::f31( // input
951  Tindex n3_V,
952  L1Station& star,
953  nsL1::vector<fvec>::TSimd& u_front_,
954  nsL1::vector<fvec>::TSimd& u_back_,
955  nsL1::vector<fvec>::TSimd& z_Pos,
956  // nsL1::vector<fvec>::TSimd& dx_,
957  // nsL1::vector<fvec>::TSimd& dy_,
958  nsL1::vector<fvec>::TSimd& dv_,
959  nsL1::vector<fvec>::TSimd& du_,
960  nsL1::vector<fvec>::TSimd& timeR,
961  nsL1::vector<fvec>::TSimd& timeER,
962  // output
963  // L1TrackPar *T_3
964  nsL1::vector<L1TrackPar>::TSimd& T_3) {
965  for (Tindex i3_V = 0; i3_V < n3_V; ++i3_V) {
966  fvec dz = z_Pos[i3_V] - T_3[i3_V].z;
967 
968 
969  L1ExtrapolateTime(T_3[i3_V], dz);
970  L1ExtrapolateLine(T_3[i3_V], z_Pos[i3_V]);
971 
972  L1UMeasurementInfo info = star.frontInfo;
973 
974  if (fUseHitErrors) info.sigma2 = du_[i3_V] * du_[i3_V];
975 
976  if (fGlobal || fmCBMmode) {
977  if ((&star - vStations) < NFStations)
978  L1Filter(T_3[i3_V], info, u_front_[i3_V]);
979  else {
980  L1FilterNoField(T_3[i3_V], info, u_front_[i3_V]);
981  }
982  } else
983  L1Filter(T_3[i3_V], info, u_front_[i3_V]);
984 
985 
986  info = star.backInfo;
987 
988  if (fUseHitErrors) info.sigma2 = dv_[i3_V] * dv_[i3_V];
989 
990  if (fGlobal || fmCBMmode) {
991  if ((&star - vStations) < NFStations)
992  L1Filter(T_3[i3_V], info, u_back_[i3_V]);
993  else
994  L1FilterNoField(T_3[i3_V], info, u_back_[i3_V]);
995  } else
996  L1Filter(T_3[i3_V], info, u_back_[i3_V]);
997 
998  FilterTime(T_3[i3_V], timeR[i3_V], timeER[i3_V]);
999  }
1000 }
1001 
1002 
1004 inline void L1Algo::f32( // input // TODO not updated after gaps introduction
1005  Tindex n3,
1006  int istal,
1007  nsL1::vector<L1TrackPar>::TSimd& T_3,
1008  vector<THitI>& hitsl_3,
1009  vector<THitI>& hitsm_3,
1010  vector<THitI>& hitsr_3,
1011  int nIterations) {
1012  const int NHits = 3; // triplets
1013 
1014  // prepare data
1015  int ista[NHits] = {istal, istal + 1, istal + 2};
1016 
1017  L1Station sta[3];
1018  for (int is = 0; is < NHits; ++is) {
1019  sta[is] = vStations[ista[is]];
1020  };
1021 
1022  for (int i3 = 0; i3 < n3; ++i3) {
1023  int i3_V = i3 / fvecLen;
1024  int i3_4 = i3 % fvecLen;
1025 
1026  L1TrackPar& T3 = T_3[i3_V];
1027  fscal& qp0 = T3.qp[i3_4];
1028 
1029  // prepare data
1030  THitI ihit[NHits] = {
1031  (*RealIHitP)[hitsl_3[i3] + StsHitsUnusedStartIndex[ista[0]]],
1032  (*RealIHitP)[hitsm_3[i3] + StsHitsUnusedStartIndex[ista[1]]],
1033  (*RealIHitP)[hitsr_3[i3] + StsHitsUnusedStartIndex[ista[2]]]};
1034 
1035  fscal u[NHits], v[NHits], x[NHits], y[NHits], z[NHits];
1036  for (int ih = 0; ih < NHits; ++ih) {
1037  const L1StsHit& hit = (*vStsHits)[ihit[ih]];
1038  u[ih] = (*vStsStrips)[hit.f];
1039  v[ih] = (*vStsStripsB)[hit.b];
1040  StripsToCoor(u[ih], v[ih], x[ih], y[ih], sta[ih]);
1041  z[ih] = (*vStsZPos)[hit.iz];
1042  };
1043 
1044  // initialize parameters
1045  L1TrackPar T;
1046 
1047  const fvec vINF = .1;
1048  T.x = x[0];
1049  T.y = y[0];
1050 
1051  fvec dz01 = 1. / (z[1] - z[0]);
1052  T.tx = (x[1] - x[0]) * dz01;
1053  T.ty = (y[1] - y[0]) * dz01;
1054 
1055  T.qp = qp0;
1056  T.z = z[0];
1057  T.chi2 = 0.;
1058  T.NDF = 2.;
1059  T.C00 = sta[0].XYInfo.C00;
1060  T.C10 = sta[0].XYInfo.C10;
1061  T.C11 = sta[0].XYInfo.C11;
1062 
1063  T.C20 = T.C21 = 0;
1064  T.C30 = T.C31 = T.C32 = 0;
1065  T.C40 = T.C41 = T.C42 = T.C43 = 0;
1066  T.C22 = T.C33 = vINF;
1067  T.C44 = 1.;
1068 
1069  // find field along the track
1070  L1FieldValue B[3] _fvecalignment;
1071  L1FieldRegion fld _fvecalignment;
1072 
1073  fvec tx[3] = {(x[1] - x[0]) / (z[1] - z[0]),
1074  (x[2] - x[0]) / (z[2] - z[0]),
1075  (x[2] - x[1]) / (z[2] - z[1])};
1076  fvec ty[3] = {(y[1] - y[0]) / (z[1] - z[0]),
1077  (y[2] - y[0]) / (z[2] - z[0]),
1078  (y[2] - y[1]) / (z[2] - z[1])};
1079  for (int ih = 0; ih < NHits; ++ih) {
1080  fvec dz = (sta[ih].z - z[ih]);
1081  sta[ih].fieldSlice.GetFieldValue(
1082  x[ih] + tx[ih] * dz, y[ih] + ty[ih] * dz, B[ih]);
1083  };
1084  fld.Set(B[0], sta[0].z, B[1], sta[1].z, B[2], sta[2].z);
1085 
1086  // fit
1087  for (int ih = 1; ih < NHits; ++ih) {
1088  L1Extrapolate(T, z[ih], T.qp, fld);
1089 #ifdef USE_RL_TABLE
1090  L1AddMaterial(T, fRadThick[ista[ih]].GetRadThick(T.x, T.y), T.qp);
1091 #else
1092  L1AddMaterial(T, sta[ih].materialInfo, T.qp);
1093 #endif
1094  if (ista[ih] == NMvdStations - 1) L1AddPipeMaterial(T, T.qp);
1095 
1096  L1Filter(T, sta[ih].frontInfo, u[ih]);
1097  L1Filter(T, sta[ih].backInfo, v[ih]);
1098  }
1099 
1100  // repeat several times in order to improve precision
1101  for (int iiter = 0; iiter < nIterations; ++iiter) {
1102  // fit backward
1103  // keep tx,ty,q/p
1104  int ih = NHits - 1;
1105  T.x = x[ih];
1106  T.y = y[ih];
1107  T.z = z[ih];
1108  T.chi2 = 0.;
1109  T.NDF = 2.;
1110  T.C00 = sta[ih].XYInfo.C00;
1111  T.C10 = sta[ih].XYInfo.C10;
1112  T.C11 = sta[ih].XYInfo.C11;
1113 
1114  T.C20 = T.C21 = 0;
1115  T.C30 = T.C31 = T.C32 = 0;
1116  T.C40 = T.C41 = T.C42 = T.C43 = 0;
1117  T.C22 = T.C33 = vINF;
1118  T.C44 = 1.;
1119 
1120  // L1Filter( T, sta[ih].frontInfo, u[ih] );
1121  // L1Filter( T, sta[ih].backInfo, v[ih] );
1122  for (ih = NHits - 2; ih >= 0; ih--) {
1123  L1Extrapolate(T, z[ih], T.qp, fld);
1124 #ifdef USE_RL_TABLE
1125  L1AddMaterial(T, fRadThick[ista[ih]].GetRadThick(T.x, T.y), T.qp);
1126 #else
1127  L1AddMaterial(T, sta[ih].materialInfo, T.qp);
1128 #endif
1129  if (ista[ih] == NMvdStations) L1AddPipeMaterial(T, T.qp);
1130 
1131  L1Filter(T, sta[ih].frontInfo, u[ih]);
1132  L1Filter(T, sta[ih].backInfo, v[ih]);
1133  }
1134  // fit forward
1135  ih = 0;
1136  T.x = x[ih];
1137  T.y = y[ih];
1138  T.z = z[ih];
1139  T.chi2 = 0.;
1140  T.NDF = 2.;
1141  T.C00 = sta[ih].XYInfo.C00;
1142  T.C10 = sta[ih].XYInfo.C10;
1143  T.C11 = sta[ih].XYInfo.C11;
1144 
1145  T.C20 = T.C21 = 0;
1146  T.C30 = T.C31 = T.C32 = 0;
1147  T.C40 = T.C41 = T.C42 = T.C43 = 0;
1148  T.C22 = T.C33 = vINF;
1149  T.C44 = 1.;
1150 
1151  // L1Filter( T, sta[ih].frontInfo, u[ih] );
1152  // L1Filter( T, sta[ih].backInfo, v[ih] );
1153  for (ih = 1; ih < NHits; ++ih) {
1154  L1Extrapolate(T, z[ih], T.qp, fld);
1155 #ifdef USE_RL_TABLE
1156  L1AddMaterial(T, fRadThick[ista[ih]].GetRadThick(T.x, T.y), T.qp);
1157 #else
1158  L1AddMaterial(T, sta[ih].materialInfo, T.qp);
1159 #endif
1160  if (ista[ih] == NMvdStations + 1) L1AddPipeMaterial(T, T.qp);
1161 
1162  L1Filter(T, sta[ih].frontInfo, u[ih]);
1163  L1Filter(T, sta[ih].backInfo, v[ih]);
1164  }
1165  } // for iiter
1166 
1167  T3.SetOneEntry(i3_4, T, i3_4);
1168  } //i3
1169 } // f32
1170 
1171 
1173 inline void L1Algo::f4( // input
1174  Tindex n3,
1175  int istal,
1176  int istam,
1177  int istar,
1178  nsL1::vector<L1TrackPar>::TSimd& T_3,
1179  vector<THitI>& hitsl_3,
1180  vector<THitI>& hitsm_3,
1181  vector<THitI>& hitsr_3,
1182  // output
1183  Tindex& nstaltriplets,
1184  vector<THitI>* /*hitsn_3*/,
1185  vector<THitI>* /*hitsr_5*/
1186 ) {
1187  THitI ihitl_priv = 0;
1188 
1189  unsigned int Station = 0;
1190  unsigned int Thread = 0;
1191  unsigned int Triplet = 0;
1192 
1193  unsigned int Location = 0;
1194 
1195  unsigned char level = 0;
1196 
1197 
1198  for (Tindex i3 = 0; i3 < n3; ++i3) {
1199  const Tindex& i3_V = i3 / fvecLen;
1200  const Tindex& i3_4 = i3 % fvecLen;
1201 
1202 
1203  L1TrackPar& T3 = T_3[i3_V];
1204 
1205  // select
1206  fscal& chi2 = T3.chi2[i3_4];
1207 
1208 
1209  Station = istal;
1210 
1211 #ifdef _OPENMP
1212  Thread = omp_get_thread_num();
1213 #else
1214  Thread = 0;
1215 #endif
1216 
1217 
1218  TripletsLocal1[Station][Thread][nTripletsThread[istal][Thread]].SetLevel(0);
1219 
1220  TripletsLocal1[Station][Thread][nTripletsThread[istal][Thread]]
1221  .SetFNeighbour(0);
1222  TripletsLocal1[Station][Thread][nTripletsThread[istal][Thread]]
1223  .SetNNeighbours(0);
1224 
1225  Triplet = nTripletsThread[istal][Thread];
1226 
1227  Location = Triplet + Station * 100000000 + Thread * 1000000;
1228 
1229  if (ihitl_priv == 0 || ihitl_priv != hitsl_3[i3]) {
1230  TripForHit[0][hitsl_3[i3] + StsHitsUnusedStartIndex[istal]] = Location;
1231 
1232  TripForHit[1][hitsl_3[i3] + StsHitsUnusedStartIndex[istal]] = Location;
1233  }
1234 
1235  ihitl_priv = hitsl_3[i3];
1236 
1237 
1238 #ifdef DO_NOT_SELECT_TRIPLETS
1239  if (isec != TRACKS_FROM_TRIPLETS_ITERATION)
1240 #endif
1241  if (!finite(chi2) || chi2 < 0 || chi2 > TRIPLET_CHI2_CUT) continue;
1242 
1243 
1244  // prepare data
1245  fscal MaxInvMomS = MaxInvMom[0];
1246  fscal qp = MaxInvMomS + T3.qp[i3_4];
1247  if (qp < 0) qp = 0;
1248  if (qp > MaxInvMomS * 2) qp = MaxInvMomS * 2;
1249  fscal Cqp = 5. * sqrt(fabs(T3.C44[i3_4]));
1250 
1251 
1252  fscal scale = 255 / (MaxInvMom[0] * 2);
1253 
1254  qp = (static_cast<unsigned int>(qp * scale)) % 256;
1255  Cqp = (static_cast<unsigned int>(Cqp * scale)) % 256;
1256  Cqp += 1;
1257 
1258  if (Cqp < 0) Cqp = 0;
1259  if (Cqp > 20) Cqp = 20;
1260  qp = static_cast<unsigned char>(qp);
1261 
1262 
1263  const THitI& ihitl = hitsl_3[i3] + StsHitsUnusedStartIndex[istal];
1264  const THitI& ihitm = hitsm_3[i3] + StsHitsUnusedStartIndex[istam];
1265  const THitI& ihitr = hitsr_3[i3] + StsHitsUnusedStartIndex[istar];
1266  L1_ASSERT(ihitl < StsHitsUnusedStopIndex[istal],
1267  ihitl << " < " << StsHitsUnusedStopIndex[istal]);
1268  L1_ASSERT(ihitm < StsHitsUnusedStopIndex[istam],
1269  ihitm << " < " << StsHitsUnusedStopIndex[istam]);
1270  L1_ASSERT(ihitr < StsHitsUnusedStopIndex[istar],
1271  ihitr << " < " << StsHitsUnusedStopIndex[istar]);
1272 
1273  fscal& time = T3.t[i3_4];
1274  // int n = T3.n[i3_4];
1275 
1276 
1277  L1Triplet& tr1 =
1278  TripletsLocal1[Station][Thread][nTripletsThread[istal][Thread]];
1279 
1280 
1281  tr1.SetLevel(0);
1282 
1283 
1284  tr1.Set(
1285  ihitl, ihitm, ihitr, istal, istam, istar, 0, qp, chi2, time, Cqp, 0);
1286 
1287  tr1.tx = sqrt(fabs(T3.tx[i3_4]));
1288  tr1.ty = sqrt(fabs(T3.ty[i3_4]));
1289  tr1.Ctx = sqrt(fabs(T3.C22[i3_4]));
1290  tr1.Cty = sqrt(fabs(T3.C33[i3_4]));
1291 
1292 
1293  ++nstaltriplets;
1294 
1295 
1296  nTripletsThread[istal][Thread]++;
1297 
1298  Triplet = nTripletsThread[istal][Thread];
1299 
1300  Location = Triplet + Station * 100000000 + Thread * 1000000;
1301 
1302 
1303  TripForHit[1][hitsl_3[i3] + StsHitsUnusedStartIndex[istal]] = Location;
1304 
1305 
1306  if (istal > (NStations - 4)) continue;
1307 
1308  unsigned int Neighbours = TripForHit[1][ihitm] - TripForHit[0][ihitm];
1309 
1310 
1311  level = 0;
1312 
1313  for (unsigned int iN = 0; iN < Neighbours; ++iN) {
1314  Location = TripForHit[0][ihitm] + iN;
1315 
1316 
1317  Station = Location / 100000000;
1318  Thread = (Location - Station * 100000000) / 1000000;
1319  Triplet = (Location - Station * 100000000 - Thread * 1000000);
1320 
1321  L1Triplet& curNeighbour = TripletsLocal1[Station][Thread][Triplet];
1322 
1323  if ((curNeighbour.GetMHit() != ihitr)) continue;
1324 
1325  if (tr1.GetFNeighbour() == 0) tr1.SetFNeighbour(Location);
1326 
1327  tr1.SetNNeighbours(Location - tr1.GetFNeighbour() + 1);
1328 
1329  const unsigned char& jlevel = curNeighbour.GetLevel();
1330 
1331  if (level <= jlevel) level = jlevel + 1;
1332  }
1333  tr1.SetLevel(level);
1334  }
1335 }
1336 
1338 inline void L1Algo::f5( // input
1339  // output
1340 
1341  int* nlevel) {
1342 #ifdef TRACKS_FROM_TRIPLETS
1343  if (isec != TRACKS_FROM_TRIPLETS_ITERATION)
1344 #endif
1345 
1346  for (int istal = NStations - 4; istal >= FIRSTCASTATION; istal--) {
1347  for (
1348  int tripType = 0; tripType < 3;
1349  tripType++) // tT = 0 - 123triplet, tT = 1 - 124triplet, tT = 2 - 134triplet
1350  {
1351  if ((((isec != kFastPrimJumpIter) && (isec != kAllPrimJumpIter)
1352  && (isec != kAllSecJumpIter))
1353  && (tripType != 0))
1354  || (((isec == kFastPrimJumpIter) || (isec == kAllPrimJumpIter)
1355  || (isec == kAllSecJumpIter))
1356  && (tripType != 0) && (istal == NStations - 4)))
1357  continue;
1358 
1359  int istam = istal + 1;
1360  int istar = istal + 2;
1361  switch (tripType) {
1362  case 1: istar++; break;
1363  case 2:
1364  istam++;
1365  istar++;
1366  break;
1367  }
1368 
1369  for (Tindex ip = 0; ip < fNThreads; ++ip) {
1370  for (Tindex itrip = 0; itrip < nTripletsThread[istal][ip]; ++itrip) {
1371  L1Triplet* trip = &(TripletsLocal1[istal][ip][itrip]);
1372  if (istam != trip->GetMSta()) continue;
1373  if (istar != trip->GetRSta()) continue;
1374 
1375  unsigned char level = 0;
1376  // float chi2 = trip->GetChi2();
1377  unsigned char qp = trip->GetQp();
1378 
1379  THitI ihitl = trip->GetLHit();
1380  THitI ihitm = trip->GetMHit();
1381  THitI ihitr = trip->GetRHit();
1382  L1_ASSERT(ihitl < StsHitsUnusedStopIndex[istal],
1383  ihitl << " < " << StsHitsUnusedStopIndex[istal]);
1384  L1_ASSERT(ihitm < StsHitsUnusedStopIndex[istam],
1385  ihitm << " < " << StsHitsUnusedStopIndex[istam]);
1386  L1_ASSERT(ihitr < StsHitsUnusedStopIndex[istar],
1387  ihitr << " < " << StsHitsUnusedStopIndex[istar]);
1388 
1389  vector<unsigned int> neighCands; // save neighbour candidates
1390  neighCands.reserve(8); // ~average is 1-2 for central, up to 5
1391 
1392  unsigned int Neighbours =
1393  TripForHit[1][ihitm] - TripForHit[0][ihitm];
1394 
1395  for (unsigned int iN = 0; iN < Neighbours; ++iN) {
1396  // for (iN = first_triplet; iN <= last_triplet; ++iN){
1397  int Location = TripForHit[0][ihitm] + iN;
1398 
1399 
1400  int Station = Location / 100000000;
1401  int Thread = (Location - Station * 100000000) / 1000000;
1402  int Triplet = (Location - Station * 100000000 - Thread * 1000000);
1403 
1404  L1Triplet& triplet = TripletsLocal1[Station][Thread][Triplet];
1405 
1406  // if (triplet.GetMSta() != istar) continue; // neighbours should have 2 common hits
1407  // if (triplet.GetMHit() != ihitr) continue; //!!!
1408 
1409  L1Triplet* tripn = &triplet;
1410 
1411  fscal qp2 = tripn->GetQp();
1412  fscal Cqp1 = trip->Cqp;
1413  fscal Cqp2 = tripn->Cqp;
1414  if (fabs(qp - qp2) > PickNeighbour * (Cqp1 + Cqp2))
1415  continue; // neighbours should have same qp
1416 
1417  // calculate level
1418  unsigned char jlevel = tripn->GetLevel();
1419  if (level <= jlevel) level = jlevel + 1;
1420  if (level == jlevel + 1) neighCands.push_back(Location);
1421  }
1422 
1423  // trip->neighbours.resize(0);
1424 
1425  // for (unsigned int in = 0; in < neighCands.size(); in++)
1426  // {
1427  // int Location = neighCands[in];
1428  //
1429  // int Station = Location/100000000;
1430  // int Thread = (Location -Station*100000000)/1000000;
1431  // int Triplet = (Location- Station*100000000-Thread*1000000);
1432 
1433  // const int nLevel = TripletsLocal1[Station][Thread][Triplet].GetLevel();
1434  // if (level == nLevel + 1) trip->neighbours.push_back(Location);
1435  // }
1436  nlevel[level]++;
1437  } // vTriplets
1438  }
1439  } // tripType
1440  } // istal
1441 }
1442 
1444 
1445 inline void L1Algo::
1446  DupletsStaPort(
1447  int istal,
1448  int istam,
1449  Tindex ip,
1450  vector<Tindex>& n_g,
1451  Tindex* portionStopIndex_,
1453  L1TrackPar*
1454  T_1,
1455  L1FieldRegion*
1456  fld_1,
1457  THitI* hitsl_1,
1458  vector<bool>& lmDuplets,
1459  Tindex& n_2,
1460  vector<THitI>& i1_2,
1461  vector<THitI>& hitsm_2) {
1462  if (istam < NStations) {
1463  L1Station& stam = vStations[istam];
1464 
1465  // prepare data
1466  L1HitPoint* vStsHits_l =
1467  &((*vStsHitPointsUnused)[0]) + StsHitsUnusedStartIndex[istal];
1468  L1HitPoint* vStsHits_m =
1469  &((*vStsHitPointsUnused)[0]) + StsHitsUnusedStartIndex[istam];
1470 
1471  fvec u_front[Portion / fvecLen], u_back[Portion / fvecLen];
1472  fvec dx0[Portion / fvecLen], dy0[Portion / fvecLen],
1473  dxy0[Portion / fvecLen];
1474  fvec dv0[Portion / fvecLen], du0[Portion / fvecLen];
1475  fvec zPos[Portion / fvecLen];
1476  fvec HitTime[Portion / fvecLen];
1477  fvec HitTimeEr[Portion / fvecLen];
1478  fvec Event[Portion / fvecLen];
1479 
1481  Tindex& n1 = n_g[ip];
1482 
1483  f10( // input
1484  (ip - portionStopIndex_[istal + 1]) * Portion,
1485  n1,
1486  vStsHits_l,
1487  // output
1488  u_front,
1489  u_back,
1490  zPos,
1491  hitsl_1,
1492  HitTime,
1493  HitTimeEr,
1494  Event,
1495  dx0,
1496  dy0,
1497  dxy0,
1498  du0,
1499  dv0);
1500 
1501  for (Tindex i = 0; i < n1; ++i)
1502  L1_ASSERT(hitsl_1[i] < StsHitsUnusedStopIndex[istal]
1503  - StsHitsUnusedStartIndex[istal],
1504  hitsl_1[i] << " < "
1505  << StsHitsUnusedStopIndex[istal]
1506  - StsHitsUnusedStartIndex[istal]);
1507 
1508  Tindex n1_V = (n1 + fvecLen - 1) / fvecLen;
1509 
1511 
1512  f11(istal,
1513  istam,
1514  n1_V,
1515 
1516  u_front,
1517  u_back,
1518  zPos,
1519  HitTime,
1520  HitTimeEr,
1521  // output
1522  T_1,
1523  fld_1,
1524  dx0,
1525  dy0,
1526  dxy0,
1527  du0,
1528  dv0);
1529 
1531 
1532 #ifdef DOUB_PERFORMANCE
1533  vector<THitI> hitsl_2;
1534 #endif // DOUB_PERFORMANCE
1535 
1536  f20( // input
1537  n1,
1538  stam,
1539  vStsHits_m,
1540  T_1,
1541  hitsl_1,
1542  // output
1543  n_2,
1544  i1_2,
1545 #ifdef DOUB_PERFORMANCE
1546  hitsl_2,
1547 #endif // DOUB_PERFORMANCE
1548  hitsm_2,
1549  Event,
1550  lmDuplets);
1551 
1552  for (Tindex i = 0; i < static_cast<Tindex>(hitsm_2.size()); ++i)
1553  L1_ASSERT(hitsm_2[i] < StsHitsUnusedStopIndex[istam]
1554  - StsHitsUnusedStartIndex[istam],
1555  hitsm_2[i] << " "
1556  << StsHitsUnusedStopIndex[istam]
1557  - StsHitsUnusedStartIndex[istam]);
1558 
1559 #ifdef DOUB_PERFORMANCE
1560  THitI* RealIHitL = &((*RealIHitP)[StsHitsUnusedStartIndex[istal]]);
1561  THitI* RealIHitM = &((*RealIHitP)[StsHitsUnusedStartIndex[istam]]);
1562  for (Tindex i = 0; i < n2; ++i) {
1563  // int i_4 = i%4;
1564  // int i_V = i/4;
1565  THitI iHits[2] = {RealIHitL[hitsl_2[i]], RealIHitM[hitsm_2[i]]};
1566  fL1Eff_doublets->AddOne(iHits);
1567  }
1568 #endif // DOUB_PERFORMANCE
1569  }
1570 }
1571 
1572 
1574 
1575 inline void
1576 L1Algo::TripletsStaPort(
1577  int istal,
1578  int istam,
1579  int istar,
1580 
1583 
1584 
1585  Tindex& nstaltriplets,
1586  L1TrackPar* T_1,
1587  L1FieldRegion* fld_1,
1588  THitI* hitsl_1,
1589 
1590  Tindex& n_2,
1591  vector<THitI>& i1_2,
1592  vector<THitI>& hitsm_2,
1593 
1594  const vector<bool>& mrDuplets
1595 
1597 
1598 ) {
1599  if (istar < NStations) {
1600  // prepare data
1601  L1Station& stam = vStations[istam];
1602  L1Station& star = vStations[istar];
1603 
1604  L1HitPoint* vStsHits_m =
1605  &((*vStsHitPointsUnused)[0]) + StsHitsUnusedStartIndex[istam];
1606 
1607  L1HitPoint* vStsHits_r = 0;
1608  vStsHits_r = &((*vStsHitPointsUnused)[0]) + StsHitsUnusedStartIndex[istar];
1609 
1610  Tindex n3 = 0, n3_V;
1611 
1613 
1614 
1615 #ifdef _OPENMP
1616  int Thread = omp_get_thread_num();
1617 #else
1618  int Thread = 0;
1619 #endif
1620 
1621  nsL1::vector<L1TrackPar>::TSimd& T_3 = fT_3[Thread];
1622  vector<THitI>& hitsl_3 = fhitsl_3[Thread];
1623  vector<THitI>& hitsm_3 = fhitsm_3[Thread];
1624  vector<THitI>& hitsr_3 = fhitsr_3[Thread];
1625  nsL1::vector<fvec>::TSimd& u_front3 = fu_front3[Thread];
1626  nsL1::vector<fvec>::TSimd& u_back3 = fu_back3[Thread];
1627  nsL1::vector<fvec>::TSimd& z_pos3 = fz_pos3[Thread];
1628  nsL1::vector<fvec>::TSimd& timeR = fTimeR[Thread];
1629  nsL1::vector<fvec>::TSimd& timeER = fTimeER[Thread];
1630  // nsL1::vector<fvec>::TSimd& dx3 = dx[Thread];
1631  // nsL1::vector<fvec>::TSimd& dy3 = dy[Thread];
1632  nsL1::vector<fvec>::TSimd& du3 = du[Thread];
1633  nsL1::vector<fvec>::TSimd& dv3 = dv[Thread];
1634 
1635  T_3.clear();
1636  hitsl_3.clear();
1637  hitsm_3.clear();
1638  hitsr_3.clear();
1639  u_front3.clear();
1640  u_back3.clear();
1641  z_pos3.clear();
1642  // dx3.clear();
1643  // dy3.clear();
1644  du3.clear();
1645  dv3.clear();
1646  timeR.clear();
1647  timeER.clear();
1648 
1650 
1651  f30( // input
1652  vStsHits_r,
1653  stam,
1654  star,
1655 
1656  istam,
1657  istar,
1658  vStsHits_m,
1659  T_1,
1660  fld_1,
1661  hitsl_1,
1662 
1663  n_2,
1664  hitsm_2,
1665  i1_2,
1666 
1667  mrDuplets,
1668  // output
1669  n3,
1670  T_3,
1671  hitsl_3,
1672  hitsm_3,
1673  hitsr_3,
1674  u_front3,
1675  u_back3,
1676  z_pos3,
1677  // dx3,
1678  // dy3,
1679  du3,
1680  dv3,
1681  timeR,
1682  timeER);
1683 
1684 
1685  n3_V = (n3 + fvecLen - 1) / fvecLen;
1686 
1687  for (Tindex i = 0; i < static_cast<Tindex>(hitsl_3.size()); ++i)
1688  L1_assert(hitsl_3[i] < StsHitsUnusedStopIndex[istal]
1689  - StsHitsUnusedStartIndex[istal]);
1690  for (Tindex i = 0; i < static_cast<Tindex>(hitsm_3.size()); ++i)
1691  L1_assert(hitsm_3[i] < StsHitsUnusedStopIndex[istam]
1692  - StsHitsUnusedStartIndex[istam]);
1693  for (Tindex i = 0; i < static_cast<Tindex>(hitsr_3.size()); ++i)
1694  L1_assert(hitsr_3[i] < StsHitsUnusedStopIndex[istar]
1695  - StsHitsUnusedStartIndex[istar]);
1696 
1697  // if (n3 >= MaxPortionTriplets) cout << "isec: " << isec << " stantion: " << istal << " portion number: " << ip << " CATrackFinder: Warning: Too many Triplets created in portion" << endl;
1698 
1700  f31( // input
1701  n3_V,
1702  star,
1703  u_front3,
1704  u_back3,
1705  z_pos3,
1706  // dx3,
1707  // dy3,
1708  du3,
1709  dv3,
1710  timeR,
1711  timeER,
1712  // output
1713  T_3);
1714 
1715 
1717  // f32( n3, istal, _RealIHit, T_3, hitsl_3, hitsm_3, hitsr_3, 0 );
1718 
1719 #ifdef TRIP_PERFORMANCE
1720  THitI* RealIHitL = &((*RealIHitP)[StsHitsUnusedStartIndex[istal]]);
1721  THitI* RealIHitM = &((*RealIHitP)[StsHitsUnusedStartIndex[istam]]);
1722  THitI* RealIHitR = &((*RealIHitP)[StsHitsUnusedStartIndex[istar]]);
1723  for (Tindex i = 0; i < n3; ++i) {
1724  Tindex i_4 = i % 4;
1725  Tindex i_V = i / 4;
1726  THitI iHits[3] = {
1727  RealIHitL[hitsl_3[i]], RealIHitM[hitsm_3[i]], RealIHitR[hitsr_3[i]]};
1728 #ifdef PULLS
1729  if (fL1Eff_triplets->AddOne(iHits))
1730  fL1Pulls->AddOne(T_3[i_V], i_4, iHits[2]);
1731 #else
1732  fL1Eff_triplets->AddOne(iHits);
1733 #endif
1734  if (T_3[i_V].chi2[i_4] < TRIPLET_CHI2_CUT)
1735  fL1Eff_triplets2->AddOne(iHits);
1736  }
1737 #endif // TRIP_PERFORMANCE
1738 
1740  f4( // input
1741  n3,
1742  istal,
1743  istam,
1744  istar,
1745  T_3,
1746  hitsl_3,
1747  hitsm_3,
1748  hitsr_3,
1749  // output
1750  nstaltriplets
1751 
1752  );
1753  }
1754 }
1755 
1756 
1757 // hitCheck::hitCheck()
1758 // {
1759 // omp_init_lock(&Occupied);
1760 // trackCandidateIndex = -1;
1761 // UsedByTrack=0;
1762 // Chi2Track = 100000000;
1763 // Length = 0;
1764 // ista = 1000;
1765 //
1766 // }
1767 // hitCheck::~hitCheck()
1768 // {
1769 // omp_destroy_lock(&Occupied);
1770 // }
1784 
1785 
1786 void L1Algo::CATrackFinder() {
1787 
1788 #ifdef _OPENMP
1789  omp_set_num_threads(fNThreads);
1790 #endif
1791 
1792 
1793 #ifdef PULLS
1794  static L1AlgoPulls* l1Pulls_ = new L1AlgoPulls();
1795  fL1Pulls = l1Pulls_;
1796  fL1Pulls->Init();
1797 #endif
1798 #ifdef TRIP_PERFORMANCE
1799  static L1AlgoEfficiencyPerformance<3>* l1Eff_triplets_ =
1800  new L1AlgoEfficiencyPerformance<3>();
1801  fL1Eff_triplets = l1Eff_triplets_;
1802  fL1Eff_triplets->Init();
1803  static L1AlgoEfficiencyPerformance<3>* l1Eff_triplets2_ =
1804  new L1AlgoEfficiencyPerformance<3>();
1805  fL1Eff_triplets2 = l1Eff_triplets2_;
1806  fL1Eff_triplets2->Init();
1807 #endif
1808 #ifdef DOUB_PERFORMANCE
1809  static L1AlgoEfficiencyPerformance<2>* l1Eff_doublets_ =
1810  new L1AlgoEfficiencyPerformance<2>();
1811  fL1Eff_doublets = l1Eff_doublets_;
1812  fL1Eff_doublets->Init();
1813 #endif
1814 
1815 #ifdef DRAW
1816  if (!draw) draw = new L1AlgoDraw;
1817  draw->InitL1Draw(this);
1818 #endif
1819 
1820  TStopwatch c_time; // for performance time
1821 #if defined(XXX) || defined(COUNTERS)
1822  static unsigned int stat_N = 0; // number of events
1823  stat_N++;
1824 #endif
1825 
1826 #ifdef XXX
1827  TStopwatch c_timerG; // global
1828  TStopwatch c_timerI; // for iterations
1829 
1830  L1CATFIterTimerInfo gti; // global
1831  gti.Add("init ");
1832  gti.Add("iterts");
1833  gti.Add("merge ");
1834 
1835  L1CATFTimerInfo ti;
1836  ti.SetNIter(fNFindIterations); // for iterations
1837  ti.Add("init ");
1838  // ti.Add("dblte1");
1839  // ti.Add("dblte2");
1840  ti.Add("tripl1");
1841 
1842 
1843  ti.Add("tracks");
1844  ti.Add("table");
1845  ti.Add("save");
1846  ti.Add("delete");
1847  ti.Add("copy");
1848  ti.Add("finish");
1849 
1850  static L1CATFIterTimerInfo stat_gti = gti;
1851  static L1CATFTimerInfo stat_ti = ti;
1852 
1853 
1854 #endif
1855 
1856 #ifdef COUNTERS
1857  static Tindex stat_nStartHits = 0;
1858  static Tindex stat_nHits[fNFindIterations] = {0};
1859 
1860  static Tindex stat_nSinglets[fNFindIterations] = {0};
1861  // static Tindex stat_nDoublets[fNFindIterations] = {0};
1862  static Tindex stat_nTriplets[fNFindIterations] = {0};
1863 
1864  static Tindex stat_nLevels[MaxNStations - 2][fNFindIterations] = {{0}};
1865  static Tindex stat_nCalls[fNFindIterations] = {0}; // n calls of CAFindTrack
1866  static Tindex stat_nTrCandidates[fNFindIterations] = {0};
1867 #endif
1868 
1869  RealIHitP = &RealIHit_v;
1870  RealIHitPBuf = &RealIHit_v_buf;
1871  vStsHitsUnused =
1872  &vStsDontUsedHits_B;
1873  vector<L1StsHit>* vStsHitsUnused_buf =
1874  &vStsDontUsedHits_A;
1875 
1876  vStsHitPointsUnused =
1877  &vStsDontUsedHitsxy_B;
1878  vector<L1HitPoint>* vStsHitPointsUnused_buf = &vStsDontUsedHitsxy_A;
1879 
1880  NHitsIsecAll = 0;
1881  NTracksIsecAll = 0;
1882 
1883  int nDontUsedHits = 0;
1884 
1885  // #pragma omp parallel for reduction(+:nDontUsedHits)
1886  for (int ista = 0; ista < NStations; ++ista) {
1887  nDontUsedHits += (StsHitsStopIndex[ista] - StsHitsStartIndex[ista]);
1888  StsHitsUnusedStartIndex[ista] = StsHitsStartIndex[ista];
1889  StsHitsUnusedStopIndex[ista] = StsHitsStopIndex[ista];
1890  }
1891 
1892  float lasttime = 0;
1893 
1894  for (int ist = 0; ist < NStations; ++ist)
1895  for (THitI ih = StsHitsStartIndex[ist]; ih < StsHitsStopIndex[ist]; ++ih) {
1896  if ((lasttime < (*vStsHits)[ih].t_reco)
1897  && (!isinf((*vStsHits)[ih].t_reco)))
1898  lasttime = (*vStsHits)[ih].t_reco;
1899  if (ist < NMvdStations) {
1900  L1StsHit& h = *(const_cast<L1StsHit*>(&((*vStsHits)[ih])));
1901  h.t_reco = 0;
1902  h.t_er = 100;
1903  }
1904  }
1905 
1906 
1907 #ifdef XXX
1908  c_time.Start();
1909  c_timerG.Start();
1910 
1911 #endif
1912 
1913  float yStep = 1.5 / sqrt(nDontUsedHits); // empirics. 0.01*sqrt(2374) ~= 0.5
1914 
1915 
1916  // float yStep = 0.5 / sqrt(nDontUsedHits); // empirics. 0.01*sqrt(2374) ~= 0.5
1917  if (yStep > 0.3) yStep = 0.3;
1918  float xStep = yStep * 3;
1919  // float xStep = yStep * 3;
1920 
1921  // yStep = 0.0078;
1922  // const float hitDensity = sqrt( nDontUsedHits );
1923 
1924  // float yStep = 0.7*4/hitDensity; // empirics. 0.01*sqrt(2374) ~= 0.5
1925  // if (yStep > 0.3)
1926  // yStep = 1.25;
1927  // xStep = 2.05;
1928 
1929  vStsHitsUnused = &vStsDontUsedHits_Buf;
1930 
1931 
1932  for (int iS = 0; iS < NStations; ++iS) {
1933  vGridTime[iS].BuildBins(
1934  -1, 1, -0.6, 0.6, 0, lasttime, xStep, yStep, (lasttime + 1) / 35);
1935  vGridTime[iS].StoreHits(
1936  StsHitsUnusedStopIndex[iS] - StsHitsUnusedStartIndex[iS],
1937  &((*vStsHits)[StsHitsUnusedStartIndex[iS]]),
1938  iS,
1939  *this,
1940  StsHitsUnusedStartIndex[iS],
1941  &(vStsDontUsedHits_Buf[StsHitsUnusedStartIndex[iS]]),
1942  &((*vStsHits)[StsHitsUnusedStartIndex[iS]]),
1943  &(RealIHit_v[StsHitsUnusedStartIndex[iS]]));
1944  }
1945 
1946 
1947  for (int ist = 0; ist < NStations; ++ist)
1948  for (THitI ih = StsHitsStartIndex[ist]; ih < StsHitsStopIndex[ist]; ++ih) {
1949  L1StsHit& h = *(const_cast<L1StsHit*>(&((*vStsHits)[ih])));
1950  SetFUnUsed(const_cast<unsigned char&>((*vSFlag)[h.f]));
1951  SetFUnUsed(const_cast<unsigned char&>((*vSFlagB)[h.b]));
1952  }
1953 
1954  for (int ista = 0; ista < NStations; ++ista) {
1955 
1956 
1957 #ifdef _OPENMP
1958 #pragma omp parallel for schedule(dynamic, 5)
1959 #endif
1960  for (THitI ih = StsHitsStartIndex[ista]; ih < StsHitsStopIndex[ista]; ++ih)
1961  CreateHitPoint(vStsDontUsedHits_Buf[ih], ista, vStsDontUsedHitsxy_B[ih]);
1962  }
1963 
1964 #ifdef COUNTERS
1965  stat_nStartHits += nDontUsedHits;
1966 #endif
1967 
1968 #ifdef XXX
1969  c_timerG.Stop();
1970  gti["init "] = c_timerG;
1971  c_timerG.Start();
1972 #endif
1973 
1974  TStopwatch c_time1;
1975  c_time1.Start();
1976 
1977  for (isec = 0; isec < fNFindIterations; ++isec) // all finder
1978  {
1979  if (fmCBMmode)
1980  if (isec != 0) continue;
1981 
1982  n_g1.assign(n_g1.size(), Portion);
1983 
1984  for (int n = 0; n < nTh; n++)
1985  for (int j = 0; j < 12; j++)
1986  nTripletsThread[j][n] = 0;
1987 
1989 #ifdef COUNTERS
1990  Tindex nSinglets = 0;
1991 #endif
1992 
1993  if (isec != 0) {
1994  vector<THitI>* RealIHitPTmp = RealIHitP;
1995  RealIHitP = RealIHitPBuf;
1996  RealIHitPBuf = RealIHitPTmp;
1997 
1998  vector<L1StsHit>* vStsHitsUnused_temp = vStsHitsUnused;
1999  vStsHitsUnused = vStsHitsUnused_buf;
2000  vStsHitsUnused_buf = vStsHitsUnused_temp;
2001 
2002  vector<L1HitPoint>* vStsHitsUnused_temp2 = vStsHitPointsUnused;
2003  vStsHitPointsUnused = vStsHitPointsUnused_buf;
2004  vStsHitPointsUnused_buf = vStsHitsUnused_temp2;
2005  }
2006 
2007  {
2008  // #pragma omp task
2009  {
2010  // --- SET PARAMETERS FOR THE ITERATION ---
2011 
2012  FIRSTCASTATION = 0;
2013  // if ( (isec == kAllSecIter) || (isec == kAllSecEIter) || (isec == kAllSecJumpIter) )
2014  // FIRSTCASTATION = 2;
2015 
2016  DOUBLET_CHI2_CUT = 11.3449 * 2.f / 3.f; // prob = 0.1
2017 
2018  TRIPLET_CHI2_CUT = 21.1075; // prob = 0.01%
2019 
2020  switch (isec) {
2021  case kFastPrimIter:
2022  TRIPLET_CHI2_CUT = 7.815 * 3; // prob = 0.05
2023  break;
2024  case kAllPrimIter:
2025  case kAllPrimEIter:
2026  TRIPLET_CHI2_CUT = 7.815 * 3; // prob = 0.05
2027  break;
2028  case kAllPrimJumpIter:
2029  TRIPLET_CHI2_CUT = 6.252 * 3; // prob = 0.1
2030  break;
2031  case kAllSecIter:
2032  case kAllSecEIter:
2033  TRIPLET_CHI2_CUT = 6.252 * 3; //2.706; // prob = 0.1
2034  break;
2035  }
2036 
2037  Pick_gather =
2038  3.0;
2039  if ((isec == kAllPrimIter) || (isec == kAllPrimEIter)
2040  || (isec == kAllPrimJumpIter) || (isec == kAllSecIter)
2041  || (isec == kAllSecEIter) || (isec == kAllSecJumpIter))
2042  Pick_gather = 4.0;
2043 
2044  PickNeighbour =
2045  1.0; // (PickNeighbour < dp/dp_error) => triplets are neighbours
2046  // if ( (isec == kFastPrimIter) )
2047  // PickNeighbour = 0.5; // TODO understand why works with 0.2
2048 
2049  MaxInvMom = 1.0 / 0.5; // max considered q/p
2050  if ((isec == kAllPrimJumpIter) || (isec == kAllSecIter)
2051  || (isec == kAllSecJumpIter))
2052  MaxInvMom = 1.0 / 0.1;
2053  if ((isec == kAllPrimIter) || (isec == kAllPrimEIter)
2054  || (isec == kAllSecEIter))
2055  MaxInvMom = 1. / 0.05;
2056 
2057  MaxSlope = 1.1;
2058  if ( // (isec == kAllPrimIter) || (isec == kAllPrimEIter) || (isec == kAllPrimJumpIter) ||
2059  (isec == kAllSecIter) || (isec == kAllSecEIter)
2060  || (isec == kAllSecJumpIter))
2061  MaxSlope = 1.5;
2062 
2063  // define the target
2064  targX = 0;
2065  targY = 0;
2066  targZ = 0; // suppose, what target will be at (0,0,0)
2067 
2068  float SigmaTargetX = 0, SigmaTargetY = 0; // target constraint [cm]
2069  if ((isec == kFastPrimIter) || (isec == kFastPrimIter2)
2070  || (isec == kFastPrimJumpIter) || (isec == kAllPrimIter)
2071  || (isec == kAllPrimEIter)
2072  || (isec == kAllPrimJumpIter)) { // target
2073  targB = vtxFieldValue;
2074  if ((isec == kFastPrimIter) || (isec == kAllPrimIter)
2075  || (isec == kAllPrimEIter))
2076  SigmaTargetX = SigmaTargetY = 1; // target
2077  else
2078  SigmaTargetX = SigmaTargetY = 5;
2079  }
2080  if (
2081  (isec == kAllSecIter) || (isec == kAllSecEIter)
2082  || (isec
2083  == kAllSecJumpIter)) { //use outer radius of the 1st station as a constraint
2084  L1Station& st = vStations[0];
2085  SigmaTargetX = SigmaTargetY = 10; //st.Rmax[0];
2086  targZ = 0.; //-1.;
2087  st.fieldSlice.GetFieldValue(0, 0, targB);
2088  }
2089 
2090  TargetXYInfo.C00 = SigmaTargetX * SigmaTargetX;
2091  TargetXYInfo.C10 = 0;
2092  TargetXYInfo.C11 = SigmaTargetY * SigmaTargetY;
2093 
2097  MaxDZ = 0;
2098  if ((isec == kAllPrimIter) || (isec == kAllPrimEIter)
2099  || (isec == kAllPrimJumpIter) || (isec == kAllSecIter)
2100  || (isec == kAllSecEIter) || (isec == kAllSecJumpIter))
2101  MaxDZ = 0.1;
2102 
2103  if (NStations > MaxNStations)
2104  cout << " CATrackFinder: Error: Too many Stantions" << endl;
2105  }
2106 
2107 #ifndef L1_NO_ASSERT
2108  for (int istal = NStations - 1; istal >= 0; istal--) {
2109  L1_ASSERT(StsHitsUnusedStopIndex[istal]
2110  >= StsHitsUnusedStartIndex[istal],
2111  StsHitsUnusedStopIndex[istal]
2112  << " >= " << StsHitsUnusedStartIndex[istal]);
2113  L1_ASSERT(StsHitsUnusedStopIndex[istal] <= (*vStsHitsUnused).size(),
2114  StsHitsUnusedStopIndex[istal] << " <= "
2115  << (*vStsHitsUnused).size());
2116  }
2117 #endif // L1_NO_ASSERT
2118  }
2119 
2120 
2121  {
2123  portionStopIndex[NStations - 1] = 0;
2124  unsigned int ip = 0; //index of curent portion
2125 
2126  for (int istal = NStations - 2; istal >= FIRSTCASTATION;
2127  istal--) { //start downstream chambers
2128  int NHits_l =
2129  StsHitsUnusedStopIndex[istal] - StsHitsUnusedStartIndex[istal];
2130 
2131  int NHits_l_P = NHits_l / Portion;
2132 
2133  for (int ipp = 0; ipp < NHits_l_P; ipp++) {
2134  // n_g1[ip++] = Portion;
2135  n_g1[ip] = (Portion);
2136 
2137  ip++;
2138  } // start_lh - portion of left hits
2139 
2140  // n_g1[ip++] = NHits_l - NHits_l_P*Portion;
2141  n_g1[ip] = (NHits_l - NHits_l_P * Portion);
2142 
2143  ip++;
2144  portionStopIndex[istal] = ip;
2145  } // lstations
2146 
2147 
2148 #ifdef COUNTERS
2149  stat_nSinglets[isec] += nSinglets;
2150 #endif
2151  }
2152 
2153  /* {
2155  portionStopIndex[NStations-1] = 0;
2156  unsigned int ip = 0; //index of curent portion
2157 
2158  for (int istal = NStations-2; istal >= FIRSTCASTATION; istal--) //start downstream chambers
2159  {
2160  int nHits = StsHitsUnusedStopIndex[istal] - StsHitsUnusedStartIndex[istal];
2161 
2162  int NHits_P = nHits/Portion;
2163 
2164  for( int ipp = 0; ipp < NHits_P; ipp++ )
2165  {
2166  n_g1[ip] = Portion;
2167  ip++;
2168  } // start_lh - portion of left hits
2169 
2170  n_g1[ip] = nHits - NHits_P*Portion;
2171 
2172  ip++;
2173  portionStopIndex[istal] = ip;
2174  }// lstations
2175  // nPortions = ip;
2176  } */
2177 
2179 
2180 #ifdef XXX
2181  TStopwatch c_timer;
2182  c_timer.Start();
2183 #endif
2184 
2185 
2186  const Tindex vPortion = Portion / fvecLen;
2187  L1TrackPar T_1[vPortion];
2188  L1FieldRegion fld_1[vPortion];
2189  THitI hitsl_1[Portion];
2190  L1TrackPar TG_1[vPortion];
2191  L1FieldRegion fldG_1[vPortion];
2192  THitI hitslG_1[Portion];
2193 
2194  vector<THitI>
2195  hitsm_2;
2196  vector<THitI>
2197  i1_2;
2198 
2199 
2200  vector<THitI>
2201  hitsmG_2;
2202  vector<THitI>
2203  i1G_2;
2204  vector<bool> lmDuplets
2205  [MaxNStations]; // is exist a doublet started from indexed by left hit
2206  vector<bool> lmDupletsG
2207  [MaxNStations]; // is exist a doublet started from indexed by left hit
2208  hitsm_2.reserve(3500);
2209  i1_2.reserve(3500);
2210 
2211  hitsmG_2.reserve(800);
2212  i1G_2.reserve(800);
2213 
2214  for (int istal = NStations - 2; istal >= FIRSTCASTATION;
2215  istal--) // //start downstream chambers
2216  {
2217 
2218 #ifdef _OPENMP
2219 #pragma omp parallel for firstprivate(T_1, \
2220  fld_1, \
2221  hitsl_1, \
2222  hitsm_2, \
2223  i1_2, \
2224  TG_1, \
2225  fldG_1, \
2226  hitslG_1, \
2227  hitsmG_2, \
2228  i1G_2) //schedule(dynamic, 2)
2229 #endif
2230  for (Tindex ip = portionStopIndex[istal + 1];
2231  ip < portionStopIndex[istal];
2232  ++ip) {
2233  Tindex n_2;
2234  int NHitsSta = StsHitsStopIndex[istal] - StsHitsStartIndex[istal];
2235  lmDuplets[istal].resize(NHitsSta);
2236  lmDupletsG[istal].resize(NHitsSta);
2237 
2238  hitsm_2.clear();
2239  i1_2.clear();
2240 
2241 
2242  DupletsStaPort(istal,
2243  istal + 1,
2244  ip,
2245  n_g1,
2246  portionStopIndex,
2247 
2248  // output
2249  T_1,
2250  fld_1,
2251  hitsl_1,
2252 
2253  lmDuplets[istal],
2254 
2255 
2256  n_2,
2257  i1_2,
2258  hitsm_2);
2259 
2260 
2261  Tindex nstaltriplets = 0;
2262 
2263  TripletsStaPort( // input
2264  istal,
2265  istal + 1,
2266  istal + 2,
2267  nstaltriplets,
2268  T_1,
2269  fld_1,
2270  hitsl_1,
2271 
2272  n_2,
2273  i1_2,
2274  hitsm_2,
2275 
2276  lmDuplets[istal + 1]
2277  // output
2278  );
2279 
2280 
2281  if ((isec == kFastPrimJumpIter) || (isec == kAllPrimJumpIter)
2282  || (isec == kAllSecJumpIter)) {
2283  Tindex nG_2;
2284  hitsmG_2.clear();
2285  i1G_2.clear();
2286 
2287  DupletsStaPort( // input
2288  istal,
2289  istal + 2,
2290  ip,
2291  n_g1,
2292  portionStopIndex,
2293  // output
2294  TG_1,
2295  fldG_1,
2296  hitslG_1,
2297 
2298  lmDupletsG[istal],
2299 
2300  nG_2,
2301  i1G_2,
2302  hitsmG_2);
2303 
2304  TripletsStaPort( // input
2305  istal,
2306  istal + 1,
2307  istal + 3,
2308  nstaltriplets,
2309  T_1,
2310  fld_1,
2311  hitsl_1,
2312 
2313  n_2,
2314  i1_2,
2315  hitsm_2,
2316  lmDupletsG[istal + 1]);
2317 
2318  TripletsStaPort( // input
2319  istal,
2320  istal + 2,
2321  istal + 3,
2322  nstaltriplets,
2323  TG_1,
2324  fldG_1,
2325  hitslG_1,
2326 
2327  nG_2,
2328  i1G_2,
2329  hitsmG_2,
2330  lmDuplets[istal + 2]
2331 
2332  );
2333  }
2334  } //
2335  }
2336 
2337  // int nlevels[MaxNStations]; // number of triplets with some number of neighbours.
2338  // for (int il = 0; il < NStations; ++il) nlevels[il] = 0;
2339  //
2340  // f5( // input
2341  // // output
2342  // 0,
2343  // nlevels
2344  // );
2345 
2346 #ifdef XXX
2347  c_timer.Stop();
2348  ti[isec]["tripl1"] = c_timer;
2349  c_timer.Start();
2350 #endif
2351 
2357 
2358  // #ifdef XXX
2359  // cout<<"---- Collect track candidates. ----"<<endl;
2360  // #endif
2361 
2362  int min_level =
2363  0; // min level for start triplet. So min track length = min_level+3.
2364  // if (isec == kFastPrimJumpIter) min_level = 1;
2365  if ((isec == kAllSecIter) || (isec == kAllSecEIter)
2366  || (isec == kAllSecJumpIter))
2367  min_level = 1; // only the long low momentum tracks
2368 
2369 #ifdef TRACKS_FROM_TRIPLETS
2370  if (isec == TRACKS_FROM_TRIPLETS_ITERATION) min_level = 0;
2371 #endif
2372 
2373  // int min_level = 1; // min level for start triplet. So min track length = min_level+3.
2374  // if (isec == kAllPrimJumpIter) min_level = 1;
2375  // if ( (isec == kAllSecIter) || (isec == kAllSecJumpIter) ) min_level = 2; // only the long low momentum tracks
2376  // // if (isec == -1) min_level = NStations-3 - 3; //only the longest tracks
2377 
2378 
2379  L1Branch curr_tr;
2380  L1Branch new_tr[MaxNStations];
2381  L1Branch best_tr;
2382  fscal curr_chi2 = 0;
2383 
2384  fscal best_chi2 = 0;
2385  unsigned char best_L = 0;
2386 
2387  unsigned char curr_L = 1;
2388  int ndf = 1;
2389 
2390  vStripToTrack.assign(StsHitsStopIndex[NStations - 1], -1);
2391  vStripToTrackB.assign(StsHitsStopIndex[NStations - 1], -1);
2392 
2393  // collect consequtive: the longest tracks, shorter, more shorter and so on
2394  for (int ilev = NStations - 3; ilev >= min_level; ilev--) {
2395  // choose length in triplets number - ilev - the maximum possible triplet level among all triplets in the searched candidate
2396  TStopwatch Time;
2397 
2398  // how many levels to check
2399  int nlevel = (NStations - 2) - ilev + 1;
2400 
2401  const unsigned char min_best_l =
2402  (ilev > min_level) ? ilev + 2 : min_level + 3; // loose maximum
2403 
2404  for (int i = 0; i < fNThreads; ++i)
2405  numberCandidateThread[i] = 0;
2406 
2407  for (int istaF = FIRSTCASTATION; istaF <= NStations - 3 - ilev; ++istaF) {
2408 
2409 #ifdef _OPENMP
2410 #pragma omp parallel for firstprivate(curr_tr, \
2411  new_tr, \
2412  best_tr, \
2413  curr_chi2, \
2414  best_chi2, \
2415  best_L, \
2416  curr_L, \
2417  ndf) // schedule(dynamic, 10)
2418 #endif
2419  for (Tindex ip = 0; ip < fNThreads; ++ip) {
2420  for (Tindex itrip = 0; itrip < nTripletsThread[istaF][ip]; ++itrip) {
2421 
2422 #ifdef _OPENMP
2423  int thread_num = omp_get_thread_num();
2424 #else
2425  int thread_num = 0;
2426 #endif
2427  L1Triplet& first_trip = (TripletsLocal1[istaF][ip][itrip]);
2428 
2429  if (GetFUsed(
2430  (*vSFlag)[(*vStsHitsUnused)[first_trip.GetLHit()].f]
2431  | (*vSFlagB)[(*vStsHitsUnused)[first_trip.GetLHit()].b]))
2432  continue;
2433 
2434 
2435  // ghost supression !!!
2436 #ifndef FIND_GAPED_TRACKS
2437  if (/*(isec == kFastPrimIter) ||*/ (isec == kAllPrimIter)
2438  || (isec == kAllPrimEIter) || (isec == kAllSecIter)
2439  || (isec == kAllSecEIter) || (isec == kAllSecJumpIter)) {
2440 #else
2441  if ((isec == kFastPrimIter) || (isec == kFastPrimIter2)
2442  || (isec == kFastPrimJumpIter) || (isec == kAllPrimIter)
2443  || (isec == kAllPrimEIter) || (isec == kAllPrimJumpIter)
2444  || (isec == kAllSecIter) || (isec == kAllSecEIter)
2445  || (isec == kAllSecJumpIter)) {
2446 #endif
2447 #ifdef TRACKS_FROM_TRIPLETS
2448  if (isec != TRACKS_FROM_TRIPLETS_ITERATION)
2449 #endif
2450  { // ghost supression !!!
2451  if (isec != kFastPrimIter && isec != kAllPrimIter
2452  && isec != kAllPrimEIter && isec != kAllSecEIter)
2453  if (first_trip.GetLevel() == 0)
2454  continue; // ghost suppression // find track with 3 hits only if it was created from a chain of triplets, but not from only one triplet
2455 
2456  if (!fmCBMmode)
2457  if ((ilev == 0)
2458  && (GetFStation((
2459  *vSFlag)[(*vStsHitsUnused)[first_trip.GetLHit()].f])
2460  != 0))
2461  continue; // ghost supression // collect only MAPS tracks-triplets CHECK!!!
2462  }
2463  if (first_trip.GetLevel() < ilev)
2464  continue; // try only triplets, which can start track with ilev+3 length. w\o it have more ghosts, but efficiency either
2465  }
2466 
2467 
2468  // curr_tr.Momentum = 0.f;
2469  curr_tr.chi2 = 0.f;
2470  // curr_tr.Lengtha = 0;
2471  curr_tr.ista = 0;
2472 
2473  (curr_tr).StsHits[0] = ((*RealIHitP)[first_trip.GetLHit()]);
2474 
2475  (curr_tr).NHits = 1;
2476 
2477  curr_L = 1;
2478  curr_chi2 = first_trip.GetChi2();
2479 
2480  best_tr = (curr_tr);
2481 
2482  best_chi2 = curr_chi2;
2483  best_L = curr_L;
2484 
2485  CAFindTrack(
2486  istaF,
2487  best_tr,
2488  best_L,
2489  best_chi2,
2490  &first_trip,
2491  (curr_tr),
2492  curr_L,
2493  curr_chi2,
2494  min_best_l,
2495  new_tr);
2496 
2497  // if ( best_L < min_best_l ) continue;
2498  if (best_L < ilev + 2) continue; // lose maximum one hit
2499 
2500  if (best_L < min_level + 3)
2501  continue; // should find all hits for min_level
2502 
2503  ndf = best_L * 2 - 5;
2504  best_chi2 = best_chi2 / ndf; //normalize
2505 
2506 #ifndef TRACKS_FROM_TRIPLETS
2507  if (fGhostSuppression) {
2508  if (best_L == 3) {
2509  // if( isec == kAllSecIter ) continue; // too /*short*/ secondary track
2510  if (((isec == kAllSecIter) || (isec == kAllSecEIter)
2511  || (isec == kAllSecJumpIter))
2512  && (istaF != 0))
2513  continue; // too /*short*/ non-MAPS track
2514  if ((isec != kAllSecIter) && (isec != kAllSecEIter)
2515  && (isec != kAllSecJumpIter) && (best_chi2 > 5.0))
2516  continue;
2517  }
2518  }
2519 #endif
2520  best_tr.Set(istaF, best_L, best_chi2, first_trip.GetQpOrig());
2521  L1Branch& tr =
2522  CandidatesTrack[thread_num][numberCandidateThread[thread_num]];
2523  tr = best_tr;
2524  tr.CandIndex =
2525  numberCandidateThread[thread_num] + thread_num * 100000;
2526 
2527  bool check = 1;
2528 
2529  for (vector<THitI>::iterator phitIt =
2530  tr.StsHits.begin();
2531  phitIt != tr.StsHits.begin() + tr.NHits;
2532  ++phitIt) {
2533  const L1StsHit& h = (*vStsHits)[*phitIt];
2534 #ifdef _OPENMP
2535  omp_set_lock(&hitToBestTrackB[h.b]);
2536 #endif
2537  int& strip1 = (vStripToTrackB)[h.b];
2538 
2539  if (strip1 != -1) {
2540  int thread = strip1 / 100000;
2541  int num = (strip1 - thread * 100000);
2542 
2543  if (L1Branch::compareCand(tr, CandidatesTrack[thread][num])) {
2544  CandidatesTrack[thread][num].CandIndex = -1;
2545  strip1 = tr.CandIndex;
2546  } else {
2547  check = 0;
2548 #ifdef _OPENMP
2549  omp_unset_lock(&hitToBestTrackB[h.b]);
2550 #endif
2551  break;
2552  }
2553  } else
2554  strip1 = tr.CandIndex;
2555 #ifdef _OPENMP
2556  omp_unset_lock(&hitToBestTrackB[h.b]);
2557 #endif
2558 
2559  if (check) {
2560 #ifdef _OPENMP
2561  omp_set_lock(&hitToBestTrackF[h.f]);
2562 #endif
2563  int& strip2 = (vStripToTrack)[h.f];
2564  if (strip2 != -1) {
2565  int thread = strip2 / 100000;
2566  int num = (strip2 - thread * 100000);
2567 
2568  if (L1Branch::compareCand(tr, CandidatesTrack[thread][num])) {
2569  CandidatesTrack[thread][num].CandIndex = -1;
2570  strip2 = tr.CandIndex;
2571 
2572  } else {
2573  check = 0;
2574 #ifdef _OPENMP
2575  omp_unset_lock(&hitToBestTrackF[h.f]);
2576 #endif
2577  break;
2578  }
2579  } else
2580  strip2 = tr.CandIndex;
2581 #ifdef _OPENMP
2582  omp_unset_lock(&hitToBestTrackF[h.f]);
2583 #endif
2584  }
2585  }
2586 
2587  if (check) numberCandidateThread[thread_num]++;
2588  } // itrip
2589  }
2590  }
2591 
2592  if (--nlevel == 0) break;
2593 
2594  for (int i = 0; i < fNThreads; ++i) {
2595  SavedCand[i] = 0;
2596  SavedHits[i] = 0;
2597  }
2598 
2599  for (int i = 0; i < fNThreads; ++i) {
2600  L1Track t;
2601 
2602 #ifdef _OPENMP
2603 #pragma omp parallel for schedule(dynamic, 10) firstprivate(t)
2604 #endif
2605  for (Tindex iCandidate = 0; iCandidate < numberCandidateThread[i];
2606  ++iCandidate) {
2607  L1Branch& tr = CandidatesTrack[i][iCandidate];
2608 
2609  bool check = 1;
2610 
2611  if (CandidatesTrack[i][iCandidate].CandIndex != -1) {
2612  for (vector<THitI>::iterator phIt =
2613  tr.StsHits.begin();
2614  phIt != tr.StsHits.begin() + tr.NHits;
2615  ++phIt) {
2616  const L1StsHit& h = (((*vStsHits))[*phIt]);
2617  if (((vStripToTrackB)[h.b] != tr.CandIndex)
2618  || ((vStripToTrack)[h.f] != tr.CandIndex)) {
2619  check = 0;
2620  break;
2621  }
2622  }
2623 
2624  if (tr.NHits < 3) check = 0;
2625 
2626  if (check) {
2627 #ifdef EXTEND_TRACKS
2628  if (!fmCBMmode)
2629  if (tr.NHits != NStations) BranchExtender(tr);
2630 #endif
2631  float sumTime = 0;
2632 
2633 
2634 #ifdef _OPENMP
2635 
2636  int num_thread = omp_get_thread_num();
2637 
2638 #else
2639 
2640  int num_thread = 0;
2641 
2642 #endif
2643 
2644  for (vector<THitI>::iterator phIt =
2645  tr.StsHits.begin();
2646  phIt != tr.StsHits.begin() + tr.NHits;
2647  ++phIt) {
2648  L1StsHit& h = *(const_cast<L1StsHit*>(&(((*vStsHits))[*phIt])));
2649 
2650 
2651  SetFUsed(const_cast<unsigned char&>((*vSFlag)[h.f]));
2652  SetFUsed(const_cast<unsigned char&>((*vSFlagB)[h.b]));
2653 
2654 
2655  vRecoHits_local[num_thread][SavedHits[num_thread]] = (*phIt);
2656 
2657  SavedHits[num_thread]++;
2658 
2659  const L1StsHit& hit = (*vStsHits)[*phIt];
2660 
2661 
2662  L1HitPoint tempPoint = CreateHitPoint(
2663  hit, 0); //TODO take number of station from hit
2664 
2665  float xcoor, ycoor = 0;
2666  L1Station stah = vStations[0];
2667  StripsToCoor(tempPoint.U(), tempPoint.V(), xcoor, ycoor, stah);
2668  float zcoor = tempPoint.Z();
2669 
2670  float timeFlight =
2671  sqrt(xcoor * xcoor + ycoor * ycoor + zcoor * zcoor)
2672  / 30.f; // c = 30[cm/ns]
2673  sumTime += (hit.t_reco - timeFlight);
2674  }
2675 
2676  t.NHits = tr.NHits;
2677  // t.Momentum = tr.Momentum;
2678  t.fTrackTime = sumTime / t.NHits;
2679 
2680  vTracks_local[num_thread][SavedCand[num_thread]] = (t);
2681  SavedCand[num_thread]++;
2682  }
2683  }
2684  }
2685  }
2686 
2687 
2688 #ifdef XXX
2689  Time.Stop();
2690  ti[isec]["table"] += Time;
2691 
2692  Time.Start();
2693 
2694 #endif
2695 
2696  vector<int> offset_tracks(nTh, NTracksIsecAll);
2697  vector<int> offset_hits(nTh, NHitsIsecAll);
2698 
2699  NTracksIsecAll += SavedCand[0];
2700  NHitsIsecAll += SavedHits[0];
2701 
2702 
2703  for (int i = 1; i < nTh; ++i) {
2704  offset_tracks[i] = offset_tracks[i - 1] + SavedCand[i - 1];
2705  offset_hits[i] = offset_hits[i - 1] + SavedHits[i - 1];
2706  NTracksIsecAll += SavedCand[i];
2707  NHitsIsecAll += SavedHits[i];
2708  }
2709 
2710 #ifdef _OPENMP
2711 #pragma omp parallel for
2712 #endif
2713  for (int i = 0; i < nTh; ++i) {
2714  for (Tindex iC = 0; iC < SavedCand[i]; ++iC)
2715  vTracks[offset_tracks[i] + iC] = (vTracks_local[i][iC]);
2716 
2717  for (Tindex iH = 0; iH < SavedHits[i]; ++iH)
2718  vRecoHits[offset_hits[i] + iH] = (vRecoHits_local[i][iH]);
2719  }
2720  } //istaf
2721 
2722 #ifdef XXX
2723  c_timer.Stop();
2724  ti[isec]["tracks"] = c_timer;
2725  c_timer.Start();
2726 #endif
2727 
2728 
2729  if (isec < (fNFindIterations - 1)) {
2730  int NHitsOnStation = 0;
2731 
2732  for (int ista = 0; ista < NStations; ++ista) {
2733  int Nelements =
2734  StsHitsUnusedStopIndex[ista] - StsHitsUnusedStartIndex[ista];
2735  int NHitsOnStationTmp = NHitsOnStation;
2736  vGridTime[ista].UpdateIterGrid(
2737  Nelements,
2738  &((*vStsHitsUnused)[StsHitsUnusedStartIndex[ista]]),
2739  RealIHitPBuf,
2740  &((*RealIHitP)[StsHitsUnusedStartIndex[ista]]),
2741  vStsHitsUnused_buf,
2742  vStsHitPointsUnused_buf,
2743  &((*vStsHitPointsUnused)[StsHitsUnusedStartIndex[ista]]),
2744  NHitsOnStation,
2745  ista,
2746  *this,
2747  vSFlag,
2748  vSFlagB);
2749  StsHitsUnusedStartIndex[ista] = NHitsOnStationTmp;
2750  StsHitsUnusedStopIndex[ista] = NHitsOnStation;
2751  }
2752 
2753 #ifdef XXX
2754  c_timer.Stop();
2755  ti[isec]["finish"] = c_timer;
2756 #endif
2757 
2758 #ifdef XXX
2759  // if( stat_max_n_trip<stat_n_trip ) stat_max_n_trip = vTriplets.size();
2760  // Tindex tsize = vTripletsP.size()*sizeof(L1Triplet);
2761  // if( stat_max_trip_size < tsize ) stat_max_trip_size = tsize;
2762 #endif
2763  // #ifdef DRAW
2764  // draw->ClearVeiw();
2765  // // draw->DrawInfo();
2766  // draw->DrawRestHits(StsHitsUnusedStartIndex, StsHitsUnusedStopIndex, RealIHit);
2767  // draw->DrawRecoTracks();
2768  // draw->SaveCanvas("Reco_"+isec+"_");
2769  // draw->DrawAsk();
2770  // #endif
2771 
2772  // #ifdef PULLS
2773  // fL1Pulls->Build(1);
2774  // #endif
2775 #ifdef COUNTERS
2776  // stat_nHits[isec] += (vStsHitsUnused*)->Size();
2777 
2778  cout << "iter = " << isec << endl;
2779  cout << " NSinglets = " << stat_nSinglets[isec] / stat_N << endl;
2780  // cout << " NDoublets = " << stat_nDoublets[isec]/stat_N << endl;
2781  cout << " NTriplets = " << stat_nTriplets[isec] / stat_N << endl;
2782  cout << " NHitsUnused = " << stat_nHits[isec] / stat_N << endl;
2783 
2784 #endif // COUNTERS
2785  }
2786  } // for (int isec
2787 
2788 #ifdef XXX
2789  c_timerG.Stop();
2790  gti["iterts"] = c_timerG;
2791  c_timerG.Start();
2792 #endif
2793 
2794 #ifdef MERGE_CLONES
2795  CAMergeClones();
2796 #endif
2797 
2798 #ifdef XXX
2799  c_timerG.Stop();
2800  gti["merge "] = c_timerG;
2801 #endif
2802 
2803  //==================================
2804 
2805  c_time.Stop();
2806 
2807  // cout << "End TrackFinder" << endl;
2808  // CATime = (double(c_time.CpuTime()));
2809  CATime = (double(c_time.RealTime()));
2810 
2811 #ifdef XXX
2812 
2813 
2814  cout << endl << " --- Timers, ms --- " << endl;
2815  ti.Calc();
2816  stat_ti += ti;
2817  L1CATFTimerInfo tmp_ti = stat_ti / 0.001 / stat_N; // ms
2818 
2819  tmp_ti.PrintReal();
2820  stat_gti += gti;
2821  L1CATFIterTimerInfo tmp_gti = stat_gti / 0.001 / stat_N; // ms
2822  tmp_gti.PrintReal(1);
2823  fstream filestr;
2824  filestr.open("speedUp.log", fstream::out | fstream::app);
2825  float tripl_speed = 1000. / (tmp_ti.GetTimerAll()["tripl1"].Real());
2826  filestr << tripl_speed << " ";
2827  filestr.close();
2828 
2829 
2830 #if 0
2831  static long int NTimes =0, NHits=0, HitSize =0, NStrips=0, StripSize =0, NStripsB=0, StripSizeB =0,
2832  NDup=0, DupSize=0, NTrip=0, TripSize=0, NBranches=0, BranchSize=0, NTracks=0, TrackSize=0 ;
2833 
2834  NTimes++;
2835  NHits += vStsHitsUnused->size();
2836  HitSize += vStsHitsUnused->size()*sizeof(L1StsHit);
2837  NStrips+= vStsStrips.size();
2838  StripSize += vStsStrips.size()*sizeof(fscal) + (*vSFlag).size()*sizeof(unsigned char);
2839  NStripsB+= (*vSFlagB).size();
2840  StripSizeB += vStsStripsB.size()*sizeof(fscal) + (*vSFlagB).size()*sizeof(unsigned char);
2841  NDup += stat_max_n_dup;
2842  DupSize += stat_max_n_dup*sizeof(/*short*/ int);
2843  NTrip += stat_max_n_trip;
2844  TripSize += stat_max_trip_size;
2845 
2846  NBranches += stat_max_n_branches;
2847  BranchSize += stat_max_BranchSize;
2848  NTracks += vTracks.size();
2849  TrackSize += sizeof(L1Track)*vTracks.size() + sizeof(THitI)*vRecoHits.size();
2850  int k = 1024*NTimes;
2851 
2852  cout<<"L1 Event size: \n"
2853  <<HitSize/k<<"kB for "<<NHits/NTimes<<" hits, \n"
2854  <<StripSize/k<<"kB for "<<NStrips/NTimes<<" strips, \n"
2855  <<StripSizeB/k<<"kB for "<<NStripsB/NTimes<<" stripsB, \n"
2856  <<DupSize/k<<"kB for "<<NDup/NTimes<<" doublets, \n"
2857  <<TripSize/k<<"kB for "<<NTrip/NTimes<<" triplets\n"
2858  <<BranchSize/k<<"kB for "<<NBranches/NTimes<<" branches, \n"
2859  <<TrackSize/k<<"kB for "<<NTracks/NTimes<<" tracks. "<<endl;
2860  cout<<" L1 total event size = "
2861  <<( HitSize + StripSize + StripSizeB + DupSize + TripSize + BranchSize + TrackSize )/k
2862  <<" Kb"<<endl;
2863 #endif // 0
2864 #endif
2865 
2866 #ifdef DRAW
2867  draw->ClearVeiw();
2868  // draw->DrawInputHits();
2869  // draw->DrawInfo();
2870  draw->DrawRecoTracks();
2871 
2872  draw->SaveCanvas("Reco_");
2873  draw->DrawAsk();
2874 #endif
2875 #ifdef PULLS
2876  static int iEvee = 0;
2877  iEvee++;
2878  if (iEvee % 1 == 0) fL1Pulls->Build(1);
2879 #endif
2880 #ifdef DOUB_PERFORMANCE
2881  fL1Eff_doublets->CalculateEff();
2882  fL1Eff_doublets->Print("Doublets performance.", 1);
2883 #endif
2884 #ifdef TRIP_PERFORMANCE
2885  fL1Eff_triplets->CalculateEff();
2886  fL1Eff_triplets->Print("Triplet performance", 1);
2887  // fL1Eff_triplets2->CalculateEff();
2888  // fL1Eff_triplets2->Print("Triplet performance. After chi2 cut");
2889 #endif
2890 }
2891 
2892 
2902 inline void L1Algo::CAFindTrack(int ista,
2903  L1Branch& best_tr,
2904  unsigned char& best_L,
2905  fscal& best_chi2,
2906  const L1Triplet* curr_trip,
2907  L1Branch& curr_tr,
2908  unsigned char& curr_L,
2909  fscal& curr_chi2,
2910  unsigned char min_best_l,
2911  L1Branch* new_tr)
2915 {
2916  if (curr_trip->GetLevel() == 0) // the end of the track -> check and store
2917  {
2918 
2919 
2920  // -- finish with current track
2921  // add rest of hits
2922  const THitI& ihitm = curr_trip->GetMHit();
2923  const THitI& ihitr = curr_trip->GetRHit();
2924 
2925 
2926  if (!GetFUsed((*vSFlag)[(*vStsHitsUnused)[ihitm].f]
2927  | (*vSFlagB)[(*vStsHitsUnused)[ihitm].b])) {
2928 
2929  // curr_tr.StsHits.push_back((*RealIHitP)[ihitm]);
2930 
2931  curr_tr.StsHits[curr_tr.NHits] = ((*RealIHitP)[ihitm]);
2932 
2933  curr_tr.NHits++;
2934 
2935  curr_L++;
2936  }
2937 
2938  if (!GetFUsed((*vSFlag)[(*vStsHitsUnused)[ihitr].f]
2939  | (*vSFlagB)[(*vStsHitsUnused)[ihitr].b])) {
2940 
2941  //curr_tr.StsHits.push_back((*RealIHitP)[ihitr]);
2942  curr_tr.StsHits[curr_tr.NHits] = ((*RealIHitP)[ihitr]);
2943 
2944  curr_tr.NHits++;
2945 
2946  curr_L++;
2947  }
2948 
2949  //if( curr_L < min_best_l - 1 ) return; // suppouse that only one hit can be added by extender
2950  if (curr_chi2 > TRACK_CHI2_CUT * (curr_L * 2 - 5)) return;
2951 
2952  if (fmCBMmode)
2953  if (curr_chi2 > TRACK_CHI2_CUT * (curr_L * 2 - 5.0)) return;
2954 
2955  // // try to find more hits
2956  // #ifdef EXTEND_TRACKS
2957  // // if( curr_L < min_best_l )
2958  // if (isec != kFastPrimJumpIter && isec != kAllPrimJumpIter && isec != kAllSecJumpIter && curr_L >= 3){
2959  // //curr_chi2 = BranchExtender(curr_tr);
2960  // BranchExtender(curr_tr);
2961  // curr_L = curr_tr.StsHits.size();
2962  // // if( 2*curr_chi2 > (2*(curr_L*2-5) + 1) * 4*4 ) return;
2963  // }
2964  // #endif // EXTEND_TRACKS
2965 
2966  // -- select the best
2967  if ((curr_L > best_L) || ((curr_L == best_L) && (curr_chi2 < best_chi2))) {
2968 
2969  best_tr = curr_tr;
2970  best_chi2 = curr_chi2;
2971  best_L = curr_L;
2972  }
2973 
2974  return;
2975  } else //MEANS level ! = 0
2976  {
2977  unsigned int Station = 0;
2978  unsigned int Thread = 0;
2979  unsigned int Triplet = 0;
2980  unsigned int Location = 0;
2981  int N_neighbour = (curr_trip->GetNNeighbours());
2982 
2983 
2984  for (Tindex in = 0; in < N_neighbour; in++) {
2985  Location = curr_trip->GetFNeighbour() + in;
2986 
2987  // Location = curr_trip->neighbours[in];
2988  // const fscal &qp2 = curr_trip->GetQp();
2989  // fscal &Cqp2 = curr_trip->Cqp;
2990  // if (( fabs(qp - qp2) > PickNeighbour * (Cqp + Cqp2) ) ) continue;
2991 
2992  Station = Location / 100000000;
2993  Thread = (Location - Station * 100000000) / 1000000;
2994  Triplet = (Location - Station * 100000000 - Thread * 1000000);
2995 
2996 
2997  const L1Triplet& new_trip = TripletsLocal1[Station][Thread][Triplet];
2998 
2999  if ((new_trip.GetMHit() != curr_trip->GetRHit())) continue;
3000  if ((new_trip.GetLHit() != curr_trip->GetMHit())) continue;
3001 
3002  const fscal& qp1 = curr_trip->GetQp();
3003  const fscal& qp2 = new_trip.GetQp();
3004  fscal dqp = fabs(qp1 - qp2);
3005  fscal Cqp = curr_trip->Cqp;
3006  Cqp += new_trip.Cqp;
3007 
3008  if (!fmCBMmode)
3009  if (dqp > PickNeighbour * Cqp)
3010  continue; // bad neighbour // CHECKME why do we need recheck it?? (it really change result)
3011 
3012  const fscal& tx1 = curr_trip->tx;
3013  const fscal& tx2 = new_trip.tx;
3014  fscal dtx = fabs(tx1 - tx2);
3015  fscal Ctx = curr_trip->Ctx;
3016  Ctx += new_trip.Ctx;
3017 
3018  const fscal& ty1 = curr_trip->ty;
3019  const fscal& ty2 = new_trip.ty;
3020  fscal dty = fabs(ty1 - ty2);
3021  fscal Cty = curr_trip->Cty;
3022  Cty += new_trip.Cty;
3023 
3024  if (fGlobal || fmCBMmode) {
3025  if (dty > PickNeighbour * Cty) continue;
3026  if (dtx > PickNeighbour * Ctx) continue;
3027  }
3028 
3029 
3030  if (GetFUsed((*vSFlag)[(*vStsHitsUnused)[new_trip.GetLHit()].f]
3031  | (*vSFlagB)[(*vStsHitsUnused)[new_trip.GetLHit()]
3032  .b])) { //hits are used
3033  // no used hits allowed -> compare and store track
3034  if ((curr_L > best_L)
3035  || ((curr_L == best_L) && (curr_chi2 < best_chi2))) {
3036  best_tr = curr_tr;
3037 
3038  best_chi2 = curr_chi2;
3039  best_L = curr_L;
3040  }
3041  } else { // if hits are used add new triplet to the current track
3042 
3043 
3044  new_tr[ista] = curr_tr;
3045 
3046 
3047  unsigned char new_L = curr_L;
3048  fscal new_chi2 = curr_chi2;
3049 
3050  // add new hit
3051  new_tr[ista].StsHits[new_tr[ista].NHits] =
3052  ((*RealIHitP)[new_trip.GetLHit()]);
3053  new_tr[ista].NHits++;
3054  new_L += 1;
3055  dqp = dqp / Cqp * 5.; // CHECKME: understand 5, why no sqrt(5)?
3056 
3057  dtx = dtx / Ctx;
3058  dty = dty / Cty;
3059 
3060  if (fGlobal || fmCBMmode) {
3061  new_chi2 += dtx * dtx;
3062  new_chi2 += dty * dty;
3063  } else
3064  new_chi2 += dqp * dqp;
3065 
3066  if (new_chi2 > TRACK_CHI2_CUT * new_L) continue;
3067 
3068 
3069  const int new_ista = ista + new_trip.GetMSta() - new_trip.GetLSta();
3070 
3071  CAFindTrack(new_ista,
3072  best_tr,
3073  best_L,
3074  best_chi2,
3075  &new_trip,
3076  new_tr[ista],
3077  new_L,
3078  new_chi2,
3079  min_best_l,
3080  new_tr);
3081  } // add triplet to track
3082  } // for neighbours
3083  } // level = 0
3084 }
3085 
3086 #ifdef DRAW
3087 void L1Algo::DrawRecoTracksTime(const vector<CbmL1Track>& tracks) {
3088  draw->DrawRecoTracksTime(tracks);
3089  draw->SaveCanvas(" ");
3090 }
3091 #endif
L1Algo::kAllPrimIter
@ kAllPrimIter
Definition: L1Algo.h:872
L1Algo::SavedCand
int SavedCand[nTh]
Definition: L1Algo.h:266
L1ExtrapolateLine
void L1ExtrapolateLine(L1TrackPar &T, fvec z_out)
Definition: L1Extrapolation.h:814
L1Algo::triplet
unsigned int & triplet
Definition: L1Algo.h:574
L1Algo::vStsHitsUnused
vector< L1StsHit > * vStsHitsUnused
Definition: L1Algo.h:421
L1TrackPar::C54
fvec C54
Definition: L1TrackPar.h:10
L1ExtrapolateJXY0
void L1ExtrapolateJXY0(fvec &tx, fvec &ty, fvec dz, L1FieldRegion &F, fvec &extrDx, fvec &extrDy, fvec &J04, fvec &J14)
Definition: L1Extrapolation.h:1022
L1FilterNoField
void L1FilterNoField(L1TrackPar &T, L1UMeasurementInfo &info, fvec u, fvec w=1.)
Definition: L1Filtration.h:147
L1Algo::fGhostSuppression
int fGhostSuppression
Definition: L1Algo.h:943
L1Algo::vStsHits
const vector< L1StsHit > * vStsHits
Definition: L1Algo.h:344
L1Portion.h
L1TrackPar::C10
fvec C10
Definition: L1TrackPar.h:9
L1Algo::TripletsLocal1
L1Vector< L1Triplet > TripletsLocal1[nSta][nTh]
Definition: L1Algo.h:259
fscal
float fscal
Definition: L1/vectors/P4_F32vec4.h:250
L1Algo::fmCBMmode
bool fmCBMmode
Definition: L1Algo.h:391
L1Algo::f4
void f4(Tindex n3, int istal, int istam, int istar, nsL1::vector< L1TrackPar >::TSimd &T_3, vector< THitI > &hitsl_3, vector< THitI > &hitsm_3, vector< THitI > &hitsr_3, Tindex &nstaltriplets, vector< THitI > *hitsn_3=0, vector< THitI > *hitsr_5=0)
Select triplets. Save them into vTriplets.
Definition: L1CATrackFinder.cxx:1173
L1TrackPar::qp
fvec qp
Definition: L1TrackPar.h:9
h
Generates beam ions for transport simulation.
Definition: CbmBeamGenerator.h:17
L1HitPoint
contain strips positions and coordinates of hit
Definition: L1HitPoint.h:6
L1HitPoint::dU
fscal dU() const
Definition: L1HitPoint.h:56
L1TrackPar::t
fvec t
Definition: L1TrackPar.h:9
L1Algo::vGridTime
L1Grid vGridTime[MaxNStations]
Definition: L1Algo.h:347
L1Algo::CAMergeClones
void CAMergeClones()
Definition: L1CAMergeClones.cxx:255
L1Algo::fNThreads
int fNThreads
Definition: L1Algo.h:389
L1Algo::fz_pos3
nsL1::vector< fvec >::TSimd fz_pos3[nTh]
Definition: L1Algo.h:443
L1Algo::vStsDontUsedHits_Buf
vector< L1StsHit > vStsDontUsedHits_Buf
Definition: L1Algo.h:368
f
float f
Definition: L1/vectors/P4_F32vec4.h:24
L1Algo.h
L1FilterChi2XYC00C10C11
void L1FilterChi2XYC00C10C11(const L1UMeasurementInfo &info, fvec &x, fvec &y, fvec &C00, fvec &C10, fvec &C11, fvec &chi2, const fvec &u)
Definition: L1Filtration.h:239
L1HitPoint::V
fscal V() const
Definition: L1HitPoint.h:60
L1Algo::dUdV_to_dY
void dUdV_to_dY(const fvec &u, const fvec &v, fvec &_y, const L1Station &sta)
Definition: L1Algo.cxx:279
L1TrackPar::C20
fvec C20
Definition: L1TrackPar.h:9
L1Branch::ista
unsigned char ista
Definition: L1Branch.h:41
L1AlgoEfficiencyPerformance.h
F32vec4
Definition: L1/vectors/P4_F32vec4.h:47
L1Triplet::Set
void Set(unsigned int iHitL, unsigned int iHitM, unsigned int iHitR, unsigned int iStaL, unsigned int iStaM, unsigned int iStaR, unsigned char Level, fscal Qp, fscal Chi2, fscal=0, fscal _Cqp=0, int=0)
Definition: L1Triplet.h:87
L1Triplet::Cty
fscal Cty
Definition: L1Triplet.h:67
L1Station::materialInfo
L1MaterialInfo materialInfo
Definition: L1Station.h:29
L1Algo::fu_back3
nsL1::vector< fvec >::TSimd fu_back3[nTh]
Definition: L1Algo.h:443
L1Algo::f20
void f20(Tindex n1, L1Station &stam, L1HitPoint *vStsHits_m, L1TrackPar *T_1, THitI *hitsl_1, Tindex &n2, vector< THitI > &i1_2, vector< THitI > &hitsm_2, fvec *Event, vector< bool > &lmDuplets)
Find the doublets. Reformat data in the portion of doublets.
Definition: L1CATrackFinder.cxx:435
L1Algo::vStripToTrackB
L1Vector< int > vStripToTrackB
Definition: L1Algo.h:387
L1XYMeasurementInfo::C11
fvec C11
Definition: L1XYMeasurementInfo.h:11
L1Algo::CandidatesTrack
L1Vector< L1Branch > CandidatesTrack[nTh]
Definition: L1Algo.h:260
L1Track::fTrackTime
float fTrackTime
Definition: L1Track.h:24
L1XYMeasurementInfo::C10
fvec C10
Definition: L1XYMeasurementInfo.h:11
L1TrackPar::C41
fvec C41
Definition: L1TrackPar.h:10
L1Algo::StripsToCoor
void StripsToCoor(const fscal &u, const fscal &v, fscal &x, fscal &y, const L1Station &sta) const
Definition: L1Algo.cxx:249
L1Triplet
Definition: L1Triplet.h:4
L1TrackPar::C51
fvec C51
Definition: L1TrackPar.h:10
L1CATFTimerInfo::GetTimerAll
L1CATFIterTimerInfo & GetTimerAll()
Definition: L1Timer.h:119
L1Station
Definition: L1Station.h:9
L1Station::fieldSlice
L1FieldSlice fieldSlice
Definition: L1Station.h:30
sqrt
friend F32vec4 sqrt(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:41
L1TrackPar::C30
fvec C30
Definition: L1TrackPar.h:9
L1AlgoPulls
Definition: L1AlgoPulls.h:81
L1TrackPar::C31
fvec C31
Definition: L1TrackPar.h:9
L1Algo::Pick_gather
float Pick_gather
parameters which are different for different iterations. Set in the begin of CAL1TrackFinder
Definition: L1Algo.h:922
L1TrackPar::C53
fvec C53
Definition: L1TrackPar.h:10
fvec
F32vec4 fvec
Definition: L1/vectors/P4_F32vec4.h:249
L1FilterChi2
void L1FilterChi2(const L1UMeasurementInfo &info, const fvec &x, const fvec &y, const fvec &C00, const fvec &C10, const fvec &C11, fvec &chi2, const fvec &u)
Definition: L1Filtration.h:217
FilterTime
void FilterTime(L1TrackPar &T, fvec t0, fvec dt0, fvec w=1.)
Definition: L1Filtration.h:13
L1Algo::StsHitsUnusedStopIndex
THitI StsHitsUnusedStopIndex[MaxNStations+1]
Definition: L1Algo.h:427
L1Extrapolate
void L1Extrapolate(L1TrackPar &T, fvec z_out, fvec qp0, const L1FieldRegion &F, fvec *w=0)
Definition: L1Extrapolation.h:314
L1Algo::f30
void f30(L1HitPoint *vStsHits_r, L1Station &stam, L1Station &star, int istam, int istar, L1HitPoint *vStsHits_m, L1TrackPar *T_1, L1FieldRegion *fld_1, THitI *hitsl_1, Tindex n2, vector< THitI > &hitsm_2, vector< THitI > &i1_2, const vector< bool > &mrDuplets, Tindex &n3, nsL1::vector< L1TrackPar >::TSimd &T_3, vector< THitI > &hitsl_3, vector< THitI > &hitsm_3, vector< THitI > &hitsr_3, nsL1::vector< fvec >::TSimd &u_front_3, nsL1::vector< fvec >::TSimd &u_back_3, nsL1::vector< fvec >::TSimd &z_Pos_3, nsL1::vector< fvec >::TSimd &du_, nsL1::vector< fvec >::TSimd &dv_, nsL1::vector< fvec >::TSimd &timeR, nsL1::vector< fvec >::TSimd &timeER)
Definition: L1CATrackFinder.cxx:609
L1Branch::compareCand
static bool compareCand(const L1Branch &a, const L1Branch &b)
Definition: L1Branch.h:87
L1Algo::vStsDontUsedHitsxy_B
vector< L1HitPoint > vStsDontUsedHitsxy_B
Definition: L1Algo.h:371
L1Algo::vStsDontUsedHits_B
vector< L1StsHit > vStsDontUsedHits_B
Definition: L1Algo.h:367
L1AlgoTBB.h
L1Algo::vRecoHits_local
L1Vector< THitI > vRecoHits_local[nTh]
Definition: L1Algo.h:373
L1AlgoEfficiencyPerformance::Init
void Init()
Definition: L1AlgoEfficiencyPerformance.h:156
L1Branch::StsHits
L1Vector< THitI > StsHits
Definition: L1Branch.h:50
L1Algo::StsHitsStartIndex
const THitI * StsHitsStartIndex
Definition: L1Algo.h:358
L1Triplet::GetFNeighbour
unsigned int GetFNeighbour() const
Definition: L1Triplet.h:124
L1HitPoint::Z
fscal Z() const
Definition: L1HitPoint.h:58
L1Algo::vStsDontUsedHitsxy_A
vector< L1HitPoint > vStsDontUsedHitsxy_A
Definition: L1Algo.h:369
L1AlgoDraw::InitL1Draw
void InitL1Draw(L1Algo *algo_)
Definition: L1AlgoDraw.h:147
L1Algo::dUdV_to_dX
void dUdV_to_dX(const fvec &u, const fvec &v, fvec &_x, const L1Station &sta)
Definition: L1Algo.cxx:287
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
L1Algo::vStsDontUsedHits_A
vector< L1StsHit > vStsDontUsedHits_A
Definition: L1Algo.h:366
L1Algo::portionStopIndex
Tindex portionStopIndex[nSta]
Definition: L1Algo.h:262
L1Algo::fTimeR
nsL1::vector< fvec >::TSimd fTimeR[nTh]
Definition: L1Algo.h:444
L1TrackPar::C42
fvec C42
Definition: L1TrackPar.h:10
L1Algo::thread
unsigned int unsigned int unsigned int & thread
Definition: L1Algo.h:576
L1Algo::MaxDZ
fvec MaxDZ
Definition: L1Algo.h:918
L1Triplet::GetRSta
int GetRSta() const
Definition: L1Triplet.h:160
L1Algo::vStsStrips
const vector< L1Strip > * vStsStrips
Definition: L1Algo.h:340
L1HitAreaTime::GetNext
bool GetNext(THitI &i)
Definition: L1HitArea.h:194
L1Algo::RealIHitPBuf
vector< THitI > * RealIHitPBuf
Definition: L1Algo.h:423
L1CATFTimerInfo::SetNIter
void SetNIter(int n)
Definition: L1Timer.h:112
L1Station::backInfo
L1UMeasurementInfo backInfo
Definition: L1Station.h:31
L1Algo::StsHitsStopIndex
const THitI * StsHitsStopIndex
Definition: L1Algo.h:359
L1ExtrapolateC10Line
void L1ExtrapolateC10Line(const L1TrackPar &T, fvec z_out, fvec &C10)
Definition: L1Extrapolation.h:855
L1TrackPar::z
fvec z
Definition: L1TrackPar.h:9
L1Algo::Portion
@ Portion
Definition: L1Algo.h:405
L1UMeasurementInfo::sigma2
fvec sigma2
Definition: L1UMeasurementInfo.h:12
L1FilterVtx
void L1FilterVtx(L1TrackPar &T, fvec x, fvec y, L1XYMeasurementInfo &info, fvec extrDx, fvec extrDy, fvec J[])
Definition: L1Filtration.h:273
L1TrackPar::SetOneEntry
void SetOneEntry(const int i0, const L1TrackPar &T1, const int i1)
Definition: L1TrackPar.h:122
L1FilterXY
void L1FilterXY(L1TrackPar &T, fvec x, fvec y, L1XYMeasurementInfo &info)
Definition: L1Filtration.h:358
L1Algo::vTracks_local
L1Vector< L1Track > vTracks_local[nTh]
Definition: L1Algo.h:372
L1CATFTimerInfo::Calc
void Calc()
Definition: L1Timer.h:136
L1StsHit
Definition: L1StsHit.h:12
L1Algo::StsHitsUnusedStartIndex
THitI StsHitsUnusedStartIndex[MaxNStations+1]
Definition: L1Algo.h:426
L1Algo::fu_front3
nsL1::vector< fvec >::TSimd fu_front3[nTh]
Definition: L1Algo.h:443
L1Timer.h
L1Triplet::GetRHit
THitI GetRHit() const
Definition: L1Triplet.h:121
L1Algo::nTripletsThread
int nTripletsThread[nSta][nTh]
Definition: L1Algo.h:271
L1Algo::TripletsStaPort
void TripletsStaPort(int istal, int istam, int istar, Tindex &nstaltriplets, L1TrackPar *T_1, L1FieldRegion *fld_1, THitI *hitsl_1, Tindex &n_2, vector< THitI > &i1_2, vector< THitI > &hitsm_2, const vector< bool > &mrDuplets)
Find triplets on station.
Definition: L1CATrackFinder.cxx:1576
L1Algo::targY
fvec targY
Definition: L1Algo.h:927
L1TrackPar::ty
fvec ty
Definition: L1TrackPar.h:9
L1CATFTimerInfo::Add
void Add(string name)
Definition: L1Timer.h:114
L1TrackPar::C32
fvec C32
Definition: L1TrackPar.h:9
L1AddMaterial
void L1AddMaterial(L1TrackPar &T, fvec radThick, fvec qp0, fvec w=1, fvec mass2=0.10565f *0.10565f)
Definition: L1AddMaterial.h:559
L1Algo::vStsStripsB
const vector< L1Strip > * vStsStripsB
Definition: L1Algo.h:341
L1HitArea.h
L1Algo::fUseHitErrors
bool fUseHitErrors
Definition: L1Algo.h:390
L1Algo::MaxInvMom
fvec MaxInvMom
Definition: L1Algo.h:925
L1Algo::NFStations
int NFStations
Definition: L1Algo.h:335
L1Track.h
L1StsHit::b
TStripI b
Definition: L1StsHit.h:15
L1Algo::SavedHits
int SavedHits[nTh]
Definition: L1Algo.h:267
L1AlgoEfficiencyPerformance
Definition: L1AlgoEfficiencyPerformance.h:123
L1Algo::BranchExtender
fscal BranchExtender(L1Branch &t)
Try to extrapolate and find additional hits on other stations.
Definition: L1TrackExtender.cxx:396
L1Algo::vStsHitPointsUnused
vector< L1HitPoint > * vStsHitPointsUnused
Definition: L1Algo.h:424
L1Algo::kFastPrimJumpIter
@ kFastPrimJumpIter
Definition: L1Algo.h:878
L1TrackPar::y
fvec y
Definition: L1TrackPar.h:9
L1Triplet::SetFNeighbour
void SetFNeighbour(unsigned int n)
Definition: L1Triplet.h:125
L1Algo::kAllPrimEIter
@ kAllPrimEIter
Definition: L1Algo.h:875
L1Triplet::Cqp
fscal Cqp
Definition: L1Triplet.h:66
L1ExtrapolateTime
void L1ExtrapolateTime(L1TrackPar &T, fvec dz)
Definition: L1Extrapolation.h:781
L1Algo::NTracksIsecAll
unsigned int NTracksIsecAll
Definition: L1Algo.h:364
L1_assert
#define L1_assert(v)
Definition: CbmL1Def.h:56
THitI
unsigned int THitI
Definition: L1StsHit.h:8
L1Triplet::GetChi2
fscal GetChi2() const
Definition: L1Triplet.h:141
L1Algo::vSFlagB
const vector< unsigned char > * vSFlagB
Definition: L1Algo.h:351
L1Algo::NStations
int NStations
Definition: L1Algo.h:333
L1Triplet::Ctx
fscal Ctx
Definition: L1Triplet.h:67
finite
T finite(T x)
Definition: CbmL1Def.h:21
L1TrackPar::C44
fvec C44
Definition: L1TrackPar.h:10
L1Branch::NHits
char NHits
Definition: L1Branch.h:44
L1Algo::targX
fvec targX
Definition: L1Algo.h:927
L1Triplet::tx
fscal tx
Definition: L1Triplet.h:67
L1Triplet::ty
fscal ty
Definition: L1Triplet.h:67
L1Branch::Set
void Set(unsigned char iStation, unsigned char Length, float Chi2, float Qp)
Definition: L1Branch.h:117
L1HitPoint::time
float time
Definition: L1HitPoint.h:101
L1FieldRegion
Definition: L1Field.h:85
L1Algo::fT_3
nsL1::vector< L1TrackPar >::TSimd fT_3[nTh]
Definition: L1Algo.h:439
fvecLen
const int fvecLen
Definition: L1/vectors/P4_F32vec4.h:251
L1Algo::kFastPrimIter2
@ kFastPrimIter2
Definition: L1Algo.h:879
L1Station::z
fvec z
Definition: L1Station.h:28
L1Algo::NHitsIsecAll
unsigned int NHitsIsecAll
Definition: L1Algo.h:363
L1AddPipeMaterial
void L1AddPipeMaterial(L1TrackPar &T, fvec qp0, fvec w=1, fvec mass2=0.10565f *0.10565f)
Definition: L1AddMaterial.h:737
L1TrackPar::chi2
fvec chi2
Definition: L1TrackPar.h:10
L1ExtrapolateXC00Line
void L1ExtrapolateXC00Line(const L1TrackPar &T, fvec z_out, fvec &x, fvec &C00)
Definition: L1Extrapolation.h:842
L1Triplet::SetLevel
void SetLevel(unsigned char Level)
Definition: L1Triplet.h:116
L1CATFIterTimerInfo
Definition: L1Timer.h:62
L1FieldRegion::Set
void Set(const L1FieldValue &B0, const fvec B0z, const L1FieldValue &B1, const fvec B1z, const L1FieldValue &B2, const fvec B2z)
Definition: L1Field.h:135
L1Algo::fTimeER
nsL1::vector< fvec >::TSimd fTimeER[nTh]
Definition: L1Algo.h:444
tracks
TClonesArray * tracks
Definition: Analyze_matching.h:17
L1Filter
void L1Filter(L1TrackPar &T, L1UMeasurementInfo &info, fvec u, fvec w=1.)
Definition: L1Filtration.h:80
L1TrackPar::NDF
fvec NDF
Definition: L1TrackPar.h:10
L1AlgoPulls.h
L1Track::NHits
unsigned char NHits
Definition: L1Track.h:22
L1Algo::dv
nsL1::vector< fvec >::TSimd dv[nTh]
Definition: L1Algo.h:444
L1Algo::numberCandidateThread
int numberCandidateThread[nTh]
Definition: L1Algo.h:269
L1Algo::PickNeighbour
float PickNeighbour
Definition: L1Algo.h:924
L1ExtrapolateYC11Line
void L1ExtrapolateYC11Line(const L1TrackPar &T, fvec z_out, fvec &y, fvec &C11)
Definition: L1Extrapolation.h:849
L1Triplet::GetLevel
unsigned char GetLevel() const
Definition: L1Triplet.h:134
L1Algo::targZ
fvec targZ
Definition: L1Algo.h:927
L1Grid::UpdateIterGrid
void UpdateIterGrid(unsigned int Nelements, L1StsHit *hits, vector< THitI > *indicesBuf, THitI *indices, vector< L1StsHit > *hits2, vector< L1HitPoint > *pointsBuf, L1HitPoint *points, int &NHitsOnStation, char iS, L1Algo &Algo, const vector< unsigned char > *vSFlag, const vector< unsigned char > *vSFlagB)
Definition: L1Grid.cxx:34
L1Branch::chi2
fscal chi2
Definition: L1Branch.h:46
L1Grid::BuildBins
void BuildBins(float yMin, float yMax, float zMin, float zMax, float tMin, float tMax, float sy, float sz, float st)
Definition: L1Grid.cxx:146
L1TrackPar::tx
fvec tx
Definition: L1TrackPar.h:9
L1HitPoint::dV
fscal dV() const
Definition: L1HitPoint.h:57
L1Algo::f10
void f10(Tindex start_lh, Tindex n1_l, L1HitPoint *StsHits_l, fvec *u_front_l, fvec *u_back_l, fvec *zPos_l, THitI *hitsl, fvec *HitTime_l, fvec *HitTimeEr, fvec *Event_l, fvec *d_x, fvec *d_y, fvec *d_xy, fvec *d_u, fvec *d_v)
Prepare the portion of left hits data.
Definition: L1CATrackFinder.cxx:73
L1TrackPar.h
L1Algo::TRACK_CHI2_CUT
float TRACK_CHI2_CUT
Definition: L1Algo.h:911
L1Algo::fhitsm_3
vector< THitI > fhitsm_3[nTh]
Definition: L1Algo.h:441
L1HitPoint.h
L1TrackPar::C50
fvec C50
Definition: L1TrackPar.h:10
L1Algo::CAFindTrack
void CAFindTrack(int ista, L1Branch &best_tr, unsigned char &best_L, fscal &best_chi2, const L1Triplet *curr_trip, L1Branch &curr_tr, unsigned char &curr_L, fscal &curr_chi2, unsigned char min_best_l, L1Branch *new_tr)
================================= FUNCTIONAL PART =================================
Definition: L1CATrackFinder.cxx:2902
L1HitAreaTime
Definition: L1HitArea.h:104
L1TrackPar::C21
fvec C21
Definition: L1TrackPar.h:9
L1Algo::vSFlag
const vector< unsigned char > * vSFlag
Definition: L1Algo.h:350
L1Algo::vTracks
L1Vector< L1Track > vTracks
Definition: L1Algo.h:355
L1Algo::kAllSecEIter
@ kAllSecEIter
Definition: L1Algo.h:876
L1Branch
Definition: L1Branch.h:27
L1Algo::vStripToTrack
L1Vector< int > vStripToTrack
Definition: L1Algo.h:386
L1TrackPar::x
fvec x
Definition: L1TrackPar.h:9
L1Algo::MaxSlope
fvec MaxSlope
Definition: L1Algo.h:926
L1Algo::RealIHit_v_buf
vector< THitI > RealIHit_v_buf
Definition: L1Algo.h:376
L1TrackPar::C00
fvec C00
Definition: L1TrackPar.h:9
L1Triplet::GetLHit
THitI GetLHit() const
Definition: L1Triplet.h:119
L1StsHit::t_reco
float t_reco
Definition: L1StsHit.h:17
L1Algo::FIRSTCASTATION
Tindex FIRSTCASTATION
================================= DATA PART =================================
Definition: L1Algo.h:865
L1FieldSlice::GetFieldValue
void GetFieldValue(const fvec &x, const fvec &y, L1FieldValue &B) const
Definition: L1Field.h:41
L1Extrapolation.h
L1Algo::fhitsl_3
vector< THitI > fhitsl_3[nTh]
Definition: L1Algo.h:441
L1Algo::TRIPLET_CHI2_CUT
float TRIPLET_CHI2_CUT
Definition: L1Algo.h:913
L1Algo::fhitsr_3
vector< THitI > fhitsr_3[nTh]
Definition: L1Algo.h:441
L1UMeasurementInfo
Definition: L1UMeasurementInfo.h:7
L1Branch::CandIndex
int CandIndex
Definition: L1Branch.h:47
L1Algo::RealIHit_v
vector< THitI > RealIHit_v
Definition: L1Algo.h:375
L1Algo::TripForHit
vector< int > TripForHit[2]
Definition: L1Algo.h:432
L1Triplet::GetMSta
int GetMSta() const
Definition: L1Triplet.h:158
L1TrackPar::C33
fvec C33
Definition: L1TrackPar.h:9
L1Track
Definition: L1Track.h:20
L1Triplet::GetQp
fscal GetQp() const
Definition: L1Triplet.h:139
L1Triplet::GetNNeighbours
char GetNNeighbours() const
Definition: L1Triplet.h:122
L1TrackPar::C52
fvec C52
Definition: L1TrackPar.h:10
L1HitsSortHelper.h
L1Algo::CATime
double CATime
Definition: L1Algo.h:353
L1Algo::DupletsStaPort
void DupletsStaPort(int istal, int istam, Tindex ip, vector< Tindex > &n_g, Tindex *portionStopIndex_, L1TrackPar *T_1, L1FieldRegion *fld_1, THitI *hitsl_1, vector< bool > &lmDuplets, Tindex &n_2, vector< THitI > &i1_2, vector< THitI > &hitsm_2)
Find doublets on station.
Definition: L1CATrackFinder.cxx:1446
L1CATFTimerInfo
Definition: L1Timer.h:109
L1Algo::vRecoHits
L1Vector< THitI > vRecoHits
Definition: L1Algo.h:356
L1Algo::RealIHitP
vector< THitI > * RealIHitP
Definition: L1Algo.h:422
L1Algo::NMvdStations
int NMvdStations
Definition: L1Algo.h:334
Tindex
int Tindex
Definition: L1Algo.h:76
v
__m128 v
Definition: L1/vectors/P4_F32vec4.h:1
L1Algo::_fvecalignment
L1Station vStations[MaxNStations] _fvecalignment
Definition: L1Algo.h:336
L1AddMaterial.h
L1CATFIterTimerInfo::Add
void Add(string name)
Definition: L1Timer.h:65
L1FieldRegion::SetOneEntry
void SetOneEntry(const int i0, const L1FieldRegion &f1, const int i1)
Definition: L1Field.h:197
L1Algo::nTh
static const int nTh
Definition: L1Algo.h:256
x
Double_t x
Definition: CbmMvdSensorDigiToHitTask.cxx:68
L1Algo::n_g1
L1Vector< Tindex > n_g1
Definition: L1Algo.h:263
L1TrackPar
Definition: L1TrackPar.h:6
fabs
friend F32vec4 fabs(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:60
L1AlgoPulls::Init
void Init()
Definition: L1AlgoPulls.h:115
L1Algo::f31
void f31(Tindex n3_V, L1Station &star, nsL1::vector< fvec >::TSimd &u_front_3, nsL1::vector< fvec >::TSimd &u_back_3, nsL1::vector< fvec >::TSimd &z_Pos_3, nsL1::vector< fvec >::TSimd &du_, nsL1::vector< fvec >::TSimd &dv_, nsL1::vector< fvec >::TSimd &timeR, nsL1::vector< fvec >::TSimd &timeER, nsL1::vector< L1TrackPar >::TSimd &T_3)
Add the right hits to parameters estimation.
Definition: L1CATrackFinder.cxx:950
nsL1::vector::TSimd
std::vector< T > TSimd
Definition: L1/vectors/PSEUDO_F32vec1.h:148
L1TrackPar::C55
fvec C55
Definition: L1TrackPar.h:10
y
Double_t y
Definition: CbmMvdSensorDigiToHitTask.cxx:68
L1TrackPar::C40
fvec C40
Definition: L1TrackPar.h:10
NS_L1TrackFitter::vINF
const fvec vINF
Definition: L1TrackFitter.cxx:34
L1TrackPar::C11
fvec C11
Definition: L1TrackPar.h:9
L1Algo::dUdV_to_dXdY
void dUdV_to_dXdY(const fvec &u, const fvec &v, fvec &_xy, const L1Station &sta)
Definition: L1Algo.cxx:295
L1Algo::kAllPrimJumpIter
@ kAllPrimJumpIter
Definition: L1Algo.h:873
L1Algo::kFastPrimIter
@ kFastPrimIter
Definition: L1Algo.h:871
L1Branch.h
L1Triplet::GetMHit
THitI GetMHit() const
Definition: L1Triplet.h:120
L1Algo::fNFindIterations
@ fNFindIterations
Definition: L1Algo.h:889
L1Algo::kAllSecIter
@ kAllSecIter
Definition: L1Algo.h:874
L1Triplet::GetQpOrig
fscal GetQpOrig()
Definition: L1Triplet.h:154
L1CATFIterTimerInfo::PrintReal
void PrintReal(int f=0)
Definition: L1Timer.h:85
L1Triplet::SetNNeighbours
void SetNNeighbours(unsigned char n)
Definition: L1Triplet.h:123
L1Algo::MaxNStations
@ MaxNStations
Definition: L1Algo.h:331
L1Algo::du
nsL1::vector< fvec >::TSimd du[nTh]
Definition: L1Algo.h:444
L1Algo::f11
void f11(int istal, int istam, Tindex n1_V, fvec *u_front_l, fvec *u_back_l, fvec *zPos_l, fvec *HitTime_l, fvec *HitTimeEr, L1TrackPar *T_1, L1FieldRegion *fld_1, fvec *d_x, fvec *d_y, fvec *d_xy, fvec *d_u, fvec *d_v)
Get the field approximation. Add the target to parameters estimation. Propagate to middle station.
Definition: L1CATrackFinder.cxx:124
L1Algo::f32
void f32(Tindex n3, int istal, nsL1::vector< L1TrackPar >::TSimd &T_3, vector< THitI > &hitsl_3, vector< THitI > &hitsm_3, vector< THitI > &hitsr_3, int nIterations=0)
Refit Triplets.
Definition: L1CATrackFinder.cxx:1004
L1AlgoDraw
Definition: L1AlgoDraw.h:29
L1Algo::isec
int isec
— data used during finding iterations
Definition: L1Algo.h:420
L1Algo::fRadThick
vector< L1Material > fRadThick
Definition: L1Algo.h:337
L1Extrapolate0
void L1Extrapolate0(L1TrackPar &T, fvec z_out, L1FieldRegion &F)
Definition: L1Extrapolation.h:681
L1HitPoint::U
fscal U() const
Definition: L1HitPoint.h:59
L1TrackPar::C43
fvec C43
Definition: L1TrackPar.h:10
L1CATFTimerInfo::PrintReal
void PrintReal()
Definition: L1Timer.h:143
L1FieldValue
Definition: L1Field.h:11
L1XYMeasurementInfo::C00
fvec C00
Definition: L1XYMeasurementInfo.h:11
L1StsHit::f
TStripI f
Definition: L1StsHit.h:15
L1StsHit::iz
TZPosI iz
Definition: L1StsHit.h:24
L1Grid::StoreHits
void StoreHits(THitI nhits, const L1StsHit *hits, char iS, L1Algo &Algo, THitI n, L1StsHit *hitsBuf1, const L1StsHit *hits1, THitI *indices1)
Definition: L1Grid.cxx:175
L1AlgoDraw.h
L1Algo::CreateHitPoint
L1HitPoint CreateHitPoint(const L1StsHit &hit, char ista)
full the hit point by hit information: takes hit as input (2 strips) and creates hit_point with all c...
Definition: L1Algo.cxx:318
L1Algo::DOUBLET_CHI2_CUT
float DOUBLET_CHI2_CUT
Definition: L1Algo.h:914
L1Triplet::GetLSta
int GetLSta() const
Definition: L1Triplet.h:156
L1Algo::fGlobal
bool fGlobal
Definition: L1Algo.h:392
L1Algo::CATrackFinder
void CATrackFinder()
The main procedure - find tracks.
Definition: L1CATrackFinder.cxx:1786
L1Filtration.h
L1_ASSERT
#define L1_ASSERT(v, msg)
Definition: CbmL1Def.h:48
L1HitPoint::timeEr
float timeEr
Definition: L1HitPoint.h:101
L1Station::frontInfo
L1UMeasurementInfo frontInfo
Definition: L1Station.h:31
L1Algo::f5
void f5(int *nlevel)
Find neighbours of triplets. Calculate level of triplets.
Definition: L1CATrackFinder.cxx:1338
L1Grid.h
L1TrackPar::C22
fvec C22
Definition: L1TrackPar.h:9
L1Algo::kAllSecJumpIter
@ kAllSecJumpIter
Definition: L1Algo.h:880
L1Station::XYInfo
L1XYMeasurementInfo XYInfo
Definition: L1Station.h:32