ATLAS Offline Software
Loading...
Searching...
No Matches
Trk::TrkDistanceFinderNeutralCharged Class Reference

#include <TrkDistanceFinderNeutralCharged.h>

Inheritance diagram for Trk::TrkDistanceFinderNeutralCharged:
Collaboration diagram for Trk::TrkDistanceFinderNeutralCharged:

Public Member Functions

 TrkDistanceFinderNeutralCharged (const std::string &t, const std::string &n, const IInterface *p)
 ~TrkDistanceFinderNeutralCharged ()
std::pair< Amg::Vector3D, double > getPointAndDistance (const Trk::NeutralTrack &, const Trk::Perigee &, double &, MagField::AtlasFieldCache &fieldCache) 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

Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
 specialization for handling Gaudi::Property<SG::VarHandleKey>

Private Attributes

double m_precision
double m_maxloopnumber
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 22 of file TrkDistanceFinderNeutralCharged.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

◆ TrkDistanceFinderNeutralCharged()

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

Definition at line 22 of file TrkDistanceFinderNeutralCharged.cxx.

22 :
23 AthAlgTool(t,n,p),
24 m_precision(1e-8),
26{
27 declareProperty("Precision",m_precision);
29 declareInterface<TrkDistanceFinderNeutralCharged>(this);
30
31}
AthAlgTool()
Default constructor:
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)

◆ ~TrkDistanceFinderNeutralCharged()

Trk::TrkDistanceFinderNeutralCharged::~TrkDistanceFinderNeutralCharged ( )
default

Member Function Documentation

◆ 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 }

◆ 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.

◆ 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

◆ getPointAndDistance()

std::pair< Amg::Vector3D, double > Trk::TrkDistanceFinderNeutralCharged::getPointAndDistance ( const Trk::NeutralTrack & neutraltrk,
const Trk::Perigee & chargedtrk,
double & distanceOnAxis,
MagField::AtlasFieldCache & fieldCache ) const

Definition at line 37 of file TrkDistanceFinderNeutralCharged.cxx.

40 {
41
42 double b_phi0=chargedtrk.parameters()[Trk::phi0];
43 double b_cosphi0=cos(b_phi0);
44 double b_sinphi0=sin(b_phi0);
45 double b_x0=chargedtrk.associatedSurface().center().x()-
46 chargedtrk.parameters()[Trk::d0]*b_sinphi0;
47 double b_y0=chargedtrk.associatedSurface().center().y()+
48 chargedtrk.parameters()[Trk::d0]*b_cosphi0;
49 double b_z0=chargedtrk.associatedSurface().center().z()+
50 chargedtrk.parameters()[Trk::z0];
51
52 Amg::Vector3D magnFieldVect;
53 fieldCache.getField(chargedtrk.associatedSurface().center().data(),magnFieldVect.data());
54
55 //Magnetic field at (x0,y0,z0)
56 double Bz=magnFieldVect.z()*299.792;//B field in Gev/mm
57
58 double b_Rt=getRadiusOfCurvature(chargedtrk,Bz);
59
60 double b_cottheta=1./tan(chargedtrk.parameters()[Trk::theta]);
61
62 //for neutrals
63 double a_x0=neutraltrk.position()[0];
64 double a_y0=neutraltrk.position()[1];
65 double a_z0=neutraltrk.position()[2];
66 double a_px=neutraltrk.momentum()[0];
67 double a_py=neutraltrk.momentum()[1];
68 double a_pz=neutraltrk.momentum()[2];
69
70
71 //more elaborate pieces for later
72 double DxNoVar=b_x0+b_Rt*b_sinphi0-a_x0;
73 double DyNoVar=b_y0-b_Rt*b_cosphi0-a_y0;
74 double DzNoVar=b_z0+b_Rt*b_cottheta*b_phi0-a_z0;
75
76 //start in both cases from the vertex
77 double b_phi=b_phi0;
78 double b_cosphi=cos(b_phi);
79 double b_sinphi=sin(b_phi);
80 double a_lambda=0.;
81
82 double deltab_phi=0.;
83 double deltaa_lambda=0.;
84
85 int loopsnumber=0;
86
87 bool isok=false;
88
89 bool secondTimeSaddleProblem=false;
90
91 while (!isok) {
92
93 // std::cout << "Loop number: " << loopsnumber << std::endl;
94 // std::cout << "phi is: " << b_phi << " lambda " << a_lambda << std::endl;
95
96 //ONLY FOR DEBUG * from here
97
98 double x1fin=b_x0+b_Rt*(sin(b_phi0)-sin(b_phi));
99 double y1fin=b_y0+b_Rt*(cos(b_phi)-cos(b_phi0));
100 double z1fin=b_z0+b_Rt*b_cottheta*(b_phi0-b_phi);
101
102 double x2fin=a_x0+a_px*a_lambda;
103 double y2fin=a_y0+a_py*a_lambda;
104 double z2fin=a_z0+a_pz*a_lambda;
105
106 ATH_MSG_VERBOSE ("x1fin " << x1fin << "y1fin " << y1fin << "z1fin " << z1fin);
107 ATH_MSG_VERBOSE ("x2fin " << x2fin << "y2fin " << y2fin << "z2fin " << z2fin);
108
109 double distance = sqrt((x1fin-x2fin)*(x1fin-x2fin)+
110 (y1fin-y2fin)*(y1fin-y2fin)+
111 (z1fin-z2fin)*(z1fin-z2fin));
112
113 ATH_MSG_VERBOSE ("distance " << distance);
114
115 //to here (please delete afterwards)
116
117 loopsnumber+=1;
118
119 double d1db_phi=((DxNoVar-a_px*a_lambda)*(-b_Rt*b_cosphi)+
120 (DyNoVar-a_py*a_lambda)*(-b_Rt*b_sinphi)+
121 (DzNoVar-b_Rt*b_phi*b_cottheta-a_pz*a_lambda)*(-b_Rt*b_cottheta));
122
123 double d1da_lambda=((DxNoVar-b_Rt*b_sinphi-a_px*a_lambda)*(-a_px)+
124 (DyNoVar+b_Rt*b_cosphi-a_py*a_lambda)*(-a_py)+
125 (DzNoVar-b_Rt*b_cottheta*b_phi-a_pz*a_lambda)*(-a_pz));
126
127 //second derivatives (d^2/d^2(a) d^2/d^2(b) d^2/d(a)d(b) )
128
129 double d2db_phi2=((DxNoVar-a_px*a_lambda)*(b_Rt*b_sinphi)+
130 (DyNoVar-a_py*a_lambda)*(-b_Rt*b_cosphi)+
131 (b_Rt*b_cottheta*b_Rt*b_cottheta));
132
133 double d2dadb_lambdaphi=a_px*b_Rt*b_cosphi+a_py*b_Rt*b_sinphi+a_pz*b_Rt*b_cottheta;
134
135 // double m_d2da_lambda2=((DxNoVar-b_Rt*b_sinphi-a_px*a_lambda)*(a_px*a_px)+
136 // (DyNoVar+b_Rt*b_cosphi-a_py*a_lambda)*(a_py*a_py)+
137 // (DzNoVar-b_Rt*b_cottheta*b_phi-a_pz*a_lambda)*(a_pz*a_pz));
138
139 //small error :-)
140
141 double d2da_lambda2=((a_px*a_px)+
142 (a_py*a_py)+
143 (a_pz*a_pz));
144
145
146 ATH_MSG_VERBOSE (" d1db_phi: " << d1db_phi << " d1da_lambda " << d1da_lambda);
147
148 ATH_MSG_VERBOSE (" d2phi2 " << d2db_phi2 << " d2lambda2 " <<d2da_lambda2 <<
149 " d2lambdaphi " << d2dadb_lambdaphi);
150
151
152 double det=d2da_lambda2*d2db_phi2-d2dadb_lambdaphi*d2dadb_lambdaphi;
153
154 if (det<0) {
155 if (!secondTimeSaddleProblem)
156 {
157 secondTimeSaddleProblem=true;
158
159 ATH_MSG_DEBUG ("Now try first to recover from the problem");
160 double denominator=1./((a_px*a_px)+(a_py*a_py));
161 double firstTerm=(DxNoVar*a_px+DyNoVar*a_py)*denominator;
162 double toSquare=(a_px*DyNoVar-a_py*DxNoVar);
163 double toSqrt=b_Rt*b_Rt/denominator-toSquare*toSquare;
164 if (toSqrt<0)
165 {
166 ATH_MSG_WARNING ("No intersection between neutral and charged track can be found. SKipping...");
167 continue;
168 }
169
170 double secondTerm=TMath::Sqrt(toSqrt)*denominator;
171 double lambda1=firstTerm+secondTerm;
172 double lambda2=firstTerm-secondTerm;
173
174 ATH_MSG_VERBOSE ("first solution: " << lambda1 << " second solution " << lambda2);
175
176 double phi1=TMath::ATan2((DxNoVar-a_px*lambda1)/b_Rt,-(DyNoVar-a_py*lambda1)/b_Rt);
177 double phi2=TMath::ATan2((DxNoVar-a_px*lambda2)/b_Rt,-(DyNoVar-a_py*lambda2)/b_Rt);
178
179 double x1case1=b_x0+b_Rt*(sin(b_phi0)-sin(phi1));
180 double y1case1=b_y0+b_Rt*(cos(phi1)-cos(b_phi0));
181 double z1case1=b_z0+b_Rt*b_cottheta*(b_phi0-phi1);
182
183 double x2case1=a_x0+a_px*lambda1;
184 double y2case1=a_y0+a_py*lambda1;
185 double z2case1=a_z0+a_pz*lambda1;
186
187 double x1case2=b_x0+b_Rt*(sin(b_phi0)-sin(phi2));
188 double y1case2=b_y0+b_Rt*(cos(phi2)-cos(b_phi0));
189 double z1case2=b_z0+b_Rt*b_cottheta*(b_phi0-phi2);
190
191 double x2case2=a_x0+a_px*lambda2;
192 double y2case2=a_y0+a_py*lambda2;
193 double z2case2=a_z0+a_pz*lambda2;
194
195 ATH_MSG_VERBOSE ("solution1 x1: " << x1case1 << " x2: " << x2case1 << " y1: " << y1case1 << " y2: " << y2case1);
196 ATH_MSG_VERBOSE ("solution2 x1: " << x1case2 << " x2: " << x2case2 << " y1: " << y1case2 << " y2: " << y2case2);
197
198 double deltaz1=fabs(z2case1-z1case1);
199 double deltaz2=fabs(z2case2-z1case2);
200
201 ATH_MSG_VERBOSE (" deltaz1: " << deltaz1 << " deltaz2: " << deltaz2);
202
203 if (deltaz1<=deltaz2)
204 {
205 b_phi=phi1;
206 a_lambda=lambda1;
207 }
208 else
209 {
210 b_phi=phi2;
211 a_lambda=lambda2;
212 }
213 //store cos and sin of phi
214 b_cosphi=cos(b_phi);
215 b_sinphi=sin(b_phi);
216 continue;
217 }
218
219
220 ATH_MSG_WARNING ("Hessian is negative: saddle point");
221 throw Error::NewtonProblem("Hessian is negative");
222
223 }
224 if (det>0&&d2da_lambda2<0) {
225 ATH_MSG_WARNING("Hessian indicates a maximum: derivative will be zero but result incorrect");
226 throw Error::NewtonProblem("Maximum point found");
227 }
228
229 //Now apply the Newton formula in more than one dimension
230 deltaa_lambda=-(d2db_phi2*d1da_lambda-d2dadb_lambdaphi*d1db_phi)/det;
231 deltab_phi=-(-d2dadb_lambdaphi*d1da_lambda+d2da_lambda2*d1db_phi)/det;
232
233 a_lambda+=deltaa_lambda;
234 b_phi+=deltab_phi;
235
236 //store cos and sin of phi
237 b_cosphi=cos(b_phi);
238 b_sinphi=sin(b_phi);
239
240 if (sqrt(std::pow(deltaa_lambda,2)+std::pow(deltab_phi,2))<m_precision ||
241 loopsnumber>m_maxloopnumber) isok=true;
242
243 }
244
245 if (loopsnumber>m_maxloopnumber) throw Error::NewtonProblem("Could not find minimum distance: max loops number reached"); //now return error, see how to do it...
246
247 double x1fin=b_x0+b_Rt*(sin(b_phi0)-sin(b_phi));
248 double y1fin=b_y0+b_Rt*(cos(b_phi)-cos(b_phi0));
249 double z1fin=b_z0+b_Rt*b_cottheta*(b_phi0-b_phi);
250
251 double x2fin=a_x0+a_px*a_lambda;
252 double y2fin=a_y0+a_py*a_lambda;
253 double z2fin=a_z0+a_pz*a_lambda;
254
255 ATH_MSG_VERBOSE ("x1fin " << x1fin << "y1fin " << y1fin << "z1fin " << z1fin);
256 ATH_MSG_VERBOSE ("x2fin " << x2fin << "y2fin " << y2fin << "z2fin " << z2fin);
257
258 double distance = sqrt((x1fin-x2fin)*(x1fin-x2fin)+
259 (y1fin-y2fin)*(y1fin-y2fin)+
260 (z1fin-z2fin)*(z1fin-z2fin));
261
262 ATH_MSG_DEBUG ("distance " << distance << std::endl
263 << "fitted dist from primary vertex " <<
264 a_lambda/fabs(a_lambda)*sqrt(std::pow(a_px*a_lambda,2)+
265 std::pow(a_py*a_lambda,2)+std::pow(a_pz*a_lambda,2)) );
266
267 distanceOnAxis=a_lambda/fabs(a_lambda)*sqrt(std::pow(a_px*a_lambda,2)+
268 std::pow(a_py*a_lambda,2)+
269 std::pow(a_pz*a_lambda,2));
270
271 Amg::Vector3D r1(x1fin,y1fin,z1fin);
272 Amg::Vector3D r2(x2fin,y2fin,z2fin);
273
274 return std::pair<Amg::Vector3D,double>(r2,distance);//give back position on neutral track... (jet axis)
275
276}
#define ATH_MSG_VERBOSE(x)
#define ATH_MSG_WARNING(x)
#define ATH_MSG_DEBUG(x)
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,...
const Amg::Vector3D & momentum() const
const Amg::Vector3D & position() const
virtual const S & associatedSurface() const override final
Access to the Surface method.
const Amg::Vector3D & center() const
Returns the center position of the Surface.
float distance(const Amg::Vector3D &p1, const Amg::Vector3D &p2)
calculates the distance between two point in 3D space
Eigen::Matrix< double, 3, 1 > Vector3D
@ phi0
Definition ParamDefs.h:65
@ theta
Definition ParamDefs.h:66
@ d0
Definition ParamDefs.h:63
@ z0
Definition ParamDefs.h:64

◆ 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 & Trk::TrkDistanceFinderNeutralCharged::interfaceID ( )
inlinestatic

Definition at line 34 of file TrkDistanceFinderNeutralCharged.h.

35 {
37 };
static const InterfaceID IID_TrkDistanceFinderNeutralCharged("TrkDistanceFinderNeutralCharged", 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.

◆ 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 }

◆ 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.

◆ 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_maxloopnumber

double Trk::TrkDistanceFinderNeutralCharged::m_maxloopnumber
private

Definition at line 43 of file TrkDistanceFinderNeutralCharged.h.

◆ m_precision

double Trk::TrkDistanceFinderNeutralCharged::m_precision
private

Definition at line 42 of file TrkDistanceFinderNeutralCharged.h.

◆ 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.


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