ATLAS Offline Software
Loading...
Searching...
No Matches
SiTrajectoryElement_xk.icc
Go to the documentation of this file.
1/*
2 Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
3*/
4
5#include <type_traits>
6
7///////////////////////////////////////////////////////////////////
8// Search clusters compatible with track
9///////////////////////////////////////////////////////////////////
10
11template<typename T>
12int InDet::SiTrajectoryElement_xk::searchClustersSub
13(Trk::PatternTrackParameters& Tp, SiClusterLink_xk* L) {
14 /// Case a): No PRD association in use
15 if (not m_useassoTool) {
16 /// a1): No stereo
17 if (not m_stereo) {
18 if (m_ndf == 2) {
19 return searchClustersWithoutStereoPIX<T>(Tp, L);
20 } else {
21 return searchClustersWithoutStereoSCT<T>(Tp, L);
22 }
23 }
24 /// a2): With stereo
25 else {
26 return searchClustersWithStereo<T>(Tp, L);
27 }
28 }
29 /// Case b) Using PRD association tool
30 else {
31 /// b1): no stereo
32 if (not m_stereo) {
33 if (m_ndf == 2) {
34 return searchClustersWithoutStereoAssPIX<T>(Tp, L,*m_prdToTrackMap);
35 } else {
36 return searchClustersWithoutStereoAssSCT<T>(Tp, L,*m_prdToTrackMap);
37 }
38 }
39 /// b2): with stereo
40 else {
41 return searchClustersWithStereoAss<T>(Tp, L, *m_prdToTrackMap);
42 }
43 }
44}
45
46///////////////////////////////////////////////////////////////////
47// Search closest cluster with stereo angle
48///////////////////////////////////////////////////////////////////
49
50template <typename T>
51int InDet::SiTrajectoryElement_xk::searchClustersWithStereo
52(Trk::PatternTrackParameters& Tp,InDet::SiClusterLink_xk* L)
53{
54 if (m_detstatus<=0) return 0;
55
56 const AmgSymMatrix(5) & Tp_cov = *Tp.covariance();
57
58 int nLinksFound = 0;
59 /// predicted local position
60 double P0 = Tp.parameters()[0];
61 double P1 = Tp.parameters()[1];
62 /// covariance for prediction
63 double PV0 = Tp_cov(0, 0); /// X,X
64 double PV1 = Tp_cov(0, 1); /// X,Y
65 double PV2 = Tp_cov(1, 1); /// Y,Y
66 /// chi² cut values
67 double Xc = m_xi2maxlink;
68 double bestX2 = m_xi2maxlink;
69 double Xm = m_xi2max ;
70 int nMaxClusters = m_tools->isITkGeometry() ? m_tools->maxclusters() : 9;
71 int nMaxClustersForBestChi2 = m_tools->isITkGeometry() ? nMaxClusters-1 : nMaxClusters;
72
73 /// holder for best cluster seen not pasing cuts
74 const InDet::SiCluster* cl = nullptr;
75
76 T* sibegin = std::any_cast<T>(&m_sibegin);
77 T* siend = std::any_cast<T>(&m_siend);
78 if (sibegin==nullptr or siend==nullptr) return 0;
79
80 /// loop over all clusters on this element
81 for (T p=*sibegin; p!=*siend; ++p) {
82 const InDet::SiCluster* c = static_cast<const InDet::SiCluster*>(*p);
83 const Amg::Vector2D& M = c->localPosition();
84 const Amg::MatrixX& V = c->localCovariance();
85
86 double MV0 = V(0,0);
87 double MV1 = V(1,0);
88 double MV2 = V(1,1);
89 /// add the two covariances
90 double v0 = MV0+PV0;
91 double v1 = MV1+PV1;
92 double v2 = MV2+PV2;
93 /// get the distance from the cluster to the predicted location
94 double r0 = M[0]-P0;
95 double r1 = M[1]-P1;
96 /// determinant of the cov matrix
97 double d = v0*v2-v1*v1;
98 if(d<=0.) continue;
99 d=1./d;
100
101 /// chi² calculation
102 double x = (r0*(r0*v2-r1*v1)+r1*(r1*v0-r0*v1))*d;
103 if(x > Xc) continue; /// if we don't satisfy the minimum, don't bother further
104
105 /// refine y estimate
106 r1 = fabs(r1+d*((PV1*v2-PV2*v1)*r0+(PV2*v0-PV1*v1)*r1));
107 x -= (r1*r1)/MV2 ;
108 r1 -= m_halflength ;
109
110 /// update the chi² estimate, and check again the cut
111 if(r1 > 0.){
112 x+=(r1*r1)/PV2;
113 if (x > Xc) continue;
114 }
115 /// if we satisfy the chi2 criterion:
116 if(x < Xm) {
117 /// create a cluster link to this cluster
118 InDet::SiClusterLink_xk l(c,x);
119
120 /// This misleadingly named method actually is a sorting
121 /// Mechanism: Our current link will take the place of the one
122 /// it replaces in ascending chi², and we will be
123 /// left with a link to the worst chi² one
124 for(int i=0; i!=nLinksFound; ++i) L[i].Comparison(l);
125 /// Now our link points not the the cluster we just processed but the worst
126 /// so far
127
128 /// If we have less than nm saved, save it in the next slot
129 if(nLinksFound<=nMaxClusters) {
130 L[nLinksFound++]=l;
131 } else {
132 /// otherwise, don't add it, and update the chi2 cut to
133 /// the worst one in our list, to get rid of further
134 /// candidates early on.
135 Xm=L[nMaxClustersForBestChi2].xi2();
136 }
137
138 /// finally, update the chi2 cut - don't bother with anything worse
139 /// than the current candidate by 6 units in chi2
140 Xc = Xm+6.;
141 } /// done with branch if we satisfy the chi2 criterion
142 /// if we have no links yet and this is the best chi2 so far
143 else if(!nLinksFound && x < bestX2) {
144 /// update best chi2
145 bestX2 = x;
146 /// update max cut
147 Xc = x+6.;
148 /// update best cluster pointer
149 cl = c;
150 }
151 } /// end cluster loop
152 /// If we didn't find any cluster satisfying the chi2 cut,
153 /// but have at least one cluster,
154 /// we write one link to the best (of the bad) candidate
155 if(cl && !nLinksFound) {L[nLinksFound++].Set(cl,bestX2);}
156 /// return number of found links
157 return nLinksFound;
158}
159
160///////////////////////////////////////////////////////////////////
161// Search closest cluster without stereo angle for pixels
162///////////////////////////////////////////////////////////////////
163
164template <typename T>
165int InDet::SiTrajectoryElement_xk::searchClustersWithoutStereoPIX
166(Trk::PatternTrackParameters& Tp,InDet::SiClusterLink_xk* L)
167{
168 if (m_detstatus<=0) return 0;
169
170 const AmgSymMatrix(5) & Tp_cov = *Tp.covariance();
171 int nLinksFound = 0;
172 /// predicted local position
173 double P0 = Tp.parameters()[0];
174 double P1 = Tp.parameters()[1];
175 /// covariance of prediction
176 double PV0 = Tp_cov(0, 0);
177 double PV1 = Tp_cov(0, 1);
178 double PV2 = Tp_cov(1, 1);
179 /// chi2 cuts
180 double Xc = m_xi2maxlink;
181 double Xm = m_xi2max ;
182 int nMaxClusters = m_tools->isITkGeometry() ? m_tools->maxclusters() : 9;
183 int nMaxClustersForBestChi2 = m_tools->isITkGeometry() ? nMaxClusters-1 : nMaxClusters;
184
185 /// best cluster seen
186 const InDet::SiCluster* cl = nullptr;
187
188 T* sibegin = std::any_cast<T>(&m_sibegin);
189 T* siend = std::any_cast<T>(&m_siend);
190 if (sibegin==nullptr or siend==nullptr) return 0;
191
192 /// loop over all hits
193 for (T p=*sibegin; p!=*siend; ++p) {
194 const InDet::SiCluster* c = static_cast<const InDet::SiCluster*>(*p);
195 const Amg::Vector2D& M = c->localPosition();
196
197 /// estimate the covariance by the pitch
198 /// Factors 1/sqrt(12) will be added later
199 double MV0 = c->width().phiR();
200 double MV2 = c->width().z ();
201
202 /// squared distance between prediction and hit location,
203 double r0 = M[0]-P0,
204 r02 = r0*r0;
205 double r1 = M[1]-P1,
206 r12 = r1*r1;
207
208 /// sum cov on X - using 1/sqrt(12) update of hit cov
209 double v0 = s_oneOverTwelve*(MV0*MV0)+PV0;
210 /// bail out if X alone sufficient to fail chi2
211 if(r02 >(Xc*v0)) continue;
212
213 /// sum cov on Y - using 1/sqrt(12) update of hit cov
214 double v2 = s_oneOverTwelve*(MV2*MV2)+PV2;
215 /// bail out if Y alone sufficient to fail chi2
216 if(r12 >(Xc*v2)) continue;
217 double v1 = PV1;
218 /// determinant of cov matrix
219 double d = v0*v2-v1*v1;
220 if( d<=0. ) continue;
221 /// chi2
222 double x = (r02*v2+r12*v0-(r0*r1)*(2.*v1))/d;
223
224 /// skip clusters failing the cut
225 if(x>Xc) continue;
226
227 /// if we satisfy the tighter requirement:
228 if(x < Xm) {
229 /// set up a cluster link
230 InDet::SiClusterLink_xk l(c,x);
231 /// This misleadingly named method actually is a sorting
232 /// Mechanism: Our current link will take the place of the one
233 /// it replaces in ascending chi², and we will be
234 /// left with a link to the worst chi² one
235 for(int i=0; i!=nLinksFound; ++i) L[i].Comparison(l);
236 /// and place the worst link (remember that l is now a different link, see above)!
237 /// If within capacity, increment counter and add
238 if(nLinksFound<=nMaxClusters)
239 L[nLinksFound++]=l;
240 /// otherwise increment cut to avoid bothering with lower qualities in the future
241 else Xm=L[nMaxClustersForBestChi2].xi2();
242 /// update looser chi2 cut to match the tighter one - we found at least one good cluster
243 Xc = Xm;
244 }
245 /// if we found no clusters yet satisfying the chi2 cut, use this one
246 else if(!nLinksFound) {
247 Xc = x; cl = c;
248 }
249 } /// end of cluster loop
250 /// if we found no cluster satisfying the chi2, we return one link to the best (least bad) we saw
251 if(cl && !nLinksFound) {
252 L[nLinksFound++].Set(cl,Xc);
253 }
254 return nLinksFound;
255}
256
257///////////////////////////////////////////////////////////////////
258// Search closest cluster without stereo angle for SCT
259///////////////////////////////////////////////////////////////////
260
261template <typename T>
262int InDet::SiTrajectoryElement_xk::searchClustersWithoutStereoSCT
263(Trk::PatternTrackParameters& Tp,InDet::SiClusterLink_xk* L)
264{
265 if (m_detstatus<=0) return 0;
266
267 const AmgSymMatrix(5) & Tp_cov = *Tp.covariance();
268 int nLinksFound = 0;
269 /// predicted local position
270 double P0 = Tp.parameters()[0];
271 double P1 = Tp.parameters()[1];
272 /// and cov of prediction
273 double PV0 = Tp_cov(0, 0);
274 double PV1 = Tp_cov(0, 1);
275 double PV2 = Tp_cov(1, 1);
276 /// chi2 cuts
277 double Xc = m_xi2maxlink;
278 double Xm = m_xi2max ;
279 int nMaxClusters = m_tools->isITkGeometry() ? m_tools->maxclusters() : 9;
280 int nMaxClustersForBestChi2 = m_tools->isITkGeometry() ? nMaxClusters-1 : nMaxClusters;
281 /// best cluster seen
282 const InDet::SiCluster* cl = nullptr;
283
284 T* sibegin = std::any_cast<T>(&m_sibegin);
285 T* siend = std::any_cast<T>(&m_siend);
286 if (sibegin==nullptr or siend==nullptr) return 0;
287
288 /// loop over all hits on this element
289 for (T p=*sibegin; p!=*siend; ++p) {
290 const InDet::SiCluster* c = static_cast<const InDet::SiCluster*>(*p);
291 const Amg::Vector2D& M = c->localPosition();
292
293 /// estimate covariance by pitch
294 double MV0 = c->width().phiR() ;
295 /// include 1/sqrt(12) factor and add to contribution from pred
296 double v0 = s_oneOverTwelve*(MV0*MV0)+PV0;
297 /// get residual in X
298 double r0 = M[0]-P0;
299 /// one over cov
300 double d = 1./v0;
301 /// chi2
302 double x = (r0*r0)*d;
303
304 /// stop if chi2 too large
305 if(x>Xc) continue;
306 /// estimate of Y residual
307 double dP1 = (P1-M[1])+PV1*d*r0;
308
309 /// if residual beyond half-length: assume hit is on very end of strip and
310 /// re-evaluate chi2 with a penalty term on the Y coordinate
311 if(std::abs(dP1) > m_halflength) {
312 double r1 = M[1]-P1;
313 if(m_tools->isITkGeometry()){
314 if(dP1 > m_halflength) r1 += m_halflength;
315 else r1 -= m_halflength;
316 }
317 else{
318 r1 = (dP1 > m_halflength ?m_halflength-P1 : -(m_halflength+P1));
319 }
320 double v1 = PV1;
321 double v2 = PV2;
322 /// determinant of 2D covariance (only contribution from hit on X,X)
323 d = v0*v2-v1*v1 ;
324 if(d<=0.) continue;
325 /// updated chi2 estimate, with penalty for being at end of strip
326 x = (r0*(r0*v2-r1*v1)+r1*(r1*v0-r0*v1))/d;
327 /// and re-apply cut
328 if(x>Xc) continue;
329 } /// end of branch for Y residual beyond half-length
330
331 /// if we satisfy the minimum criterion:
332 if(x < Xm) {
333 /// create a cluster link for our firend
334 InDet::SiClusterLink_xk l(c,x);
335 /// This misleadingly named method actually is a sorting
336 /// Mechanism: Our current link will take the place of the one
337 /// it replaces in ascending chi², and we will be
338 /// left with a link to the worst chi² one
339 for(int i=0; i!=nLinksFound; ++i) L[i].Comparison(l);
340 /// and place the worst link so far at the end if space.
341 if(nLinksFound<=nMaxClusters)
342 L[nLinksFound++]=l;
343 /// otherwise, update chi2 cut to not bother with worse clusters anymore
344 else Xm=L[nMaxClustersForBestChi2].xi2();
345 /// update looser chi2 to match the tighter one - we have at least one good link now
346 Xc = Xm;
347 }
348 /// if we didn't find any good links yet, keep this one as best so far
349 /// update looser chi2 cut to only keep better candidates in the future
350 else if(!nLinksFound) {Xc = x; cl = c;}
351 }
352 /// if we found no clusters satisfying the chi2 cut, we return
353 /// a single link to the best (least bad) candidate
354 if(cl && !nLinksFound) {L[nLinksFound++].Set(cl,Xc);}
355 return nLinksFound;
356}
357
358///////////////////////////////////////////////////////////////////
359// Search closest cluster with stereo angle
360///////////////////////////////////////////////////////////////////
361
362template <typename T>
363int InDet::SiTrajectoryElement_xk::searchClustersWithStereoAss
364(Trk::PatternTrackParameters& Tp,InDet::SiClusterLink_xk* L, const Trk::PRDtoTrackMap &prd_to_track_map)
365{
366
367 /// for detailed documentation, please see searchClustersWithStereo
368
369 if (m_detstatus<=0) return 0;
370
371 const AmgSymMatrix(5) & Tp_cov = *Tp.covariance();
372 int nLinksFound = 0;
373 double P0 = Tp.parameters()[0];
374 double P1 = Tp.parameters()[1];
375 double PV0 = Tp_cov(0, 0);
376 double PV1 = Tp_cov(0, 1);
377 double PV2 = Tp_cov(1, 1);
378 double Xc = m_xi2maxlink;
379 double Xl = m_xi2maxlink;
380 double Xm = m_xi2max ;
381 int nMaxClusters = m_tools->isITkGeometry() ? m_tools->maxclusters() : 9;
382 int nMaxClustersForBestChi2 = m_tools->isITkGeometry() ? nMaxClusters-1 : nMaxClusters;
383
384 const InDet::SiCluster* cl = nullptr;
385
386 T* sibegin = std::any_cast<T>(&m_sibegin);
387 T* siend = std::any_cast<T>(&m_siend);
388 if (sibegin==nullptr or siend==nullptr) return 0;
389 for (T p=*sibegin; p!=*siend; ++p) {
390 const InDet::SiCluster* c = static_cast<const InDet::SiCluster*>(*p);
391
392 /// check the PRD association map, skip clusters already in use
393 if (prd_to_track_map.isUsed(*c)) continue;
394
395 const Amg::Vector2D& M = c->localPosition();
396 const Amg::MatrixX& V = c->localCovariance();
397
398 double MV0 = V(0,0);
399 double MV1 = V(1,0);
400 double MV2 = V(1,1);
401 double v0 = MV0+PV0;
402 double v1 = MV1+PV1;
403 double v2 = MV2+PV2;
404 double r0 = M[0]-P0;
405 double r1 = M[1]-P1;
406 double d = v0*v2-v1*v1;
407 if(d<=0.) continue;
408 d=1./d;
409 double x = (r0*(r0*v2-r1*v1)+r1*(r1*v0-r0*v1))*d;
410 if(x > Xc) continue;
411
412 r1 = fabs(r1+d*((PV1*v2-PV2*v1)*r0+(PV2*v0-PV1*v1)*r1));
413 x -= (r1*r1)/MV2 ;
414 r1 -= m_halflength ;
415
416 if(r1 > 0. && (x+=((r1*r1)/PV2)) > Xc) continue;
417
418 if(x < Xm) {
419 InDet::SiClusterLink_xk l(c,x);
420 /// This misleadingly named method actually is a sorting
421 /// Mechanism: Our current link will take the place of the one
422 /// it replaces in ascending chi², and we will be
423 /// left with a link to the worst chi² one
424 for(int i=0; i!=nLinksFound; ++i) L[i].Comparison(l);
425 if(nLinksFound<=nMaxClusters)
426 L[nLinksFound++]=l;
427 else Xm=L[nMaxClustersForBestChi2].xi2();
428 Xc = Xm+6.;
429 }
430 else if(!nLinksFound && x < Xl) {Xl = x; Xc = x+6.; cl = c;}
431 }
432 if(cl && !nLinksFound) {L[nLinksFound++].Set(cl,Xl);}
433 return nLinksFound;
434}
435
436///////////////////////////////////////////////////////////////////
437// Search closest cluster without stereo angle for pixels
438///////////////////////////////////////////////////////////////////
439
440template <typename T>
441int InDet::SiTrajectoryElement_xk::searchClustersWithoutStereoAssPIX
442(Trk::PatternTrackParameters& Tp,InDet::SiClusterLink_xk* L, const Trk::PRDtoTrackMap &prd_to_track_map)
443{
444 /// for detailed documentation, please see searchClustersWithoutStereoPIX
445 if (m_detstatus<=0) return 0;
446
447 const AmgSymMatrix(5) & Tp_cov = *Tp.covariance();
448 int nLinksFound = 0;
449 double P0 = Tp.parameters()[0];
450 double P1 = Tp.parameters()[1];
451 double PV0 = Tp_cov(0, 0);
452 double PV1 = Tp_cov(0, 1);
453 double PV2 = Tp_cov(1, 1);
454 double Xc = m_xi2maxlink;
455 double Xm = m_xi2max ;
456 int nMaxClusters = m_tools->isITkGeometry() ? m_tools->maxclusters() : 9;
457 int nMaxClustersForBestChi2 = m_tools->isITkGeometry() ? nMaxClusters-1 : nMaxClusters;
458
459 const InDet::SiCluster* cl = nullptr;
460
461 T* sibegin = std::any_cast<T>(&m_sibegin);
462 T* siend = std::any_cast<T>(&m_siend);
463 if (sibegin==nullptr or siend==nullptr) return 0;
464
465 for (T p=*sibegin; p!=*siend; ++p) {
466 const InDet::SiCluster* c = static_cast<const InDet::SiCluster*>(*p);
467 if (prd_to_track_map.isUsed(*c)) continue;
468 const Amg::Vector2D& M = c->localPosition();
469
470 double MV0 = c->width().phiR();
471 double MV2 = c->width().z ();
472
473 double r0 = M[0]-P0, r02 = r0*r0;
474 double r1 = M[1]-P1, r12 = r1*r1;
475
476 double v0 = s_oneOverTwelve*(MV0*MV0)+PV0; if(r02 >(Xc*v0)) continue;
477 double v2 = s_oneOverTwelve*(MV2*MV2)+PV2; if(r12 >(Xc*v2)) continue;
478 double v1 = PV1;
479 double d = v0*v2-v1*v1; if( d<=0. ) continue;
480 double x = (r02*v2+r12*v0-(r0*r1)*(2.*v1))/d;
481
482 if(x>Xc) continue;
483
484 if(x < Xm) {
485 InDet::SiClusterLink_xk l(c,x);
486 /// This misleadingly named method actually is a sorting
487 /// Mechanism: Our current link will take the place of the one
488 /// it replaces in ascending chi², and we will be
489 /// left with a link to the worst chi² one
490 for(int i=0; i!=nLinksFound; ++i) L[i].Comparison(l);
491 if(nLinksFound<=nMaxClusters)
492 L[nLinksFound++]=l;
493 else Xm=L[nMaxClustersForBestChi2].xi2();
494 Xc = Xm;
495 }
496 else if(!nLinksFound) {Xc = x; cl = c;}
497 }
498 if(cl && !nLinksFound) {L[nLinksFound++].Set(cl,Xc);}
499 return nLinksFound;
500}
501
502///////////////////////////////////////////////////////////////////
503// Search closest cluster without stereo angle for SCT
504///////////////////////////////////////////////////////////////////
505
506template <typename T>
507int InDet::SiTrajectoryElement_xk::searchClustersWithoutStereoAssSCT
508(Trk::PatternTrackParameters& Tp,InDet::SiClusterLink_xk* L, const Trk::PRDtoTrackMap &prd_to_track_map)
509{
510 /// for detailed documentation, please see searchClustersWithoutStereoSCT
511 if (m_detstatus<=0) return 0;
512
513 const AmgSymMatrix(5) & Tp_cov = *Tp.covariance();
514 int nLinksFound = 0;
515 double P0 = Tp.parameters()[0];
516 double P1 = Tp.parameters()[1];
517 double PV0 = Tp_cov(0, 0);
518 double PV1 = Tp_cov(0, 1);
519 double PV2 = Tp_cov(1, 1);
520 double Xc = m_xi2maxlink;
521 double Xm = m_xi2max ;
522 int nMaxClusters = m_tools->isITkGeometry() ? m_tools->maxclusters() : 9;
523 int nMaxClustersForBestChi2 = m_tools->isITkGeometry() ? nMaxClusters-1 : nMaxClusters;
524
525 const InDet::SiCluster* cl = nullptr;
526
527 T* sibegin = std::any_cast<T>(&m_sibegin);
528 T* siend = std::any_cast<T>(&m_siend);
529 if (sibegin==nullptr or siend==nullptr) return 0;
530
531 for (T p=*sibegin; p!=*siend; ++p) {
532 const InDet::SiCluster* c = static_cast<const InDet::SiCluster*>(*p);
533 if (prd_to_track_map.isUsed(*c)) continue;
534 const Amg::Vector2D& M = c->localPosition();
535
536 double MV0 = c->width().phiR() ;
537 double v0 = s_oneOverTwelve*(MV0*MV0)+PV0;
538 double r0 = M[0]-P0;
539 double d = 1./v0;
540 double x = (r0*r0)*d;
541
542 if(x>Xc) continue;
543
544 double dP1 = (P1-M[1])+PV1*d*r0;
545
546 if(fabs(dP1) > m_halflength) {
547
548 double r1 = M[1]-P1;
549 if(m_tools->isITkGeometry()){
550 if(dP1 > m_halflength) r1 += m_halflength;
551 else r1 -= m_halflength;
552 }
553 else{
554 r1 = (dP1 > m_halflength ?m_halflength-P1 : -(m_halflength+P1));
555 }
556
557 double v1 = PV1;
558 double v2 = PV2;
559 d = v0*v2-v1*v1 ; if(d<=0.) continue;
560 x = (r0*(r0*v2-r1*v1)+r1*(r1*v0-r0*v1))/d;
561 if(x>Xc) continue;
562 }
563
564 if(x < Xm) {
565 InDet::SiClusterLink_xk l(c,x);
566 /// This misleadingly named method actually is a sorting
567 /// Mechanism: Our current link will take the place of the one
568 /// it replaces in ascending chi², and we will be
569 /// left with a link to the worst chi² one
570 for(int i=0; i!=nLinksFound; ++i) L[i].Comparison(l);
571 if(nLinksFound<=nMaxClusters)
572 L[nLinksFound++]=l;
573 else Xm=L[nMaxClustersForBestChi2].xi2();
574 Xc = Xm;
575 }
576 else if(!nLinksFound) {Xc = x; cl = c;}
577 }
578 if(cl && !nLinksFound) {L[nLinksFound++].Set(cl,Xc);}
579 return nLinksFound;
580}
581
582///////////////////////////////////////////////////////////////////
583// Set trajectory element
584///////////////////////////////////////////////////////////////////
585
586// T = InDet::SiClusterCollection::const_iterator or
587// InDet::PixelClusterCollection::const_iterator or
588// InDet::SCT_ClusterCollection::const_iterator
589template <typename T>
590bool InDet::SiTrajectoryElement_xk::set
591(int st,
592 const InDet::SiDetElementBoundaryLink_xk*& dl,
593 const T& sb,
594 const T& se,
595 const InDet::SiCluster* si,
596 [[maybe_unused]] const EventContext& ctx)
597{
598 if (std::is_same<T, InDet::SiClusterCollection::const_iterator>::value) m_itType = SiClusterColl;
599 else if (std::is_same<T, InDet::PixelClusterCollection::const_iterator>::value) m_itType = PixelClusterColl;
600 else if (std::is_same<T, InDet::SCT_ClusterCollection::const_iterator>::value) m_itType = SCT_ClusterColl;
601
602 m_fieldMode = false;
603 if(m_tools->fieldTool().magneticFieldMode()!=0) m_fieldMode = true;
604 m_status = 0 ;
605 m_detstatus = st ;
606 m_nMissing = 0 ;
607 m_nlinksForward = 0 ;
608 m_nlinksBackward = 0 ;
609 m_nholesForward = 0 ;
610 m_nholesBackward = 0 ;
611 m_dholesForward = 0 ;
612 m_dholesBackward = 0 ;
613 m_nclustersForward = 0 ;
614 m_nclustersBackward = 0 ;
615 m_npixelsBackward = 0 ;
616 m_ndfForward = 0 ;
617 m_ndfBackward = 0 ;
618 m_ntsos = 0 ;
619 m_detelement = dl->detElement() ;
620 m_detlink = dl ;
621 m_surface = &m_detelement->surface(); if(!m_surface && m_tools->isITkGeometry()) return false;
622 m_sibegin = sb ;
623 m_siend = se ;
624 m_cluster = si ;
625 m_clusterOld = si ;
626 m_clusterNoAdd = 0 ;
627 m_stereo = false ;
628 m_xi2Forward = 10000. ;
629 m_xi2Backward = 10000. ;
630 m_xi2totalForward = 0. ;
631 m_xi2totalBackward = 0. ;
632 m_tools->electron() ? m_xi2max = m_tools->xi2maxBrem() : m_xi2max = m_tools->xi2max();
633 m_halflength = 0. ;
634 m_detelement->isSCT() ? m_ndf=1 : m_ndf=2;
635
636 if(m_tools->heavyion()) {
637 if(m_ndf==2) {m_xi2max = 13.; m_xi2maxNoAdd = 13.;}
638 else {m_xi2max = 4.; m_xi2maxNoAdd = 8.;}
639 }
640
641 noiseInitiate() ;
642 (m_detelement->isSCT() && (m_detelement->design().shape()==InDetDD::Trapezoid || m_detelement->design().shape()==InDetDD::Annulus) ) ?
643 m_stereo = true : m_stereo = false;
644
645 if(m_detstatus && m_ndf == 1){
646 m_halflength = si ? si->width().z()*.5 : (*sb)->width().z()*.5;
647 }
648
649 if(!m_detstatus) {
650
651 IdentifierHash idHash = m_detelement->identifyHash();
652 if (m_ndf==2) {
653 if (m_tools->pixelStatus() ? !m_tools->pixelStatus()->at(idHash) : !m_tools->pixcond()->isGood(idHash, ctx) ) {
654 m_detstatus = -1;
655 }
656 }
657 else {
658 if (m_tools->sctStatus() ? !m_tools->sctStatus()->at(idHash) : !m_tools->sctcond()->isGood(idHash, ctx) ) {
659 m_detstatus = -1;
660 }
661 }
662 VALIDATE_STATUS_ARRAY(m_tools->pixelStatus() && m_tools->sctStatus(), ((m_ndf==2) ? m_tools->pixelStatus() : m_tools->sctStatus())->at(idHash), ( (m_ndf==2) ? m_tools->pixcond() : m_tools->sctcond())->isGood(idHash, ctx));
663
664 }
665
666 const Amg::Transform3D& tr = m_surface->transform();
667 //Initialise material corrections
668 m_radlength = 0.03;
669 if(m_tools->isITkGeometry()) m_radlength = 0.04;
670 /// add material for the forward pixels
671 if(m_ndf == 2 && std::abs(tr(2,2)) >= 1.) m_radlength = .07;
672 if(m_tools->isITkGeometry()) m_radlengthN=m_radlength;
673
674 m_localTransform[ 0] = tr(0,0); m_localTransform[ 1]=tr(1,0); m_localTransform[ 2]=tr(2,0);
675 m_localTransform[ 3] = tr(0,1); m_localTransform[ 4]=tr(1,1); m_localTransform[ 5]=tr(2,1);
676 m_localTransform[ 6] = tr(0,2); m_localTransform[ 7]=tr(1,2); m_localTransform[ 8]=tr(2,2);
677 m_localTransform[ 9] = tr(0,3); m_localTransform[10]=tr(1,3); m_localTransform[11]=tr(2,3);
678 m_localTransform[12] = m_localTransform[ 9]*m_localTransform[ 6]+m_localTransform[10]*m_localTransform[ 7]+m_localTransform[11]*m_localTransform[ 8];
679 m_localDir[0] = 1.;
680 m_localDir[1] = 0.;
681 m_localDir[2] = 0.;
682 return true;
683}
684
685///////////////////////////////////////////////////////////////////
686// Forward propagation for search closest cluster
687///////////////////////////////////////////////////////////////////
688
689// T = InDet::SiClusterCollection::const_iterator or
690// InDet::PixelClusterCollection::const_iterator or
691// InDet::SCT_ClusterCollection::const_iterator
692template <typename T>
693bool InDet::SiTrajectoryElement_xk::ForwardPropagationForClusterSeach
694(int n,
695 const Trk::TrackParameters& Tpa,
696 const InDet::SiDetElementBoundaryLink_xk*& dl,
697 const T& sb,
698 const T& se,
699 const EventContext& ctx)
700{
701 // remove case if you have trajectory element without actual detector element
702 // this happens if you have added a dead cylinder
703 if(!m_detelement) {
704 return false;
705 }
706
707 if (std::is_same<T, InDet::SiClusterCollection::const_iterator>::value) m_itType = SiClusterColl;
708 else if (std::is_same<T, InDet::PixelClusterCollection::const_iterator>::value) m_itType = PixelClusterColl;
709 else if (std::is_same<T, InDet::SCT_ClusterCollection::const_iterator>::value) m_itType = SCT_ClusterColl;
710
711 m_detstatus = 1 ;
712 m_sibegin = sb ;
713 m_siend = se ;
714 m_detelement = dl->detElement() ;
715 m_detlink = dl ;
716 m_surface = &m_detelement->surface();
717 m_detelement->isSCT() ? m_ndf=1 : m_ndf=2;
718 m_halflength = 0. ;
719 m_stereo = false ;
720
721 (m_detelement->isSCT() && m_detelement->design().shape()==InDetDD::Trapezoid) ?
722 m_stereo = true : m_stereo = false;
723
724 if(m_detstatus && m_ndf == 1) m_halflength = (*sb)->width().z()*.5;
725
726 if(!n) {
727 Trk::PatternTrackParameters Tp; if(!Tp.production(&Tpa)) return false;
728 if(!propagateParameters(Tp,m_parametersPredForward,m_step)) return false;
729
730 if(!m_parametersPredForward.iscovariance()) {
731
732 double cv[15]={ .02 ,
733 0., .02,
734 0., 0.,.000001,
735 0., 0., 0.,.000001,
736 0., 0., 0., 0.,.000000001};
737
738 m_parametersPredForward.setCovariance(cv);
739 }
740 }
741 else {
742 if(!propagate(m_parametersPredForward,m_parametersPredForward,m_step,ctx)) return false;
743 }
744 m_nlinksForward=searchClusters(m_parametersPredForward,m_linkForward);
745 return true;
746}
747
748///////////////////////////////////////////////////////////////////
749// Forward propagation for search closest cluster
750///////////////////////////////////////////////////////////////////
751
752// T = InDet::SiClusterCollection::const_iterator or
753// InDet::PixelClusterCollection::const_iterator or
754// InDet::SCT_ClusterCollection::const_iterator
755template <typename T>
756void InDet::SiTrajectoryElement_xk::CloseClusterSeach
757(Trk::PatternTrackParameters& Tpa,
758 const InDet::SiDetElementBoundaryLink_xk*& dl,
759 const T& sb,
760 const T& se)
761{
762 if (std::is_same<T, InDet::SiClusterCollection::const_iterator>::value) m_itType = SiClusterColl;
763 else if (std::is_same<T, InDet::PixelClusterCollection::const_iterator>::value) m_itType = PixelClusterColl;
764 else if (std::is_same<T, InDet::SCT_ClusterCollection::const_iterator>::value) m_itType = SCT_ClusterColl;
765
766 m_detstatus = 1 ;
767 m_cluster = 0 ;
768 m_sibegin = sb ;
769 m_siend = se ;
770 m_detelement = dl->detElement() ;
771 m_detlink = dl ;
772 m_surface = &m_detelement->surface();
773 m_detelement->isSCT() ? m_ndf=1 : m_ndf=2;
774 m_halflength = 0. ;
775 m_stereo = false ;
776
777 (m_detelement->isSCT() && m_detelement->design().shape()==InDetDD::Trapezoid) ?
778 m_stereo = true : m_stereo = false;
779 if(m_detstatus && m_ndf == 1) m_halflength = (*sb)->width().z()*.5;
780 /// in this case, we do not wish to overwrite the m_inside member.
781 /// So store the old value and restore after the check.
782 int inside_originalVal = m_inside;
783 checkBoundaries(Tpa);
784 if(m_inside > 0 || !searchClusters(Tpa,m_linkForward)) {
785 m_inside = inside_originalVal;
786 return;
787 }
788 m_inside = inside_originalVal;
789 m_cluster = m_linkForward[0].cluster();
790 m_xi2Forward = m_linkForward[0].xi2 ();
791}
792
793
794inline const Trk::PatternTrackParameters *
795InDet::SiTrajectoryElement_xk::parameters() const
796{
797 // logic from SiTrajectoryElement_xk::trackParameters(bool, int Q) with Q==1
798
799 if (m_status == 1) {
800 if (m_cluster) {
801 return &m_parametersUpdatedForward;
802 } else {
803 return &m_parametersPredForward;
804 }
805 } else if (m_status == 2) {
806 if (m_cluster) {
807 return &m_parametersUpdatedBackward;
808 } else {
809 return &m_parametersPredBackward;
810 }
811 } else if (m_status == 3) {
812 if (m_cluster) {
813 return &m_parametersUpdatedBackward;
814 }
815 }
816 return nullptr;
817}
818
Generated on for ATLAS Offline Software by doxygen 1.14.0