ATLAS Offline Software
Loading...
Searching...
No Matches
TrigInDetTrackFollowingTool Class Reference

#include <TrigInDetTrackFollowingTool.h>

Inheritance diagram for TrigInDetTrackFollowingTool:
Collaboration diagram for TrigInDetTrackFollowingTool:

Public Member Functions

 TrigInDetTrackFollowingTool (const std::string &, const std::string &, const IInterface *)
virtual StatusCode initialize ()
virtual StatusCode finalize ()
virtual Trk::TrackgetTrack (const std::vector< const Trk::SpacePoint * > &, const std::vector< const InDetDD::SiDetectorElement * > &, const EventContext &) const
ServiceHandle< StoreGateSvc > & evtStore ()
 The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.
const ServiceHandle< StoreGateSvc > & detStore () const
 The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.
virtual StatusCode sysInitialize () override
 Perform system initialization for an algorithm.
virtual StatusCode sysStart () override
 Handle START transition.
virtual std::vector< Gaudi::DataHandle * > inputHandles () const override
 Return this algorithm's input handles.
virtual std::vector< Gaudi::DataHandle * > outputHandles () const override
 Return this algorithm's output handles.
Gaudi::Details::PropertyBase & declareProperty (Gaudi::Property< T, V, H > &t)
void updateVHKA (Gaudi::Details::PropertyBase &)
MsgStream & msg () const
bool msgLvl (const MSG::Level lvl) const

Static Public Member Functions

static const InterfaceID & interfaceID ()

Protected Member Functions

void renounceArray (SG::VarHandleKeyArray &handlesArray)
 remove all handles from I/O resolution
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce (T &h)
void extraDeps_update_handler (Gaudi::Details::PropertyBase &ExtraDeps)
 Add StoreName to extra input/output deps as needed.

Private Types

typedef ServiceHandle< StoreGateSvcStoreGateSvc_t

Private Member Functions

std::unique_ptr< TrigFTF_ExtendedTrackStatefitTheSeed (const std::vector< const Trk::SpacePoint * > &, MagField::AtlasFieldCache &) const
const Trk::PrepRawDataupdateTrackState (const InDet::PixelCluster *, TrigFTF_ExtendedTrackState &, bool) const
const Trk::PrepRawDataupdateTrackState (const InDet::SCT_Cluster *, int, TrigFTF_ExtendedTrackState &) const
int extrapolateTrackState (TrigFTF_ExtendedTrackState &, const Trk::PlaneSurface *, MagField::AtlasFieldCache &) const
int RungeKutta34 (double *, double *, const Trk::PlaneSurface *, MagField::AtlasFieldCache &, bool) const
bool tentativeExtrapolation (double const *, double *, const Trk::PlaneSurface *, const Trk::PlaneSurface *, MagField::AtlasFieldCache &) const
double processHit (const InDet::PixelCluster *, double *, double *, const TrigFTF_ExtendedTrackState &) const
double processHit (const InDet::SCT_Cluster *, int, double &, double &, double *, const TrigFTF_ExtendedTrackState &) const
void findNearestHit (int, const InDet::PixelClusterCollection *, const double *, std::vector< std::tuple< double, const Trk::PrepRawData *, int > > &) const
void findNearestHit (int, const InDet::SCT_ClusterCollection *, int, const double *, std::vector< std::tuple< double, const Trk::PrepRawData *, int > > &) const
void crossProduct (double const *, double const *, double *) const
double estimateRK_Step (const Trk::PlaneSurface *, double const *) const
Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
 specialization for handling Gaudi::Property<SG::VarHandleKey>

Private Attributes

ToolHandle< ITrigL2LayerNumberToolm_layerNumberTool {this, "LayerNumberTool", "TrigL2LayerNumberToolITk"}
Gaudi::Property< int > m_nClustersMin {this, "nClustersMin", 7, "Minimum number of clusters on track"}
Gaudi::Property< int > m_nHolesMax {this, "nHolesMax", 100, "Maximum number of holes on track"}
Gaudi::Property< double > m_maxChi2Dist_Pixels {this, "Chi2MaxPixels", 25.0, "the Pixel hit chi2 cut"}
Gaudi::Property< double > m_maxChi2Dist_Strips {this, "Chi2MaxStrips", 12.0, "the Strip hit chi2 cut"}
Gaudi::Property< double > m_winX_Pixels {this, "XSearchWindowPixels", 3.0, "x-size of hit search window for Pixels"}
Gaudi::Property< double > m_winY_Pixels {this, "YSearchWindowPixels", 3.0, "y-size of hit search window for Pixels"}
Gaudi::Property< double > m_winX_Strips {this, "XSearchWindowStrips", 3.0, "x-size of hit search window for Strips"}
Gaudi::Property< bool > m_useHitErrors {this, "UseHitErrors", false, "use PrepRawData errors"}
Gaudi::Property< bool > m_useDetectorThickness {this, "UseDetectorThickness", false, "get Si-modules thickness from InDet Geometry"}
Gaudi::Property< double > m_nominalRadLength {this, "ModuleRadLength", 0.04, "fixed radiation thickness of the detector modules"}
SG::ReadCondHandleKey< AtlasFieldCacheCondObjm_fieldCondObjInputKey {this, "AtlasFieldCacheCondObj", "fieldCondObj", "Name of the Magnetic Field conditions object key"}
SG::ReadHandleKey< InDet::PixelClusterContainer > m_pixcontainerkey {this, "PixelClusterContainer","ITkPixelClusters"}
SG::ReadHandleKey< InDet::SCT_ClusterContainer > m_sctcontainerkey {this, "SCT_ClusterContainer", "ITkStripClusters"}
StoreGateSvc_t m_evtStore
 Pointer to StoreGate (event store by default)
StoreGateSvc_t m_detStore
 Pointer to StoreGate (detector store by default)
std::vector< SG::VarHandleKeyArray * > m_vhka
bool m_varHandleArraysDeclared

Detailed Description

Definition at line 83 of file TrigInDetTrackFollowingTool.h.

Member Typedef Documentation

◆ StoreGateSvc_t

typedef ServiceHandle<StoreGateSvc> AthCommonDataStore< AthCommonMsg< AlgTool > >::StoreGateSvc_t
privateinherited

Definition at line 388 of file AthCommonDataStore.h.

Constructor & Destructor Documentation

◆ TrigInDetTrackFollowingTool()

TrigInDetTrackFollowingTool::TrigInDetTrackFollowingTool ( const std::string & t,
const std::string & n,
const IInterface * p )

Definition at line 166 of file TrigInDetTrackFollowingTool.cxx.

168 : AthAlgTool(t,n,p)
169{
170 declareInterface< ITrigInDetTrackFollowingTool >( this );
171}
AthAlgTool()
Default constructor:

Member Function Documentation

◆ crossProduct()

void TrigInDetTrackFollowingTool::crossProduct ( double const * B,
double const * V,
double * A ) const
inlineprivate

Definition at line 1665 of file TrigInDetTrackFollowingTool.cxx.

1665 {
1666 A[0] = -B[1]*V[2] + B[2]*V[1];
1667 A[1] = B[0]*V[2] - B[2]*V[0];
1668 A[2] = -B[0]*V[1] + B[1]*V[0];
1669}

◆ declareGaudiProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T, V, H > & hndl,
const SG::VarHandleKeyType &  )
inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKey>

Definition at line 156 of file AthCommonDataStore.h.

158 {
160 hndl.value(),
161 hndl.documentation());
162
163 }
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)

◆ declareProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( Gaudi::Property< T, V, H > & t)
inlineinherited

Definition at line 145 of file AthCommonDataStore.h.

