CMS 3D CMS Logo

PixelTriplet.h
Go to the documentation of this file.
1 #ifndef RecoTracker_LSTCore_src_alpaka_PixelTriplet_h
2 #define RecoTracker_LSTCore_src_alpaka_PixelTriplet_h
3 
12 
14 
15  template <typename TAcc>
16  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool runTripletDefaultAlgoPPBB(TAcc const& acc,
20  SegmentsConst segments,
21  SegmentsPixelConst segmentsPixel,
22  uint16_t pixelModuleIndex,
23  uint16_t outerInnerLowerModuleIndex,
24  uint16_t outerOuterLowerModuleIndex,
25  unsigned int innerSegmentIndex,
26  unsigned int outerSegmentIndex,
27  unsigned int firstMDIndex,
28  unsigned int secondMDIndex,
29  unsigned int thirdMDIndex,
30  unsigned int fourthMDIndex);
31  template <typename TAcc>
32  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool runTripletDefaultAlgoPPEE(TAcc const& acc,
36  SegmentsConst segments,
37  SegmentsPixelConst segmentsPixel,
38  uint16_t pixelModuleIndex,
39  uint16_t outerInnerLowerModuleIndex,
40  uint16_t outerOuterLowerModuleIndex,
41  unsigned int innerSegmentIndex,
42  unsigned int outerSegmentIndex,
43  unsigned int firstMDIndex,
44  unsigned int secondMDIndex,
45  unsigned int thirdMDIndex,
46  unsigned int fourthMDIndex);
47 
48  ALPAKA_FN_ACC ALPAKA_FN_INLINE void addPixelTripletToMemory(MiniDoubletsConst mds,
49  SegmentsConst segments,
50  TripletsConst triplets,
51  PixelTriplets pixelTriplets,
52  unsigned int pixelSegmentIndex,
53  unsigned int tripletIndex,
54  float pixelRadius,
55  float tripletRadius,
56  float centerX,
57  float centerY,
58  float rPhiChiSquared,
59  float rPhiChiSquaredInwards,
60  float rzChiSquared,
61  unsigned int pixelTripletIndex,
62  float pt,
63  float eta,
64  float phi,
65  float eta_pix,
66  float phi_pix,
67  float score) {
68  pixelTriplets.pixelSegmentIndices()[pixelTripletIndex] = pixelSegmentIndex;
69  pixelTriplets.tripletIndices()[pixelTripletIndex] = tripletIndex;
70  pixelTriplets.pixelRadius()[pixelTripletIndex] = __F2H(pixelRadius);
71  pixelTriplets.tripletRadius()[pixelTripletIndex] = __F2H(tripletRadius);
72  pixelTriplets.pt()[pixelTripletIndex] = __F2H(pt);
73  pixelTriplets.eta()[pixelTripletIndex] = __F2H(eta);
74  pixelTriplets.phi()[pixelTripletIndex] = __F2H(phi);
75  pixelTriplets.eta_pix()[pixelTripletIndex] = __F2H(eta_pix);
76  pixelTriplets.phi_pix()[pixelTripletIndex] = __F2H(phi_pix);
77  pixelTriplets.isDup()[pixelTripletIndex] = false;
78  pixelTriplets.score()[pixelTripletIndex] = __F2H(score);
79 
80  pixelTriplets.centerX()[pixelTripletIndex] = __F2H(centerX);
81  pixelTriplets.centerY()[pixelTripletIndex] = __F2H(centerY);
82  pixelTriplets.logicalLayers()[pixelTripletIndex][0] = 0;
83  pixelTriplets.logicalLayers()[pixelTripletIndex][1] = 0;
84  pixelTriplets.logicalLayers()[pixelTripletIndex][2] = triplets.logicalLayers()[tripletIndex][0];
85  pixelTriplets.logicalLayers()[pixelTripletIndex][3] = triplets.logicalLayers()[tripletIndex][1];
86  pixelTriplets.logicalLayers()[pixelTripletIndex][4] = triplets.logicalLayers()[tripletIndex][2];
87 
88  pixelTriplets.lowerModuleIndices()[pixelTripletIndex][0] = segments.innerLowerModuleIndices()[pixelSegmentIndex];
89  pixelTriplets.lowerModuleIndices()[pixelTripletIndex][1] = segments.outerLowerModuleIndices()[pixelSegmentIndex];
90  pixelTriplets.lowerModuleIndices()[pixelTripletIndex][2] = triplets.lowerModuleIndices()[tripletIndex][0];
91  pixelTriplets.lowerModuleIndices()[pixelTripletIndex][3] = triplets.lowerModuleIndices()[tripletIndex][1];
92  pixelTriplets.lowerModuleIndices()[pixelTripletIndex][4] = triplets.lowerModuleIndices()[tripletIndex][2];
93 
94  unsigned int pixelInnerMD = segments.mdIndices()[pixelSegmentIndex][0];
95  unsigned int pixelOuterMD = segments.mdIndices()[pixelSegmentIndex][1];
96 
97  pixelTriplets.hitIndices()[pixelTripletIndex][0] = mds.anchorHitIndices()[pixelInnerMD];
98  pixelTriplets.hitIndices()[pixelTripletIndex][1] = mds.outerHitIndices()[pixelInnerMD];
99  pixelTriplets.hitIndices()[pixelTripletIndex][2] = mds.anchorHitIndices()[pixelOuterMD];
100  pixelTriplets.hitIndices()[pixelTripletIndex][3] = mds.outerHitIndices()[pixelOuterMD];
101 
102  pixelTriplets.hitIndices()[pixelTripletIndex][4] = triplets.hitIndices()[tripletIndex][0];
103  pixelTriplets.hitIndices()[pixelTripletIndex][5] = triplets.hitIndices()[tripletIndex][1];
104  pixelTriplets.hitIndices()[pixelTripletIndex][6] = triplets.hitIndices()[tripletIndex][2];
105  pixelTriplets.hitIndices()[pixelTripletIndex][7] = triplets.hitIndices()[tripletIndex][3];
106  pixelTriplets.hitIndices()[pixelTripletIndex][8] = triplets.hitIndices()[tripletIndex][4];
107  pixelTriplets.hitIndices()[pixelTripletIndex][9] = triplets.hitIndices()[tripletIndex][5];
108  pixelTriplets.rPhiChiSquared()[pixelTripletIndex] = rPhiChiSquared;
109  pixelTriplets.rPhiChiSquaredInwards()[pixelTripletIndex] = rPhiChiSquaredInwards;
110  pixelTriplets.rzChiSquared()[pixelTripletIndex] = rzChiSquared;
111  }
112 
113  template <typename TAcc>
114  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool runPixelTrackletDefaultAlgopT3(TAcc const& acc,
117  MiniDoubletsConst mds,
118  SegmentsConst segments,
119  SegmentsPixelConst segmentsPixel,
120  uint16_t pixelLowerModuleIndex,
121  uint16_t outerInnerLowerModuleIndex,
122  uint16_t outerOuterLowerModuleIndex,
123  unsigned int innerSegmentIndex,
124  unsigned int outerSegmentIndex) {
125  short outerInnerLowerModuleSubdet = modules.subdets()[outerInnerLowerModuleIndex];
126  short outerOuterLowerModuleSubdet = modules.subdets()[outerOuterLowerModuleIndex];
127 
128  unsigned int firstMDIndex = segments.mdIndices()[innerSegmentIndex][0];
129  unsigned int secondMDIndex = segments.mdIndices()[innerSegmentIndex][1];
130 
131  unsigned int thirdMDIndex = segments.mdIndices()[outerSegmentIndex][0];
132  unsigned int fourthMDIndex = segments.mdIndices()[outerSegmentIndex][1];
133 
134  if (outerInnerLowerModuleSubdet == Barrel and
135  (outerOuterLowerModuleSubdet == Barrel or outerOuterLowerModuleSubdet == Endcap)) {
136  return runTripletDefaultAlgoPPBB(acc,
137  modules,
138  ranges,
139  mds,
140  segments,
141  segmentsPixel,
142  pixelLowerModuleIndex,
143  outerInnerLowerModuleIndex,
144  outerOuterLowerModuleIndex,
145  innerSegmentIndex,
146  outerSegmentIndex,
147  firstMDIndex,
148  secondMDIndex,
149  thirdMDIndex,
150  fourthMDIndex);
151  } else if (outerInnerLowerModuleSubdet == Endcap and outerOuterLowerModuleSubdet == Endcap) {
152  return runTripletDefaultAlgoPPEE(acc,
153  modules,
154  ranges,
155  mds,
156  segments,
157  segmentsPixel,
158  pixelLowerModuleIndex,
159  outerInnerLowerModuleIndex,
160  outerOuterLowerModuleIndex,
161  innerSegmentIndex,
162  outerSegmentIndex,
163  firstMDIndex,
164  secondMDIndex,
165  thirdMDIndex,
166  fourthMDIndex);
167  }
168  return false;
169  }
170 
171  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passPT3RZChiSquaredCuts(ModulesConst modules,
172  uint16_t lowerModuleIndex1,
173  uint16_t lowerModuleIndex2,
174  uint16_t lowerModuleIndex3,
175  float rzChiSquared) {
176  const int layer1 =
177  modules.layers()[lowerModuleIndex1] + 6 * (modules.subdets()[lowerModuleIndex1] == Endcap) +
178  5 * (modules.subdets()[lowerModuleIndex1] == Endcap and modules.moduleType()[lowerModuleIndex1] == TwoS);
179  const int layer2 =
180  modules.layers()[lowerModuleIndex2] + 6 * (modules.subdets()[lowerModuleIndex2] == Endcap) +
181  5 * (modules.subdets()[lowerModuleIndex2] == Endcap and modules.moduleType()[lowerModuleIndex2] == TwoS);
182  const int layer3 =
183  modules.layers()[lowerModuleIndex3] + 6 * (modules.subdets()[lowerModuleIndex3] == Endcap) +
184  5 * (modules.subdets()[lowerModuleIndex3] == Endcap and modules.moduleType()[lowerModuleIndex3] == TwoS);
185 
186  if (layer1 == 8 and layer2 == 9 and layer3 == 10) {
187  return rzChiSquared < 13.6067f;
188  } else if (layer1 == 8 and layer2 == 9 and layer3 == 15) {
189  return rzChiSquared < 5.5953f;
190  } else if (layer1 == 7 and layer2 == 8 and layer3 == 9) {
191  return rzChiSquared < 3.9263f;
192  }
193  /*
194  else if(layer1 == 7 and layer2 == 8 and layer3 == 14)
195  {
196  // PS+PS+2S in endcap layers 1+2+3, which is not really feasible in the current geometry,
197  // without skipping barrel layers 1 and 2 (not allowed by algorithm logic).
198  }
199  */
200  else if (layer1 == 1 and layer2 == 2 and layer3 == 3) {
201  return rzChiSquared < 9.4377f;
202  } else if (layer1 == 1 and layer2 == 2 and layer3 == 7) {
203  return rzChiSquared < 9.9975f;
204  } else if (layer1 == 1 and layer2 == 7 and layer3 == 8) {
205  return rzChiSquared < 8.6369f;
206  } else if (layer1 == 2 and layer2 == 3 and layer3 == 7) {
207  return rzChiSquared < 37.945f;
208  } else if (layer1 == 2 and layer2 == 3 and layer3 == 12) {
209  return rzChiSquared < 43.0167f;
210  } else if (layer1 == 2 and layer2 == 3 and layer3 == 4) {
211  return rzChiSquared < 8.6923f;
212  } else if (layer1 == 2 and layer2 == 7 and layer3 == 8) {
213  return rzChiSquared < 11.9672f;
214  } else if (layer1 == 2 and layer2 == 7 and layer3 == 13) {
215  return rzChiSquared < 16.2133f;
216  }
217 
218  //default - category not found!
219  return true;
220  }
221 
222  template <typename TAcc>
223  ALPAKA_FN_ACC ALPAKA_FN_INLINE float computeChiSquaredpT3(TAcc const& acc,
224  unsigned int nPoints,
225  float* xs,
226  float* ys,
227  float* delta1,
228  float* delta2,
229  float* slopes,
230  bool* isFlat,
231  float g,
232  float f,
233  float radius) {
234  //given values of (g, f, radius) and a set of points (and its uncertainties)
235  //compute chi squared
236  float c = g * g + f * f - radius * radius;
237  float chiSquared = 0.f;
238  float absArctanSlope, angleM, xPrime, yPrime, sigma2;
239  for (size_t i = 0; i < nPoints; i++) {
240  absArctanSlope = ((slopes[i] != kVerticalModuleSlope) ? alpaka::math::abs(acc, alpaka::math::atan(acc, slopes[i]))
241  : kPi / 2.f);
242  if (xs[i] > 0 and ys[i] > 0) {
243  angleM = kPi / 2.f - absArctanSlope;
244  } else if (xs[i] < 0 and ys[i] > 0) {
245  angleM = absArctanSlope + kPi / 2.f;
246  } else if (xs[i] < 0 and ys[i] < 0) {
247  angleM = -(absArctanSlope + kPi / 2.f);
248  } else if (xs[i] > 0 and ys[i] < 0) {
249  angleM = -(kPi / 2.f - absArctanSlope);
250  } else {
251  angleM = 0;
252  }
253 
254  if (not isFlat[i]) {
255  xPrime = xs[i] * alpaka::math::cos(acc, angleM) + ys[i] * alpaka::math::sin(acc, angleM);
256  yPrime = ys[i] * alpaka::math::cos(acc, angleM) - xs[i] * alpaka::math::sin(acc, angleM);
257  } else {
258  xPrime = xs[i];
259  yPrime = ys[i];
260  }
261  sigma2 = 4 * ((xPrime * delta1[i]) * (xPrime * delta1[i]) + (yPrime * delta2[i]) * (yPrime * delta2[i]));
262  chiSquared += (xs[i] * xs[i] + ys[i] * ys[i] - 2 * g * xs[i] - 2 * f * ys[i] + c) *
263  (xs[i] * xs[i] + ys[i] * ys[i] - 2 * g * xs[i] - 2 * f * ys[i] + c) / sigma2;
264  }
265  return chiSquared;
266  }
267 
268  //TODO: merge this one and the pT5 function later into a single function
269  template <typename TAcc>
270  ALPAKA_FN_ACC ALPAKA_FN_INLINE float computePT3RPhiChiSquared(TAcc const& acc,
272  uint16_t* lowerModuleIndices,
273  float g,
274  float f,
275  float radius,
276  float* xs,
277  float* ys) {
278  float delta1[3]{}, delta2[3]{}, slopes[3]{};
279  bool isFlat[3]{};
280  float chiSquared = 0;
281  float inv1 = kWidthPS / kWidth2S;
282  float inv2 = kPixelPSZpitch / kWidth2S;
283  for (size_t i = 0; i < 3; i++) {
284  ModuleType moduleType = modules.moduleType()[lowerModuleIndices[i]];
285  short moduleSubdet = modules.subdets()[lowerModuleIndices[i]];
286  short moduleSide = modules.sides()[lowerModuleIndices[i]];
287  float drdz = modules.drdzs()[lowerModuleIndices[i]];
288  slopes[i] = modules.dxdys()[lowerModuleIndices[i]];
289  //category 1 - barrel PS flat
290  if (moduleSubdet == Barrel and moduleType == PS and moduleSide == Center) {
291  delta1[i] = inv1;
292  delta2[i] = inv1;
293  slopes[i] = -999;
294  isFlat[i] = true;
295  }
296  //category 2 - barrel 2S
297  else if (moduleSubdet == Barrel and moduleType == TwoS) {
298  delta1[i] = 1;
299  delta2[i] = 1;
300  slopes[i] = -999;
301  isFlat[i] = true;
302  }
303  //category 3 - barrel PS tilted
304  else if (moduleSubdet == Barrel and moduleType == PS and moduleSide != Center) {
305  delta1[i] = inv1;
306  isFlat[i] = false;
307  delta2[i] = (inv2 * drdz / alpaka::math::sqrt(acc, 1 + drdz * drdz));
308  }
309  //category 4 - endcap PS
310  else if (moduleSubdet == Endcap and moduleType == PS) {
311  delta1[i] = inv1;
312  isFlat[i] = false;
313 
314  /*
315  despite the type of the module layer of the lower module index, all anchor
316  hits are on the pixel side and all non-anchor hits are on the strip side!
317  */
318  delta2[i] = inv2;
319  }
320  //category 5 - endcap 2S
321  else if (moduleSubdet == Endcap and moduleType == TwoS) {
322  delta1[i] = 1;
323  delta2[i] = 500 * inv1;
324  isFlat[i] = false;
325  }
326 #ifdef WARNINGS
327  else {
328  printf("ERROR!!!!! I SHOULDN'T BE HERE!!!! subdet = %d, type = %d, side = %d\n",
329  moduleSubdet,
330  moduleType,
331  moduleSide);
332  }
333 #endif
334  }
335  chiSquared = computeChiSquaredpT3(acc, 3, xs, ys, delta1, delta2, slopes, isFlat, g, f, radius);
336 
337  return chiSquared;
338  }
339 
340  ALPAKA_FN_ACC ALPAKA_FN_INLINE float computePT3RPhiChiSquaredInwards(
341  float g, float f, float r, float* xPix, float* yPix) {
342  float residual = (xPix[0] - g) * (xPix[0] - g) + (yPix[0] - f) * (yPix[0] - f) - r * r;
343  float chiSquared = residual * residual;
344  residual = (xPix[1] - g) * (xPix[1] - g) + (yPix[1] - f) * (yPix[1] - f) - r * r;
345  chiSquared += residual * residual;
346 
347  chiSquared *= 0.5f;
348  return chiSquared;
349  }
350 
351  //90pc threshold
352  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passPT3RPhiChiSquaredCuts(ModulesConst modules,
353  uint16_t lowerModuleIndex1,
354  uint16_t lowerModuleIndex2,
355  uint16_t lowerModuleIndex3,
356  float chiSquared) {
357  const int layer1 =
358  modules.layers()[lowerModuleIndex1] + 6 * (modules.subdets()[lowerModuleIndex1] == Endcap) +
359  5 * (modules.subdets()[lowerModuleIndex1] == Endcap and modules.moduleType()[lowerModuleIndex1] == TwoS);
360  const int layer2 =
361  modules.layers()[lowerModuleIndex2] + 6 * (modules.subdets()[lowerModuleIndex2] == Endcap) +
362  5 * (modules.subdets()[lowerModuleIndex2] == Endcap and modules.moduleType()[lowerModuleIndex2] == TwoS);
363  const int layer3 =
364  modules.layers()[lowerModuleIndex3] + 6 * (modules.subdets()[lowerModuleIndex3] == Endcap) +
365  5 * (modules.subdets()[lowerModuleIndex3] == Endcap and modules.moduleType()[lowerModuleIndex3] == TwoS);
366 
367  if (layer1 == 8 and layer2 == 9 and layer3 == 10) {
368  return chiSquared < 7.003f;
369  } else if (layer1 == 8 and layer2 == 9 and layer3 == 15) {
370  return chiSquared < 0.5f;
371  } else if (layer1 == 7 and layer2 == 8 and layer3 == 9) {
372  return chiSquared < 8.046f;
373  } else if (layer1 == 7 and layer2 == 8 and layer3 == 14) {
374  return chiSquared < 0.575f;
375  } else if (layer1 == 1 and layer2 == 2 and layer3 == 7) {
376  return chiSquared < 5.304f;
377  } else if (layer1 == 1 and layer2 == 2 and layer3 == 3) {
378  return chiSquared < 10.6211f;
379  } else if (layer1 == 1 and layer2 == 7 and layer3 == 8) {
380  return chiSquared < 4.617f;
381  } else if (layer1 == 2 and layer2 == 7 and layer3 == 8) {
382  return chiSquared < 8.046f;
383  } else if (layer1 == 2 and layer2 == 7 and layer3 == 13) {
384  return chiSquared < 0.435f;
385  } else if (layer1 == 2 and layer2 == 3 and layer3 == 7) {
386  return chiSquared < 9.244f;
387  } else if (layer1 == 2 and layer2 == 3 and layer3 == 12) {
388  return chiSquared < 0.287f;
389  } else if (layer1 == 2 and layer2 == 3 and layer3 == 4) {
390  return chiSquared < 18.509f;
391  }
392 
393  return true;
394  }
395 
396  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passPT3RPhiChiSquaredInwardsCuts(ModulesConst modules,
397  uint16_t lowerModuleIndex1,
398  uint16_t lowerModuleIndex2,
399  uint16_t lowerModuleIndex3,
400  float chiSquared) {
401  const int layer1 =
402  modules.layers()[lowerModuleIndex1] + 6 * (modules.subdets()[lowerModuleIndex1] == Endcap) +
403  5 * (modules.subdets()[lowerModuleIndex1] == Endcap and modules.moduleType()[lowerModuleIndex1] == TwoS);
404  const int layer2 =
405  modules.layers()[lowerModuleIndex2] + 6 * (modules.subdets()[lowerModuleIndex2] == Endcap) +
406  5 * (modules.subdets()[lowerModuleIndex2] == Endcap and modules.moduleType()[lowerModuleIndex2] == TwoS);
407  const int layer3 =
408  modules.layers()[lowerModuleIndex3] + 6 * (modules.subdets()[lowerModuleIndex3] == Endcap) +
409  5 * (modules.subdets()[lowerModuleIndex3] == Endcap and modules.moduleType()[lowerModuleIndex3] == TwoS);
410 
411  if (layer1 == 7 and layer2 == 8 and layer3 == 9) // endcap layer 1,2,3, ps
412  {
413  return chiSquared < 22016.8055f;
414  } else if (layer1 == 7 and layer2 == 8 and layer3 == 14) // endcap layer 1,2,3 layer3->2s
415  {
416  return chiSquared < 935179.56807f;
417  } else if (layer1 == 8 and layer2 == 9 and layer3 == 10) // endcap layer 2,3,4
418  {
419  return chiSquared < 29064.12959f;
420  } else if (layer1 == 8 and layer2 == 9 and layer3 == 15) // endcap layer 2,3,4, layer3->2s
421  {
422  return chiSquared < 935179.5681f;
423  } else if (layer1 == 1 and layer2 == 2 and layer3 == 3) // barrel 1,2,3
424  {
425  return chiSquared < 1370.0113195101474f;
426  } else if (layer1 == 1 and layer2 == 2 and layer3 == 7) // barrel 1,2 endcap 1
427  {
428  return chiSquared < 5492.110048314815f;
429  } else if (layer1 == 2 and layer2 == 3 and layer3 == 4) // barrel 2,3,4
430  {
431  return chiSquared < 4160.410806470067f;
432  } else if (layer1 == 1 and layer2 == 7 and layer3 == 8) // barrel 1, endcap 1,2
433  {
434  return chiSquared < 29064.129591225726f;
435  } else if (layer1 == 2 and layer2 == 3 and layer3 == 7) // barrel 2,3 endcap 1
436  {
437  return chiSquared < 12634.215376250893f;
438  } else if (layer1 == 2 and layer2 == 3 and layer3 == 12) // barrel 2,3, endcap 1->2s
439  {
440  return chiSquared < 353821.69361145404f;
441  } else if (layer1 == 2 and layer2 == 7 and layer3 == 8) // barrel2, endcap 1,2
442  {
443  return chiSquared < 33393.26076341235f;
444  } else if (layer1 == 2 and layer2 == 7 and layer3 == 13) //barrel 2, endcap 1, endcap2->2s
445  {
446  return chiSquared < 935179.5680742573f;
447  }
448 
449  return true;
450  }
451 
452  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool checkIntervalOverlappT3(float firstMin,
453  float firstMax,
454  float secondMin,
455  float secondMax) {
456  return ((firstMin <= secondMin) && (secondMin < firstMax)) || ((secondMin < firstMin) && (firstMin < secondMax));
457  }
458 
459  /*bounds for high Pt taken from : http://uaf-10.t2.ucsd.edu/~bsathian/SDL/T5_efficiency/efficiencies/new_efficiencies/efficiencies_20210513_T5_recovering_high_Pt_efficiencies/highE_radius_matching/highE_bounds.txt */
460  template <typename TAcc>
461  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passRadiusCriterionBBB(TAcc const& acc,
462  float pixelRadius,
463  float pixelRadiusError,
464  float tripletRadius) {
465  float tripletInvRadiusErrorBound = 0.15624f;
466  float pixelInvRadiusErrorBound = 0.17235f;
467 
468  if (pixelRadius > 2.0f * kR1GeVf) {
469  pixelInvRadiusErrorBound = 0.6375f;
470  tripletInvRadiusErrorBound = 0.6588f;
471  }
472 
473  float tripletRadiusInvMax = (1 + tripletInvRadiusErrorBound) / tripletRadius;
474  float tripletRadiusInvMin = alpaka::math::max(acc, (1 - tripletInvRadiusErrorBound) / tripletRadius, 0.0f);
475 
476  float pixelRadiusInvMax =
477  alpaka::math::max(acc, (1 + pixelInvRadiusErrorBound) / pixelRadius, 1.f / (pixelRadius - pixelRadiusError));
478  float pixelRadiusInvMin =
479  alpaka::math::min(acc, (1 - pixelInvRadiusErrorBound) / pixelRadius, 1.f / (pixelRadius + pixelRadiusError));
480 
481  return checkIntervalOverlappT3(tripletRadiusInvMin, tripletRadiusInvMax, pixelRadiusInvMin, pixelRadiusInvMax);
482  }
483 
484  template <typename TAcc>
485  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passRadiusCriterionBBE(TAcc const& acc,
486  float pixelRadius,
487  float pixelRadiusError,
488  float tripletRadius) {
489  float tripletInvRadiusErrorBound = 0.45972f;
490  float pixelInvRadiusErrorBound = 0.19644f;
491 
492  if (pixelRadius > 2.0f * kR1GeVf) {
493  pixelInvRadiusErrorBound = 0.6805f;
494  tripletInvRadiusErrorBound = 0.8557f;
495  }
496 
497  float tripletRadiusInvMax = (1 + tripletInvRadiusErrorBound) / tripletRadius;
498  float tripletRadiusInvMin = alpaka::math::max(acc, (1 - tripletInvRadiusErrorBound) / tripletRadius, 0.0f);
499 
500  float pixelRadiusInvMax =
501  alpaka::math::max(acc, (1 + pixelInvRadiusErrorBound) / pixelRadius, 1.f / (pixelRadius - pixelRadiusError));
502  float pixelRadiusInvMin =
503  alpaka::math::min(acc, (1 - pixelInvRadiusErrorBound) / pixelRadius, 1.f / (pixelRadius + pixelRadiusError));
504 
505  return checkIntervalOverlappT3(tripletRadiusInvMin, tripletRadiusInvMax, pixelRadiusInvMin, pixelRadiusInvMax);
506  }
507 
508  template <typename TAcc>
509  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passRadiusCriterionBEE(TAcc const& acc,
510  float pixelRadius,
511  float pixelRadiusError,
512  float tripletRadius) {
513  float tripletInvRadiusErrorBound = 1.59294f;
514  float pixelInvRadiusErrorBound = 0.255181f;
515 
516  if (pixelRadius > 2.0f * kR1GeVf) //as good as not having selections
517  {
518  pixelInvRadiusErrorBound = 2.2091f;
519  tripletInvRadiusErrorBound = 2.3548f;
520  }
521 
522  float tripletRadiusInvMax = (1 + tripletInvRadiusErrorBound) / tripletRadius;
523  float tripletRadiusInvMin = alpaka::math::max(acc, (1 - tripletInvRadiusErrorBound) / tripletRadius, 0.0f);
524 
525  float pixelRadiusInvMax =
526  alpaka::math::max(acc, (1 + pixelInvRadiusErrorBound) / pixelRadius, 1.f / (pixelRadius - pixelRadiusError));
527  float pixelRadiusInvMin =
528  alpaka::math::min(acc, (1 - pixelInvRadiusErrorBound) / pixelRadius, 1.f / (pixelRadius + pixelRadiusError));
529  pixelRadiusInvMin = alpaka::math::max(acc, pixelRadiusInvMin, 0.0f);
530 
531  return checkIntervalOverlappT3(tripletRadiusInvMin, tripletRadiusInvMax, pixelRadiusInvMin, pixelRadiusInvMax);
532  }
533 
534  template <typename TAcc>
535  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passRadiusCriterionEEE(TAcc const& acc,
536  float pixelRadius,
537  float pixelRadiusError,
538  float tripletRadius) {
539  float tripletInvRadiusErrorBound = 1.7006f;
540  float pixelInvRadiusErrorBound = 0.26367f;
541 
542  if (pixelRadius > 2.0f * kR1GeVf) //as good as not having selections
543  {
544  pixelInvRadiusErrorBound = 2.286f;
545  tripletInvRadiusErrorBound = 2.436f;
546  }
547 
548  float tripletRadiusInvMax = (1 + tripletInvRadiusErrorBound) / tripletRadius;
549  float tripletRadiusInvMin = alpaka::math::max(acc, (1 - tripletInvRadiusErrorBound) / tripletRadius, 0.0f);
550 
551  float pixelRadiusInvMax =
552  alpaka::math::max(acc, (1 + pixelInvRadiusErrorBound) / pixelRadius, 1.f / (pixelRadius - pixelRadiusError));
553  float pixelRadiusInvMin =
554  alpaka::math::min(acc, (1 - pixelInvRadiusErrorBound) / pixelRadius, 1.f / (pixelRadius + pixelRadiusError));
555  pixelRadiusInvMin = alpaka::math::max(acc, 0.0f, pixelRadiusInvMin);
556 
557  return checkIntervalOverlappT3(tripletRadiusInvMin, tripletRadiusInvMax, pixelRadiusInvMin, pixelRadiusInvMax);
558  }
559 
560  template <typename TAcc>
561  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passRadiusCriterion(TAcc const& acc,
563  float pixelRadius,
564  float pixelRadiusError,
565  float tripletRadius,
566  int16_t lowerModuleIndex,
567  uint16_t middleModuleIndex,
568  uint16_t upperModuleIndex) {
569  if (modules.subdets()[lowerModuleIndex] == Endcap) {
570  return passRadiusCriterionEEE(acc, pixelRadius, pixelRadiusError, tripletRadius);
571  } else if (modules.subdets()[middleModuleIndex] == Endcap) {
572  return passRadiusCriterionBEE(acc, pixelRadius, pixelRadiusError, tripletRadius);
573  } else if (modules.subdets()[upperModuleIndex] == Endcap) {
574  return passRadiusCriterionBBE(acc, pixelRadius, pixelRadiusError, tripletRadius);
575  } else {
576  return passRadiusCriterionBBB(acc, pixelRadius, pixelRadiusError, tripletRadius);
577  }
578  }
579 
580  template <typename TAcc>
581  ALPAKA_FN_ACC ALPAKA_FN_INLINE float computePT3RZChiSquared(TAcc const& acc,
583  const uint16_t* lowerModuleIndices,
584  const float* rtPix,
585  const float* xPix,
586  const float* yPix,
587  const float* zPix,
588  const float* rts,
589  const float* xs,
590  const float* ys,
591  const float* zs,
592  float pixelSegmentPt,
593  float pixelSegmentPx,
594  float pixelSegmentPy,
595  float pixelSegmentPz,
596  int pixelSegmentCharge) {
597  float residual = 0;
598  float error2 = 0;
599  float RMSE = 0;
600 
601  float Px = pixelSegmentPx, Py = pixelSegmentPy, Pz = pixelSegmentPz;
602  int charge = pixelSegmentCharge;
603  float x1 = xPix[1] / 100;
604  float y1 = yPix[1] / 100;
605  float z1 = zPix[1] / 100;
606  float r1 = rtPix[1] / 100;
607 
608  float a = -2.f * k2Rinv1GeVf * 100 * charge; // multiply by 100 to make the correct length units
609 
610  for (size_t i = 0; i < Params_T3::kLayers; i++) {
611  float zsi = zs[i] / 100;
612  float rtsi = rts[i] / 100;
613  uint16_t lowerModuleIndex = lowerModuleIndices[i];
614  const int moduleType = modules.moduleType()[lowerModuleIndex];
615  const int moduleSide = modules.sides()[lowerModuleIndex];
616  const int moduleSubdet = modules.subdets()[lowerModuleIndex];
617 
618  // calculation is detailed documented here https://indico.cern.ch/event/1185895/contributions/4982756/attachments/2526561/4345805/helix%20pT3%20summarize.pdf
619  float diffr, diffz;
620  float p = alpaka::math::sqrt(acc, Px * Px + Py * Py + Pz * Pz);
621 
622  float rou = a / p;
623  if (moduleSubdet == Endcap) {
624  float s = (zsi - z1) * p / Pz;
625  float x = x1 + Px / a * alpaka::math::sin(acc, rou * s) - Py / a * (1 - alpaka::math::cos(acc, rou * s));
626  float y = y1 + Py / a * alpaka::math::sin(acc, rou * s) + Px / a * (1 - alpaka::math::cos(acc, rou * s));
627  diffr = alpaka::math::abs(acc, rtsi - alpaka::math::sqrt(acc, x * x + y * y)) * 100;
628  }
629 
630  if (moduleSubdet == Barrel) {
631  float paraA = r1 * r1 + 2 * (Px * Px + Py * Py) / (a * a) + 2 * (y1 * Px - x1 * Py) / a - rtsi * rtsi;
632  float paraB = 2 * (x1 * Px + y1 * Py) / a;
633  float paraC = 2 * (y1 * Px - x1 * Py) / a + 2 * (Px * Px + Py * Py) / (a * a);
634  float A = paraB * paraB + paraC * paraC;
635  float B = 2 * paraA * paraB;
636  float C = paraA * paraA - paraC * paraC;
637  float sol1 = (-B + alpaka::math::sqrt(acc, B * B - 4 * A * C)) / (2 * A);
638  float sol2 = (-B - alpaka::math::sqrt(acc, B * B - 4 * A * C)) / (2 * A);
639  float solz1 = alpaka::math::asin(acc, sol1) / rou * Pz / p + z1;
640  float solz2 = alpaka::math::asin(acc, sol2) / rou * Pz / p + z1;
641  float diffz1 = alpaka::math::abs(acc, solz1 - zsi) * 100;
642  float diffz2 = alpaka::math::abs(acc, solz2 - zsi) * 100;
643  diffz = alpaka::math::min(acc, diffz1, diffz2);
644  }
645 
646  residual = moduleSubdet == Barrel ? diffz : diffr;
647 
648  //PS Modules
649  if (moduleType == 0) {
650  error2 = kPixelPSZpitch * kPixelPSZpitch;
651  } else //2S modules
652  {
653  error2 = kStrip2SZpitch * kStrip2SZpitch;
654  }
655 
656  //special dispensation to tilted PS modules!
657  if (moduleType == 0 and moduleSubdet == Barrel and moduleSide != Center) {
658  float drdz = modules.drdzs()[lowerModuleIndex];
659  error2 /= (1 + drdz * drdz);
660  }
661  RMSE += (residual * residual) / error2;
662  }
663 
664  RMSE = alpaka::math::sqrt(acc, 0.2f * RMSE); // Divided by the degree of freedom 5.
665 
666  return RMSE;
667  }
668 
669  template <typename TAcc>
670  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool runPixelTripletDefaultAlgo(TAcc const& acc,
673  MiniDoubletsConst mds,
674  SegmentsConst segments,
675  SegmentsPixelConst segmentsPixel,
676  TripletsConst triplets,
677  unsigned int pixelSegmentIndex,
678  unsigned int tripletIndex,
679  float& pixelRadius,
680  float& tripletRadius,
681  float& centerX,
682  float& centerY,
683  float& rzChiSquared,
684  float& rPhiChiSquared,
685  float& rPhiChiSquaredInwards,
686  bool runChiSquaredCuts = true) {
687  //run pT4 compatibility between the pixel segment and inner segment, and between the pixel and outer segment of the triplet
688  uint16_t pixelModuleIndex = segments.innerLowerModuleIndices()[pixelSegmentIndex];
689 
690  uint16_t lowerModuleIndex = triplets.lowerModuleIndices()[tripletIndex][0];
691  uint16_t middleModuleIndex = triplets.lowerModuleIndices()[tripletIndex][1];
692  uint16_t upperModuleIndex = triplets.lowerModuleIndices()[tripletIndex][2];
693 
694  {
695  // pixel segment vs inner segment of the triplet
696  if (not runPixelTrackletDefaultAlgopT3(acc,
697  modules,
698  ranges,
699  mds,
700  segments,
701  segmentsPixel,
702  pixelModuleIndex,
703  lowerModuleIndex,
704  middleModuleIndex,
705  pixelSegmentIndex,
706  triplets.segmentIndices()[tripletIndex][0]))
707  return false;
708 
709  //pixel segment vs outer segment of triplet
710  if (not runPixelTrackletDefaultAlgopT3(acc,
711  modules,
712  ranges,
713  mds,
714  segments,
715  segmentsPixel,
716  pixelModuleIndex,
717  middleModuleIndex,
718  upperModuleIndex,
719  pixelSegmentIndex,
720  triplets.segmentIndices()[tripletIndex][1]))
721  return false;
722  }
723 
724  //pt matching between the pixel ptin and the triplet circle pt
725  unsigned int pixelSegmentArrayIndex = pixelSegmentIndex - ranges.segmentModuleIndices()[pixelModuleIndex];
726  float pixelSegmentPt = segmentsPixel.ptIn()[pixelSegmentArrayIndex];
727  float pixelSegmentPtError = segmentsPixel.ptErr()[pixelSegmentArrayIndex];
728  float pixelSegmentPx = segmentsPixel.px()[pixelSegmentArrayIndex];
729  float pixelSegmentPy = segmentsPixel.py()[pixelSegmentArrayIndex];
730  float pixelSegmentPz = segmentsPixel.pz()[pixelSegmentArrayIndex];
731  int pixelSegmentCharge = segmentsPixel.charge()[pixelSegmentArrayIndex];
732 
733  float pixelG = segmentsPixel.circleCenterX()[pixelSegmentArrayIndex];
734  float pixelF = segmentsPixel.circleCenterY()[pixelSegmentArrayIndex];
735  float pixelRadiusPCA = segmentsPixel.circleRadius()[pixelSegmentArrayIndex];
736 
737  unsigned int pixelInnerMDIndex = segments.mdIndices()[pixelSegmentIndex][0];
738  unsigned int pixelOuterMDIndex = segments.mdIndices()[pixelSegmentIndex][1];
739 
740  pixelRadius = pixelSegmentPt * kR1GeVf;
741  float pixelRadiusError = pixelSegmentPtError * kR1GeVf;
742  unsigned int tripletInnerSegmentIndex = triplets.segmentIndices()[tripletIndex][0];
743  unsigned int tripletOuterSegmentIndex = triplets.segmentIndices()[tripletIndex][1];
744 
745  unsigned int firstMDIndex = segments.mdIndices()[tripletInnerSegmentIndex][0];
746  unsigned int secondMDIndex = segments.mdIndices()[tripletInnerSegmentIndex][1];
747  unsigned int thirdMDIndex = segments.mdIndices()[tripletOuterSegmentIndex][1];
748 
749  float xs[Params_T3::kLayers] = {
750  mds.anchorX()[firstMDIndex], mds.anchorX()[secondMDIndex], mds.anchorX()[thirdMDIndex]};
751  float ys[Params_T3::kLayers] = {
752  mds.anchorY()[firstMDIndex], mds.anchorY()[secondMDIndex], mds.anchorY()[thirdMDIndex]};
753 
754  float g, f;
755  tripletRadius = triplets.radius()[tripletIndex];
756  g = triplets.centerX()[tripletIndex];
757  f = triplets.centerY()[tripletIndex];
758 
759  if (not passRadiusCriterion(acc,
760  modules,
761  pixelRadius,
762  pixelRadiusError,
763  tripletRadius,
764  lowerModuleIndex,
765  middleModuleIndex,
766  upperModuleIndex))
767  return false;
768 
769  uint16_t lowerModuleIndices[Params_T3::kLayers] = {lowerModuleIndex, middleModuleIndex, upperModuleIndex};
770 
771  if (runChiSquaredCuts and pixelSegmentPt < 5.0f) {
772  float rts[Params_T3::kLayers] = {
773  mds.anchorRt()[firstMDIndex], mds.anchorRt()[secondMDIndex], mds.anchorRt()[thirdMDIndex]};
774  float zs[Params_T3::kLayers] = {
775  mds.anchorZ()[firstMDIndex], mds.anchorZ()[secondMDIndex], mds.anchorZ()[thirdMDIndex]};
776  float rtPix[Params_pLS::kLayers] = {mds.anchorRt()[pixelInnerMDIndex], mds.anchorRt()[pixelOuterMDIndex]};
777  float xPix[Params_pLS::kLayers] = {mds.anchorX()[pixelInnerMDIndex], mds.anchorX()[pixelOuterMDIndex]};
778  float yPix[Params_pLS::kLayers] = {mds.anchorY()[pixelInnerMDIndex], mds.anchorY()[pixelOuterMDIndex]};
779  float zPix[Params_pLS::kLayers] = {mds.anchorZ()[pixelInnerMDIndex], mds.anchorZ()[pixelOuterMDIndex]};
780 
781  rzChiSquared = computePT3RZChiSquared(acc,
782  modules,
783  lowerModuleIndices,
784  rtPix,
785  xPix,
786  yPix,
787  zPix,
788  rts,
789  xs,
790  ys,
791  zs,
792  pixelSegmentPt,
793  pixelSegmentPx,
794  pixelSegmentPy,
795  pixelSegmentPz,
796  pixelSegmentCharge);
797  if (not passPT3RZChiSquaredCuts(modules, lowerModuleIndex, middleModuleIndex, upperModuleIndex, rzChiSquared))
798  return false;
799  } else {
800  rzChiSquared = -1;
801  }
802 
803  rPhiChiSquared = computePT3RPhiChiSquared(acc, modules, lowerModuleIndices, pixelG, pixelF, pixelRadiusPCA, xs, ys);
804 
805  if (runChiSquaredCuts and pixelSegmentPt < 5.0f) {
806  if (not passPT3RPhiChiSquaredCuts(modules, lowerModuleIndex, middleModuleIndex, upperModuleIndex, rPhiChiSquared))
807  return false;
808  }
809 
810  float xPix[Params_pLS::kLayers] = {mds.anchorX()[pixelInnerMDIndex], mds.anchorX()[pixelOuterMDIndex]};
811  float yPix[Params_pLS::kLayers] = {mds.anchorY()[pixelInnerMDIndex], mds.anchorY()[pixelOuterMDIndex]};
812  rPhiChiSquaredInwards = computePT3RPhiChiSquaredInwards(g, f, tripletRadius, xPix, yPix);
813 
814  if (runChiSquaredCuts and pixelSegmentPt < 5.0f) {
816  modules, lowerModuleIndex, middleModuleIndex, upperModuleIndex, rPhiChiSquaredInwards))
817  return false;
818  }
819  centerX = 0;
820  centerY = 0;
821  return true;
822  }
823 
825  template <typename TAcc>
826  ALPAKA_FN_ACC void operator()(TAcc const& acc,
828  ModulesPixelConst modulesPixel,
830  MiniDoubletsConst mds,
831  SegmentsConst segments,
832  SegmentsPixelConst segmentsPixel,
833  Triplets triplets,
834  TripletsOccupancyConst tripletsOccupancy,
835  PixelTriplets pixelTriplets,
836  unsigned int* connectedPixelSize,
837  unsigned int* connectedPixelIndex,
838  unsigned int nPixelSegments) const {
839  auto const globalBlockIdx = alpaka::getIdx<alpaka::Grid, alpaka::Blocks>(acc);
840  auto const globalThreadIdx = alpaka::getIdx<alpaka::Grid, alpaka::Threads>(acc);
841  auto const gridBlockExtent = alpaka::getWorkDiv<alpaka::Grid, alpaka::Blocks>(acc);
842  auto const gridThreadExtent = alpaka::getWorkDiv<alpaka::Grid, alpaka::Threads>(acc);
843 
844  for (unsigned int i_pLS = globalThreadIdx[1]; i_pLS < nPixelSegments; i_pLS += gridThreadExtent[1]) {
845  auto iLSModule_max = connectedPixelIndex[i_pLS] + connectedPixelSize[i_pLS];
846 
847  for (unsigned int iLSModule = connectedPixelIndex[i_pLS] + globalBlockIdx[0]; iLSModule < iLSModule_max;
848  iLSModule += gridBlockExtent[0]) {
849  uint16_t tripletLowerModuleIndex =
850  modulesPixel.connectedPixels()
851  [iLSModule]; //connected pixels will have the appropriate lower module index by default!
852 #ifdef WARNINGS
853  if (tripletLowerModuleIndex >= modules.nLowerModules()) {
854  printf("tripletLowerModuleIndex %d >= modules.nLowerModules %d \n",
855  tripletLowerModuleIndex,
856  modules.nLowerModules());
857  continue; //sanity check
858  }
859 #endif
860  //Removes 2S-2S :FIXME: filter these out in the pixel map
861  if (modules.moduleType()[tripletLowerModuleIndex] == TwoS)
862  continue;
863 
864  uint16_t pixelModuleIndex = modules.nLowerModules();
865  unsigned int nOuterTriplets = tripletsOccupancy.nTriplets()[tripletLowerModuleIndex];
866  if (nOuterTriplets == 0)
867  continue;
868 
869  unsigned int pixelSegmentIndex = ranges.segmentModuleIndices()[pixelModuleIndex] + i_pLS;
870 
871  if (segmentsPixel.isDup()[i_pLS])
872  continue;
873  if (segmentsPixel.partOfPT5()[i_pLS])
874  continue; //don't make pT3s for those pixels that are part of pT5
875 
876  short layer2_adjustment;
877  if (modules.layers()[tripletLowerModuleIndex] == 1) {
878  layer2_adjustment = 1;
879  } //get upper segment to be in second layer
880  else if (modules.layers()[tripletLowerModuleIndex] == 2) {
881  layer2_adjustment = 0;
882  } // get lower segment to be in second layer
883  else {
884  continue;
885  }
886 
887  //fetch the triplet
888  for (unsigned int outerTripletArrayIndex = globalThreadIdx[2]; outerTripletArrayIndex < nOuterTriplets;
889  outerTripletArrayIndex += gridThreadExtent[2]) {
890  unsigned int outerTripletIndex =
891  ranges.tripletModuleIndices()[tripletLowerModuleIndex] + outerTripletArrayIndex;
892  if (modules.moduleType()[triplets.lowerModuleIndices()[outerTripletIndex][1]] == TwoS)
893  continue; //REMOVES PS-2S
894 
895  if (triplets.partOfPT5()[outerTripletIndex])
896  continue; //don't create pT3s for T3s accounted in pT5s
897 
898  float pixelRadius, tripletRadius, rPhiChiSquared, rzChiSquared, rPhiChiSquaredInwards, centerX, centerY;
900  modules,
901  ranges,
902  mds,
903  segments,
904  segmentsPixel,
905  triplets,
906  pixelSegmentIndex,
907  outerTripletIndex,
908  pixelRadius,
909  tripletRadius,
910  centerX,
911  centerY,
912  rzChiSquared,
913  rPhiChiSquared,
914  rPhiChiSquaredInwards);
915 
916  if (success) {
917  float phi =
918  mds.anchorPhi()[segments
919  .mdIndices()[triplets.segmentIndices()[outerTripletIndex][0]][layer2_adjustment]];
920  float eta =
921  mds.anchorEta()[segments
922  .mdIndices()[triplets.segmentIndices()[outerTripletIndex][0]][layer2_adjustment]];
923  float eta_pix = segmentsPixel.eta()[i_pLS];
924  float phi_pix = segmentsPixel.phi()[i_pLS];
925  float pt = segmentsPixel.ptIn()[i_pLS];
926  float score = rPhiChiSquared + rPhiChiSquaredInwards;
927  unsigned int totOccupancyPixelTriplets =
928  alpaka::atomicAdd(acc, &pixelTriplets.totOccupancyPixelTriplets(), 1u, alpaka::hierarchy::Threads{});
929  if (totOccupancyPixelTriplets >= n_max_pixel_triplets) {
930 #ifdef WARNINGS
931  printf("Pixel Triplet excess alert!\n");
932 #endif
933  } else {
934  unsigned int pixelTripletIndex =
935  alpaka::atomicAdd(acc, &pixelTriplets.nPixelTriplets(), 1u, alpaka::hierarchy::Threads{});
937  segments,
938  triplets,
939  pixelTriplets,
940  pixelSegmentIndex,
941  outerTripletIndex,
942  pixelRadius,
943  tripletRadius,
944  centerX,
945  centerY,
946  rPhiChiSquared,
947  rPhiChiSquaredInwards,
948  rzChiSquared,
949  pixelTripletIndex,
950  pt,
951  eta,
952  phi,
953  eta_pix,
954  phi_pix,
955  score);
956  triplets.partOfPT3()[outerTripletIndex] = true;
957  }
958  }
959  } // for outerTripletArrayIndex
960  } // for iLSModule < iLSModule_max
961  } // for i_pLS
962  }
963  };
964 
965  template <typename TAcc>
966  ALPAKA_FN_ACC ALPAKA_FN_INLINE void runDeltaBetaIterationspT3(TAcc const& acc,
967  float& betaIn,
968  float& betaOut,
969  float betaAv,
970  float& pt_beta,
971  float sdIn_dr,
972  float sdOut_dr,
973  float dr,
974  float lIn) {
975  if (lIn == 0) {
976  betaOut += alpaka::math::copysign(
977  acc,
978  alpaka::math::asin(
979  acc, alpaka::math::min(acc, sdOut_dr * k2Rinv1GeVf / alpaka::math::abs(acc, pt_beta), kSinAlphaMax)),
980  betaOut);
981  return;
982  }
983 
984  if (betaIn * betaOut > 0.f and
985  (alpaka::math::abs(acc, pt_beta) < 4.f * kPt_betaMax or
986  (lIn >= 11 and alpaka::math::abs(acc, pt_beta) <
987  8.f * kPt_betaMax))) //and the pt_beta is well-defined; less strict for endcap-endcap
988  {
989  const float betaInUpd =
990  betaIn +
991  alpaka::math::copysign(
992  acc,
993  alpaka::math::asin(
994  acc, alpaka::math::min(acc, sdIn_dr * k2Rinv1GeVf / alpaka::math::abs(acc, pt_beta), kSinAlphaMax)),
995  betaIn); //FIXME: need a faster version
996  const float betaOutUpd =
997  betaOut +
998  alpaka::math::copysign(
999  acc,
1000  alpaka::math::asin(
1001  acc, alpaka::math::min(acc, sdOut_dr * k2Rinv1GeVf / alpaka::math::abs(acc, pt_beta), kSinAlphaMax)),
1002  betaOut); //FIXME: need a faster version
1003  betaAv = 0.5f * (betaInUpd + betaOutUpd);
1004 
1005  //1st update
1006  const float pt_beta_inv =
1007  1.f / alpaka::math::abs(acc, dr * k2Rinv1GeVf / alpaka::math::sin(acc, betaAv)); //get a better pt estimate
1008 
1009  betaIn += alpaka::math::copysign(
1010  acc,
1011  alpaka::math::asin(acc, alpaka::math::min(acc, sdIn_dr * k2Rinv1GeVf * pt_beta_inv, kSinAlphaMax)),
1012  betaIn); //FIXME: need a faster version
1013  betaOut += alpaka::math::copysign(
1014  acc,
1015  alpaka::math::asin(acc, alpaka::math::min(acc, sdOut_dr * k2Rinv1GeVf * pt_beta_inv, kSinAlphaMax)),
1016  betaOut); //FIXME: need a faster version
1017  //update the av and pt
1018  betaAv = 0.5f * (betaIn + betaOut);
1019  //2nd update
1020  pt_beta = dr * k2Rinv1GeVf / alpaka::math::sin(acc, betaAv); //get a better pt estimate
1021  } else if (lIn < 11 && alpaka::math::abs(acc, betaOut) < 0.2f * alpaka::math::abs(acc, betaIn) &&
1022  alpaka::math::abs(acc, pt_beta) < 12.f * kPt_betaMax) //use betaIn sign as ref
1023  {
1024  const float pt_betaIn = dr * k2Rinv1GeVf / alpaka::math::sin(acc, betaIn);
1025 
1026  const float betaInUpd =
1027  betaIn +
1028  alpaka::math::copysign(
1029  acc,
1030  alpaka::math::asin(
1031  acc, alpaka::math::min(acc, sdIn_dr * k2Rinv1GeVf / alpaka::math::abs(acc, pt_betaIn), kSinAlphaMax)),
1032  betaIn); //FIXME: need a faster version
1033  const float betaOutUpd =
1034  betaOut +
1035  alpaka::math::copysign(
1036  acc,
1037  alpaka::math::asin(
1038  acc,
1039  alpaka::math::min(acc, sdOut_dr * k2Rinv1GeVf / alpaka::math::abs(acc, pt_betaIn), kSinAlphaMax)),
1040  betaIn); //FIXME: need a faster version
1041  betaAv = (alpaka::math::abs(acc, betaOut) > 0.2f * alpaka::math::abs(acc, betaIn))
1042  ? (0.5f * (betaInUpd + betaOutUpd))
1043  : betaInUpd;
1044 
1045  //1st update
1046  pt_beta = dr * k2Rinv1GeVf / alpaka::math::sin(acc, betaAv); //get a better pt estimate
1047  betaIn += alpaka::math::copysign(
1048  acc,
1049  alpaka::math::asin(
1050  acc, alpaka::math::min(acc, sdIn_dr * k2Rinv1GeVf / alpaka::math::abs(acc, pt_beta), kSinAlphaMax)),
1051  betaIn); //FIXME: need a faster version
1052  betaOut += alpaka::math::copysign(
1053  acc,
1054  alpaka::math::asin(
1055  acc, alpaka::math::min(acc, sdOut_dr * k2Rinv1GeVf / alpaka::math::abs(acc, pt_beta), kSinAlphaMax)),
1056  betaIn); //FIXME: need a faster version
1057  //update the av and pt
1058  betaAv = 0.5f * (betaIn + betaOut);
1059  //2nd update
1060  pt_beta = dr * k2Rinv1GeVf / alpaka::math::sin(acc, betaAv); //get a better pt estimate
1061  }
1062  }
1063 
1064  template <typename TAcc>
1065  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool runTripletDefaultAlgoPPBB(TAcc const& acc,
1068  MiniDoubletsConst mds,
1069  SegmentsConst segments,
1070  SegmentsPixelConst segmentsPixel,
1071  uint16_t pixelModuleIndex,
1072  uint16_t outerInnerLowerModuleIndex,
1073  uint16_t outerOuterLowerModuleIndex,
1074  unsigned int innerSegmentIndex,
1075  unsigned int outerSegmentIndex,
1076  unsigned int firstMDIndex,
1077  unsigned int secondMDIndex,
1078  unsigned int thirdMDIndex,
1079  unsigned int fourthMDIndex) {
1080  float dPhi, betaIn, betaOut, pt_beta, zLo, zHi, zLoPointed, zHiPointed, dPhiCut, betaOutCut;
1081 
1082  bool isPS_OutLo = (modules.moduleType()[outerInnerLowerModuleIndex] == PS);
1083 
1084  float rt_InLo = mds.anchorRt()[firstMDIndex];
1085  float rt_InUp = mds.anchorRt()[secondMDIndex];
1086  float rt_OutLo = mds.anchorRt()[thirdMDIndex];
1087  float rt_OutUp = mds.anchorRt()[fourthMDIndex];
1088 
1089  float z_InUp = mds.anchorZ()[secondMDIndex];
1090  float z_OutLo = mds.anchorZ()[thirdMDIndex];
1091 
1092  float x_InLo = mds.anchorX()[firstMDIndex];
1093  float x_InUp = mds.anchorX()[secondMDIndex];
1094  float x_OutLo = mds.anchorX()[thirdMDIndex];
1095  float x_OutUp = mds.anchorX()[fourthMDIndex];
1096 
1097  float y_InLo = mds.anchorY()[firstMDIndex];
1098  float y_InUp = mds.anchorY()[secondMDIndex];
1099  float y_OutLo = mds.anchorY()[thirdMDIndex];
1100  float y_OutUp = mds.anchorY()[fourthMDIndex];
1101 
1102  float rt_InOut = rt_InUp;
1103 
1104  if (alpaka::math::abs(acc, deltaPhi(acc, x_InUp, y_InUp, x_OutLo, y_OutLo)) > kPi / 2.f)
1105  return false;
1106 
1107  unsigned int pixelSegmentArrayIndex = innerSegmentIndex - ranges.segmentModuleIndices()[pixelModuleIndex];
1108  float ptIn = segmentsPixel.ptIn()[pixelSegmentArrayIndex];
1109  float ptSLo = ptIn;
1110  float px = segmentsPixel.px()[pixelSegmentArrayIndex];
1111  float py = segmentsPixel.py()[pixelSegmentArrayIndex];
1112  float pz = segmentsPixel.pz()[pixelSegmentArrayIndex];
1113  float ptErr = segmentsPixel.ptErr()[pixelSegmentArrayIndex];
1114  float etaErr = segmentsPixel.etaErr()[pixelSegmentArrayIndex];
1115  ptSLo = alpaka::math::max(acc, ptCut, ptSLo - 10.0f * alpaka::math::max(acc, ptErr, 0.005f * ptSLo));
1116  ptSLo = alpaka::math::min(acc, 10.0f, ptSLo);
1117 
1118  float alpha1GeV_OutLo =
1119  alpaka::math::asin(acc, alpaka::math::min(acc, rt_OutLo * k2Rinv1GeVf / ptCut, kSinAlphaMax));
1120  const float rtRatio_OutLoInOut =
1121  rt_OutLo / rt_InOut; // Outer segment beginning rt divided by inner segment beginning rt;
1122 
1123  float dzDrtScale =
1124  alpaka::math::tan(acc, alpha1GeV_OutLo) / alpha1GeV_OutLo; // The track can bend in r-z plane slightly
1125  const float zpitch_InLo = 0.05f;
1126  const float zpitch_InOut = 0.05f;
1127  float zpitch_OutLo = (isPS_OutLo ? kPixelPSZpitch : kStrip2SZpitch);
1128  float zGeom = zpitch_InLo + zpitch_OutLo;
1129  zHi = z_InUp + (z_InUp + kDeltaZLum) * (rtRatio_OutLoInOut - 1.f) * (z_InUp < 0.f ? 1.f : dzDrtScale) +
1130  (zpitch_InOut + zpitch_OutLo);
1131  zLo = z_InUp + (z_InUp - kDeltaZLum) * (rtRatio_OutLoInOut - 1.f) * (z_InUp > 0.f ? 1.f : dzDrtScale) -
1132  (zpitch_InOut + zpitch_OutLo); //slope-correction only on outer end
1133 
1134  if ((z_OutLo < zLo) || (z_OutLo > zHi))
1135  return false;
1136 
1137  const float cosh2Eta = 1.f + (pz * pz) / (ptIn * ptIn);
1138 
1139  const float drt_OutLo_InUp = (rt_OutLo - rt_InUp);
1140 
1141  const float r3_InUp = alpaka::math::sqrt(acc, z_InUp * z_InUp + rt_InUp * rt_InUp);
1142 
1143  float drt_InSeg = rt_InOut - rt_InLo;
1144 
1145  const float thetaMuls2 =
1146  (kMulsInGeV * kMulsInGeV) * (0.1f + 0.2f * (rt_OutLo - rt_InUp) / 50.f) * (r3_InUp / rt_InUp);
1147  const float muls2 = thetaMuls2 * 9.f / (ptCut * ptCut) * 16.f;
1148 
1149  float dzErr = (drt_OutLo_InUp * drt_OutLo_InUp) * (etaErr * etaErr) * cosh2Eta;
1150  dzErr += 0.03f * 0.03f; // Approximately account for IT module size
1151  dzErr *= 9.f; // 3 sigma
1152  dzErr += muls2 * (drt_OutLo_InUp * drt_OutLo_InUp) / 3.f * cosh2Eta;
1153  dzErr += zGeom * zGeom;
1154  dzErr = alpaka::math::sqrt(acc, dzErr);
1155 
1156  const float dzDrIn = pz / ptIn;
1157  const float zWindow = dzErr / drt_InSeg * drt_OutLo_InUp + zGeom;
1158  const float dzMean = dzDrIn * drt_OutLo_InUp *
1159  (1.f + drt_OutLo_InUp * drt_OutLo_InUp * 4 * k2Rinv1GeVf * k2Rinv1GeVf / ptIn / ptIn /
1160  24.f); // with curved path correction
1161  // Constructing upper and lower bound
1162  zLoPointed = z_InUp + dzMean - zWindow;
1163  zHiPointed = z_InUp + dzMean + zWindow;
1164 
1165  if ((z_OutLo < zLoPointed) || (z_OutLo > zHiPointed))
1166  return false;
1167 
1168  const float pvOffset = 0.1f / rt_OutLo;
1169  dPhiCut = alpha1GeV_OutLo + alpaka::math::sqrt(acc, muls2 + pvOffset * pvOffset);
1170 
1171  //no dphipos cut
1172  float midPointX = 0.5f * (x_InLo + x_OutLo);
1173  float midPointY = 0.5f * (y_InLo + y_OutLo);
1174 
1175  float diffX = x_OutLo - x_InLo;
1176  float diffY = y_OutLo - y_InLo;
1177 
1178  dPhi = deltaPhi(acc, midPointX, midPointY, diffX, diffY);
1179 
1180  if (alpaka::math::abs(acc, dPhi) > dPhiCut)
1181  return false;
1182 
1183  //lots of array accesses below this...
1184 
1185  float alpha_InLo = __H2F(segments.dPhiChanges()[innerSegmentIndex]);
1186  float alpha_OutLo = __H2F(segments.dPhiChanges()[outerSegmentIndex]);
1187 
1188  bool isEC_lastLayer = modules.subdets()[outerOuterLowerModuleIndex] == Endcap and
1189  modules.moduleType()[outerOuterLowerModuleIndex] == TwoS;
1190 
1191  float alpha_OutUp, alpha_OutUp_highEdge, alpha_OutUp_lowEdge;
1192  alpha_OutUp = deltaPhi(acc, x_OutUp, y_OutUp, x_OutUp - x_OutLo, y_OutUp - y_OutLo);
1193 
1194  alpha_OutUp_highEdge = alpha_OutUp;
1195  alpha_OutUp_lowEdge = alpha_OutUp;
1196 
1197  float tl_axis_x = x_OutUp - x_InUp;
1198  float tl_axis_y = y_OutUp - y_InUp;
1199 
1200  float tl_axis_highEdge_x = tl_axis_x;
1201  float tl_axis_highEdge_y = tl_axis_y;
1202 
1203  float tl_axis_lowEdge_x = tl_axis_x;
1204  float tl_axis_lowEdge_y = tl_axis_y;
1205 
1206  betaIn = -deltaPhi(acc, px, py, tl_axis_x, tl_axis_y);
1207  float betaInRHmin = betaIn;
1208  float betaInRHmax = betaIn;
1209 
1210  betaOut = -alpha_OutUp + deltaPhi(acc, x_OutUp, y_OutUp, tl_axis_x, tl_axis_y);
1211 
1212  float betaOutRHmin = betaOut;
1213  float betaOutRHmax = betaOut;
1214 
1215  if (isEC_lastLayer) {
1216  alpha_OutUp_highEdge = deltaPhi(acc,
1217  mds.anchorHighEdgeX()[fourthMDIndex],
1218  mds.anchorHighEdgeY()[fourthMDIndex],
1219  mds.anchorHighEdgeX()[fourthMDIndex] - x_OutLo,
1220  mds.anchorHighEdgeY()[fourthMDIndex] - y_OutLo);
1221  alpha_OutUp_lowEdge = deltaPhi(acc,
1222  mds.anchorLowEdgeX()[fourthMDIndex],
1223  mds.anchorLowEdgeY()[fourthMDIndex],
1224  mds.anchorLowEdgeX()[fourthMDIndex] - x_OutLo,
1225  mds.anchorLowEdgeY()[fourthMDIndex] - y_OutLo);
1226 
1227  tl_axis_highEdge_x = mds.anchorHighEdgeX()[fourthMDIndex] - x_InUp;
1228  tl_axis_highEdge_y = mds.anchorHighEdgeY()[fourthMDIndex] - y_InUp;
1229  tl_axis_lowEdge_x = mds.anchorLowEdgeX()[fourthMDIndex] - x_InUp;
1230  tl_axis_lowEdge_y = mds.anchorLowEdgeY()[fourthMDIndex] - y_InUp;
1231 
1232  betaOutRHmin = -alpha_OutUp_highEdge + deltaPhi(acc,
1233  mds.anchorHighEdgeX()[fourthMDIndex],
1234  mds.anchorHighEdgeY()[fourthMDIndex],
1235  tl_axis_highEdge_x,
1236  tl_axis_highEdge_y);
1237  betaOutRHmax = -alpha_OutUp_lowEdge + deltaPhi(acc,
1238  mds.anchorLowEdgeX()[fourthMDIndex],
1239  mds.anchorLowEdgeY()[fourthMDIndex],
1240  tl_axis_lowEdge_x,
1241  tl_axis_lowEdge_y);
1242  }
1243 
1244  //beta computation
1245  float drt_tl_axis = alpaka::math::sqrt(acc, tl_axis_x * tl_axis_x + tl_axis_y * tl_axis_y);
1246 
1247  //innerOuterAnchor - innerInnerAnchor
1248  const float rt_InSeg =
1249  alpaka::math::sqrt(acc, (x_InUp - x_InLo) * (x_InUp - x_InLo) + (y_InUp - y_InLo) * (y_InUp - y_InLo));
1250 
1251  //no betaIn cut for the pixels
1252  float betaAv = 0.5f * (betaIn + betaOut);
1253  pt_beta = ptIn;
1254 
1255  int lIn = 0;
1256  int lOut = isEC_lastLayer ? 11 : 5;
1257  float sdOut_dr =
1258  alpaka::math::sqrt(acc, (x_OutUp - x_OutLo) * (x_OutUp - x_OutLo) + (y_OutUp - y_OutLo) * (y_OutUp - y_OutLo));
1259  float sdOut_d = rt_OutUp - rt_OutLo;
1260 
1261  runDeltaBetaIterationspT3(acc, betaIn, betaOut, betaAv, pt_beta, rt_InSeg, sdOut_dr, drt_tl_axis, lIn);
1262 
1263  const float betaInMMSF = (alpaka::math::abs(acc, betaInRHmin + betaInRHmax) > 0)
1264  ? (2.f * betaIn / alpaka::math::abs(acc, betaInRHmin + betaInRHmax))
1265  : 0.; //mean value of min,max is the old betaIn
1266  const float betaOutMMSF = (alpaka::math::abs(acc, betaOutRHmin + betaOutRHmax) > 0)
1267  ? (2.f * betaOut / alpaka::math::abs(acc, betaOutRHmin + betaOutRHmax))
1268  : 0.;
1269  betaInRHmin *= betaInMMSF;
1270  betaInRHmax *= betaInMMSF;
1271  betaOutRHmin *= betaOutMMSF;
1272  betaOutRHmax *= betaOutMMSF;
1273 
1274  float min_ptBeta_ptBetaMax = alpaka::math::min(
1275  acc, alpaka::math::abs(acc, pt_beta), kPt_betaMax); //need to confirm the range-out value of 7 GeV
1276  const float dBetaMuls2 = thetaMuls2 * 16.f / (min_ptBeta_ptBetaMax * min_ptBeta_ptBetaMax);
1277  const float alphaInAbsReg =
1278  alpaka::math::max(acc,
1279  alpaka::math::abs(acc, alpha_InLo),
1280  alpaka::math::asin(acc, alpaka::math::min(acc, rt_InUp * k2Rinv1GeVf / 3.0f, kSinAlphaMax)));
1281  const float alphaOutAbsReg =
1282  alpaka::math::max(acc,
1283  alpaka::math::abs(acc, alpha_OutLo),
1284  alpaka::math::asin(acc, alpaka::math::min(acc, rt_OutLo * k2Rinv1GeVf / 3.0f, kSinAlphaMax)));
1285  const float dBetaInLum = lIn < 11 ? 0.0f : alpaka::math::abs(acc, alphaInAbsReg * kDeltaZLum / z_InUp);
1286  const float dBetaOutLum = lOut < 11 ? 0.0f : alpaka::math::abs(acc, alphaOutAbsReg * kDeltaZLum / z_OutLo);
1287  const float dBetaLum2 = (dBetaInLum + dBetaOutLum) * (dBetaInLum + dBetaOutLum);
1288 
1289  const float sinDPhi = alpaka::math::sin(acc, dPhi);
1290  const float dBetaRIn2 = 0; // TODO-RH
1291 
1292  float dBetaROut = 0;
1293  if (isEC_lastLayer) {
1294  dBetaROut = (alpaka::math::sqrt(acc,
1295  mds.anchorHighEdgeX()[fourthMDIndex] * mds.anchorHighEdgeX()[fourthMDIndex] +
1296  mds.anchorHighEdgeY()[fourthMDIndex] * mds.anchorHighEdgeY()[fourthMDIndex]) -
1297  alpaka::math::sqrt(acc,
1298  mds.anchorLowEdgeX()[fourthMDIndex] * mds.anchorLowEdgeX()[fourthMDIndex] +
1299  mds.anchorLowEdgeY()[fourthMDIndex] * mds.anchorLowEdgeY()[fourthMDIndex])) *
1300  sinDPhi / drt_tl_axis;
1301  }
1302 
1303  const float dBetaROut2 = dBetaROut * dBetaROut;
1304 
1305  //FIXME: need faster version
1306  betaOutCut = alpaka::math::asin(acc, alpaka::math::min(acc, drt_tl_axis * k2Rinv1GeVf / ptCut, kSinAlphaMax)) +
1307  (0.02f / sdOut_d) + alpaka::math::sqrt(acc, dBetaLum2 + dBetaMuls2);
1308 
1309  //Cut #6: The real beta cut
1310  if (alpaka::math::abs(acc, betaOut) >= betaOutCut)
1311  return false;
1312  const float dBetaRes = 0.02f / alpaka::math::min(acc, sdOut_d, drt_InSeg);
1313  const float dBetaCut2 =
1314  (dBetaRes * dBetaRes * 2.0f + dBetaMuls2 + dBetaLum2 + dBetaRIn2 + dBetaROut2 +
1315  0.25f *
1316  (alpaka::math::abs(acc, betaInRHmin - betaInRHmax) + alpaka::math::abs(acc, betaOutRHmin - betaOutRHmax)) *
1317  (alpaka::math::abs(acc, betaInRHmin - betaInRHmax) + alpaka::math::abs(acc, betaOutRHmin - betaOutRHmax)));
1318  float dBeta = betaIn - betaOut;
1319  return dBeta * dBeta <= dBetaCut2;
1320  }
1321 
1322  template <typename TAcc>
1323  ALPAKA_FN_ACC ALPAKA_FN_INLINE bool runTripletDefaultAlgoPPEE(TAcc const& acc,
1326  MiniDoubletsConst mds,
1327  SegmentsConst segments,
1328  SegmentsPixelConst segmentsPixel,
1329  uint16_t pixelModuleIndex,
1330  uint16_t outerInnerLowerModuleIndex,
1331  uint16_t outerOuterLowerModuleIndex,
1332  unsigned int innerSegmentIndex,
1333  unsigned int outerSegmentIndex,
1334  unsigned int firstMDIndex,
1335  unsigned int secondMDIndex,
1336  unsigned int thirdMDIndex,
1337  unsigned int fourthMDIndex) {
1338  float dPhi, betaIn, betaOut, pt_beta, rtLo, rtHi, dPhiCut, betaOutCut;
1339 
1340  bool isPS_OutLo = (modules.moduleType()[outerInnerLowerModuleIndex] == PS);
1341 
1342  float z_InUp = mds.anchorZ()[secondMDIndex];
1343  float z_OutLo = mds.anchorZ()[thirdMDIndex];
1344 
1345  if (z_InUp * z_OutLo <= 0)
1346  return false;
1347 
1348  float rt_InLo = mds.anchorRt()[firstMDIndex];
1349  float rt_InUp = mds.anchorRt()[secondMDIndex];
1350  float rt_OutLo = mds.anchorRt()[thirdMDIndex];
1351  float rt_OutUp = mds.anchorRt()[fourthMDIndex];
1352 
1353  float x_InLo = mds.anchorX()[firstMDIndex];
1354  float x_InUp = mds.anchorX()[secondMDIndex];
1355  float x_OutLo = mds.anchorX()[thirdMDIndex];
1356  float x_OutUp = mds.anchorX()[fourthMDIndex];
1357 
1358  float y_InLo = mds.anchorY()[firstMDIndex];
1359  float y_InUp = mds.anchorY()[secondMDIndex];
1360  float y_OutLo = mds.anchorY()[thirdMDIndex];
1361  float y_OutUp = mds.anchorY()[fourthMDIndex];
1362 
1363  unsigned int pixelSegmentArrayIndex = innerSegmentIndex - ranges.segmentModuleIndices()[pixelModuleIndex];
1364 
1365  float ptIn = segmentsPixel.ptIn()[pixelSegmentArrayIndex];
1366  float ptSLo = ptIn;
1367  float px = segmentsPixel.px()[pixelSegmentArrayIndex];
1368  float py = segmentsPixel.py()[pixelSegmentArrayIndex];
1369  float pz = segmentsPixel.pz()[pixelSegmentArrayIndex];
1370  float ptErr = segmentsPixel.ptErr()[pixelSegmentArrayIndex];
1371  float etaErr = segmentsPixel.etaErr()[pixelSegmentArrayIndex];
1372 
1373  ptSLo = alpaka::math::max(acc, ptCut, ptSLo - 10.0f * alpaka::math::max(acc, ptErr, 0.005f * ptSLo));
1374  ptSLo = alpaka::math::min(acc, 10.0f, ptSLo);
1375 
1376  const float zpitch_InLo = 0.05f;
1377  float zpitch_OutLo = (isPS_OutLo ? kPixelPSZpitch : kStrip2SZpitch);
1378  float zGeom = zpitch_InLo + zpitch_OutLo;
1379 
1380  const float slope = alpaka::math::asin(acc, alpaka::math::min(acc, rt_OutLo * k2Rinv1GeVf / ptCut, kSinAlphaMax));
1381  const float dzDrtScale = alpaka::math::tan(acc, slope) / slope; //FIXME: need approximate value
1382 
1383  const float dLum = alpaka::math::copysign(acc, kDeltaZLum, z_InUp);
1384  bool isOutSgInnerMDPS = modules.moduleType()[outerInnerLowerModuleIndex] == PS;
1385 
1386  const float rtGeom1 = isOutSgInnerMDPS
1387  ? kPixelPSZpitch
1388  : kStrip2SZpitch; //FIXME: make this chosen by configuration for lay11,12 full PS
1389  const float zGeom1 = alpaka::math::copysign(acc, zGeom, z_InUp); //used in B-E region
1390  rtLo = rt_InUp * (1.f + (z_OutLo - z_InUp - zGeom1) / (z_InUp + zGeom1 + dLum) / dzDrtScale) -
1391  rtGeom1; //slope correction only on the lower end
1392 
1393  float zInForHi = z_InUp - zGeom1 - dLum;
1394  if (zInForHi * z_InUp < 0)
1395  zInForHi = alpaka::math::copysign(acc, 0.1f, z_InUp);
1396  rtHi = rt_InUp * (1.f + (z_OutLo - z_InUp + zGeom1) / zInForHi) + rtGeom1;
1397 
1398  // Cut #2: rt condition
1399  if ((rt_OutLo < rtLo) || (rt_OutLo > rtHi))
1400  return false;
1401 
1402  const float dzOutInAbs = alpaka::math::abs(acc, z_OutLo - z_InUp);
1403  const float cosh2Eta = 1.f + (pz * pz) / (ptIn * ptIn);
1404  const float multDzDr2 = (dzOutInAbs * dzOutInAbs) * cosh2Eta / ((cosh2Eta - 1.f) * (cosh2Eta - 1.f));
1405  const float r3_InUp = alpaka::math::sqrt(acc, z_InUp * z_InUp + rt_InUp * rt_InUp);
1406  const float thetaMuls2 =
1407  (kMulsInGeV * kMulsInGeV) * (0.1f + 0.2f * (rt_OutLo - rt_InUp) / 50.f) * (r3_InUp / rt_InUp);
1408  const float muls2 = thetaMuls2 * 9.f / (ptCut * ptCut) * 16.f;
1409 
1410  float drtErr = (etaErr * etaErr) * multDzDr2;
1411  drtErr += 0.03f * 0.03f; // Approximately account for IT module size
1412  drtErr *= 9.f; // 3 sigma
1413  drtErr += muls2 * multDzDr2 / 3.f * cosh2Eta;
1414  drtErr = alpaka::math::sqrt(acc, drtErr);
1415  const float drtDzIn = alpaka::math::abs(acc, ptIn / pz);
1416 
1417  const float drt_OutLo_InUp = (rt_OutLo - rt_InUp); // drOutIn
1418 
1419  const float rtWindow = drtErr + rtGeom1;
1420  const float drtMean = drtDzIn * dzOutInAbs *
1421  (1.f - drt_OutLo_InUp * drt_OutLo_InUp * 4 * k2Rinv1GeVf * k2Rinv1GeVf / ptIn / ptIn /
1422  24.f); // with curved path correction
1423  const float rtLo_point = rt_InUp + drtMean - rtWindow;
1424  const float rtHi_point = rt_InUp + drtMean + rtWindow;
1425 
1426  // Cut #3: rt-z pointed
1427  if ((rt_OutLo < rtLo_point) || (rt_OutLo > rtHi_point))
1428  return false;
1429 
1430  const float alpha1GeV_OutLo =
1431  alpaka::math::asin(acc, alpaka::math::min(acc, rt_OutLo * k2Rinv1GeVf / ptCut, kSinAlphaMax));
1432  const float pvOffset = 0.1f / rt_OutLo;
1433  dPhiCut = alpha1GeV_OutLo + alpaka::math::sqrt(acc, muls2 + pvOffset * pvOffset);
1434 
1435  float midPointX = 0.5f * (x_InLo + x_OutLo);
1436  float midPointY = 0.5f * (y_InLo + y_OutLo);
1437 
1438  float diffX = x_OutLo - x_InLo;
1439  float diffY = y_OutLo - y_InLo;
1440 
1441  dPhi = deltaPhi(acc, midPointX, midPointY, diffX, diffY);
1442 
1443  // Cut #5: deltaPhiChange
1444  if (alpaka::math::abs(acc, dPhi) > dPhiCut)
1445  return false;
1446 
1447  float alpha_InLo = __H2F(segments.dPhiChanges()[innerSegmentIndex]);
1448  float alpha_OutLo = __H2F(segments.dPhiChanges()[outerSegmentIndex]);
1449 
1450  bool isEC_lastLayer = modules.subdets()[outerOuterLowerModuleIndex] == Endcap and
1451  modules.moduleType()[outerOuterLowerModuleIndex] == TwoS;
1452 
1453  float alpha_OutUp, alpha_OutUp_highEdge, alpha_OutUp_lowEdge;
1454 
1455  alpha_OutUp = deltaPhi(acc, x_OutUp, y_OutUp, x_OutUp - x_OutLo, y_OutUp - y_OutLo);
1456  alpha_OutUp_highEdge = alpha_OutUp;
1457  alpha_OutUp_lowEdge = alpha_OutUp;
1458 
1459  float tl_axis_x = x_OutUp - x_InUp;
1460  float tl_axis_y = y_OutUp - y_InUp;
1461 
1462  float tl_axis_highEdge_x = tl_axis_x;
1463  float tl_axis_highEdge_y = tl_axis_y;
1464 
1465  float tl_axis_lowEdge_x = tl_axis_x;
1466  float tl_axis_lowEdge_y = tl_axis_y;
1467 
1468  betaIn = -deltaPhi(acc, px, py, tl_axis_x, tl_axis_y);
1469  float betaInRHmin = betaIn;
1470  float betaInRHmax = betaIn;
1471 
1472  betaOut = -alpha_OutUp + deltaPhi(acc, x_OutUp, y_OutUp, tl_axis_x, tl_axis_y);
1473  float betaOutRHmin = betaOut;
1474  float betaOutRHmax = betaOut;
1475 
1476  if (isEC_lastLayer) {
1477  alpha_OutUp_highEdge = deltaPhi(acc,
1478  mds.anchorHighEdgeX()[fourthMDIndex],
1479  mds.anchorHighEdgeY()[fourthMDIndex],
1480  mds.anchorHighEdgeX()[fourthMDIndex] - x_OutLo,
1481  mds.anchorHighEdgeY()[fourthMDIndex] - y_OutLo);
1482  alpha_OutUp_lowEdge = deltaPhi(acc,
1483  mds.anchorLowEdgeX()[fourthMDIndex],
1484  mds.anchorLowEdgeY()[fourthMDIndex],
1485  mds.anchorLowEdgeX()[fourthMDIndex] - x_OutLo,
1486  mds.anchorLowEdgeY()[fourthMDIndex] - y_OutLo);
1487 
1488  tl_axis_highEdge_x = mds.anchorHighEdgeX()[fourthMDIndex] - x_InUp;
1489  tl_axis_highEdge_y = mds.anchorHighEdgeY()[fourthMDIndex] - y_InUp;
1490  tl_axis_lowEdge_x = mds.anchorLowEdgeX()[fourthMDIndex] - x_InUp;
1491  tl_axis_lowEdge_y = mds.anchorLowEdgeY()[fourthMDIndex] - y_InUp;
1492 
1493  betaOutRHmin = -alpha_OutUp_highEdge + deltaPhi(acc,
1494  mds.anchorHighEdgeX()[fourthMDIndex],
1495  mds.anchorHighEdgeY()[fourthMDIndex],
1496  tl_axis_highEdge_x,
1497  tl_axis_highEdge_y);
1498  betaOutRHmax = -alpha_OutUp_lowEdge + deltaPhi(acc,
1499  mds.anchorLowEdgeX()[fourthMDIndex],
1500  mds.anchorLowEdgeY()[fourthMDIndex],
1501  tl_axis_lowEdge_x,
1502  tl_axis_lowEdge_y);
1503  }
1504 
1505  //beta computation
1506  float drt_tl_axis = alpaka::math::sqrt(acc, tl_axis_x * tl_axis_x + tl_axis_y * tl_axis_y);
1507  //no betaIn cut for the pixels
1508  const float rt_InSeg =
1509  alpaka::math::sqrt(acc, (x_InUp - x_InLo) * (x_InUp - x_InLo) + (y_InUp - y_InLo) * (y_InUp - y_InLo));
1510 
1511  float betaAv = 0.5f * (betaIn + betaOut);
1512  pt_beta = ptIn;
1513 
1514  int lIn = 0;
1515  int lOut = isEC_lastLayer ? 11 : 5;
1516  float sdOut_dr =
1517  alpaka::math::sqrt(acc, (x_OutUp - x_OutLo) * (x_OutUp - x_OutLo) + (y_OutUp - y_OutLo) * (y_OutUp - y_OutLo));
1518  float sdOut_d = rt_OutUp - rt_OutLo;
1519 
1520  runDeltaBetaIterationspT3(acc, betaIn, betaOut, betaAv, pt_beta, rt_InSeg, sdOut_dr, drt_tl_axis, lIn);
1521 
1522  const float betaInMMSF = (alpaka::math::abs(acc, betaInRHmin + betaInRHmax) > 0)
1523  ? (2.f * betaIn / alpaka::math::abs(acc, betaInRHmin + betaInRHmax))
1524  : 0.; //mean value of min,max is the old betaIn
1525  const float betaOutMMSF = (alpaka::math::abs(acc, betaOutRHmin + betaOutRHmax) > 0)
1526  ? (2.f * betaOut / alpaka::math::abs(acc, betaOutRHmin + betaOutRHmax))
1527  : 0.;
1528  betaInRHmin *= betaInMMSF;
1529  betaInRHmax *= betaInMMSF;
1530  betaOutRHmin *= betaOutMMSF;
1531  betaOutRHmax *= betaOutMMSF;
1532 
1533  float min_ptBeta_ptBetaMax = alpaka::math::min(
1534  acc, alpaka::math::abs(acc, pt_beta), kPt_betaMax); //need to confirm the range-out value of 7 GeV
1535  const float dBetaMuls2 = thetaMuls2 * 16.f / (min_ptBeta_ptBetaMax * min_ptBeta_ptBetaMax);
1536 
1537  const float alphaInAbsReg =
1538  alpaka::math::max(acc,
1539  alpaka::math::abs(acc, alpha_InLo),
1540  alpaka::math::asin(acc, alpaka::math::min(acc, rt_InUp * k2Rinv1GeVf / 3.0f, kSinAlphaMax)));
1541  const float alphaOutAbsReg =
1542  alpaka::math::max(acc,
1543  alpaka::math::abs(acc, alpha_OutLo),
1544  alpaka::math::asin(acc, alpaka::math::min(acc, rt_OutLo * k2Rinv1GeVf / 3.0f, kSinAlphaMax)));
1545  const float dBetaInLum = lIn < 11 ? 0.0f : alpaka::math::abs(acc, alphaInAbsReg * kDeltaZLum / z_InUp);
1546  const float dBetaOutLum = lOut < 11 ? 0.0f : alpaka::math::abs(acc, alphaOutAbsReg * kDeltaZLum / z_OutLo);
1547  const float dBetaLum2 = (dBetaInLum + dBetaOutLum) * (dBetaInLum + dBetaOutLum);
1548 
1549  const float sinDPhi = alpaka::math::sin(acc, dPhi);
1550  const float dBetaRIn2 = 0; // TODO-RH
1551 
1552  float dBetaROut = 0;
1553  if (isEC_lastLayer) {
1554  dBetaROut = (alpaka::math::sqrt(acc,
1555  mds.anchorHighEdgeX()[fourthMDIndex] * mds.anchorHighEdgeX()[fourthMDIndex] +
1556  mds.anchorHighEdgeY()[fourthMDIndex] * mds.anchorHighEdgeY()[fourthMDIndex]) -
1557  alpaka::math::sqrt(acc,
1558  mds.anchorLowEdgeX()[fourthMDIndex] * mds.anchorLowEdgeX()[fourthMDIndex] +
1559  mds.anchorLowEdgeY()[fourthMDIndex] * mds.anchorLowEdgeY()[fourthMDIndex])) *
1560  sinDPhi / drt_tl_axis;
1561  }
1562 
1563  const float dBetaROut2 = dBetaROut * dBetaROut;
1564 
1565  betaOutCut =
1566  alpaka::math::asin(
1567  acc, alpaka::math::min(acc, drt_tl_axis * k2Rinv1GeVf / ptCut, kSinAlphaMax)) //FIXME: need faster version
1568  + (0.02f / sdOut_d) + alpaka::math::sqrt(acc, dBetaLum2 + dBetaMuls2);
1569 
1570  //Cut #6: The real beta cut
1571  if (alpaka::math::abs(acc, betaOut) >= betaOutCut)
1572  return false;
1573 
1574  float drt_InSeg = rt_InUp - rt_InLo;
1575 
1576  const float dBetaRes = 0.02f / alpaka::math::min(acc, sdOut_d, drt_InSeg);
1577  const float dBetaCut2 =
1578  (dBetaRes * dBetaRes * 2.0f + dBetaMuls2 + dBetaLum2 + dBetaRIn2 + dBetaROut2 +
1579  0.25f *
1580  (alpaka::math::abs(acc, betaInRHmin - betaInRHmax) + alpaka::math::abs(acc, betaOutRHmin - betaOutRHmax)) *
1581  (alpaka::math::abs(acc, betaInRHmin - betaInRHmax) + alpaka::math::abs(acc, betaOutRHmin - betaOutRHmax)));
1582  float dBeta = betaIn - betaOut;
1583  return dBeta * dBeta <= dBetaCut2;
1584  }
1585 
1586 } // namespace ALPAKA_ACCELERATOR_NAMESPACE::lst
1587 #endif
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passRadiusCriterionBBE(TAcc const &acc, float pixelRadius, float pixelRadiusError, float tripletRadius)
Definition: PixelTriplet.h:485
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passRadiusCriterionEEE(TAcc const &acc, float pixelRadius, float pixelRadiusError, float tripletRadius)
Definition: PixelTriplet.h:535
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float kR1GeVf
Definition: Common.h:50
Definition: APVGainStruct.h:7
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passRadiusCriterionBEE(TAcc const &acc, float pixelRadius, float pixelRadiusError, float tripletRadius)
Definition: PixelTriplet.h:509
ALPAKA_FN_ACC ALPAKA_FN_INLINE float computeChiSquaredpT3(TAcc const &acc, unsigned int nPoints, float *xs, float *ys, float *delta1, float *delta2, float *slopes, bool *isFlat, float g, float f, float radius)
Definition: PixelTriplet.h:223
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float k2Rinv1GeVf
Definition: Common.h:49
ALPAKA_FN_ACC ALPAKA_FN_INLINE float computePT3RZChiSquared(TAcc const &acc, ModulesConst modules, const uint16_t *lowerModuleIndices, const float *rtPix, const float *xPix, const float *yPix, const float *zPix, const float *rts, const float *xs, const float *ys, const float *zs, float pixelSegmentPt, float pixelSegmentPx, float pixelSegmentPy, float pixelSegmentPz, int pixelSegmentCharge)
Definition: PixelTriplet.h:581
ALPAKA_FN_ACC ALPAKA_FN_INLINE float computePT3RPhiChiSquared(TAcc const &acc, ModulesConst modules, uint16_t *lowerModuleIndices, float g, float f, float radius, float *xs, float *ys)
Definition: PixelTriplet.h:270
static const double slope[3]
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float kWidthPS
Definition: Common.h:58
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool runTripletDefaultAlgoPPEE(TAcc const &acc, ModulesConst modules, ObjectRangesConst ranges, MiniDoubletsConst mds, SegmentsConst segments, SegmentsPixelConst segmentsPixel, uint16_t pixelModuleIndex, uint16_t outerInnerLowerModuleIndex, uint16_t outerOuterLowerModuleIndex, unsigned int innerSegmentIndex, unsigned int outerSegmentIndex, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex)
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float kWidth2S
Definition: Common.h:57
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float kVerticalModuleSlope
Definition: Common.h:61
constexpr unsigned int n_max_pixel_triplets
Definition: Common.h:35
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e g
Definition: Activities.doc:4
ALPAKA_FN_ACC ALPAKA_FN_INLINE float computePT3RPhiChiSquaredInwards(float g, float f, float r, float *xPix, float *yPix)
Definition: PixelTriplet.h:340
#define __H2F
Definition: Common.h:50
#define __F2H
Definition: Common.h:49
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float ptCut
Definition: Common.h:52
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool runTripletDefaultAlgoPPBB(TAcc const &acc, ModulesConst modules, ObjectRangesConst ranges, MiniDoubletsConst mds, SegmentsConst segments, SegmentsPixelConst segmentsPixel, uint16_t pixelModuleIndex, uint16_t outerInnerLowerModuleIndex, uint16_t outerOuterLowerModuleIndex, unsigned int innerSegmentIndex, unsigned int outerSegmentIndex, unsigned int firstMDIndex, unsigned int secondMDIndex, unsigned int thirdMDIndex, unsigned int fourthMDIndex)
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool runPixelTripletDefaultAlgo(TAcc const &acc, ModulesConst modules, ObjectRangesConst ranges, MiniDoubletsConst mds, SegmentsConst segments, SegmentsPixelConst segmentsPixel, TripletsConst triplets, unsigned int pixelSegmentIndex, unsigned int tripletIndex, float &pixelRadius, float &tripletRadius, float &centerX, float &centerY, float &rzChiSquared, float &rPhiChiSquared, float &rPhiChiSquaredInwards, bool runChiSquaredCuts=true)
Definition: PixelTriplet.h:670
TripletsSoA::ConstView TripletsConst
Definition: TripletsSoA.h:31
ALPAKA_FN_ACC void operator()(TAcc const &acc, ModulesConst modules, ModulesPixelConst modulesPixel, ObjectRangesConst ranges, MiniDoubletsConst mds, SegmentsConst segments, SegmentsPixelConst segmentsPixel, Triplets triplets, TripletsOccupancyConst tripletsOccupancy, PixelTriplets pixelTriplets, unsigned int *connectedPixelSize, unsigned int *connectedPixelIndex, unsigned int nPixelSegments) const
Definition: PixelTriplet.h:826
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passPT3RPhiChiSquaredCuts(ModulesConst modules, uint16_t lowerModuleIndex1, uint16_t lowerModuleIndex2, uint16_t lowerModuleIndex3, float chiSquared)
Definition: PixelTriplet.h:352
ALPAKA_FN_HOST_ACC ALPAKA_FN_INLINE float phi(TAcc const &acc, float x, float y)
Definition: Hit.h:29
T sqrt(T t)
Definition: SSEVec.h:23
ModulesPixelSoA::ConstView ModulesPixelConst
Definition: ModulesSoA.h:49
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventID const &, edm::Timestamp const & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
ModuleType
Definition: DQModule.h:24
double f[11][100]
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passPT3RZChiSquaredCuts(ModulesConst modules, uint16_t lowerModuleIndex1, uint16_t lowerModuleIndex2, uint16_t lowerModuleIndex3, float rzChiSquared)
Definition: PixelTriplet.h:171
MiniDoubletsSoA::ConstView MiniDoubletsConst
ModulesSoA::ConstView ModulesConst
Definition: ModulesSoA.h:47
PixelTripletsSoA::View PixelTriplets
string ranges
Definition: diffTwoXMLs.py:79
SegmentsPixelSoA::ConstView SegmentsPixelConst
Definition: SegmentsSoA.h:53
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float kMulsInGeV
Definition: Common.h:41
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float kPixelPSZpitch
Definition: Common.h:54
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool checkIntervalOverlappT3(float firstMin, float firstMax, float secondMin, float secondMax)
Definition: PixelTriplet.h:452
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float kPt_betaMax
Definition: Common.h:59
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passRadiusCriterion(TAcc const &acc, ModulesConst modules, float pixelRadius, float pixelRadiusError, float tripletRadius, int16_t lowerModuleIndex, uint16_t middleModuleIndex, uint16_t upperModuleIndex)
Definition: PixelTriplet.h:561
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float kPi
Definition: Common.h:39
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float kStrip2SZpitch
Definition: Common.h:56
TripletsSoA::View Triplets
Definition: TripletsSoA.h:30
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passPT3RPhiChiSquaredInwardsCuts(ModulesConst modules, uint16_t lowerModuleIndex1, uint16_t lowerModuleIndex2, uint16_t lowerModuleIndex3, float chiSquared)
Definition: PixelTriplet.h:396
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool passRadiusCriterionBBB(TAcc const &acc, float pixelRadius, float pixelRadiusError, float tripletRadius)
Definition: PixelTriplet.h:461
TripletsOccupancySoA::ConstView TripletsOccupancyConst
Definition: TripletsSoA.h:39
ObjectRangesSoA::ConstView ObjectRangesConst
double a
Definition: hdecay.h:121
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float kSinAlphaMax
Definition: Common.h:51
SegmentsSoA::ConstView SegmentsConst
Definition: SegmentsSoA.h:49
float x
Definition: APVGainStruct.h:7
ALPAKA_STATIC_ACC_MEM_GLOBAL constexpr float kDeltaZLum
Definition: Common.h:53
ALPAKA_FN_ACC ALPAKA_FN_INLINE bool runPixelTrackletDefaultAlgopT3(TAcc const &acc, ModulesConst modules, ObjectRangesConst ranges, MiniDoubletsConst mds, SegmentsConst segments, SegmentsPixelConst segmentsPixel, uint16_t pixelLowerModuleIndex, uint16_t outerInnerLowerModuleIndex, uint16_t outerOuterLowerModuleIndex, unsigned int innerSegmentIndex, unsigned int outerSegmentIndex)
Definition: PixelTriplet.h:114
ALPAKA_FN_HOST_ACC ALPAKA_FN_INLINE float deltaPhi(TAcc const &acc, float x1, float y1, float x2, float y2)
Definition: Hit.h:34
T1 atomicAdd(T1 *a, T2 b)
Definition: cudaCompat.h:61
ALPAKA_FN_ACC ALPAKA_FN_INLINE void addPixelTripletToMemory(MiniDoubletsConst mds, SegmentsConst segments, TripletsConst triplets, PixelTriplets pixelTriplets, unsigned int pixelSegmentIndex, unsigned int tripletIndex, float pixelRadius, float tripletRadius, float centerX, float centerY, float rPhiChiSquared, float rPhiChiSquaredInwards, float rzChiSquared, unsigned int pixelTripletIndex, float pt, float eta, float phi, float eta_pix, float phi_pix, float score)
Definition: PixelTriplet.h:48
ALPAKA_FN_ACC ALPAKA_FN_INLINE void runDeltaBetaIterationspT3(TAcc const &acc, float &betaIn, float &betaOut, float betaAv, float &pt_beta, float sdIn_dr, float sdOut_dr, float dr, float lIn)
Definition: PixelTriplet.h:966
ALPAKA_FN_HOST_ACC ALPAKA_FN_INLINE float eta(TAcc const &acc, float x, float y, float z)
Definition: Hit.h:11