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

#include <MuonSegmentMomentumFromField.h>

Inheritance diagram for MuonSegmentMomentumFromField:
Collaboration diagram for MuonSegmentMomentumFromField:

Public Member Functions

 MuonSegmentMomentumFromField (const std::string &, const std::string &, const IInterface *)
 constructor
MuonSegmentMomentumFromFieldoperator= (const MuonSegmentMomentumFromField &)=delete
virtual StatusCode initialize ()
 to initiate private members
virtual void fitMomentum2Segments (const EventContext &ctx, const Muon::MuonSegment *segment1, const Muon::MuonSegment *segment2, double &signedMomentum) const
 fits a momentum to 2 segments
virtual void fitMomentum2Segments_old (const EventContext &ctx, const Muon::MuonSegment *segment1, const Muon::MuonSegment *segment2, double &signedMomentum) const
virtual void fitMomentumVectorSegments (const EventContext &, const std::vector< const Muon::MuonSegment * > &, double &signedMomentum) const
 fits a momentum to a vector of segments
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

double fieldIntegralEstimate (const EventContext &ctx, const Muon::MuonSegment *segment1, const Muon::MuonSegment *segment2) const
double fieldIntegralEstimate_old (const EventContext &ctx, const Muon::MuonSegment *segment1, const Muon::MuonSegment *segment2) const
Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
 specialization for handling Gaudi::Property<SG::VarHandleKey>

Private Attributes

ToolHandle< Trk::IPropagatorm_propagator
ToolHandle< Trk::INavigatorm_navigator
ServiceHandle< Muon::IMuonIdHelperSvcm_idHelperSvc {this, "MuonIdHelperSvc", "Muon::MuonIdHelperSvc/MuonIdHelperSvc"}
SG::ReadCondHandleKey< AtlasFieldCacheCondObjm_fieldCondObjInputKey {this, "AtlasFieldCacheCondObj", "fieldCondObj"}
Gaudi::Property< bool > m_doOld {this, "DoOld", false, "Use old fitMomentum2Segments"}
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 23 of file MuonSegmentMomentumFromField.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

◆ MuonSegmentMomentumFromField()

MuonSegmentMomentumFromField::MuonSegmentMomentumFromField ( const std::string & type,
const std::string & name,
const IInterface * parent )

constructor

Definition at line 20 of file MuonSegmentMomentumFromField.cxx.

21 :AthAlgTool(type,name,parent)
22{
23 declareInterface<IMuonSegmentMomentumEstimator>(this);
24}
AthAlgTool()
Default constructor:

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

◆ 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

◆ fieldIntegralEstimate()

double MuonSegmentMomentumFromField::fieldIntegralEstimate ( const EventContext & ctx,
const Muon::MuonSegment * segment1,
const Muon::MuonSegment * segment2 ) const
private

Definition at line 66 of file MuonSegmentMomentumFromField.cxx.

67{
68
69 Amg::Vector3D pos1=segment1->globalPosition();
70 Amg::Vector3D pos2=segment2->globalPosition();
71 Amg::Vector3D posdiff=pos2-pos1;
72 if (posdiff.dot(segment1->globalDirection())<0){
73 posdiff=-posdiff;
74 pos1=segment2->globalPosition();
75 pos2=segment1->globalPosition();
76 }
77 Amg::Vector3D point1=pos1+.25*posdiff;
78 Amg::Vector3D point2=pos1+.5*posdiff;
79 Amg::Vector3D point3=pos1+.75*posdiff;
80 Amg::Vector3D field1,field2,field3;
81
82 MagField::AtlasFieldCache fieldCache;
83 SG::ReadCondHandle<AtlasFieldCacheCondObj> readHandle{m_fieldCondObjInputKey, ctx};
84 const AtlasFieldCacheCondObj* fieldCondObj{*readHandle};
85 if (!fieldCondObj) {
86 throw std::runtime_error(Form("File: %s, Line: %d\nMuonSegmentMomentumFromField::fieldIntegralEstimate() - Failed to retrieve AtlasFieldCacheCondObj with key %s", __FILE__, __LINE__, (m_fieldCondObjInputKey.key()).c_str()));
87 }
88 fieldCondObj->getInitializedCache(fieldCache);
89
90 fieldCache.getField(point1.data(),field1.data()); // field in kiloTesla
91 fieldCache.getField(point2.data(),field2.data());
92 fieldCache.getField(point3.data(),field3.data());
93 ATH_MSG_DEBUG("Mid Point " << Amg::toString(point2) << " field " << Amg::toString(field2) );
94 field1[2]=field2[2]=field3[2]=0;
95 Amg::Vector3D fieldsum=field1+field2+field3;
96 Amg::Vector3D crossvec=posdiff.unit().cross(field1+field2+field3);
97 Amg::Vector2D rphidir(-posdiff.y(),posdiff.x());
98 double averagelcrossB=crossvec.mag()/3;
99 if (rphidir.x()*fieldsum.x()+rphidir.y()*fieldsum.y()<0) averagelcrossB=-averagelcrossB;
100 ATH_MSG_DEBUG("field integral " << averagelcrossB << "dist " << posdiff.mag() << " tot " << averagelcrossB*posdiff.mag());
101 return averagelcrossB*posdiff.mag();
102}
#define ATH_MSG_DEBUG(x)
void getInitializedCache(MagField::AtlasFieldCache &cache) const
get B field cache for evaluation as a function of 2-d or 3-d position.
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,...
SG::ReadCondHandleKey< AtlasFieldCacheCondObj > m_fieldCondObjInputKey
virtual const Amg::Vector3D & globalPosition() const override final
global position
std::string toString(const Translation3D &translation, int precision=4)
GeoPrimitvesToStringConverter.
Eigen::Matrix< double, 2, 1 > Vector2D
Eigen::Matrix< double, 3, 1 > Vector3D

◆ fieldIntegralEstimate_old()

double MuonSegmentMomentumFromField::fieldIntegralEstimate_old ( const EventContext & ctx,
const Muon::MuonSegment * segment1,
const Muon::MuonSegment * segment2 ) const
private

Definition at line 244 of file MuonSegmentMomentumFromField.cxx.

245{
246
247 Amg::Vector3D pos1=segment1->globalPosition();
248 Amg::Vector3D pos2=segment2->globalPosition();
249 Amg::Vector3D posdiff=pos2-pos1;
250 if (posdiff.dot(segment1->globalDirection())<0){
251 posdiff=-posdiff;
252 pos1=segment2->globalPosition();
253 pos2=segment1->globalPosition();
254 }
255 Amg::Vector3D point1=pos1+.25*posdiff;
256 Amg::Vector3D point2=pos1+.5*posdiff;
257 Amg::Vector3D point3=pos1+.75*posdiff;
258 Amg::Vector3D field1,field2,field3;
259
260 MagField::AtlasFieldCache fieldCache;
261 SG::ReadCondHandle<AtlasFieldCacheCondObj> readHandle{m_fieldCondObjInputKey, ctx};
262 const AtlasFieldCacheCondObj* fieldCondObj{*readHandle};
263 if (!fieldCondObj) {
264 throw std::runtime_error(Form("File: %s, Line: %d\nMuonSegmentMomentumFromField::fieldIntegralEstimate_old() - Failed to retrieve AtlasFieldCacheCondObj with key %s", __FILE__, __LINE__, (m_fieldCondObjInputKey.key()).c_str()));
265 }
266 fieldCondObj->getInitializedCache(fieldCache);
267
268 fieldCache.getField(point1.data(),field1.data()); // field in kiloTesla
269 fieldCache.getField(point2.data(),field2.data());
270 fieldCache.getField(point3.data(),field3.data());
271 ATH_MSG_DEBUG("Mid Point " << Amg::toString(point2) << " field " << Amg::toString(field2) );
272 double averageBcrossl=(field1.cross(posdiff.unit()).mag()+field2.cross(posdiff.unit()).mag()+field3.cross(posdiff.unit()).mag())/3;
273 ATH_MSG_DEBUG("field integral " << averageBcrossl << "dist " << posdiff.mag() << " tot " << averageBcrossl*posdiff.mag());
274 return averageBcrossl*posdiff.mag();
275}

◆ fitMomentum2Segments()

void MuonSegmentMomentumFromField::fitMomentum2Segments ( const EventContext & ctx,
const Muon::MuonSegment * segment1,
const Muon::MuonSegment * segment2,
double & signedMomentum ) const
virtual

fits a momentum to 2 segments

Estimate signed momentum for two segments by fitting 2 segments to one approximate track model

Implements Muon::IMuonSegmentMomentumEstimator.

Definition at line 104 of file MuonSegmentMomentumFromField.cxx.

105{
106
109
110 ATH_MSG_VERBOSE(" Executing MuonSegmentMomentumTool fitMomentum2Segments ");
111
112 const Muon::MuonSegment* myseg1=segment1,*myseg2=segment2;
113 if (myseg1->globalDirection().dot(myseg2->globalPosition()-myseg1->globalPosition())<0) {
114 myseg1=segment2;
115 myseg2=segment1;
116 }
117 double fieldintegral=fieldIntegralEstimate(ctx, myseg1,myseg2);
118 double theta1=myseg1->globalDirection().theta();
119 double theta2=myseg2->globalDirection().theta();
120 double phi1=myseg1->globalDirection().phi();
121 double phi2=myseg2->globalDirection().phi();
122 double deltaphi=std::abs(phi2-phi1);
123 double deltatheta=99999999;
124 if (std::abs(deltaphi-2*M_PI)<deltaphi) deltaphi=std::abs(deltaphi-2*M_PI);
125 if (deltaphi<M_PI/2) deltatheta=theta2-theta1;
126 else {
127 if (theta1>M_PI/2) deltatheta=2*M_PI-theta2-theta1;
128 else deltatheta=-theta1-theta2;
129 }
130 double dist1=-1,dist2=-1;
131 const Trk::RIO_OnTrack *firstphi1=nullptr,*lastphi1=nullptr,*firstphi2=nullptr,*lastphi2=nullptr;
132
133 const Muon::MuonSegment *bestseg=myseg1,*worstseg=myseg2;
134 for (int i=0;i<(int)myseg1->numberOfMeasurementBases();i++){
135 const Trk::RIO_OnTrack *rot=dynamic_cast<const Trk::RIO_OnTrack *>(myseg1->measurement(i));
136 if (!rot){
137 const Trk::CompetingRIOsOnTrack *crot=dynamic_cast<const Trk::CompetingRIOsOnTrack *>(myseg1->measurement(i));
138 if (crot) rot=&crot->rioOnTrack(0);
139 }
140 if (!rot) continue;
141 Identifier id=rot->identify();
142
143 if ((m_idHelperSvc->isRpc(id) && m_idHelperSvc->rpcIdHelper().measuresPhi(id)) || (m_idHelperSvc->isCsc(id) && m_idHelperSvc->cscIdHelper().measuresPhi(id)) || (m_idHelperSvc->isTgc(id) && m_idHelperSvc->tgcIdHelper().isStrip(id))
144 || (m_idHelperSvc->issTgc(id) && m_idHelperSvc->stgcIdHelper().measuresPhi(id) ) ){
145 if (!firstphi1) firstphi1=rot;
146 lastphi1=rot;
147 }
148 }
149 for (int i=0;i<(int)myseg2->numberOfMeasurementBases();i++){
150 const Trk::RIO_OnTrack *rot=dynamic_cast<const Trk::RIO_OnTrack *>(myseg2->measurement(i));
151 if (!rot){
152 const Trk::CompetingRIOsOnTrack *crot=dynamic_cast<const Trk::CompetingRIOsOnTrack *>(myseg2->measurement(i));
153 if (crot) rot=&crot->rioOnTrack(0);
154 }
155 if (!rot) continue;
156 Identifier id=rot->identify();
157 if ((m_idHelperSvc->isRpc(id) && m_idHelperSvc->rpcIdHelper().measuresPhi(id)) || (m_idHelperSvc->isCsc(id) && m_idHelperSvc->cscIdHelper().measuresPhi(id)) || (m_idHelperSvc->isTgc(id) && m_idHelperSvc->tgcIdHelper().isStrip(id))
158 || (m_idHelperSvc->issTgc(id) && m_idHelperSvc->stgcIdHelper().measuresPhi(id) ) ){
159 if (!firstphi2) firstphi2=rot;
160 lastphi2=rot;
161 }
162 }
163 bool flip = false;
164 if (firstphi1) dist1=std::abs((firstphi1->globalPosition()-lastphi1->globalPosition()).dot(myseg1->globalDirection()));
165 if (firstphi2) dist2=std::abs((firstphi2->globalPosition()-lastphi2->globalPosition()).dot(myseg2->globalDirection()));
166 if (dist2>dist1) {
167 flip = true;
168 bestseg=myseg2;
169 worstseg=myseg1;
170 }
171
172 for (int itry=0;itry<2;itry++){
173
174 if(itry==1) {
175// retry with other as best segmt segment
176 if(!flip) {
177 bestseg = myseg2;
178 worstseg=myseg1;
179 } else {
180 bestseg = myseg1;
181 worstseg=myseg2;
182 }
183 }
184
185 signedMomentum =-.3e3*fieldintegral/deltatheta;
186 if(std::abs(signedMomentum)<1000.) signedMomentum = 1e6;
187 ATH_MSG_DEBUG("integral: " << fieldintegral << " deltatheta: " << deltatheta << " signedmomentum : " << signedMomentum);
188 double residual = 9999.;
189 double resi[4],qoverp[4];
190 for (int i=0;i<4;i++){
191 std::optional<Trk::TransportJacobian> jac{};
192 Trk::AtaPlane startpar(bestseg->globalPosition(),bestseg->globalDirection().phi(),
193 bestseg->globalDirection().theta(),1/signedMomentum,bestseg->associatedSurface());
194 auto par=m_propagator->propagateParameters(
195 ctx, startpar,worstseg->associatedSurface(),
196 (bestseg==myseg1) ? Trk::alongMomentum : Trk::oppositeMomentum,false,
197 Trk::MagneticFieldProperties(Trk::FullField),jac,Trk::nonInteracting);
198 ATH_MSG_DEBUG("par: " << par.get() << " jac: " << *jac );
199 if (par && jac && (*jac)(1,4)!=0){
200 residual = worstseg->localParameters()[Trk::locY] - par->parameters()[Trk::locY];
201 resi[i] = residual;
202 qoverp[i] = 1/signedMomentum;
203// update qoverp
204 double delta_qoverp=residual/(*jac)(1,4);
205 double der_simple = -10.*(bestseg->globalPosition()-worstseg->globalPosition()).mag()/(.3*fieldintegral);
206 if(bestseg->globalDirection().z()<0) der_simple = - der_simple;
207 if(i>0) {
208 if(qoverp[i]!=qoverp[i-1]) {
209 double derivative = -(resi[i]-resi[i-1])/(qoverp[i]-qoverp[i-1]);
210 ATH_MSG_DEBUG(" numerical derivative " << derivative << " derivative from track " << (*jac)(1,4) << " der_simple " << der_simple);
211 if(std::abs(derivative)>std::abs((*jac)(1,4))) {
212 ATH_MSG_DEBUG(" use numerical derivative " << derivative << " derivative from track " << (*jac)(1,4));
213 delta_qoverp = residual/derivative;
214 }
215 }
216 } else {
217 if(std::abs(der_simple)>std::abs((*jac)(1,4))) {
218 ATH_MSG_DEBUG(" use simple numerical derivative " << der_simple << " derivative from track " << (*jac)(1,4));
219 delta_qoverp = residual/der_simple;
220 }
221 }
222 ATH_MSG_DEBUG("residual: " << residual << " jac " << (*jac)(1,4) << " signedmomentum: " << signedMomentum << " delta_qoverp " << delta_qoverp);
223 double signedMomentum_updated = signedMomentum/(1+signedMomentum*delta_qoverp);
224 if(std::abs(signedMomentum_updated)<1000.) {
225 ATH_MSG_DEBUG("Too low signed momentum " << signedMomentum_updated );
226// protect against too low momenta as propagation will fail
227 signedMomentum = signedMomentum_updated>0? 1000.:-1000;
228 } else {
229// do not change momentum beyond last iteration
230 if(i<3) signedMomentum = signedMomentum_updated;
231 }
232 //delete par;
233 if (std::abs(residual)<1) break;
234 }
235 }
236 if(std::abs(residual)>10.) {
237 ATH_MSG_DEBUG("NOT converged residual: " << residual << " itry " << itry);
238 if(itry==1) ATH_MSG_DEBUG("NOT converged residual after two trials ");
239 } else break;
240 }
241}
#define M_PI
#define ATH_MSG_VERBOSE(x)
double fieldIntegralEstimate(const EventContext &ctx, const Muon::MuonSegment *segment1, const Muon::MuonSegment *segment2) const
ToolHandle< Trk::IPropagator > m_propagator
ServiceHandle< Muon::IMuonIdHelperSvc > m_idHelperSvc
virtual const Trk::PlaneSurface & associatedSurface() const override final
returns the surface for the local to global transformation
virtual const RIO_OnTrack & rioOnTrack(unsigned int) const =0
returns the RIO_OnTrack (also known as ROT) objects depending on the integer.
Identifier identify() const
return the identifier -extends MeasurementBase
virtual const Amg::Vector3D & globalPosition() const override=0
Interface method to get the global Position.
const MeasurementBase * measurement(unsigned int) const
returns the Trk::MeasurementBase objects depending on the integer
unsigned int numberOfMeasurementBases() const
Return the number of contained Trk::MeasurementBase (s)
constexpr double derivative(const double x)
Evaluates the d-th derivative of the n-th Legendre polynomial at x.
@ oppositeMomentum
@ alongMomentum
@ FullField
Field is set to be realistic, but within a given Volume.
@ locY
local cartesian
Definition ParamDefs.h:38
ParametersT< TrackParametersDim, Charged, PlaneSurface > AtaPlane

◆ fitMomentum2Segments_old()

void MuonSegmentMomentumFromField::fitMomentum2Segments_old ( const EventContext & ctx,
const Muon::MuonSegment * segment1,
const Muon::MuonSegment * segment2,
double & signedMomentum ) const
virtual

Estimate signed momentum for two segments by fitting 2 segments to one approximate track model

Definition at line 277 of file MuonSegmentMomentumFromField.cxx.

278{
279
282
283 ATH_MSG_VERBOSE(" Executing MuonSegmentMomentumTool fitMomentum2Segments ");
284
285 const Muon::MuonSegment* myseg1=segment1,*myseg2=segment2;
286 if (myseg1->globalDirection().dot(myseg2->globalPosition()-myseg1->globalPosition())<0) {
287 myseg1=segment2;
288 myseg2=segment1;
289 }
290 double fieldintegral=fieldIntegralEstimate(ctx, myseg1,myseg2);
291 double theta1=myseg1->globalDirection().theta();
292 double theta2=myseg2->globalDirection().theta();
293 double phi1=myseg1->globalDirection().phi();
294 double phi2=myseg2->globalDirection().phi();
295 double deltaphi=std::abs(phi2-phi1);
296 double deltatheta=99999999;
297 if (std::abs(deltaphi-2*M_PI)<deltaphi) deltaphi=std::abs(deltaphi-2*M_PI);
298 if (deltaphi<M_PI/2) deltatheta=theta2-theta1;
299 else {
300 if (theta1>M_PI/2) deltatheta=2*M_PI-theta2-theta1;
301 else deltatheta=-theta1-theta2;
302 }
303 double dist1=-1,dist2=-1;
304 const Trk::RIO_OnTrack *firstphi1=nullptr,*lastphi1=nullptr,*firstphi2=nullptr,*lastphi2=nullptr;
305
306 const Muon::MuonSegment *bestseg=myseg1,*worstseg=myseg2;
307 for (int i=0;i<(int)myseg1->numberOfMeasurementBases();i++){
308 const Trk::RIO_OnTrack *rot=dynamic_cast<const Trk::RIO_OnTrack *>(myseg1->measurement(i));
309 if (!rot){
310 const Trk::CompetingRIOsOnTrack *crot=dynamic_cast<const Trk::CompetingRIOsOnTrack *>(myseg1->measurement(i));
311 if (crot) rot=&crot->rioOnTrack(0);
312 }
313 if (!rot) continue;
314 Identifier id=rot->identify();
315
316 if ((m_idHelperSvc->isRpc(id) && m_idHelperSvc->rpcIdHelper().measuresPhi(id)) || (m_idHelperSvc->isCsc(id) && m_idHelperSvc->cscIdHelper().measuresPhi(id)) ||
317 (m_idHelperSvc->isTgc(id) && m_idHelperSvc->tgcIdHelper().isStrip(id))){
318 if (!firstphi1) firstphi1=rot;
319 lastphi1=rot;
320 }
321 }
322 for (int i=0;i<(int)myseg2->numberOfMeasurementBases();i++){
323 const Trk::RIO_OnTrack *rot=dynamic_cast<const Trk::RIO_OnTrack *>(myseg2->measurement(i));
324 if (!rot){
325 const Trk::CompetingRIOsOnTrack *crot=dynamic_cast<const Trk::CompetingRIOsOnTrack *>(myseg2->measurement(i));
326 if (crot) rot=&crot->rioOnTrack(0);
327 }
328 if (!rot) continue;
329 Identifier id=rot->identify();
330 if ((m_idHelperSvc->isRpc(id) && m_idHelperSvc->rpcIdHelper().measuresPhi(id)) || (m_idHelperSvc->isCsc(id) && m_idHelperSvc->cscIdHelper().measuresPhi(id)) ||
331 (m_idHelperSvc->isTgc(id) && m_idHelperSvc->tgcIdHelper().isStrip(id))){
332 if (!firstphi2) firstphi2=rot;
333 lastphi2=rot;
334 }
335 }
336
337 if (firstphi1) dist1=std::abs((firstphi1->globalPosition()-lastphi1->globalPosition()).dot(myseg1->globalDirection()));
338 if (firstphi2) dist2=std::abs((firstphi2->globalPosition()-lastphi2->globalPosition()).dot(myseg2->globalDirection()));
339 if (dist2>dist1) {
340 bestseg=myseg2;
341 worstseg=myseg1;
342 }
343 signedMomentum =-.3e3*fieldintegral/deltatheta;
344 ATH_MSG_DEBUG("integral: " << fieldintegral << " deltatheta: " << deltatheta << " signedmomentum : " << signedMomentum);
345 for (int i=0;i<3;i++){
346 Trk::AtaPlane startpar(bestseg->globalPosition(),bestseg->globalDirection().phi(),bestseg->globalDirection().theta(),
347 1/signedMomentum,bestseg->associatedSurface());
348 std::optional<Trk::TransportJacobian> jac{};
349 auto par=m_propagator->propagateParameters(
350 ctx,startpar,worstseg->associatedSurface(),
351 (bestseg==myseg1) ? Trk::alongMomentum : Trk::oppositeMomentum,false,
352 Trk::MagneticFieldProperties(Trk::FullField),jac,Trk::nonInteracting);
353 if (par && jac && (*jac)(1,4)!=0){
354 double residual=worstseg->localParameters()[Trk::locY] - par->parameters()[Trk::locY];
355 double delta_qoverp=residual/(*jac)(1,4);
356 signedMomentum=1/(1/signedMomentum+delta_qoverp);
357 ATH_MSG_DEBUG("residual: " << residual << " jac " << (*jac)(1,4) << " dp " << delta_qoverp << " signedmomentum: " << signedMomentum);
358 // delete par;
359 if (std::abs(residual)<1) break;
360 }
361 }
362}

◆ fitMomentumVectorSegments()

void MuonSegmentMomentumFromField::fitMomentumVectorSegments ( const EventContext & ctx,
const std::vector< const Muon::MuonSegment * > & segments,
double & signedMomentum ) const
virtual

fits a momentum to a vector of segments

Estimate signed momentum from vector of MDT/CSC segments using fit to pairs of segments

Implements Muon::IMuonSegmentMomentumEstimator.

Definition at line 38 of file MuonSegmentMomentumFromField.cxx.

39{
40
43
44 ATH_MSG_VERBOSE(" Executing MuonSegmentMomentumTool fitMomentumVectorSegments ");
45 ATH_MSG_DEBUG(" fitMomentumVectorSegments " << segments.size() << " segments " );
46
47 std::vector<const Muon::MuonSegment*>::const_iterator it = segments.begin();
48 std::vector<const Muon::MuonSegment*>::const_iterator it2 = segments.begin();
49 std::vector<const Muon::MuonSegment*>::const_iterator it_end = segments.end();
50
51 ++it2;
52 double maxintegral=0;
53 while (it2!=it_end){
54 double integral= m_doOld ? fieldIntegralEstimate_old(ctx,*it,*it2) : fieldIntegralEstimate(ctx,*it,*it2);
55 if (std::abs(integral)>maxintegral){
56 maxintegral=std::abs(integral);
57 if( m_doOld ) fitMomentum2Segments_old(ctx, *it, *it2, signedMomentum);
58 else fitMomentum2Segments(ctx, *it, *it2, signedMomentum);
59 }
60 ++it;
61 ++it2;
62 }
63 ATH_MSG_DEBUG( " Estimated signed momentum " << signedMomentum );
64}
double fieldIntegralEstimate_old(const EventContext &ctx, const Muon::MuonSegment *segment1, const Muon::MuonSegment *segment2) const
virtual void fitMomentum2Segments_old(const EventContext &ctx, const Muon::MuonSegment *segment1, const Muon::MuonSegment *segment2, double &signedMomentum) const
virtual void fitMomentum2Segments(const EventContext &ctx, const Muon::MuonSegment *segment1, const Muon::MuonSegment *segment2, double &signedMomentum) const
fits a momentum to 2 segments
double integral(TH1 *h)
Definition computils.cxx:59

◆ initialize()

StatusCode MuonSegmentMomentumFromField::initialize ( )
virtual

to initiate private members

Definition at line 26 of file MuonSegmentMomentumFromField.cxx.

27{
28 ATH_MSG_VERBOSE(" MuonSegmentMomentumFromField::Initializing ");
29 ATH_CHECK(AthAlgTool::initialize());
30 ATH_CHECK(m_propagator.retrieve());
31 ATH_CHECK(m_navigator.retrieve());
32 ATH_CHECK(m_idHelperSvc.retrieve());
34 ATH_MSG_VERBOSE("End of Initializing");
35 return StatusCode::SUCCESS;
36}
#define ATH_CHECK
Evaluate an expression and check for errors.
ToolHandle< Trk::INavigator > m_navigator

◆ 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 & Muon::IMuonSegmentMomentumEstimator::interfaceID ( )
inlinestaticinherited

Definition at line 29 of file IMuonSegmentMomentumEstimator.h.

static const InterfaceID IID_IMuonSegmentMomentumEstimator("Muon::IMuonSegmentMomentumEstimator", 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 }

◆ operator=()

MuonSegmentMomentumFromField & MuonSegmentMomentumFromField::operator= ( const MuonSegmentMomentumFromField & )
delete

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

Gaudi::Property<bool> MuonSegmentMomentumFromField::m_doOld {this, "DoOld", false, "Use old fitMomentum2Segments"}
private

Definition at line 50 of file MuonSegmentMomentumFromField.h.

50{this, "DoOld", false, "Use old fitMomentum2Segments"};

◆ 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> MuonSegmentMomentumFromField::m_fieldCondObjInputKey {this, "AtlasFieldCacheCondObj", "fieldCondObj"}
private

Definition at line 48 of file MuonSegmentMomentumFromField.h.

48{this, "AtlasFieldCacheCondObj", "fieldCondObj"};

◆ m_idHelperSvc

ServiceHandle<Muon::IMuonIdHelperSvc> MuonSegmentMomentumFromField::m_idHelperSvc {this, "MuonIdHelperSvc", "Muon::MuonIdHelperSvc/MuonIdHelperSvc"}
private

Definition at line 47 of file MuonSegmentMomentumFromField.h.

47{this, "MuonIdHelperSvc", "Muon::MuonIdHelperSvc/MuonIdHelperSvc"};

◆ m_navigator

ToolHandle<Trk::INavigator> MuonSegmentMomentumFromField::m_navigator
private
Initial value:
{this, "NavigatorTool",
"Trk::Navigator/MuonNavigator"}

Definition at line 45 of file MuonSegmentMomentumFromField.h.

45 {this, "NavigatorTool",
46 "Trk::Navigator/MuonNavigator"};

◆ m_propagator

ToolHandle<Trk::IPropagator> MuonSegmentMomentumFromField::m_propagator
private
Initial value:
{this, "PropagatorTool",
"Trk::STEP_Propagator/MuonPropagator"}

Definition at line 43 of file MuonSegmentMomentumFromField.h.

43 {this, "PropagatorTool",
44 "Trk::STEP_Propagator/MuonPropagator"};

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