145 {
146 typedef typename SG::HandleClassifier<T>::type htype;
148 }
Gaudi::Details::PropertyBase & declareGaudiProperty(Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
specialization for handling Gaudi::Property<SG::VarHandleKey>

◆ detStore()

const ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< AlgTool > >::detStore ( ) const
inlineinherited

The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 95 of file AthCommonDataStore.h.

◆ estimateRK_Step()

double TrigInDetTrackFollowingTool::estimateRK_Step ( const Trk::PlaneSurface * pN,
double const * P ) const
private

Definition at line 1353 of file TrigInDetTrackFollowingTool.cxx.

1353 {
1354
1355 double Step = 1e8;
1356
1357 if(pN == nullptr) { //step to perigee "surface" assuming the global c.s. origin at (0,0)
1358
1359 Step = -(P[0]*P[3] + P[1]*P[4])/(1 - P[5]*P[5]);
1360 return Step;
1361 }
1362
1363 const Amg::Vector3D& normal = pN->normal();
1364 const Amg::Vector3D& center = pN->center();
1365
1366 double D = 0.0;// -(r0,n)
1367
1368 for(int i=0;i<3;i++) D += -normal[i]*center[i];
1369
1370 double Sum = D;
1371 double a = 0.0;
1372
1373 for(int i=0;i<3;i++) {
1374 a += normal[i]*P[i+3];
1375 Sum += normal[i]*P[i];
1376 }
1377 if(a==0.0) return Step;
1378
1379 Step = -Sum/a;
1380
1381 return Step;
1382}
static Double_t a
static Double_t P(Double_t *tt, Double_t *par)
virtual const Amg::Vector3D & normal() const
Returns the normal vector of the Surface (i.e.
const Amg::Vector3D & center() const
Returns the center position of the Surface.
Eigen::Matrix< double, 3, 1 > Vector3D

◆ evtStore()

ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< AlgTool > >::evtStore ( )
inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 85 of file AthCommonDataStore.h.

◆ extraDeps_update_handler()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::extraDeps_update_handler ( Gaudi::Details::PropertyBase & ExtraDeps)
protectedinherited

Add StoreName to extra input/output deps as needed.

use the logic of the VarHandleKey to parse the DataObjID keys supplied via the ExtraInputs and ExtraOuputs Properties to add the StoreName if it's not explicitly given

◆ extrapolateTrackState()

int TrigInDetTrackFollowingTool::extrapolateTrackState ( TrigFTF_ExtendedTrackState & ETS,
const Trk::PlaneSurface * pN,
MagField::AtlasFieldCache & fieldCache ) const
private

Definition at line 1031 of file TrigInDetTrackFollowingTool.cxx.

1031 {
1032
1033 double Re[5];
1034
1035 memcpy(&Re[0], &ETS.m_Xk[0], sizeof(Re));
1036
1037 double P[8]; //parameters + path
1038 double Jm[40];//Jacobian
1039
1040 memset(&P[0],0,sizeof(P));
1041
1042 const Amg::Transform3D& Trf = ETS.m_pS->transform();
1043 double Ax[3] = {Trf(0,0),Trf(1,0),Trf(2,0)}; //loc1-axis in the global frame
1044 double Ay[3] = {Trf(0,1),Trf(1,1),Trf(2,1)}; //loc2-axis in the global frame
1045
1046 double gP[3];
1047 gP[0] = Trf(0,3) + Ax[0]*Re[0] + Ay[0]*Re[1];
1048 gP[1] = Trf(1,3) + Ax[1]*Re[0] + Ay[1]*Re[1];
1049 gP[2] = Trf(2,3) + Ax[2]*Re[0] + Ay[2]*Re[1];
1050
1051 double sinf, cosf;
1052 double sint, cost;
1053
1054 sincos(Re[2], &sinf, &cosf);
1055 sincos(Re[3], &sint, &cost);
1056
1057 double gV[3] = {cosf*sint, sinf*sint, cost};
1058
1059 memset(&Jm[0],0,sizeof(Jm));
1060
1061 //Track state and Jacobian initialization
1062
1063 P[6] = Re[4];
1064 P[7] = 0.0;
1065
1066 for(int i=0;i<3;i++) {
1067 P[i] = gP[i];
1068 P[i+3] = gV[i];
1069 Jm[i] = Ax[i];
1070 Jm[7+i] = Ay[i];
1071 }
1072
1073 Jm[17] =-P[4];//17
1074 Jm[18] = P[3];//18
1075 Jm[24] = cosf*cost;//24
1076 Jm[25] = sinf*cost;//25
1077 Jm[26] =-sint;//26
1078 Jm[34] = 1.0;//34
1079
1080 int code = RungeKutta34(P, Jm, pN, fieldCache, true);
1081
1082 if(code!=0) return code;
1083
1084 double J[21];
1085 memset(&J[0],0,sizeof(J));
1086
1087 double BigP[45];
1088 memset(&BigP[0],0,sizeof(BigP));
1089
1090 for(int i=0;i<7;i++) BigP[i] = P[i];
1091 for(int i=0;i<35;i++) BigP[i+7] = Jm[i];
1092
1093 if(pN == nullptr) {//special case: to perigee
1094 Trk::PerigeeSurface perSurf;
1095 Trk::RungeKuttaUtils::transformGlobalToLocal(&perSurf, true, BigP, Re, J);
1096 }
1097 else {//from plane to plane
1098 Trk::RungeKuttaUtils::transformGlobalToLocal(pN, true, BigP, Re, J);
1099 }
1100
1101 const double* Az = Trf.matrix().col(2).data();
1102
1103 // z-component of track direction vector in the local c.s.
1104
1105 double lV = Az[0]*gV[0] + Az[1]*gV[1] + Az[2]*gV[2];
1106
1107 double xOverX0 = m_nominalRadLength;
1108
1109 const InDetDD::SiDetectorElement* pDE = dynamic_cast<const InDetDD::SiDetectorElement*>(ETS.m_pS->associatedDetectorElement());
1110 if(pDE!=nullptr) {
1112 xOverX0 = pDE->design().thickness()/93.7;//Radiation length of silicon according to PDG
1113 }
1114 else {
1115 if(pDE->isPixel() && std::abs(Trf(2,2)) >= 1.0) xOverX0 = 0.05;//increase for endcap Pixel modules
1116 }
1117 }
1118
1119 double lenCorr = 1/std::fabs(lV);
1120
1121 double radLength = xOverX0*lenCorr;
1122
1123 double qpCorr = Re[4]*(1.0 + 0.038 * std::log(radLength));
1124 double sigmaMS2 = 185.0 * radLength * qpCorr*qpCorr; //Highland formula
1125
1126 //multiple scattering
1127
1128 ETS.m_Gk[2][2] += sigmaMS2/(sint*sint);
1129 ETS.m_Gk[3][3] += sigmaMS2;
1130
1131 //Sym. product J*C*J^T
1132
1133 double Be[5][5];//upper off-diagonal block
1134
1135 memset(&Be[0][0],0,sizeof(Be));
1136
1137 for(int i=0;i<4;i++) {
1138 for(int j=0;j<5;j++) {
1139 for(int k=0;k<5;k++) Be[i][j] += J[k+i*5]*ETS.m_Gk[k][j+5];
1140 }
1141 }
1142 for(int j=0;j<5;j++) {
1143 Be[4][j] = ETS.m_Gk[4][j+5];
1144 }
1145
1146 double JC[5][5];
1147 double Ce[5][5];//"running" diagonal block
1148
1149 memset(&JC[0][0],0,sizeof(JC));
1150 memset(&Ce[0][0],0,sizeof(Ce));
1151
1152 for(int i=0;i<4;i++) {
1153 for(int j=0;j<5;j++) {
1154 for(int k=0;k<5;k++) JC[i][j] += J[k+i*5]*ETS.m_Gk[k][j];
1155 }
1156 }
1157 for(int j=0;j<5;j++) {
1158 JC[4][j] = ETS.m_Gk[4][j];
1159 }
1160
1161 for(int i=0;i<5;i++) {
1162 for(int j=0;j<=i;j++) {
1163 if(j<4) {
1164 for(int k=0;k<5;k++) Ce[i][j] += JC[i][k]*J[k+j*5];
1165 Ce[j][i] = Ce[i][j];
1166 }
1167 else {
1168 Ce[i][4] = Ce[4][i] = JC[i][4];
1169 }
1170 }
1171 }
1172
1173 //copy back to the state
1174
1175 ETS.m_pS = pN;//moving forward
1176
1177 for(int i=0;i<5;i++) {
1178 ETS.m_Xk[i] = Re[i];
1179 for(int j=0;j<5;j++) {
1180 ETS.m_Gk[i][j] = Ce[i][j];
1181 ETS.m_Gk[i][j+5] = Be[i][j];
1182 ETS.m_Gk[j+5][i] = Be[i][j];
1183 }
1184 }
1185 return 0;
1186}
double thickness() const
Method which returns thickness of the silicon wafer.
virtual const SiDetectorDesign & design() const override final
access to the local description (inline):
int RungeKutta34(double *, double *, const Trk::PlaneSurface *, MagField::AtlasFieldCache &, bool) const
Gaudi::Property< double > m_nominalRadLength
Gaudi::Property< bool > m_useDetectorThickness
const TrkDetElementBase * associatedDetectorElement() const
return associated Detector Element
const Amg::Transform3D & transform() const
Returns HepGeom::Transform3D by reference.
int cost(std::vector< std::string > &files, node &n, const std::string &directory="", bool deleteref=false, bool relocate=false)
Definition hcg.cxx:922
Eigen::Affine3d Transform3D
void transformGlobalToLocal(double *ATH_RESTRICT, double *ATH_RESTRICT)

◆ finalize()

StatusCode TrigInDetTrackFollowingTool::finalize ( )
virtual

Definition at line 190 of file TrigInDetTrackFollowingTool.cxx.

190 {
191 return StatusCode::SUCCESS;
192}

◆ findNearestHit() [1/2]

void TrigInDetTrackFollowingTool::findNearestHit ( int moduleIdx,
const InDet::PixelClusterCollection * pColl,
const double * TP,
std::vector< std::tuple< double, const Trk::PrepRawData *, int > > & hitLinks ) const
inlineprivate

Definition at line 195 of file TrigInDetTrackFollowingTool.cxx.

195 {
196
197 const InDet::PixelCluster* bestHit = nullptr;
198
199 double bestDist = 1e8;
200
201 for(const auto pPRD : *pColl) {
202
203 double rx = std::fabs(pPRD->localPosition().x() - TP[0]);
204 if(rx > m_winX_Pixels) continue;
205
206 double ry = std::fabs(pPRD->localPosition().y() - TP[1]);
207 if(ry > m_winY_Pixels) continue;
208
209 double dist = std::pow(5*rx,2) + std::pow(ry,2);//x-residual is given more weight
210
211 if(dist < bestDist) {
212 bestDist = dist;
213 bestHit = pPRD;
214 }
215
216 }
217 if(bestHit != nullptr) {
218 hitLinks.emplace_back(std::make_tuple(bestDist, bestHit, moduleIdx));
219 }
220}
Gaudi::Property< double > m_winX_Pixels
Gaudi::Property< double > m_winY_Pixels

◆ findNearestHit() [2/2]

void TrigInDetTrackFollowingTool::findNearestHit ( int moduleIdx,
const InDet::SCT_ClusterCollection * pColl,
int shape,
const double * TP,
std::vector< std::tuple< double, const Trk::PrepRawData *, int > > & hitLinks ) const
inlineprivate

Definition at line 222 of file TrigInDetTrackFollowingTool.cxx.

222 {
223
224 const InDet::SCT_Cluster* bestHit = nullptr;
225 float bestDist = 1e8;
226
227 for(const auto pPRD : *pColl) {
228
229 double rx = 0.0;
230
231 if(shape == InDetDD::Box) {
232 rx = std::fabs(pPRD->localPosition().x() - TP[0]);
233 }
234 else {
235
236 double meas_x = pPRD->localPosition().x();
237 double meas_y = pPRD->localPosition().y();
238
239 double e00 = pPRD->localCovariance()(0, 0);
240 double e01 = pPRD->localCovariance()(0, 1);
241 double e11 = pPRD->localCovariance()(1, 1);
242
243 double beta = 0.5*std::atan(2*e01/(e00-e11));
244 double sinB, cosB;
245 sincos(beta, &sinB, &cosB);
246
247 rx = std::fabs((meas_x - TP[0])*cosB + (meas_y - TP[1])*sinB);
248 }
249
250 if(rx > m_winX_Strips) continue;
251
252 if(rx < bestDist) {
253 bestHit = pPRD;
254 bestDist = rx;
255 }
256 }
257 if(bestHit != nullptr) {
258 hitLinks.emplace_back(std::make_tuple(bestDist, bestHit, moduleIdx));
259 }
260}
Gaudi::Property< double > m_winX_Strips

◆ fitTheSeed()

std::unique_ptr< TrigFTF_ExtendedTrackState > TrigInDetTrackFollowingTool::fitTheSeed ( const std::vector< const Trk::SpacePoint * > & seed,
MagField::AtlasFieldCache & fieldCache ) const
private

Definition at line 447 of file TrigInDetTrackFollowingTool.cxx.

447 {
448
449 //track parameters are estimated at the last point of the seed (inside-out approach)
450
451 const Trk::SpacePoint& SPb = *seed.at(0);
452 const Trk::SpacePoint& SPm = *seed.at(1);
453 const Trk::SpacePoint& SPe = *seed.at(seed.size()-1);
454
455 double pb[3] = {SPb.globalPosition().x(), SPb.globalPosition().y(), SPb.globalPosition().z()};
456 double pm[3] = {SPm.globalPosition().x(), SPm.globalPosition().y(), SPm.globalPosition().z()};
457 double pe[3] = {SPe.globalPosition().x(), SPe.globalPosition().y(), SPe.globalPosition().z()};
458
459 double dsp[2] = {pe[0]-pb[0], pe[1]-pb[1]};
460 double Lsp = std::sqrt(dsp[0]*dsp[0] + dsp[1]*dsp[1]);
461 double cosA = dsp[0]/Lsp;
462 double sinA = dsp[1]/Lsp;
463 double tau0 = (pe[2]-pb[2])/(SPe.globalPosition().perp()-SPb.globalPosition().perp());
464
465 double dxm = pm[0] - pb[0];
466 double dym = pm[1] - pb[1];
467
468 double x1 = dxm*cosA + dym*sinA;
469 double m1 = -dxm*sinA + dym*cosA;
470
471 double a_parab = -2*m1/(x1*(Lsp-x1));
472 double b_parab = -0.5*a_parab*Lsp;
473
474 double Rx[3] = {0,b_parab,a_parab};
475 double Ry[2] = {pb[2],tau0};
476
477 double Cx[3][3];
478 double Cy[2][2];
479
480 memset(&Cx[0][0],0,sizeof(Cx));
481 Cx[0][0] = 0.1;
482 Cx[1][1] = 0.001;
483 Cx[2][2] = 0.001;
484
485 memset(&Cy[0][0],0,sizeof(Cy));
486 Cy[0][0] = 0.25;
487 Cy[1][1] = 0.001;
488
489 double path = -Lsp;
490 path = 0.0;
491 double radius= SPb.globalPosition().perp();
492
493 double Fx[3][3];
494 double Fy[2][2];
495
496 memset(&Fx[0][0],0,sizeof(Fx));
497 Fx[0][0] = 1.0;
498 Fx[1][1] = 1.0;
499 Fx[2][2] = 1.0;
500
501 memset(&Fy[0][0],0,sizeof(Fy));
502 Fy[0][0] = 1.0;
503 Fy[1][1] = 1.0;
504
505 double sigma_x2 = std::pow(0.08,2);
506 double sigma_y2 = std::pow(0.3,2);
507
508 for(const auto& sp : seed) {
509
510 //1. extrapolate
511
512 double pk[3] = {sp->globalPosition().x(), sp->globalPosition().y(), sp->globalPosition().z()};
513
514 double dx = pk[0] - pb[0];
515 double dy = pk[1] - pb[1];
516
517 double dist = dx*cosA + dy*sinA;
518 double measx =-dx*sinA + dy*cosA;
519 double measy = pk[2];
520
521 double ds = dist - path;
522
523 path = dist;
524
525 double rk = sp->globalPosition().perp();
526 double dr = rk - radius;
527 radius = rk;
528
529 //update Jacobians
530 Fx[0][1] = ds;
531 Fx[0][2] = 0.5*ds*ds;
532 Fx[1][2] = ds;
533
534 Fy[0][1] = dr;
535
536 Cx[1][1] += 1e-7;
537 Cy[1][1] += 1e-7;
538
539 double Rex[3] = {Rx[0], Rx[1], Rx[2]};
540
541 Rex[0] += Fx[0][1]*Rx[1] + Fx[0][2]*Rx[2];
542 Rex[1] += Fx[1][2]*Rx[2];
543
544 double Cex[3][3], CFT[3][3];
545
546 for(int i=0;i<3;i++) {
547 for(int j=0;j<3;j++) {
548 CFT[i][j] = 0;
549 for(int m=0;m<3;m++) CFT[i][j] += Cx[i][m]*Fx[j][m];
550 }
551 }
552 for(int i=0;i<3;i++) {
553 for(int j=0;j<3;j++) {
554 Cex[i][j] = 0;
555 for(int m=0;m<3;m++) Cex[i][j] += Fx[i][m]*CFT[m][j];
556 }
557 }
558
559 double Rey[2] = {Ry[0], Ry[1]};
560
561 Rey[0] += Fy[0][1]*Ry[1];
562
563 double Cey[3][3];
564
565 for(int i=0;i<2;i++) {
566 for(int j=0;j<2;j++) {
567 CFT[i][j] = 0;
568 for(int m=0;m<2;m++) CFT[i][j] += Cy[i][m]*Fy[j][m];
569 }
570 }
571 for(int i=0;i<2;i++) {
572 for(int j=0;j<2;j++) {
573 Cey[i][j] = 0;
574 for(int m=0;m<2;m++) Cey[i][j] += Fy[i][m]*CFT[m][j];
575 }
576 }
577
578 //2. update
579
580 double CHTx[3] = {Cex[0][0], Cex[0][1], Cex[0][2]};
581 double Dx = 1/(Cex[0][0] + sigma_x2);
582
583 double resid = measx - Rex[0];
584
585 double Kx[3] = {Dx*CHTx[0], Dx*CHTx[1], Dx*CHTx[2]};
586
587 for(int i=0;i<3;i++) Rx[i] = Rex[i] + Kx[i]*resid;
588 for(int i=0;i<3;i++) {
589 for(int j=0;j<3;j++) {
590 Cx[i][j] = Cex[i][j] - Kx[i]*CHTx[j];
591 }
592 }
593
594 double CHTy[2] = {Cey[0][0], Cey[0][1]};
595 double Dy = 1/(Cey[0][0] + sigma_y2);
596
597 resid = measy - Rey[0];
598
599 double Ky[2] = {Dy*CHTy[0], Dy*CHTy[1]};
600
601 for(int i=0;i<2;i++) Ry[i] = Rey[i] + Ky[i]*resid;
602 for(int i=0;i<2;i++) {
603 for(int j=0;j<2;j++) {
604 Cy[i][j] = Cey[i][j] - Ky[i]*CHTy[j];
605 }
606 }
607 }
608
609 //create initial track state at the last spacepoint of the seed
610
611 double B0[3];
612
613 fieldCache.getField(pe, B0);
614
615 double P0[6];
616
617 const Trk::PrepRawData* cle = SPe.clusterList().first;
618
619 const Trk::PlaneSurface* thePlane = static_cast<const Trk::PlaneSurface*>(&cle->detectorElement()->surface());
620
621 if(thePlane == nullptr) return nullptr;
622
623 P0[0] = cle->localPosition()[0];
624 P0[1] = cle->localPosition()[1];
625 P0[2] = std::atan2(sinA + Rx[1]*cosA, cosA - Rx[1]*sinA);//phi in the global c.s.
626 P0[3] = std::atan2(1, Ry[1]);//theta in the global c.s.
627 double coeff = 1.0/(300.0*B0[2]*std::sqrt(1+Ry[1]*Ry[1]));
628 P0[4] = -Rx[2]*coeff;//qOverP estimate
629 P0[5] = Cx[2][2]*coeff*coeff;//qOverP covariance
630
631 return std::make_unique<TrigFTF_ExtendedTrackState>(P0, thePlane);
632
633}
static const int B0
Definition AtlasPID.h:122
static Double_t sp
static Double_t tau0
void getField(const double *ATH_RESTRICT xyz, double *ATH_RESTRICT bxyz, double *ATH_RESTRICT deriv=nullptr)
get B field value at given position xyz[3] is in mm, bxyz[3] is in kT if deriv[9] is given,...
virtual const TrkDetElementBase * detectorElement() const =0
return the detector element corresponding to this PRD The pointer will be zero if the det el is not d...
const Amg::Vector2D & localPosition() const
return the local position reference
virtual const Amg::Vector3D & globalPosition() const override final
Interface method to get the global Position.
const std::pair< const PrepRawData *, const PrepRawData * > & clusterList() const
return the pair of cluster pointers by reference
GeoGenfun::FunctionNoop Fy(double r, GeoGenfun::GENFUNCTION G, const double Cenx[], const double Ceny[])
GeoGenfun::FunctionNoop Fx(double r, GeoGenfun::GENFUNCTION G, const double Cenx[], const double Ceny[])
constexpr double coeff(const unsigned l, const unsigned k)
Calculates the n-th coefficient of the legendre polynomial series.
path
python interpreter configuration --------------------------------------—
Definition athena.py:128

◆ getTrack()

Trk::Track * TrigInDetTrackFollowingTool::getTrack ( const std::vector< const Trk::SpacePoint * > & seed,
const std::vector< const InDetDD::SiDetectorElement * > & road,
const EventContext & ctx ) const
virtual

Implements ITrigInDetTrackFollowingTool.

Definition at line 635 of file TrigInDetTrackFollowingTool.cxx.

635 {
636
637 //1. get magnetic field
638
639 MagField::AtlasFieldCache fieldCache;
640
641 SG::ReadCondHandle<AtlasFieldCacheCondObj> fieldCondObj{m_fieldCondObjInputKey, ctx};
642 if (!fieldCondObj.isValid()) {
643 ATH_MSG_ERROR("Failed to retrieve AtlasFieldCacheCondObj with key " << m_fieldCondObjInputKey.key());
644 return nullptr;
645 }
646
647 fieldCondObj->getInitializedCache (fieldCache);
648
649 //2. get hits
650
651 SG::ReadHandle<InDet::PixelClusterContainer> pixcontainer(m_pixcontainerkey, ctx);
652
653 if (!pixcontainer.isValid()) {
654 ATH_MSG_ERROR("Failed to retrieve Pixel Cluster Container with key " << m_pixcontainerkey.key());
655 return nullptr;
656 }
657
658 const InDet::PixelClusterContainer* p_pixcontainer = pixcontainer.ptr();
659
660 SG::ReadHandle<InDet::SCT_ClusterContainer> sctcontainer(m_sctcontainerkey, ctx);
661
662 if (!sctcontainer.isValid()) {
663 ATH_MSG_ERROR("Failed to retrieve Strip Cluster Container with key " << m_sctcontainerkey.key());
664 return nullptr;
665 }
666
667 const InDet::SCT_ClusterContainer* p_sctcontainer = sctcontainer.ptr();
668
669 //3. prepare the initial track state
670
671 int nModules = road.size();
672
673 std::vector<const Trk::PrepRawData*> assignedHits(nModules, nullptr);//assuming maximum one assigned hit per detector element
674
675 std::vector<int> moduleStatus(nModules, 0);//initial status: unchecked
676
677 std::vector<Identifier> seedIdents;
678 std::vector<const Trk::PrepRawData*> seedHits;
679
680 for(const auto& sp : seed) {
681 const Trk::PrepRawData* prd = sp->clusterList().first;
682 seedIdents.push_back(prd->detectorElement()->identify());
683 seedHits.push_back(prd);
684 prd = sp->clusterList().second;
685 if(prd == nullptr) continue;
686 seedIdents.push_back(prd->detectorElement()->identify());//the second cluster of a strip SP
687 seedHits.push_back(prd);
688 }
689
690 //pre-assigning the hits contained in the input track seed
691
692 unsigned int seedSize = seedIdents.size();
693 int nUnassigned = seedSize;
694
695 int startModuleIdx = -1;
696
697 for(int moduleIdx = 0;moduleIdx<nModules;moduleIdx++) {
698
699 Identifier ident = road.at(moduleIdx)->identify();
700
701 for(unsigned int clIdx=0;clIdx<seedSize;clIdx++) {
702
703 if(seedIdents[clIdx] != ident) continue;
704
705 assignedHits[moduleIdx] = seedHits[clIdx];
706 moduleStatus[moduleIdx] = 1;//seed hit assigned
707
708 startModuleIdx = moduleIdx;
709 --nUnassigned;
710 break;
711 }
712
713 if(nUnassigned == 0) break;
714 }
715
716 if(nUnassigned > 0) return nullptr;//bad road or bad seed
717
718 std::unique_ptr<TrigFTF_ExtendedTrackState> initialState = fitTheSeed(seed, fieldCache);
719
720 if(initialState == nullptr) return nullptr;
721
722 initialState->correctAngles();
723
724 TrigFTF_ExtendedTrackState& theState = *initialState;
725
726 //3. create layer sequences for forward and backward passes
727
728 std::vector<int> layerSequence[2]; //the order in which layers will be explored to find hits
729 std::unordered_map<int, std::vector<int> > layerMap[2]; //layer-to-modules map
730
731 const std::vector<short>& vPixelL = *(m_layerNumberTool->pixelLayers());
732 const std::vector<short>& vStripL = *(m_layerNumberTool->sctLayers());
733
734 for(int passIdx=0;passIdx<2;passIdx++) {//0: forward pass, 1: backward pass
735
736 int start = passIdx==0 ? startModuleIdx + 1 : startModuleIdx;
737 int end = passIdx==0 ? nModules : -1;
738 int step = passIdx==0 ? 1 : -1;
739
740 for(int moduleIdx = start;moduleIdx!=end;moduleIdx+=step) {//iterating over modules in the road
741 const InDetDD::SiDetectorElement* de = road.at(moduleIdx);
742 int h = de->identifyHash();
743 int l = de->isPixel() ? vPixelL.at(h) : vStripL.at(h);
744 if(!de->isPixel()) {//Strips
745 l *= 1000;
746 if(h % 2) {//odd hash means the other side of the double Strip layer
747 l += 1;
748 }
749 }
750 auto it = layerMap[passIdx].find(l);
751 if(it != layerMap[passIdx].end()) (*it).second.push_back(moduleIdx);
752 else {
753 layerSequence[passIdx].push_back(l);
754 std::vector<int> v = {moduleIdx};
755 layerMap[passIdx].insert(std::make_pair(l,v));
756 }
757 }
758 }
759
760 //4. The layer-based track following loop with two passes
761
762 for(int passIdx=0;passIdx<2;passIdx++) {
763
764 for(const auto& lkey : layerSequence[passIdx]) {
765
766 if(theState.m_nHoles > m_nHolesMax) { //bailing-out early to save CPU time
767 return nullptr;
768 }
769
770 const auto& lp = (*layerMap[passIdx].find(lkey));
771
772 //4a. collect hits from modules on the layer
773
774 std::vector<std::tuple<double, const Trk::PrepRawData*, int> > hitLinks;
775
776 for(const auto& moduleIdx : lp.second) {
777
778 if(moduleIdx == startModuleIdx) {
779 theState.SwapTheEnds();
780 }
781
782 if(moduleStatus[moduleIdx] < 0) continue;//checked and rejected
783
784 if (assignedHits[moduleIdx] != nullptr) {//we have a pre-assigned hit, so keep it
785 hitLinks.emplace_back(std::make_tuple(-1.0,assignedHits[moduleIdx],moduleIdx));// negative distance forces accept
786 continue;
787 }
788
789 //no pre-assigned hit
790
791 const InDetDD::SiDetectorElement* de = road.at(moduleIdx);// module in the road
792
793 //5a. skip empty, dead, bad, etc. modules
794
795 unsigned int moduleHash = de->identifyHash();
796
797 bool noHits = false;
798
799 if(de->isPixel()) {//Pixel module
800 noHits = (p_pixcontainer == nullptr) || ((*p_pixcontainer).indexFindPtr(moduleHash) == nullptr);
801 } else {//Strip module
802 noHits = (p_sctcontainer == nullptr) || ((*p_sctcontainer).indexFindPtr(moduleHash) == nullptr);
803 }
804
805 if (noHits) {
806 moduleStatus[moduleIdx] = -3;//skip the module
807 continue;
808 }
809
810 //5b. tentative extrapolation to check if the track can cross the module at all
811
812 const Trk::PlaneSurface* plane = static_cast<const Trk::PlaneSurface*>(&de->surface());
813
814 double trackParams[2];
815
816 if(!tentativeExtrapolation(theState.m_Xk, trackParams, theState.m_pS, plane, fieldCache)) {
817 moduleStatus[moduleIdx] = -2;//miss due to extrapolation failure
818 continue;
819 }
820
821 const double bound_tol = 0.2;
822
823 InDetDD::SiIntersect intersection = de->inDetector(Amg::Vector2D(trackParams[0], trackParams[1]), bound_tol, bound_tol);
824
825 if (intersection.out()) {
826 moduleStatus[moduleIdx] = -2;//miss
827 continue;
828 }
829
830 //5c. search for the nearest hit on the module
831
832 moduleStatus[moduleIdx] = -4;//the default: all hits are rejected
833
834 if(de->isPixel()) {
835 findNearestHit(moduleIdx, (*p_pixcontainer).indexFindPtr(moduleHash), trackParams, hitLinks);
836 }
837 else {
838 findNearestHit(moduleIdx, (*p_sctcontainer).indexFindPtr(moduleHash), de->design().shape(), trackParams, hitLinks);
839 }
840 }
841
842 //4b. check what we've found
843
844 if(hitLinks.empty()) {//add a hole and go to the next layer
845 theState.AddHole();
846 continue;
847 }
848
849 //4c. sorting by distance
850
851 std::sort(hitLinks.begin(), hitLinks.end());
852
853 //4d. update the track state using selected hits
854
855 for(const auto& hl : hitLinks) {
856
857 int moduleIdx = get<2>(hl);
858
859 const Trk::PrepRawData* pPRD = get<1>(hl);
860
861 const InDetDD::SiDetectorElement* de = road.at(moduleIdx);
862
863 const Trk::PlaneSurface* plane = static_cast<const Trk::PlaneSurface*>(&de->surface());
864
865 //precise extrapolation with covariance, material effects, etc.
866
867 int rkCode = extrapolateTrackState(theState, plane, fieldCache);
868
869 if(rkCode!=0) {
870 moduleStatus[moduleIdx] = -2;//extrapolation failure
871 continue;
872 }
873
874 const Trk::PrepRawData* acceptedHit = nullptr;
875
876 if(de->isPixel()) {//Pixel module
877 const InDet::PixelCluster* pPixelHit = dynamic_cast<const InDet::PixelCluster*>(pPRD);
878 acceptedHit = updateTrackState(pPixelHit, theState, (get<0>(hl) < 0.0));
879 }
880 else {
881 const InDet::SCT_Cluster* pStripHit = dynamic_cast<const InDet::SCT_Cluster*>(pPRD);
882 acceptedHit = updateTrackState(pStripHit, de->design().shape(), theState);
883 }
884
885 if(acceptedHit == nullptr) {//the corresponding max dchi2 exceeded
886 //break;
887 }
888 else {
889 assignedHits[moduleIdx] = acceptedHit;
890 moduleStatus[moduleIdx] = de->isPixel() ? 2 : 3;//a new Pixel/Strip hit assigned
891 }
892 }//end of the update cycle
893 }//end of layer processing cycle
894 } //end of the track following loop
895
896 //6. create output track
897
898 int nClusters = theState.m_nClusters;
899 int nHoles = theState.m_nHoles;
900 double chi2tot = theState.m_chi2;
901
902 if(nClusters < m_nClustersMin) {
903 return nullptr;
904 }
905
906 if(nHoles > m_nHolesMax) {
907 return nullptr;
908 }
909
910 int rkCode = extrapolateTrackState(theState, nullptr, fieldCache);//to perigee
911
912 if(rkCode !=0) {
913 return nullptr;
914 }
915
916 int ndoftot = theState.m_ndof;
917 double qOverP = theState.m_Xk[4];
918 double pt = std::sin(theState.m_Xk[3])/qOverP;
919 double phi0 = theState.m_Xk[2];
920 double theta = theState.m_Xk[3];
921
922 double z0 = theState.m_Xk[1];
923 double d0 = theState.m_Xk[0];
924
925 bool bad_cov = false;
926
927 auto cov = AmgSymMatrix(5){};
928
929 for(int i=0;i<5;i++) {
930 for(int j=0;j<=i;j++) {
931 double c = theState.m_Gk[i][j];
932 if (i == j && c < 0) {
933 bad_cov = true;//Diagonal elements must be positive
934 ATH_MSG_DEBUG("REGTEST: cov(" << i << "," << i << ") =" << c << " < 0, reject track");
935 break;
936 }
937 cov.fillSymmetric(i,j, c);
938 }
939 }
940
941 if((ndoftot<0) || (fabs(pt)<100.0) || (std::isnan(pt)) || bad_cov) {
942 ATH_MSG_DEBUG("Track following failed - possibly floating point problem");
943 return nullptr;
944 }
945
946 Trk::PerigeeSurface perigeeSurface;
947
948 std::unique_ptr<Trk::TrackParameters> pPP = perigeeSurface.createUniqueTrackParameters(d0, z0, phi0, theta, qOverP, cov);
949
950 std::bitset<Trk::TrackStateOnSurface::NumberOfTrackStateOnSurfaceTypes> perType;
952
953 auto pParVec = std::make_unique<Trk::TrackStates>();
954
955 pParVec->reserve(50);
956 pParVec->push_back(new Trk::TrackStateOnSurface(nullptr, std::move(pPP), nullptr, perType));
957
958 std::bitset<Trk::TrackStateOnSurface::NumberOfTrackStateOnSurfaceTypes> rioType(0);
961
962 for(const auto& ha : theState.m_track) {//loop over hit assignments
963
964 const Trk::PrepRawData* pPRD = ha.m_pPRD;
965
966 if(pPRD == nullptr) continue;//skip holes
967
968 //create track states on surface
969
970 int ndof = ha.m_ndof;
971
972 const Trk::PlaneSurface* pPS = dynamic_cast<const Trk::PlaneSurface*>(&pPRD->detectorElement()->surface());
973
974 if(pPS==nullptr) continue;
975
976 const InDet::SiCluster* pCL = dynamic_cast<const InDet::SiCluster*>(pPRD);
977
978 Trk::LocalParameters locPos = Trk::LocalParameters(pCL->localPosition());
979
980 const Amg::MatrixX& cov = pCL->localCovariance();
981
982 IdentifierHash hash = pPRD->detectorElement()->identifyHash();
983
984 std::unique_ptr<Trk::MeasurementBase> pRIO{};
985
986 if(ndof == 2) {
987 const InDet::PixelCluster* pPixel = static_cast<const InDet::PixelCluster*>(pCL);
988 if(pPixel) {
989 pRIO = std::make_unique<InDet::PixelClusterOnTrack>(pPixel, std::move(locPos), Amg::MatrixX(cov),
990 hash, pPixel->globalPosition(), pPixel->gangedPixel());
991 }
992 }
993 else {
994 const InDet::SCT_Cluster* pStrip = static_cast<const InDet::SCT_Cluster*>(pCL);
995 if(pStrip) {
996 pRIO = std::make_unique<InDet::SCT_ClusterOnTrack>(pStrip, std::move(locPos), Amg::MatrixX(cov),
997 hash, pStrip->globalPosition());
998 }
999 }
1000
1001 auto pM = AmgSymMatrix(5){};
1002
1003 int idx = 0;
1004 for(int i=0;i<5;i++) {
1005 for(int j=0;j<=i;j++) {
1006 pM.fillSymmetric(i,j,ha.m_Ck[idx++]);
1007 }
1008 }
1009
1010 std::unique_ptr<Trk::TrackParameters> pTP = pPS->createUniqueTrackParameters(ha.m_Xk[0], ha.m_Xk[1], ha.m_Xk[2], ha.m_Xk[3], ha.m_Xk[4], pM);
1011
1012 Trk::FitQualityOnSurface FQ = Trk::FitQualityOnSurface(ha.m_dchi2, ndof);
1013
1014 Trk::TrackStateOnSurface* pTSS = new Trk::TrackStateOnSurface(FQ, std::move(pRIO), std::move(pTP), nullptr, rioType);
1015
1016 pParVec->push_back(pTSS);
1017 }
1018
1019 auto pFQ = std::make_unique<Trk::FitQuality>(chi2tot, ndoftot);
1020
1021 Trk::TrackInfo info{};
1022
1023 info.setParticleHypothesis(Trk::pion);
1024 info.setPatternRecognitionInfo(Trk::TrackInfo::strategyC);
1025
1026 Trk::Track* foundTrack = new Trk::Track(info, std::move(pParVec), std::move(pFQ));
1027
1028 return foundTrack;
1029}
Scalar theta() const
theta method
#define ATH_MSG_ERROR(x)
#define ATH_MSG_DEBUG(x)
#define AmgSymMatrix(dim)
virtual DetectorShape shape() const
Shape of element.
virtual IdentifierHash identifyHash() const override final
identifier hash (inline)
SiIntersect inDetector(const Amg::Vector2D &localPosition, double phiTol, double etaTol) const
Test that it is in the active region.
Trk::Surface & surface()
Element Surface.
const Amg::Vector3D & globalPosition() const
return global position reference
bool gangedPixel() const
return the flag of this cluster containing a gangedPixel
Trk::PrepRawDataContainer< SCT_ClusterCollection > SCT_ClusterContainer
Trk::PrepRawDataContainer< PixelClusterCollection > PixelClusterContainer
int extrapolateTrackState(TrigFTF_ExtendedTrackState &, const Trk::PlaneSurface *, MagField::AtlasFieldCache &) const
std::unique_ptr< TrigFTF_ExtendedTrackState > fitTheSeed(const std::vector< const Trk::SpacePoint * > &, MagField::AtlasFieldCache &) const
void findNearestHit(int, const InDet::PixelClusterCollection *, const double *, std::vector< std::tuple< double, const Trk::PrepRawData *, int > > &) const
SG::ReadHandleKey< InDet::SCT_ClusterContainer > m_sctcontainerkey
SG::ReadHandleKey< InDet::PixelClusterContainer > m_pixcontainerkey
const Trk::PrepRawData * updateTrackState(const InDet::PixelCluster *, TrigFTF_ExtendedTrackState &, bool) const
ToolHandle< ITrigL2LayerNumberTool > m_layerNumberTool
SG::ReadCondHandleKey< AtlasFieldCacheCondObj > m_fieldCondObjInputKey
bool tentativeExtrapolation(double const *, double *, const Trk::PlaneSurface *, const Trk::PlaneSurface *, MagField::AtlasFieldCache &) const
virtual Surface::ChargedTrackParametersUniquePtr createUniqueTrackParameters(double l1, double l2, double phi, double theta, double qop, std::optional< AmgSymMatrix(5)> cov=std::nullopt) const override final
Use the Surface as a ParametersBase constructor, from local parameters - charged.
virtual Surface::ChargedTrackParametersUniquePtr createUniqueTrackParameters(double l1, double l2, double phi, double theta, double qop, std::optional< AmgSymMatrix(5)> cov=std::nullopt) const override final
Use the Surface as a ParametersBase constructor, from local parameters - charged.
const Amg::MatrixX & localCovariance() const
return const ref to the error matrix
@ Measurement
This is a measurement, and will at least contain a Trk::MeasurementBase.
@ Perigee
This represents a perigee, and so will contain a Perigee object only.
@ Scatterer
This represents a scattering point on the track, and so will contain TrackParameters and MaterialEffe...
virtual IdentifierHash identifyHash() const =0
Identifier hash.
virtual Identifier identify() const =0
Identifier.
std::vector< std::string > intersection(std::vector< std::string > &v1, std::vector< std::string > &v2)
T * get(TKey *tobj)
get a TObject* from a TKey* (why can't a TObject be a TKey?)
Definition hcg.cxx:130
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > MatrixX
Dynamic Matrix - dynamic allocation.
Eigen::Matrix< double, 2, 1 > Vector2D
float ndof(const U &p)
@ qOverP
perigee
l
Printing final latex table to .tex output file.
@ ident
Definition HitInfo.h:77
void sort(typename DataModel_detail::iterator< DVL > beg, typename DataModel_detail::iterator< DVL > end)
Specialization of sort for DataVector/List.
std::list< TrigFTF_HitAssignment > m_track

◆ initialize()

StatusCode TrigInDetTrackFollowingTool::initialize ( )
virtual

Definition at line 173 of file TrigInDetTrackFollowingTool.cxx.

173 {
174
175 StatusCode sc = m_layerNumberTool.retrieve();
176 if(sc.isFailure()) {
177 ATH_MSG_ERROR("Could not retrieve "<<m_layerNumberTool);
178 return sc;
179 } else {
180 ATH_MSG_INFO("Detector layer structure has "<<m_layerNumberTool->maxNumberOfUniqueLayers()<<" unique layers");
181 }
182
183 ATH_CHECK( m_fieldCondObjInputKey.initialize());
186
187 return StatusCode::SUCCESS;
188}
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_INFO(x)
static Double_t sc
::StatusCode StatusCode
StatusCode definition for legacy code.

◆ inputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< AlgTool > >::inputHandles ( ) const
overridevirtualinherited

Return this algorithm's input handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ interfaceID()

const InterfaceID & ITrigInDetTrackFollowingTool::interfaceID ( )
inlinestaticinherited

Definition at line 26 of file ITrigInDetTrackFollowingTool.h.

26 {
28 }
static const InterfaceID IID_ITrigInDetTrackFollowingTool("ITrigInDetTrackFollowingTool", 1, 0)

◆ msg()

MsgStream & AthCommonMsg< AlgTool >::msg ( ) const
inlineinherited

Definition at line 24 of file AthCommonMsg.h.

24 {
25 return this->msgStream();
26 }

◆ msgLvl()

bool AthCommonMsg< AlgTool >::msgLvl ( const MSG::Level lvl) const
inlineinherited

Definition at line 30 of file AthCommonMsg.h.

30 {
31 return this->msgLevel(lvl);
32 }

◆ outputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< AlgTool > >::outputHandles ( ) const
overridevirtualinherited

Return this algorithm's output handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ processHit() [1/2]

double TrigInDetTrackFollowingTool::processHit ( const InDet::PixelCluster * pPRD,
double * resid,
double * invcov,
const TrigFTF_ExtendedTrackState & ets ) const
inlineprivate

Definition at line 262 of file TrigInDetTrackFollowingTool.cxx.

262 {
263
264 double cluster_cov[2];
265
266 if(m_useHitErrors) {
267 cluster_cov[0] = pPRD->localCovariance()(0,0);
268 cluster_cov[1] = pPRD->localCovariance()(1,1);
269 }
270 else {
271 cluster_cov[0] = pPRD->width().phiR();
272 cluster_cov[1] = pPRD->width().z();
273 for(int i=0;i<2;i++) cluster_cov[i] *= cluster_cov[i]/12.0;
274 }
275
276 double covsum[2][2] = {{ets.m_Gk[0][0] + cluster_cov[0], ets.m_Gk[0][1]}, {ets.m_Gk[0][1], ets.m_Gk[1][1] + cluster_cov[1]}};
277
278 double detr = 1/(covsum[0][0]*covsum[1][1] - covsum[0][1]*covsum[1][0]);
279
280 invcov[0] = detr*covsum[1][1];
281 invcov[1] = -detr*covsum[1][0];
282 invcov[2] = detr*covsum[0][0];
283
284 resid[0] = pPRD->localPosition().x() - ets.m_Xk[0];
285 resid[1] = pPRD->localPosition().y() - ets.m_Xk[1];
286
287 return (resid[0]*(resid[0]*invcov[0] + resid[1]*invcov[1]) + resid[1]*(resid[0]*invcov[1] + resid[1]*invcov[2]));
288
289}
const InDet::SiWidth & width() const
return width class reference
double z() const
Definition SiWidth.h:131
double phiR() const
Definition SiWidth.h:126

◆ processHit() [2/2]

double TrigInDetTrackFollowingTool::processHit ( const InDet::SCT_Cluster * pPRD,
int shape,
double & resid,
double & invcov,
double * H,
const TrigFTF_ExtendedTrackState & ets ) const
inlineprivate

Definition at line 291 of file TrigInDetTrackFollowingTool.cxx.

291 {
292
293 if(shape==InDetDD::Box) {
294
295 resid = pPRD->localPosition().x() - ets.m_Xk[0];
296
297 double covX = pPRD->localCovariance()(0, 0);
298
299 if(!m_useHitErrors) {
300 covX = pPRD->width().phiR();
301 covX *= (covX/12);
302 }
303
304 invcov = 1/(ets.m_Gk[0][0] + covX);
305 H[0] = 1;
306 H[1] = 0;
307 }
308
309 else {
310
311 double meas_x = pPRD->localPosition().x();
312 double meas_y = pPRD->localPosition().y();
313
314 double e00 = pPRD->localCovariance()(0, 0);
315 double e01 = pPRD->localCovariance()(0, 1);
316 double e11 = pPRD->localCovariance()(1, 1);
317
318 double beta = 0.5*std::atan(2*e01/(e00-e11));
319 double sinB, cosB;
320 sincos(beta, &sinB, &cosB);
321
322 resid = (meas_x - ets.m_Xk[0])*cosB + (meas_y - ets.m_Xk[1])*sinB;
323
324 H[0] = cosB;
325 H[1] = sinB;
326
327 double track_cov = cosB*cosB*(ets.m_Gk[0][0]+e00) + 2*cosB*sinB*(ets.m_Gk[0][1]+e01) + sinB*sinB*(ets.m_Gk[1][1]+e11);
328
329 invcov = 1/track_cov;
330 }
331
332 return (resid * resid * invcov);
333
334}
#define H(x, y, z)
Definition MD5.cxx:114

◆ renounce()

std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > AthCommonDataStore< AthCommonMsg< AlgTool > >::renounce ( T & h)
inlineprotectedinherited

Definition at line 380 of file AthCommonDataStore.h.

381 {
382 h.renounce();
384 }
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce(T &h)

◆ renounceArray()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::renounceArray ( SG::VarHandleKeyArray & handlesArray)
inlineprotectedinherited

remove all handles from I/O resolution

Definition at line 364 of file AthCommonDataStore.h.

364 {
366 }

◆ RungeKutta34()

int TrigInDetTrackFollowingTool::RungeKutta34 ( double * P,
double * J,
const Trk::PlaneSurface * pN,
MagField::AtlasFieldCache & fieldCache,
bool withJacobian ) const
private

Definition at line 1384 of file TrigInDetTrackFollowingTool.cxx.

1384 {
1385
1386 const double coeff = 299.7;
1387 const double min_step = 3.0;
1388 const double const_field_step = 30.0;
1389 const double maxPath = 3000.0;
1390 const double minQp = 0.01; //100 MeV cut
1391 const double minRad = 300.0;
1392
1393 const int maxIter = 10;
1394
1395 double exStep = 0.0;
1396
1397 if(std::fabs(P[6]) > minQp) return -1;
1398
1399 double Step = estimateRK_Step(pN, P);
1400
1401 if(Step > 1e7) return -2;
1402
1403 double absStep = fabs(Step);
1404
1405 if(absStep <= min_step) {
1406 for(int i=0;i<3;i++) P[i] += Step*P[i+3];
1407 P[7] += Step;
1408 return 0;
1409 }
1410
1411 if(fabs(P[6]*Step) > minRad) {
1412 Step = (Step > 0.0 ? minRad : -minRad)/fabs(P[6]);
1413 }
1414
1415 int nFlips = 0;
1416
1417 double mom = P[6];
1418
1419 double Y[6];
1420
1421 memcpy(&Y[0], P, sizeof(Y));
1422
1423 double B[3];
1424
1425 fieldCache.getField(Y, B);
1426
1427 for(int i=0;i<3;i++) B[i] *= coeff;
1428
1429 if(exStep != 0) {
1430 if(absStep > fabs(exStep))
1431 Step = exStep;
1432 }
1433
1434 for(int iter=0;iter<maxIter;iter++) {//solving single-value boundary problem
1435
1436 bool useConstField = fabs(Step) < const_field_step;
1437
1438 if(!useConstField) {
1439 fieldCache.getField(Y, B);
1440 for(int i=0;i<3;i++) B[i] *= coeff;
1441 }
1442
1443 double B2[3], B3[3];
1444
1445 double H = Step;
1446 double H3 = H/3;
1447 double H23 = 2*H3;
1448 double H4 = 0.25*H;
1449 double H34 = 3*H4;
1450
1451 double H3mom = mom*H3;
1452 double H23mom = mom*H23;
1453 double H4mom = mom*H4;
1454 double H34mom = mom*H34;
1455
1456 double YB[3];
1457
1458 crossProduct(B, Y+3, YB);
1459
1460 //second point
1461
1462 double Y2[6];
1463
1464 for(int i=0;i<3;i++) Y2[i] = Y[i] + H3*Y[i+3];
1465 for(int i=3;i<6;i++) Y2[i] = Y[i] + H3mom*YB[i-3];
1466
1467 double YB2[3];
1468
1469 if(useConstField) {
1470 crossProduct(B, Y2+3, YB2);
1471 }
1472 else {
1473 fieldCache.getField(Y2, B2);
1474 for(int i=0;i<3;i++) B2[i] *= coeff;
1475 crossProduct(B2, Y2+3, YB2);
1476 }
1477
1478 //last point
1479
1480 double Y3[6];
1481
1482 for(int i=0;i<3;i++) Y3[i] = Y[i] + H23*Y2[i+3];
1483 for(int i=3;i<6;i++) Y3[i] = Y[i] + H23mom*YB2[i-3];
1484
1485 double YB3[3];
1486
1487 if(useConstField) {
1488 crossProduct(B, Y3+3, YB3);
1489 }
1490 else {
1491 fieldCache.getField(Y3, B3);
1492 for(int i=0;i<3;i++) B3[i] *= coeff;
1493 crossProduct(B3, Y3+3, YB3);
1494 }
1495
1496 double Y1[6];
1497
1498 for(int i=3;i<6;i++) Y1[i-3] = Y[i-3] + H4*(Y[i] + 3*Y3[i]);
1499 for(int i=0;i<3;i++) Y1[i+3] = Y[i+3] + H4mom*(YB[i] + 3*YB3[i]);
1500
1501 if(fabs(Y1[5])>1) return -10;
1502
1503 //Jacobian calculations go here ...
1504
1505 if(withJacobian) {
1506
1507 double J1C[9], L2C[9], J2C[9], L3C[9], J3C[9];
1508
1509 double CqB3H34[3];
1510 double CqB2H23[3];
1511 double CqBH3[3];
1512
1513 if(!useConstField) {
1514 for(int i=0;i<3;i++) CqBH3[i] = H3mom*B[i];
1515 for(int i=0;i<3;i++) CqB2H23[i] = H23mom*B2[i];
1516 for(int i=0;i<3;i++) CqB3H34[i] = H34mom*B3[i];
1517 }
1518 else {
1519 for(int i=0;i<3;i++) CqBH3[i] = H3mom*B[i];
1520 for(int i=0;i<3;i++) CqB2H23[i] = H23mom*B[i];
1521 for(int i=0;i<3;i++) CqB3H34[i] = H34mom*B[i];
1522 }
1523
1524 crossProduct(CqBH3, J+17, J1C);
1525 crossProduct(CqBH3, J+24, J1C+3);
1526 crossProduct(CqBH3, J+31, J1C+6);
1527
1528 J1C[6] += H3*YB[0];
1529 J1C[7] += H3*YB[1];
1530 J1C[8] += H3*YB[2];
1531
1532 L2C[0] = J[17] + J1C[0];
1533 L2C[1] = J[18] + J1C[1];
1534 L2C[2] = J[19] + J1C[2];
1535
1536 L2C[3] = J[24] + J1C[3];
1537 L2C[4] = J[25] + J1C[4];
1538 L2C[5] = J[26] + J1C[5];
1539
1540 L2C[6] = J[31] + J1C[6];
1541 L2C[7] = J[32] + J1C[7];
1542 L2C[8] = J[33] + J1C[8];
1543
1544 crossProduct(CqB2H23, L2C, J2C);
1545 crossProduct(CqB2H23, L2C+3, J2C+3);
1546 crossProduct(CqB2H23, L2C+6, J2C+6);
1547
1548 J2C[6] += H23*YB2[0];
1549 J2C[7] += H23*YB2[1];
1550 J2C[8] += H23*YB2[2];
1551
1552 L3C[0] = J[17] + J2C[0];
1553 L3C[1] = J[18] + J2C[1];
1554 L3C[2] = J[19] + J2C[2];
1555
1556 L3C[3] = J[24] + J2C[3];
1557 L3C[4] = J[25] + J2C[4];
1558 L3C[5] = J[26] + J2C[5];
1559
1560 L3C[6] = J[31] + J2C[6];
1561 L3C[7] = J[32] + J2C[7];
1562 L3C[8] = J[33] + J2C[8];
1563
1564 crossProduct(CqB3H34, L3C, J3C);
1565 crossProduct(CqB3H34, L3C+3, J3C+3);
1566 crossProduct(CqB3H34, L3C+6, J3C+6);
1567
1568 J3C[6] += H34*YB3[0];
1569 J3C[7] += H34*YB3[1];
1570 J3C[8] += H34*YB3[2];
1571
1572 for(int i=0;i<9;i++) J1C[i] = 0.75*J1C[i] + J3C[i];
1573
1574 for(int i=0;i<9;i++) J2C[i] *= H34;
1575
1576 J[14] += H*J[17];
1577 J[15] += H*J[18];
1578 J[16] += H*J[19];
1579
1580 J[21] += H*J[24];
1581 J[22] += H*J[25];
1582 J[23] += H*J[26];
1583
1584 J[28] += H*J[31];
1585 J[29] += H*J[32];
1586 J[30] += H*J[33];
1587
1588 J[14] += J2C[0];
1589 J[15] += J2C[1];
1590 J[16] += J2C[2];
1591
1592 J[21] += J2C[3];
1593 J[22] += J2C[4];
1594 J[23] += J2C[5];
1595
1596 J[28] += J2C[6];
1597 J[29] += J2C[7];
1598 J[30] += J2C[8];
1599
1600 J[17] += J1C[0];
1601 J[18] += J1C[1];
1602 J[19] += J1C[2];
1603
1604 J[24] += J1C[3];
1605 J[25] += J1C[4];
1606 J[26] += J1C[5];
1607
1608 J[31] += J1C[6];
1609 J[32] += J1C[7];
1610 J[33] += J1C[8];
1611
1612 }
1613
1614 P[7] += Step;
1615
1616 if(fabs(P[7]) > maxPath) return -3;
1617
1618 double norm = 1/std::sqrt(Y1[3]*Y1[3]+Y1[4]*Y1[4]+Y1[5]*Y1[5]);
1619
1620 Y1[3] *= norm; Y1[4] *= norm; Y1[5] *= norm;
1621
1622 double newStep = estimateRK_Step(pN, Y1);
1623
1624 if(newStep > 1e7) return -4;
1625
1626 double absNewStep = fabs(newStep);
1627
1628 if(absNewStep <= min_step) {//the boundary is too close, using straight line
1629
1630 if(withJacobian) {
1631 if(!useConstField) {
1632 crossProduct(B3, Y1+3, J+35);
1633 }
1634 else {
1635 crossProduct(B, Y1+3, J+35);
1636 }
1637 }
1638
1639 for(int i=0;i<3;i++) {
1640 P[i+3] = Y1[i+3];
1641 P[i] = Y1[i] + newStep*Y1[i+3];
1642 }
1643 P[7] += newStep;
1644
1645 return 0;
1646 }
1647
1648 double absStep = fabs(Step);
1649
1650 if(Step*newStep < 0.0) {//the boundary is overshot
1651 if(++nFlips > 2) return -5;//oscillations
1652 Step = absNewStep < absStep ? newStep : -Step;
1653 }
1654 else {
1655 if(absNewStep < absStep) Step = newStep;
1656 }
1657
1658 for(int i=0;i<6;i++) Y[i] = Y1[i];
1659 }
1660
1661 return -11;//max. number of iteration reached
1662
1663}
double estimateRK_Step(const Trk::PlaneSurface *, double const *) const
void crossProduct(double const *, double const *, double *) const

◆ sysInitialize()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysInitialize ( )
overridevirtualinherited

Perform system initialization for an algorithm.

We override this to declare all the elements of handle key arrays at the end of initialization. See comments on updateVHKA.

Reimplemented in asg::AsgMetadataTool, AthCheckedComponent< AthAlgTool >, AthCheckedComponent<::AthAlgTool >, and DerivationFramework::CfAthAlgTool.

◆ sysStart()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysStart ( )
overridevirtualinherited

Handle START transition.

We override this in order to make sure that conditions handle keys can cache a pointer to the conditions container.

◆ tentativeExtrapolation()

bool TrigInDetTrackFollowingTool::tentativeExtrapolation ( double const * Rk,
double * TP,
const Trk::PlaneSurface * pS,
const Trk::PlaneSurface * pN,
MagField::AtlasFieldCache & fieldCache ) const
private

Definition at line 1188 of file TrigInDetTrackFollowingTool.cxx.

1188 {
1189
1190 const double C = 299.9975;
1191 const double minStep = 100.0;
1192
1193 double Re[5];
1194
1195 memcpy(&Re[0], &Rk[0], sizeof(Re));
1196
1197 const Amg::Transform3D& Trf = pS->transform();
1198
1199 double Ax[3] = {Trf(0,0),Trf(1,0),Trf(2,0)}; //loc1-axis in the global frame
1200
1201 double Ay[3] = {Trf(0,1),Trf(1,1),Trf(2,1)}; //loc2-axis in the global frame
1202
1203 double gP[3];
1204
1205 gP[0] = Trf(0,3) + Ax[0]*Re[0] + Ay[0]*Re[1];
1206 gP[1] = Trf(1,3) + Ax[1]*Re[0] + Ay[1]*Re[1];
1207 gP[2] = Trf(2,3) + Ax[2]*Re[0] + Ay[2]*Re[1];
1208
1209 double sinf, cosf;
1210 double sint, cost;
1211
1212 sincos(Re[2], &sinf, &cosf);
1213 sincos(Re[3], &sint, &cost);
1214
1215 double gV[3] = {cosf*sint, sinf*sint, cost};
1216
1217 const Amg::Vector3D& normal = pN->normal();
1218 const Amg::Vector3D& center = pN->center();
1219
1220 double D[4] = {normal[0], normal[1], normal[2], 0.0};
1221
1222 for(int i=0;i<3;i++) D[3] += -normal[i]*center[i];
1223
1224 double CQ = C*Re[4];
1225
1226 double gB[3];
1227
1228 fieldCache.getField(gP, gB);
1229
1230 double c = D[0]*gP[0] + D[1]*gP[1] + D[2]*gP[2] + D[3];
1231 double b = D[0]*gV[0] + D[1]*gV[1] + D[2]*gV[2];
1232 double a = 0.5*CQ*(gB[0]*(D[1]*gV[2]-D[2]*gV[1]) + gB[1]*(D[2]*gV[0]-D[0]*gV[2]) + gB[2]*(D[0]*gV[1]-D[1]*gV[0]));
1233
1234 double ratio = 4*a*c/(b*b);
1235
1236 bool useExpansion = std::fabs(ratio)<0.1;
1237
1238 double sl = 0.0;
1239
1240 if(useExpansion) {
1241 sl = -c/b;
1242 sl = sl*(1-a*sl/b);
1243 }
1244 else {
1245
1246 double discr = b*b-4.0*a*c;
1247
1248 if(discr<0.0) {
1249 return false;
1250 }
1251
1252 int signb = (b<0.0)?-1:1;
1253 sl = (-b + signb*std::sqrt(discr))/(2*a);
1254 }
1255
1256 int nStepMax = 1;
1257
1258 if(std::fabs(sl)>=minStep) {
1259 nStepMax = (int)(std::fabs(sl)/minStep)+1;
1260 }
1261
1262 if((nStepMax<0)||(nStepMax>100)) {
1263 return false;
1264 }
1265
1266 double Av = sl*CQ;
1267 double Ac = 0.5*sl*Av;
1268 double DVx = gV[1]*gB[2] - gV[2]*gB[1];
1269 double DVy = gV[2]*gB[0] - gV[0]*gB[2];
1270 double DVz = gV[0]*gB[1] - gV[1]*gB[0];
1271
1272 double E[3] = {gP[0]+gV[0]*sl+Ac*DVx, gP[1]+gV[1]*sl+Ac*DVy, gP[2]+gV[2]*sl+Ac*DVz};
1273
1274 if(nStepMax == 1) {
1275 for(int i=0;i<3;i++) gP[i] = E[i];
1276 }
1277 else {
1278 double gBe[3];
1279
1280 fieldCache.getField(E, gBe);
1281
1282 double inv_step = 1/sl;
1283
1284 double dBds[3] = {inv_step*(gBe[0]-gB[0]),inv_step*(gBe[1]-gB[1]),inv_step*(gBe[2]-gB[2])};
1285
1286 int nStep = nStepMax;
1287
1288 while(nStep > 0) {
1289
1290 c = D[0]*gP[0] + D[1]*gP[1] + D[2]*gP[2] + D[3];
1291 b = D[0]*gV[0] + D[1]*gV[1] + D[2]*gV[2];
1292 a = 0.5*CQ*(gB[0]*(D[1]*gV[2]-D[2]*gV[1])+gB[1]*(D[2]*gV[0]-D[0]*gV[2])+gB[2]*(D[0]*gV[1]-D[1]*gV[0]));
1293
1294 ratio = 4*a*c/(b*b);
1295 useExpansion = std::fabs(ratio) < 0.1;
1296
1297 if(useExpansion) {
1298 sl = -c/b;
1299 sl = sl*(1-a*sl/b);
1300 }
1301 else {
1302 double discr=b*b-4.0*a*c;
1303 if(discr<0.0) {
1304 return false;
1305 }
1306 int signb = (b<0.0)?-1:1;
1307 sl = (-b+signb*std::sqrt(discr))/(2*a);
1308 }
1309
1310 double ds = sl/nStep;
1311 Av = ds*CQ;
1312 Ac = 0.5*ds*Av;
1313
1314 DVx = gV[1]*gB[2] - gV[2]*gB[1];
1315 DVy = gV[2]*gB[0] - gV[0]*gB[2];
1316 DVz = gV[0]*gB[1] - gV[1]*gB[0];
1317
1318 E[0] = gP[0] + gV[0]*ds + Ac*DVx;
1319 E[1] = gP[1] + gV[1]*ds + Ac*DVy;
1320 E[2] = gP[2] + gV[2]*ds + Ac*DVz;
1321
1322 double V[3];
1323
1324 V[0] = gV[0] + Av*DVx;
1325 V[1] = gV[1] + Av*DVy;
1326 V[2] = gV[2] + Av*DVz;
1327
1328 for(int i=0;i<3;i++) {
1329 gV[i] = V[i];gP[i] = E[i];
1330 }
1331
1332 for(int i=0;i<3;i++) gB[i] += dBds[i]*ds;//field interpolation
1333
1334 nStep--;
1335 }
1336 }
1337
1338 //global to local
1339
1340 const Amg::Transform3D& Trf2 = pN->transform();
1341
1342 const double* Ax2 = Trf2.matrix().col(0).data();
1343 const double* Ay2 = Trf2.matrix().col(1).data();
1344
1345 const double d[3] = { gP[0] - Trf2(0, 3), gP[1] - Trf2(1, 3), gP[2] - Trf2(2, 3) };
1346
1347 TP[0] = d[0] * Ax2[0] + d[1] * Ax2[1] + d[2] * Ax2[2];
1348 TP[1] = d[0] * Ay2[0] + d[1] * Ay2[1] + d[2] * Ay2[2];
1349
1350 return true;
1351}
struct color C

◆ updateTrackState() [1/2]

const Trk::PrepRawData * TrigInDetTrackFollowingTool::updateTrackState ( const InDet::PixelCluster * pInputHit,
TrigFTF_ExtendedTrackState & ets,
bool forceAccept = false ) const
private

Definition at line 337 of file TrigInDetTrackFollowingTool.cxx.

337 {
338
339 double resid[2];
340 double invcov[3];
341
342 double dchi2 = processHit(pInputHit, resid, invcov, ets);
343
344 if(!forceAccept) {
345 if(dchi2 > m_maxChi2Dist_Pixels) return nullptr;//hit rejected
346 }
347
348 double CHT[10][2];
349
350 for(int i=0;i<10;i++) {
351 CHT[i][0] = ets.m_Gk[i][0];
352 CHT[i][1] = ets.m_Gk[i][1];
353 }
354
355 double Gain[10][2];
356 double V[2][2] = {{invcov[0], invcov[1]}, {invcov[1], invcov[2]}};
357
358 for(int i=0;i<10;i++)
359 for(int j=0;j<2;j++) Gain[i][j] = CHT[i][0]*V[0][j] + CHT[i][1]*V[1][j];
360
361 for(int i=0;i<10;i++) {
362 ets.m_Xk[i] += Gain[i][0]*resid[0] + Gain[i][1]*resid[1];
363
364 for(int j=0;j<=i;j++) {
365 ets.m_Gk[i][j] = ets.m_Gk[i][j] - (Gain[i][0]*CHT[j][0] + Gain[i][1]*CHT[j][1]);
366 ets.m_Gk[j][i] = ets.m_Gk[i][j];
367 }
368 }
369
370 ets.AddHit(pInputHit, dchi2, 2);
371
372 return pInputHit;
373}
double processHit(const InDet::PixelCluster *, double *, double *, const TrigFTF_ExtendedTrackState &) const
Gaudi::Property< double > m_maxChi2Dist_Pixels
void AddHit(const Trk::PrepRawData *, double, int)

◆ updateTrackState() [2/2]

const Trk::PrepRawData * TrigInDetTrackFollowingTool::updateTrackState ( const InDet::SCT_Cluster * pInputHit,
int shape,
TrigFTF_ExtendedTrackState & ets ) const
private

Definition at line 376 of file TrigInDetTrackFollowingTool.cxx.

376 {
377
378 double resid, invcov;
379 double H[2];//linearized measurement matrix
380
381 double dchi2 = 0.0;
382
383 if(shape == InDetDD::Box) {
384
385 dchi2 = processHit(pInputHit, shape, resid, invcov, H, ets);
386
387 if(dchi2 > m_maxChi2Dist_Strips) return nullptr;
388
389 double CHT[10], Gain[10];
390
391 for(int i=0;i<10;i++) {
392 CHT[i] = ets.m_Gk[i][0];
393 Gain[i] = CHT[i] * invcov;
394 }
395
396 for(int i=0;i<10;i++) {
397 ets.m_Xk[i] += Gain[i]*resid;
398 for(int j=0;j<=i;j++) {
399 ets.m_Gk[i][j] = ets.m_Gk[i][j] - Gain[i]*CHT[j];
400 ets.m_Gk[j][i] = ets.m_Gk[i][j];
401 }
402 }
403
404 //boundary check
405
406 double covY = pInputHit->localCovariance()(1, 1);
407 double stripCentre = pInputHit->localPosition().y();
408 double stripHalfLength = std::sqrt(3*covY);
409
410 double dY = ets.m_Xk[1] - stripCentre;
411
412 if(dY > stripHalfLength) {
413 ets.m_Xk[1] = stripCentre + stripHalfLength;
414 }
415
416 if(dY < -stripHalfLength) {
417 ets.m_Xk[1] = stripCentre - stripHalfLength;
418 }
419 }
420 else { //Annulus
421
422 dchi2 = processHit(pInputHit, shape, resid, invcov, H, ets);
423
424 if(dchi2 > m_maxChi2Dist_Strips) return nullptr;
425
426 double CHT[10], Gain[10];
427
428 for(int i=0;i<10;i++) {
429 CHT[i] = H[0]*ets.m_Gk[i][0] + H[1]*ets.m_Gk[i][1];
430 Gain[i] = CHT[i] * invcov;
431 }
432
433 for(int i=0;i<10;i++) {
434 ets.m_Xk[i] += Gain[i]*resid;
435 for(int j=0;j<=i;j++) {
436 ets.m_Gk[i][j] = ets.m_Gk[i][j] - Gain[i]*CHT[j];
437 ets.m_Gk[j][i] = ets.m_Gk[i][j];
438 }
439 }
440 }
441
442 ets.AddHit(pInputHit, dchi2, 1);
443
444 return pInputHit;
445}
Gaudi::Property< double > m_maxChi2Dist_Strips

◆ updateVHKA()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::updateVHKA ( Gaudi::Details::PropertyBase & )
inlineinherited

Definition at line 308 of file AthCommonDataStore.h.

308 {
309 // debug() << "updateVHKA for property " << p.name() << " " << p.toString()
310 // << " size: " << m_vhka.size() << endmsg;
311 for (auto &a : m_vhka) {
313 for (auto k : keys) {
314 k->setOwner(this);
315 }
316 }
317 }
std::vector< SG::VarHandleKeyArray * > m_vhka

Member Data Documentation

◆ m_detStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_detStore
privateinherited

Pointer to StoreGate (detector store by default)

Definition at line 393 of file AthCommonDataStore.h.

◆ m_evtStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_evtStore
privateinherited

Pointer to StoreGate (event store by default)

Definition at line 390 of file AthCommonDataStore.h.

◆ m_fieldCondObjInputKey

SG::ReadCondHandleKey<AtlasFieldCacheCondObj> TrigInDetTrackFollowingTool::m_fieldCondObjInputKey {this, "AtlasFieldCacheCondObj", "fieldCondObj", "Name of the Magnetic Field conditions object key"}
private

Definition at line 130 of file TrigInDetTrackFollowingTool.h.

130{this, "AtlasFieldCacheCondObj", "fieldCondObj", "Name of the Magnetic Field conditions object key"};

◆ m_layerNumberTool

ToolHandle<ITrigL2LayerNumberTool> TrigInDetTrackFollowingTool::m_layerNumberTool {this, "LayerNumberTool", "TrigL2LayerNumberToolITk"}
private

Definition at line 94 of file TrigInDetTrackFollowingTool.h.

94{this, "LayerNumberTool", "TrigL2LayerNumberToolITk"};

◆ m_maxChi2Dist_Pixels

Gaudi::Property<double> TrigInDetTrackFollowingTool::m_maxChi2Dist_Pixels {this, "Chi2MaxPixels", 25.0, "the Pixel hit chi2 cut"}
private

Definition at line 119 of file TrigInDetTrackFollowingTool.h.

119{this, "Chi2MaxPixels", 25.0, "the Pixel hit chi2 cut"};

◆ m_maxChi2Dist_Strips

Gaudi::Property<double> TrigInDetTrackFollowingTool::m_maxChi2Dist_Strips {this, "Chi2MaxStrips", 12.0, "the Strip hit chi2 cut"}
private

Definition at line 120 of file TrigInDetTrackFollowingTool.h.

120{this, "Chi2MaxStrips", 12.0, "the Strip hit chi2 cut"};

◆ m_nClustersMin

Gaudi::Property<int> TrigInDetTrackFollowingTool::m_nClustersMin {this, "nClustersMin", 7, "Minimum number of clusters on track"}
private

Definition at line 117 of file TrigInDetTrackFollowingTool.h.

117{this, "nClustersMin", 7, "Minimum number of clusters on track"};

◆ m_nHolesMax

Gaudi::Property<int> TrigInDetTrackFollowingTool::m_nHolesMax {this, "nHolesMax", 100, "Maximum number of holes on track"}
private

Definition at line 118 of file TrigInDetTrackFollowingTool.h.

118{this, "nHolesMax", 100, "Maximum number of holes on track"};

◆ m_nominalRadLength

Gaudi::Property<double> TrigInDetTrackFollowingTool::m_nominalRadLength {this, "ModuleRadLength", 0.04, "fixed radiation thickness of the detector modules"}
private

Definition at line 128 of file TrigInDetTrackFollowingTool.h.

128{this, "ModuleRadLength", 0.04, "fixed radiation thickness of the detector modules"};

◆ m_pixcontainerkey

SG::ReadHandleKey<InDet::PixelClusterContainer> TrigInDetTrackFollowingTool::m_pixcontainerkey {this, "PixelClusterContainer","ITkPixelClusters"}
private

Definition at line 131 of file TrigInDetTrackFollowingTool.h.

131{this, "PixelClusterContainer","ITkPixelClusters"};

◆ m_sctcontainerkey

SG::ReadHandleKey<InDet::SCT_ClusterContainer> TrigInDetTrackFollowingTool::m_sctcontainerkey {this, "SCT_ClusterContainer", "ITkStripClusters"}
private

Definition at line 132 of file TrigInDetTrackFollowingTool.h.

132{this, "SCT_ClusterContainer", "ITkStripClusters"};

◆ m_useDetectorThickness

Gaudi::Property<bool> TrigInDetTrackFollowingTool::m_useDetectorThickness {this, "UseDetectorThickness", false, "get Si-modules thickness from InDet Geometry"}
private

Definition at line 127 of file TrigInDetTrackFollowingTool.h.

127{this, "UseDetectorThickness", false, "get Si-modules thickness from InDet Geometry"};

◆ m_useHitErrors

Gaudi::Property<bool> TrigInDetTrackFollowingTool::m_useHitErrors {this, "UseHitErrors", false, "use PrepRawData errors"}
private

Definition at line 126 of file TrigInDetTrackFollowingTool.h.

126{this, "UseHitErrors", false, "use PrepRawData errors"};

◆ m_varHandleArraysDeclared

bool AthCommonDataStore< AthCommonMsg< AlgTool > >::m_varHandleArraysDeclared
privateinherited

Definition at line 399 of file AthCommonDataStore.h.

◆ m_vhka

std::vector<SG::VarHandleKeyArray*> AthCommonDataStore< AthCommonMsg< AlgTool > >::m_vhka
privateinherited

Definition at line 398 of file AthCommonDataStore.h.

◆ m_winX_Pixels

Gaudi::Property<double> TrigInDetTrackFollowingTool::m_winX_Pixels {this, "XSearchWindowPixels", 3.0, "x-size of hit search window for Pixels"}
private

Definition at line 122 of file TrigInDetTrackFollowingTool.h.

122{this, "XSearchWindowPixels", 3.0, "x-size of hit search window for Pixels"};

◆ m_winX_Strips

Gaudi::Property<double> TrigInDetTrackFollowingTool::m_winX_Strips {this, "XSearchWindowStrips", 3.0, "x-size of hit search window for Strips"}
private

Definition at line 124 of file TrigInDetTrackFollowingTool.h.

124{this, "XSearchWindowStrips", 3.0, "x-size of hit search window for Strips"};

◆ m_winY_Pixels

Gaudi::Property<double> TrigInDetTrackFollowingTool::m_winY_Pixels {this, "YSearchWindowPixels", 3.0, "y-size of hit search window for Pixels"}
private

Definition at line 123 of file TrigInDetTrackFollowingTool.h.

123{this, "YSearchWindowPixels", 3.0, "y-size of hit search window for Pixels"};

The documentation for this class was generated from the following files: