ATLAS Offline Software
Classes | Public Types | Public Member Functions | Static Public Member Functions | Private Types | Private Member Functions | Private Attributes | List of all members
MuonR4::MdtSegmentFitter Class Reference

#include <MdtSegmentFitter.h>

Inheritance diagram for MuonR4::MdtSegmentFitter:
Collaboration diagram for MuonR4::MdtSegmentFitter:

Classes

struct  Config
 
struct  LineWithPartials
 Store the partial derivative of the line w.r.t. More...
 
struct  ResidualAuxillaries
 Helper struct carrying the space for all auxillary variables needed to calculate the residual from wire measurements. More...
 
struct  ResidualWithPartials
 Helper struct carrying the residual with its derivatives. More...
 

Public Types

using ParamDefs = SegmentFit::ParamDefs
 
using Parameters = SegmentFit::Parameters
 
using HitType = std::unique_ptr< CalibratedSpacePoint >
 
using HitVec = std::vector< HitType >
 

Public Member Functions

 MdtSegmentFitter (const std::string &name, Config &&config)
 Standard constructor. More...
 
SegmentFitResult fitSegment (const EventContext &ctx, HitVec &&calibHits, const Parameters &startPars, const Amg::Transform3D &localToGlobal) const
 
void calculateWireResiduals (const LineWithPartials &line, const CalibratedSpacePoint &spacePoint, ResidualWithPartials &residual) const
 Calculates the residuals together with the corresponding derivatives for a drift-circle measurement. More...
 
void calculateStripResiduals (const LineWithPartials &line, const CalibratedSpacePoint &spacePoint, ResidualWithPartials &residual) const
 Calculates the residual together with hte correspdonding derivatives for strip measurements. More...
 
bool msgLvl (const MSG::Level lvl) const
 Test the output level. More...
 
MsgStream & msg () const
 The standard message stream. More...
 
MsgStream & msg (const MSG::Level lvl) const
 The standard message stream. More...
 
void setLevel (MSG::Level lvl)
 Change the current logging level. More...
 

Static Public Member Functions

static void updateLinePartials (const Parameters &fitPars, LineWithPartials &line)
 Updates the line parameters together with its first & second order partial derivatives. More...
 

Private Types

enum  UpdateStatus { UpdateStatus::allOkay = 0, UpdateStatus::outOfBounds = 1, UpdateStatus::noChange = 2 }
 
using MeasCov_t = CalibratedSpacePoint::Covariance_t
 Updates the chi2, its Gradient & Hessian from the measurement residual. More...
 

Private Member Functions

bool recalibrate (const EventContext &ctx, SegmentFitResult &fitResult) const
 Recalibrate the measurements participating in the fit & shift them into the centre-of gravity frame. More...
 
bool updateHitSummary (SegmentFitResult &fitResult) const
 Brief updates the hit summary from the contributing hits. More...
 
void updateDriftSigns (const Amg::Vector3D &segPos, const Amg::Vector3D &segDir, SegmentFitResult &fitRes) const
 Update the signs of the measurement. More...
 
void updateDerivatives (const ResidualWithPartials &fitMeas, const MeasCov_t &measCovariance, AmgVector(5)&gradient, AmgSymMatrix(5)&hessian, double &chi2, int startPar) const
 
template<unsigned int nDim>
UpdateStatus updateParameters (Parameters &currentPars, Parameters &previousPars, Parameters &currGrad, Parameters &prevGrad, const AmgSymMatrix(5)&hessian) const
 Update step of the segment parameters using the Hessian and the gradient. More...
 
template<unsigned int nDim>
void blockCovariance (const AmgSymMatrix(5)&hessian, SegmentFit::Covariance &covariance) const
 
void initMessaging () const
 Initialize our message level and MessageSvc. More...
 

Private Attributes

Config m_cfg {}
 
std::string m_nm
 Message source name. More...
 
boost::thread_specific_ptr< MsgStream > m_msg_tls
 MsgStream instance (a std::cout like with print-out levels) More...
 
std::atomic< IMessageSvc * > m_imsg { nullptr }
 MessageSvc pointer. More...
 
std::atomic< MSG::Level > m_lvl { MSG::NIL }
 Current logging level. More...
 
std::atomic_flag m_initialized ATLAS_THREAD_SAFE = ATOMIC_FLAG_INIT
 Messaging initialized (initMessaging) More...
 

Detailed Description

Definition at line 19 of file MdtSegmentFitter.h.

Member Typedef Documentation

◆ HitType

Definition at line 23 of file MdtSegmentFitter.h.

◆ HitVec

Definition at line 24 of file MdtSegmentFitter.h.

◆ MeasCov_t

Updates the chi2, its Gradient & Hessian from the measurement residual.

Depending on whether the time needs to be taken into account or not the template parameter nDim is either 3 or 2, respectively.

Parameters
fitMeasResidual of the prediction from the measurement together with its derivatives
measCovarianceCovariance matrix of the measurement
gradientReference to the gradient of the chi2 which is updated by the function
hessianReference to the Hessian matrix of the chi2 which is updated by the function
chi2Reference to the chi2 value itself which is updated by the function
startParNumber of the maximum fit parameter to consider

Definition at line 163 of file MdtSegmentFitter.h.

◆ ParamDefs

Definition at line 21 of file MdtSegmentFitter.h.

◆ Parameters

Definition at line 22 of file MdtSegmentFitter.h.

Member Enumeration Documentation

◆ UpdateStatus

Enumerator
allOkay 
outOfBounds 
noChange 

Definition at line 172 of file MdtSegmentFitter.h.

172  {
173  allOkay = 0,
174  outOfBounds = 1,
175  noChange = 2,
176  };

Constructor & Destructor Documentation

◆ MdtSegmentFitter()

MuonR4::MdtSegmentFitter::MdtSegmentFitter ( const std::string &  name,
Config &&  config 
)

Standard constructor.

Parameters
nameName to be printed in the messaging
configFit configuration parameters

Definition at line 51 of file MdtSegmentFitter.cxx.

51  :
53  m_cfg{std::move(config)}{}

Member Function Documentation

◆ blockCovariance()

template<unsigned int nDim>
void MuonR4::MdtSegmentFitter::blockCovariance ( const AmgSymMatrix(5)&  hessian,
SegmentFit::Covariance covariance 
) const
private

Definition at line 606 of file MdtSegmentFitter.cxx.

607  {
608 
609  covariance.setIdentity();
610  AmgSymMatrix(nDim) miniHessian = hessian.block<nDim, nDim>(0,0);
611  if (std::abs(miniHessian.determinant()) <= detCutOff) {
612  ATH_MSG_VERBOSE("Boeser mini hessian ("<<miniHessian.determinant()<<")\n"<<miniHessian<<"\n\n"<<hessian);
613  return;
614  }
615  ATH_MSG_VERBOSE("Hessian matrix: \n"<<hessian<<",\nblock Hessian:\n"<<miniHessian<<",\n determinant: "<<miniHessian.determinant());
616  covariance.block<nDim,nDim>(0,0) = miniHessian.inverse();
617  ATH_MSG_VERBOSE("covariance: \n"<<covariance);
618  }

◆ calculateStripResiduals()

void MuonR4::MdtSegmentFitter::calculateStripResiduals ( const LineWithPartials line,
const CalibratedSpacePoint spacePoint,
ResidualWithPartials residual 
) const

Calculates the residual together with hte correspdonding derivatives for strip measurements.

The residual is calculated by a straight line propagation of the segment onto the plane & then taking the absolute disances in the [precision] & [non-precision direction]

Parameters
lineCurrent segment line together with the first & second order derivatives
spacePointReference to the measurement to calculate the residual to
residualOutput residual object to save the actual residual and the derivatives

Update the residual accordingly

Definition at line 297 of file MdtSegmentFitter.cxx.

299  {
300 
301 
302  const Amg::Vector3D& hitPos{hit.positionInChamber()};
303  const Amg::Vector3D normal = hit.type() != xAOD::UncalibMeasType::Other
304  ? hit.spacePoint()->planeNormal() : Amg::Vector3D::UnitZ();
305  const double planeOffSet = normal.dot(hitPos);
306 
307  const double normDot = normal.dot(line.dir);
308  if (std::abs(normDot) < std::numeric_limits<double>::epsilon()){
309  residual.residual.setZero();
310  for (Amg::Vector3D& deriv: residual.gradient) {
311  deriv.setZero();
312  }
313  if (m_cfg.useSecOrderDeriv) {
314  for (Amg::Vector3D& hessian : residual.hessian){
315  hessian.setZero();
316  }
317  }
318  ATH_MSG_WARNING("The hit parallel with the segment line "<<Amg::toString(line.dir));
319  return;
320  }
321  const double travelledDist = (planeOffSet - line.pos.dot(normal)) / normDot;
322 
323  const Amg::Vector3D planeIsect = line.pos + travelledDist * line.dir;
325  residual.residual.block<2,1>(0,0) = (planeIsect - hitPos).block<2,1>(0,0);
326 
327  for (unsigned fitPar = 0; fitPar < line.gradient.size(); ++fitPar) {
328  switch (fitPar) {
329  case toInt(ParamDefs::phi):
330  case toInt(ParamDefs::theta):{
331  const double partialDist = - travelledDist / normDot * normal.dot(line.gradient[fitPar]);
332  residual.gradient[fitPar].block<2,1>(0,0) = (travelledDist * line.gradient[fitPar] + partialDist * line.dir).block<2,1>(0,0);
333  break;
334  } case toInt(ParamDefs::y0):
335  case toInt(ParamDefs::x0): {
336  residual.gradient[fitPar].block<2,1>(0,0) = (line.gradient[fitPar] - line.gradient[fitPar].dot(normal) / normDot * line.dir).block<2,1>(0,0);
337  break;
338  }
339  default: {
340  break;
341  }
342  }
343  }
344  if (!m_cfg.useSecOrderDeriv) {
345  return;
346  }
347  constexpr unsigned nLinePars = LineWithPartials::nPars;
348  for (int param = toInt(ParamDefs::phi); param >=0; --param){
349  for (int param1 = param; param1>=0; --param1) {
350  const int lineIdx = vecIdxFromSymMat<nLinePars>(param, param1);
351  const int resIdx = vecIdxFromSymMat<toInt(ParamDefs::nPars)>(param, param1);
352  if ( (param == toInt(ParamDefs::theta) || param == toInt(ParamDefs::phi)) &&
353  (param1 == toInt(ParamDefs::theta) || param1 == toInt(ParamDefs::phi))) {
354  residual.hessian[resIdx].block<2,1>(0,0) =(
355  travelledDist * line.hessian[lineIdx] - travelledDist *(normal.dot(line.hessian[lineIdx])) / normDot * line.dir
356  - normal.dot(line.gradient[param1]) / normDot * residual.gradient[param]
357  - normal.dot(line.gradient[param]) / normDot * residual.gradient[param1]).block<2,1>(0,0);
358  } else if (param == toInt(ParamDefs::theta) || param == toInt(ParamDefs::phi)) {
359  const double gradientDisplace = normal.dot(line.gradient[param1]);
360  if (gradientDisplace > std::numeric_limits<float>::epsilon()){
361  residual.hessian[resIdx].block<2,1>(0,0) = gradientDisplace*( normal.dot(line.gradient[param]) / std::pow(normDot,2)*line.dir
362  - line.gradient[param] / normDot ).block<2,1>(0,0);
363  } else {
364  residual.hessian[resIdx].setZero();
365  }
366  }
367  }
368  }
369  }

◆ calculateWireResiduals()

void MuonR4::MdtSegmentFitter::calculateWireResiduals ( const LineWithPartials line,
const CalibratedSpacePoint spacePoint,
ResidualWithPartials residual 
) const

Calculates the residuals together with the corresponding derivatives for a drift-circle measurement.

The residual is calculated by projecting the hit into the plane perpendicular to the wire. If the hit is a twin-hit the distance along the wire is taken into account as well.

Parameters
lineCurrent segment line together with the first & second order derivatives
spacePointReference to the measurement to calculate the residual to
residualOutput residual object to save the actual residual and the derivatives

Fetch the hit position & direction

Cache the scalar product & the lengths as they're appearing in the formulas below more often

Project the segment onto the plane

Calculate the distance from the two reference points

Distance from the segment line to the tube wire

If the tube is a twin-tube, the hit position is no longer arbitrary along the wire. Calculate the distance along the wire towards the point of closest approach.

Calculate the first derivative of the residual

No variation of the track parameters w.r.t. the time

Loop to include the second order derivatvies

Pure angular derivatives of the residual

Angular & Spatial mixed terms

Definition at line 167 of file MdtSegmentFitter.cxx.

169  {
171  const Amg::Vector3D& hitPos{hit.positionInChamber()};
172  const Amg::Vector3D& hitDir{hit.directionInChamber()};
174  resObj.projIntoWirePlane = line.dir.dot(hitDir);
175  resObj.projDirLenSq = 1. - resObj.projIntoWirePlane*resObj.projIntoWirePlane;
176  if (resObj.projDirLenSq < std::numeric_limits<float>::epsilon()) {
177  resObj.residual.setZero();
178  for (Amg::Vector3D& grad : resObj.gradient){
179  grad.setZero();
180  }
181  if (m_cfg.useSecOrderDeriv) {
182  for (Amg::Vector3D& hess : resObj.hessian) {
183  hess.setZero();
184  }
185  }
186  ATH_MSG_WARNING("Segment is parallel along the wire");
187  return;
188  }
189  resObj.invProjLenSq = 1. / resObj.projDirLenSq;
190  resObj.invProjLen = std::sqrt(resObj.invProjLenSq);
192  resObj.projDir = (line.dir - resObj.projIntoWirePlane*hitDir) * resObj.invProjLen;
194  const Amg::Vector3D hitMinSeg = hitPos - line.pos ;
196  const double lineDist = resObj.projDir.cross(hitDir).dot(hitMinSeg);
197  const double resVal = (lineDist - hit.driftRadius());
198  resObj.residual = resVal * Amg::Vector3D::Unit(toInt(AxisDefs::eta));
199  ATH_MSG_VERBOSE("Mdt drift radius: "<<hit.driftRadius()<<" distance: "<<lineDist<<";"
200  <<Amg::signedDistance(hitPos, hitDir, line.pos, line.dir)<<", residual: "<<resVal);
203  if (hit.dimension() == 2) {
204  resObj.residual[toInt(AxisDefs::phi)] = (hitMinSeg.dot(line.dir)*resObj.projIntoWirePlane - hitMinSeg.dot(hitDir)) * resObj.invProjLenSq;
205  }
206  constexpr unsigned nLinePars = LineWithPartials::nPars;
208  for (const int param : {toInt(ParamDefs::y0), toInt(ParamDefs::x0),
210  if (!resObj.evalPhiPars && (param == toInt(ParamDefs::x0) || param == toInt(ParamDefs::phi))){
211  continue;
212  }
213  switch (param) {
214  case toInt(ParamDefs::theta):
215  case toInt(ParamDefs::phi): {
216  resObj.partWirePlaneProj[param] = line.gradient[param].dot(hitDir);
217  resObj.partProjDir[param] = (line.gradient[param] - resObj.partWirePlaneProj[param]*hitDir) * resObj.invProjLen
218  + resObj.partWirePlaneProj[param]*resObj.projIntoWirePlane* resObj.projDir * resObj.invProjLenSq;
219  const double partialDist = resObj.partProjDir[param].cross(hitDir).dot(hitMinSeg);
220 
221  resObj.gradient[param] = partialDist * Amg::Vector3D::Unit(toInt(AxisDefs::eta));
222  if (hit.dimension() == 2) {
223  resObj.gradient[param][toInt(AxisDefs::phi)] =
224  ( hitMinSeg.dot(line.gradient[param]) * resObj.projIntoWirePlane +
225  hitMinSeg.dot(line.dir)*resObj.partWirePlaneProj[param]) * resObj.invProjLenSq
226  +2.* resObj.residual[toInt(AxisDefs::phi)]*( resObj.projIntoWirePlane * resObj.partWirePlaneProj[param]) * resObj.invProjLenSq;
227  }
228  break;
229  } case toInt(ParamDefs::y0):
230  case toInt(ParamDefs::x0): {
231  const double partialDist = - resObj.projDir.cross(hitDir).dot(line.gradient[param]);
232  resObj.gradient[param] = partialDist * Amg::Vector3D::Unit(toInt(AxisDefs::eta));
233  if (hit.dimension() == 2) {
234  resObj.gradient[param][toInt(AxisDefs::phi)] = -(line.gradient[param].dot(line.dir) * resObj.projIntoWirePlane -
235  line.gradient[param].dot(hitDir)) * resObj.invProjLenSq;
236  }
237  break;
238  }
240  default:
241  break;
242  }
243  }
244  if (!m_cfg.useSecOrderDeriv) {
245  return;
246  }
248  for (int param = toInt(ParamDefs::phi); param >=0; --param){
249  if (!resObj.evalPhiPars && (param == toInt(ParamDefs::x0) || param == toInt(ParamDefs::phi))){
250  continue;
251  }
252  for (int param1 = param; param1>=0; --param1) {
253  if (!resObj.evalPhiPars && (param1 == toInt(ParamDefs::x0) || param1 == toInt(ParamDefs::phi))){
254  continue;
255  }
256  const int lineIdx = vecIdxFromSymMat<nLinePars>(param, param1);
257  const int resIdx = vecIdxFromSymMat<toInt(ParamDefs::nPars)>(param, param1);
259  if ( (param == toInt(ParamDefs::theta) || param == toInt(ParamDefs::phi)) &&
260  (param1 == toInt(ParamDefs::theta) || param1 == toInt(ParamDefs::phi))) {
261 
262  const double partSqLineProject = line.hessian[lineIdx].dot(hitDir);
263  const Amg::Vector3D projDirPartSq = (line.hessian[lineIdx] - partSqLineProject * hitDir) * resObj.invProjLen
264  + (resObj.partWirePlaneProj[param1] * resObj.projIntoWirePlane) * resObj.invProjLenSq * resObj.partProjDir[param]
265  + (resObj.partWirePlaneProj[param] * resObj.projIntoWirePlane) * resObj.invProjLenSq * resObj.partProjDir[param1]
266  + (partSqLineProject*resObj.projIntoWirePlane) * resObj.invProjLenSq * resObj.projDir
267  + (resObj.partWirePlaneProj[param1] * resObj.partWirePlaneProj[param]) * std::pow(resObj.invProjLenSq, 2) * resObj.projDir;
268 
269  const double partialSqDist = projDirPartSq.cross(hitDir).dot(hitMinSeg);
270  resObj.hessian[resIdx] = partialSqDist * Amg::Vector3D::Unit(toInt(AxisDefs::eta));
271  if (hit.dimension() == 2) {
272  const double partialSqAlongWire =2.*resObj.residual[toInt(AxisDefs::phi)]*resObj.projIntoWirePlane*partSqLineProject * resObj.invProjLenSq
273  +2.*resObj.residual[toInt(AxisDefs::phi)]*resObj.partWirePlaneProj[param]*resObj.partWirePlaneProj[param1]*resObj.invProjLenSq
274  +2.*resObj.gradient[param1][toInt(AxisDefs::phi)]*resObj.projIntoWirePlane*resObj.partWirePlaneProj[param]*resObj.invProjLenSq
275  +2.*resObj.gradient[param][toInt(AxisDefs::phi)]*resObj.projIntoWirePlane*resObj.partWirePlaneProj[param1]*resObj.invProjLenSq
276  + hitMinSeg.dot(line.hessian[lineIdx]) *resObj.projIntoWirePlane * resObj.invProjLenSq
277  + hitMinSeg.dot(line.dir)*partSqLineProject * resObj.invProjLenSq
278  + hitMinSeg.dot(line.gradient[param])*resObj.partWirePlaneProj[param1]*resObj.invProjLenSq
279  + hitMinSeg.dot(line.gradient[param1])*resObj.partWirePlaneProj[param]*resObj.invProjLenSq;
280  resObj.hessian[resIdx][toInt(AxisDefs::phi)] = partialSqAlongWire;
281  }
282  }
284  else if (param == toInt(ParamDefs::theta) || param == toInt(ParamDefs::phi)){
285  const double partialSqDist = - resObj.partProjDir[param].cross(hitDir).dot(line.gradient[param1]);
286  resObj.hessian[resIdx] = partialSqDist * Amg::Vector3D::Unit(toInt(AxisDefs::eta));
287  if (hit.dimension() == 2) {
288  const double partialSqAlongWire = -(line.gradient[param1].dot(line.gradient[param])*resObj.projIntoWirePlane +
289  line.gradient[param1].dot(line.dir)*resObj.partWirePlaneProj[param]) * resObj.invProjLenSq
290  + 2.* resObj.gradient[param1][toInt(AxisDefs::phi)]*( resObj.projIntoWirePlane * resObj.partWirePlaneProj[param]) * resObj.invProjLenSq;
291  resObj.hessian[resIdx][toInt(AxisDefs::phi)] = partialSqAlongWire;
292  }
293  }
294  }
295  }
296  }

◆ fitSegment()

SegmentFitResult MuonR4::MdtSegmentFitter::fitSegment ( const EventContext &  ctx,
HitVec &&  calibHits,
const Parameters startPars,
const Amg::Transform3D localToGlobal 
) const

Cache of the partial derivatives of the line parameters w.r.t the fit parameters

Partials of the residual w.r.t. the fit parameters

Update the partial derivatives of the direction vector

First step calibrate the hits

Reset chi2

Loop over the hits to calculate the partial derivatives

If the time is fitted the drift radii need to be updated. That's properly only possible with recalibration

Time residual calculation

Loop over hits is done. Symmetrise the Hessian

Check whether the gradient is already sufficiently small

Pure eta segment fit

In the case that the time is fit & that there are no phi measurements -> compress matrix by swaping the time column with whaever the second column is... it's zero

Subtract 1 degree of freedom to take the time into account

Calculate the chi2 per measurement

Sort the measurements by ascending z

Update the covariance

Definition at line 371 of file MdtSegmentFitter.cxx.

374  {
375 
376  const Muon::IMuonIdHelperSvc* idHelperSvc{calibHits[0]->spacePoint()->msSector()->idHelperSvc()};
378 
379  if (msgLvl(MSG::VERBOSE)) {
380  std::stringstream hitStream{};
381  const auto [startPos, startDir] = makeLine(startPars);
382  for (const HitType& hit : calibHits) {
383  hitStream<<" **** "<<(hit->type() != xAOD::UncalibMeasType::Other ? idHelperSvc->toString(hit->spacePoint()->identify()): "beamspot" )
384  <<" position: "<<Amg::toString(hit->positionInChamber());
385  if (hit->type() == xAOD::UncalibMeasType::MdtDriftCircleType) {
386  hitStream<<", driftRadius: "<<SegmentFitHelpers::driftSign(startPos,startDir, *hit, msg())*hit->driftRadius() ;
387  }
388  hitStream<<", channel dir: "<<Amg::toString(hit->directionInChamber())<<std::endl;
389  }
390  ATH_MSG_VERBOSE("Start segment fit with parameters "<<toString(startPars)
391  <<", plane location: "<<Amg::toString(localToGlobal)<<std::endl<<hitStream.str());
392 
393  }
394 
395  SegmentFitResult fitResult{};
396  fitResult.segmentPars = startPars;
397  fitResult.timeFit = m_cfg.doTimeFit;
398  fitResult.calibMeasurements = std::move(calibHits);
399  if (!updateHitSummary(fitResult)) {
400  ATH_MSG_WARNING("No valid segment seed parsed from the beginning");
401  return fitResult;
402  }
403  if (!m_cfg.reCalibrate) {
404  const auto [segPos, segDir] = makeLine(startPars);
405  updateDriftSigns(segPos, segDir, fitResult);
406  }
407  Parameters gradient{AmgVector(5)::Zero()}, prevGrad{AmgVector(5)::Zero()}, prevPars{AmgVector(5)::Zero()};
408  AmgSymMatrix(5) hessian{AmgSymMatrix(5)::Zero()};
409 
411  LineWithPartials segmentLine{};
413  ResidualWithPartials residual{};
414 
415  unsigned int noChangeIter{0};
416  while (fitResult.nIter++ < m_cfg.nMaxCalls) {
417  ATH_MSG_DEBUG("Iteration: "<<fitResult.nIter<<" parameters: "<<toString(fitResult.segmentPars)<<", chi2: "<<fitResult.chi2);
419  updateLinePartials(fitResult.segmentPars, segmentLine);
420 
422  if (m_cfg.reCalibrate && !recalibrate(ctx, fitResult)) {
423  break;
424  }
426  fitResult.chi2 = 0;
427  hessian.setZero();
428  gradient.setZero();
429 
431  for (HitType& hit : fitResult.calibMeasurements) {
432  if (hit->fitState() != State::Valid) {
433  continue;
434  }
435  const Amg::Vector3D& hitPos{hit->positionInChamber()};
436  const Amg::Vector3D& hitDir{hit->directionInChamber()};
437  ATH_MSG_DEBUG("Update chi2 from measurement "<<(hit->spacePoint() ? idHelperSvc->toString(hit->spacePoint()->identify())
438  : "pseudo meas")<<" position: "<<Amg::toString(hitPos)<<" + "<<Amg::toString(hitDir));
439 
440  const int start = toInt(fitResult.nPhiMeas ? ParamDefs::phi : ParamDefs::theta);
441  switch (hit->type()) {
445  const double dSign = (hit->driftRadius() > 0 ? 1. : -1.);
446  if (fitResult.timeFit && !m_cfg.reCalibrate && fitResult.nIter > 1) {
447  hit = m_cfg.calibrator->calibrate(ctx, *hit, segmentLine.pos, segmentLine.dir,
448  fitResult.segmentPars[toInt(ParamDefs::time)]);
449  hit->setDriftRadius(dSign*hit->driftRadius());
450  }
451  calculateWireResiduals(segmentLine, *hit, residual);
452 
453  if (!fitResult.timeFit) {
454  break;
455  }
457  residual.gradient[toInt(ParamDefs::time)] = dSign*m_cfg.calibrator->driftVelocity(ctx, *hit) * Amg::Vector3D::Unit(toInt(AxisDefs::eta));
458  if (m_cfg.useSecOrderDeriv) {
459  constexpr unsigned hessIdx = vecIdxFromSymMat<ResidualWithPartials::nPars>(toInt(ParamDefs::time), toInt(ParamDefs::time));
460  residual.hessian[hessIdx] = dSign*m_cfg.calibrator->driftAcceleration(ctx, *hit) * Amg::Vector3D::Unit(toInt(AxisDefs::eta));
461  }
462  break;
463  }
468  calculateStripResiduals(segmentLine, *hit, residual);
469  if (!fitResult.timeFit || !hit->measuresTime()) {
470  break;
471  }
472  constexpr ParamDefs par = ParamDefs::time;
473  residual.gradient[toInt(par)] = -Amg::Vector3D::Unit(toInt(AxisDefs::t0));
474  const Amg::Vector3D planeIsect = residual.residual + hitPos;
475  const double totFlightDist = (localToGlobal * planeIsect).mag() * c_inv;
476  residual.residual[toInt(AxisDefs::t0)] = hit->time() - totFlightDist * c_inv - fitResult.segmentPars[toInt(ParamDefs::time)];
477 
478  for (int p = start; p >= 0; --p) {
479  residual.gradient[p][toInt(AxisDefs::t0)] = -residual.gradient[p].perp() * c_inv;
480  ATH_MSG_VERBOSE("Partial derivative of "<<idHelperSvc->toString(hit->spacePoint()->identify())
481  <<" residual "<<Amg::toString(residual.residual)<<" "<<Amg::toString(multiply(inverse(hit->covariance()), residual.residual))
482  <<" "<<std::endl<<toString(hit->covariance())<<std::endl<<" w.r.t "
483  <<toString(static_cast<ParamDefs>(p))<<"="<<Amg::toString(residual.gradient[p]));
484  }
485  break;
487  calculateStripResiduals(segmentLine, *hit, residual);
488  break;
489  } default: {
490  const auto &measurement = *hit->spacePoint()->primaryMeasurement();
491  ATH_MSG_WARNING("MdtSegmentFitter() - Unsupported measurment type" <<typeid(measurement).name());
492  }
493  }
494 
495  ATH_MSG_VERBOSE("Update derivatives for hit "<< (hit->spacePoint() ? idHelperSvc->toString(hit->spacePoint()->identify()) : "beamspot"));
496  updateDerivatives(residual, hit->covariance(), gradient, hessian, fitResult.chi2,
497  fitResult.timeFit && hit->measuresTime() ? toInt(ParamDefs::time) : start);
498  }
499 
501  for (int k =1; k < 5 - (!fitResult.timeFit); ++k){
502  for (int l = 0; l< k; ++l){
503  hessian(l,k) = hessian(k,l);
504  }
505  }
507  if (gradient.mag() < m_cfg.tolerance) {
508  fitResult.converged = true;
509  ATH_MSG_VERBOSE("Fit converged after "<<fitResult.nIter<<" iterations with "<<fitResult.chi2);
510  break;
511  }
512  ATH_MSG_VERBOSE("Chi2: "<<fitResult.chi2<<", gradient: "<<Amg::toString(gradient)<<"hessian: "<<std::endl<<hessian);
515  if (!fitResult.nPhiMeas && !fitResult.timeFit) {
516  paramUpdate = updateParameters<2>(fitResult.segmentPars, prevPars, gradient, prevGrad, hessian);
517  } else if (!fitResult.nPhiMeas && fitResult.timeFit) {
520  hessian.col(2).swap(hessian.col(toInt(ParamDefs::time)));
521  hessian.row(2).swap(hessian.row(toInt(ParamDefs::time)));
522  std::swap(gradient[2], gradient[toInt(ParamDefs::time)]);
523  std::swap(fitResult.segmentPars[2], fitResult.segmentPars[toInt(ParamDefs::time)]);
524 
525  paramUpdate = updateParameters<3>(fitResult.segmentPars, prevPars, gradient, prevGrad, hessian);
526 
527  std::swap(fitResult.segmentPars[2], fitResult.segmentPars[toInt(ParamDefs::time)]);
528  hessian.col(2).swap(hessian.col(toInt(ParamDefs::time)));
529  hessian.row(2).swap(hessian.row(toInt(ParamDefs::time)));
530  std::swap(gradient[2], gradient[toInt(ParamDefs::time)]);
531  } else if (fitResult.nPhiMeas && !fitResult.timeFit) {
532  paramUpdate = updateParameters<4>(fitResult.segmentPars, prevPars, gradient, prevGrad, hessian);
533  } else if (fitResult.nPhiMeas && fitResult.timeFit) {
534  paramUpdate = updateParameters<5>(fitResult.segmentPars, prevPars, gradient, prevGrad, hessian);
535  }
536  if (paramUpdate == UpdateStatus::noChange) {
537  if ((++noChangeIter) >= m_cfg.noMoveIter){
538  fitResult.converged = true;
539  break;
540  }
541  } else if (paramUpdate == UpdateStatus::allOkay) {
542  noChangeIter = 0;
543  } else if (paramUpdate == UpdateStatus::outOfBounds){
544  return fitResult;
545  }
546  }
547 
549  fitResult.nDoF-=fitResult.timeFit;
550 
552  const auto [segPos, segDir] = fitResult.makeLine();
553  fitResult.chi2 =0.;
554 
555  std::optional<double> toF = fitResult.timeFit ? std::make_optional<double>((localToGlobal * segPos).mag() * c_inv) : std::nullopt;
557  std::ranges::stable_sort(fitResult.calibMeasurements, [](const HitType&a, const HitType& b){
558  return a->positionInChamber().z() < b->positionInChamber().z();
559  });
560 
561  /*** Remove the drift sign again */
562  fitResult.chi2 =0.;
563  for (const HitType& hit : fitResult.calibMeasurements) {
564  hit->setDriftRadius(std::abs(hit->driftRadius()));
565  fitResult.chi2 +=SegmentFitHelpers::chiSqTerm(segPos, segDir, fitResult.segmentPars[toInt(ParamDefs::time)], toF, *hit, msg());
566  }
568  if (!fitResult.nPhiMeas&& !fitResult.timeFit) {
569  blockCovariance<2>(std::move(hessian), /* std::move(gradient), std::move(prevGrad), */ fitResult.segmentParErrs);
570  }else if (!fitResult.nPhiMeas && fitResult.timeFit) {
571  hessian.col(2).swap(hessian.col(toInt(ParamDefs::time)));
572  hessian.row(2).swap(hessian.row(toInt(ParamDefs::time)));
573  blockCovariance<3>(std::move(hessian),/* std::move(gradient), std::move(prevGrad), */ fitResult.segmentParErrs);
574  fitResult.segmentParErrs.col(2).swap(fitResult.segmentParErrs.col(toInt(ParamDefs::time)));
575  fitResult.segmentParErrs.row(2).swap(fitResult.segmentParErrs.row(toInt(ParamDefs::time)));
576  } else if (fitResult.nPhiMeas) {
577  blockCovariance<4>(std::move(hessian), /* std::move(gradient), std::move(prevGrad), */ fitResult.segmentParErrs);
578  } else if (fitResult.nPhiMeas && fitResult.timeFit) {
579  blockCovariance<5>(std::move(hessian),/* std::move(gradient), std::move(prevGrad), */ fitResult.segmentParErrs);
580  }
581  return fitResult;
582  }

◆ initMessaging()

void AthMessaging::initMessaging ( ) const
privateinherited

Initialize our message level and MessageSvc.

This method should only be called once.

Definition at line 39 of file AthMessaging.cxx.

40 {
42  m_lvl = m_imsg ?
43  static_cast<MSG::Level>( m_imsg.load()->outputLevel(m_nm) ) :
44  MSG::INFO;
45 }

◆ msg() [1/2]

MsgStream & AthMessaging::msg ( ) const
inlineinherited

The standard message stream.

Returns a reference to the default message stream May not be invoked before sysInitialize() has been invoked.

Definition at line 164 of file AthMessaging.h.

165 {
166  MsgStream* ms = m_msg_tls.get();
167  if (!ms) {
168  if (!m_initialized.test_and_set()) initMessaging();
169  ms = new MsgStream(m_imsg,m_nm);
170  m_msg_tls.reset( ms );
171  }
172 
173  ms->setLevel (m_lvl);
174  return *ms;
175 }

◆ msg() [2/2]

MsgStream & AthMessaging::msg ( const MSG::Level  lvl) const
inlineinherited

The standard message stream.

Returns a reference to the default message stream May not be invoked before sysInitialize() has been invoked.

Definition at line 179 of file AthMessaging.h.

180 { return msg() << lvl; }

◆ msgLvl()

bool AthMessaging::msgLvl ( const MSG::Level  lvl) const
inlineinherited

Test the output level.

Parameters
lvlThe message level to test against
Returns
boolean Indicating if messages at given level will be printed
Return values
trueMessages at level "lvl" will be printed

Definition at line 151 of file AthMessaging.h.

152 {
153  if (!m_initialized.test_and_set()) initMessaging();
154  if (m_lvl <= lvl) {
155  msg() << lvl;
156  return true;
157  } else {
158  return false;
159  }
160 }

◆ recalibrate()

bool MuonR4::MdtSegmentFitter::recalibrate ( const EventContext &  ctx,
SegmentFitResult fitResult 
) const
inlineprivate

Recalibrate the measurements participating in the fit & shift them into the centre-of gravity frame.

The parsed centerOfGravity position is updated to be the average measurement position weighted by the inverse covariance. Returns true if there're enough degrees of freedom left to fit

Parameters
ctxEventContext to fetch the calibration constants
fitResultCurrent fit state
centerOfGravityVector to save the center of gravity position later on.
invCovVector to cache the inverse covariances of the problem.

Switch off the time fit if too little degrees of freedom are left

Recalibrate the measurements

Definition at line 62 of file MdtSegmentFitter.cxx.

63  {
64 
65  const auto [segPos, segDir] = makeLine(fitResult.segmentPars);
66 
67  fitResult.calibMeasurements = m_cfg.calibrator->calibrate(ctx, std::move(fitResult.calibMeasurements), segPos, segDir,
68  fitResult.segmentPars[toInt(ParamDefs::time)]);
69  if (!updateHitSummary(fitResult)){
70  return false;
71  }
72  updateDriftSigns(segPos, segDir, fitResult);
74  if (fitResult.timeFit && fitResult.nDoF <= 1) {
75  fitResult.timeFit = false;
76  ATH_MSG_DEBUG("Switch of the time fit because nDoF: "<<fitResult.nDoF);
77  fitResult.segmentPars[toInt(ParamDefs::time)] = 0.;
79  fitResult.calibMeasurements = m_cfg.calibrator->calibrate(ctx, std::move(fitResult.calibMeasurements), segPos, segDir,
80  fitResult.segmentPars[toInt(ParamDefs::time)]);
82  } else if (!fitResult.timeFit && m_cfg.doTimeFit) {
83  ATH_MSG_DEBUG("Somehow a measurement is on the narrow ridge of validity. Let's try if the time can be fitted now ");
84  fitResult.timeFit = true;
85  }
86 
87  return true;
88  }

◆ setLevel()

void AthMessaging::setLevel ( MSG::Level  lvl)
inherited

Change the current logging level.

Use this rather than msg().setLevel() for proper operation with MT.

Definition at line 28 of file AthMessaging.cxx.

29 {
30  m_lvl = lvl;
31 }

◆ updateDerivatives()

void MuonR4::MdtSegmentFitter::updateDerivatives ( const ResidualWithPartials fitMeas,
const MeasCov_t measCovariance,
AmgVector(5)&  gradient,
AmgSymMatrix(5)&  hessian,
double &  chi2,
int  startPar 
) const
private

Definition at line 584 of file MdtSegmentFitter.cxx.

587  {
588  const MeasCov_t invCov{inverse(covariance)};
589  const Amg::Vector3D covRes = multiply(invCov, fitMeas.residual);
590  chi2 += covRes.dot(fitMeas.residual);
591  for (int p = startPar; p >=0 ; --p) {
592  gradient[p] +=2.*covRes.dot(fitMeas.gradient[p]);
593  for (int k=p; k>=0; --k) {
594  hessian(p,k)+= 2.*contract(invCov, fitMeas.gradient[p], fitMeas.gradient[k]);
595  if (m_cfg.useSecOrderDeriv) {
596  const int symMatIdx = vecIdxFromSymMat<ResidualWithPartials::nPars>(p,k);
597  hessian(p,k)+=contract(invCov, covRes, fitMeas.hessian[symMatIdx]);
598  }
599  }
600  }
601  ATH_MSG_VERBOSE("After derivative update --- chi2: "<<chi2<<"("<<covRes.dot(fitMeas.residual)<<"), gradient: "
602  <<toString(gradient)<<", Hessian:\n"<<hessian<<", measurement covariance\n"<<toString(invCov));
603  }

◆ updateDriftSigns()

void MuonR4::MdtSegmentFitter::updateDriftSigns ( const Amg::Vector3D segPos,
const Amg::Vector3D segDir,
SegmentFitResult fitRes 
) const
inlineprivate

Update the signs of the measurement.

Definition at line 55 of file MdtSegmentFitter.cxx.

56  {
57  for (std::unique_ptr<CalibratedSpacePoint>& meas : fitResult.calibMeasurements) {
58  meas->setDriftRadius(SegmentFitHelpers::driftSign(segPos,segDir, *meas, msg())*meas->driftRadius());
59  }
60  }

◆ updateHitSummary()

bool MuonR4::MdtSegmentFitter::updateHitSummary ( SegmentFitResult fitResult) const
inlineprivate

Brief updates the hit summary from the contributing hits.

Returns fals if there're too little valid hits.

Parameters
fitResultReference to the container containing all the measurements

Count the phi & time measurements measurements

Mdts are already counted in the measures eta category. Don't count them twice

Definition at line 89 of file MdtSegmentFitter.cxx.

89  {
92  fitResult.nPhiMeas = fitResult.nDoF = fitResult.nTimeMeas = 0;
93  for (const HitType& hit : fitResult.calibMeasurements) {
94  if (hit->fitState() != State::Valid){
95  continue;
96  }
97 
98  fitResult.nPhiMeas+= hit->measuresPhi();
99  fitResult.nDoF+= hit->measuresPhi();
100  fitResult.nDoF+= hit->measuresEta();
101  fitResult.nPrecMeas+= (hit->type() == xAOD::UncalibMeasType::MdtDriftCircleType);
103  fitResult.nDoF += (m_cfg.doTimeFit && hit->type() != xAOD::UncalibMeasType::MdtDriftCircleType && hit->measuresTime());
104  fitResult.nTimeMeas+=hit->measuresTime();
105  }
106  if (!fitResult.nDoF) {
107  ATH_MSG_VERBOSE("Measurements rejected.");
108  return false;
109  }
110  if (!fitResult.nPhiMeas) {
111  ATH_MSG_VERBOSE("No phi measurements are left.");
112  fitResult.segmentPars[toInt(ParamDefs::phi)] = 90. * Gaudi::Units::deg;
113  const double nHits = fitResult.calibMeasurements.size();
114  double avgX{0.};
115  std::ranges::for_each(fitResult.calibMeasurements,[&avgX, nHits](const HitType& hit) {
116  return avgX += hit->positionInChamber().x() / nHits;
117  });
118  fitResult.segmentPars[toInt(ParamDefs::x0)] = avgX;
119  }
120 
121  fitResult.nDoF = fitResult.nDoF - 2 - (fitResult.nPhiMeas > 0 ? 2 : 0);
122 
123  return true;
124  }

◆ updateLinePartials()

void MuonR4::MdtSegmentFitter::updateLinePartials ( const Parameters fitPars,
LineWithPartials line 
)
static

Updates the line parameters together with its first & second order partial derivatives.

Parameters
fitParsSet of segment parameters in the current iteration
lineReference to the line object to be updated

x_{0} cos (phi) sin (theta) segPos = y_{0} , segDir = sin (phi) sin (theta) 0 cos theta

d segDir cos (phi) cos(theta) dSegDir -sin(phi) sin (theta) -------—= sin (phi) cos(theta) -----— = cos(phi) sin (theta) dTheta - sin (theta) dPhi 0

Definition at line 125 of file MdtSegmentFitter.cxx.

126  {
127 
132  line.gradient[toInt(ParamDefs::x0)] = Amg::Vector3D::UnitX();
133  line.gradient[toInt(ParamDefs::y0)] = Amg::Vector3D::UnitY();
146  line.gradient[toInt(ParamDefs::theta)] = Amg::Vector3D{phi.cs*theta.cs, phi.sn*theta.cs, -theta.sn};
147  line.gradient[toInt(ParamDefs::phi)] = Amg::Vector3D{-theta.sn *phi.sn, theta.sn*phi.cs, 0};
148  /*********************************************************************************
149  * Non-vanishing second order derivatives
150  *
151  * d^{2} segDir d^{2} segDir cos phi
152  * ------------- = - segDir , ------------ = - sin(theta) sin phi
153  * d^{2} theta d^{2} phi 0
154  *
155  * d^{2} segDir -sin phi
156  * ------------- = cos(theta) cos phi
157  * d theta dPhi 0
158  ************************************************************************************/
159  constexpr unsigned nPars = LineWithPartials::nPars;
160  constexpr int idxThetaSq = vecIdxFromSymMat<nPars>(toInt(ParamDefs::theta), toInt(ParamDefs::theta));
161  constexpr int idxPhiSq = vecIdxFromSymMat<nPars>(toInt(ParamDefs::phi), toInt(ParamDefs::phi));
162  constexpr int idxPhiTheta = vecIdxFromSymMat<nPars>(toInt(ParamDefs::theta), toInt(ParamDefs::phi));
163  line.hessian[idxThetaSq] = - Amg::Vector3D{phi.cs*theta.sn,phi.sn*theta.sn, theta.cs};
164  line.hessian[idxPhiSq] = -theta.sn * Amg::Vector3D{phi.cs, phi.sn, 0.};
165  line.hessian[idxPhiTheta] = theta.cs * Amg::Vector3D{-phi.sn, phi.cs, 0.};
166  }

◆ updateParameters()

template<unsigned int nDim>
MdtSegmentFitter::UpdateStatus MuonR4::MdtSegmentFitter::updateParameters ( Parameters currentPars,
Parameters previousPars,
Parameters currGrad,
Parameters prevGrad,
const AmgSymMatrix(5)&  hessian 
) const
private

Update step of the segment parameters using the Hessian and the gradient.

If the Hessian is definite, the currentParameters are updated according to x_{n+1} = x_{n} - H_{n}^{1} * grad(x_{n}) Otherwise, the method of steepest descent is attempted.

Parameters
currentParsBest segment estimator parameters
previousParsSegment estimator parameters from the last iteration
currGradGradient of the chi2 from this iteration
prevGradGradient of the chi2 from the previous iteration
hessianHessian estimator

Check that all parameters remain within the parameter boundary window

Recall that for 3x3 the dimensions are [y0, theta, time]

Definition at line 622 of file MdtSegmentFitter.cxx.

624  {
625 
626  AmgSymMatrix(nDim) miniHessian = currentHessian.block<nDim, nDim>(0,0);
627  ATH_MSG_DEBUG("Parameter update -- \ncurrenPars: "<<toString(currPars)<<", \ngradient: "<<toString(currGrad)
628  <<", hessian ("<<miniHessian.determinant()<<")"<<std::endl<<miniHessian);
629 
630  double changeMag{0.};
631  if (std::abs(miniHessian.determinant()) > detCutOff) {
632  prevPars.block<nDim,1>(0,0) = currPars.block<nDim,1>(0,0);
633  // Update the parameters accrodingly to the hessian
634  const AmgVector(nDim) updateMe = miniHessian.inverse()* currGrad.block<nDim, 1>(0,0);
635  changeMag = std::sqrt(updateMe.dot(updateMe));
636  currPars.block<nDim,1>(0,0) -= updateMe;
637  prevGrad.block<nDim,1>(0,0) = currGrad.block<nDim,1>(0,0);
638  ATH_MSG_DEBUG("Hessian inverse:\n"<<miniHessian.inverse()<<"\n\nUpdate the parameters by -"
639  <<Amg::toString(updateMe));
640  } else {
641  const AmgVector(nDim) gradDiff = (currGrad - prevGrad).block<nDim,1>(0,0);
642  const double gradDiffMag = gradDiff.mag2();
643  double denom = (gradDiffMag > std::numeric_limits<float>::epsilon() ? gradDiffMag : 1.);
644  const double gamma = std::abs((currPars - prevPars).block<nDim,1>(0,0).dot(gradDiff))
645  / denom;
646  ATH_MSG_VERBOSE("Hessian determinant invalid. Try deepest descent - \nprev parameters: "
647  <<toString(prevPars)<<",\nprevious gradient: "<<toString(prevGrad)<<", gamma: "<<gamma);
648  prevPars.block<nDim, 1>(0,0) = currPars.block<nDim, 1>(0,0);
649  currPars.block<nDim, 1>(0,0) -= gamma* currGrad.block<nDim, 1>(0,0);
650  prevGrad.block<nDim, 1>(0,0) = currGrad.block<nDim,1>(0,0);
651  changeMag = std::abs(gamma) * currGrad.block<nDim, 1>(0,0).mag();
652  }
654  unsigned int nOutOfBound{0};
655  for (unsigned int p = 0; p< nDim; ++p) {
656  double& parValue{currPars[p]};
657  if constexpr(nDim == 3) {
659  if (p == 2) {
661  }
662  }
663  if (m_cfg.ranges[p][0] > parValue || m_cfg.ranges[p][1] < parValue) {
664  ATH_MSG_VERBOSE("The "<<p<<"-th parameter "<<toString(static_cast<ParamDefs>(p))<<" is out of range "<<parValue
665  <<" allowed ["<<m_cfg.ranges[p][0]<<"-"<<m_cfg.ranges[p][1]<<"]");
666  ++nOutOfBound;
667  parValue = std::clamp(parValue, m_cfg.ranges[p][0], m_cfg.ranges[p][1]);
668  }
669 
670  }
671  if (nOutOfBound > m_cfg.nParsOutOfBounds){
673  }
674  if (changeMag <= m_cfg.tolerance) {
675  return UpdateStatus::noChange;
676  }
677  return UpdateStatus::allOkay;
678  }

Member Data Documentation

◆ ATLAS_THREAD_SAFE

std::atomic_flag m_initialized AthMessaging::ATLAS_THREAD_SAFE = ATOMIC_FLAG_INIT
mutableprivateinherited

Messaging initialized (initMessaging)

Definition at line 141 of file AthMessaging.h.

◆ m_cfg

Config MuonR4::MdtSegmentFitter::m_cfg {}
private

Definition at line 111 of file MdtSegmentFitter.h.

◆ m_imsg

std::atomic<IMessageSvc*> AthMessaging::m_imsg { nullptr }
mutableprivateinherited

MessageSvc pointer.

Definition at line 135 of file AthMessaging.h.

◆ m_lvl

std::atomic<MSG::Level> AthMessaging::m_lvl { MSG::NIL }
mutableprivateinherited

Current logging level.

Definition at line 138 of file AthMessaging.h.

◆ m_msg_tls

boost::thread_specific_ptr<MsgStream> AthMessaging::m_msg_tls
mutableprivateinherited

MsgStream instance (a std::cout like with print-out levels)

Definition at line 132 of file AthMessaging.h.

◆ m_nm

std::string AthMessaging::m_nm
privateinherited

Message source name.

Definition at line 129 of file AthMessaging.h.


The documentation for this class was generated from the following files:
AthMessaging::m_lvl
std::atomic< MSG::Level > m_lvl
Current logging level.
Definition: AthMessaging.h:138
MuonR4::MdtSegmentFitter::updateDriftSigns
void updateDriftSigns(const Amg::Vector3D &segPos, const Amg::Vector3D &segDir, SegmentFitResult &fitRes) const
Update the signs of the measurement.
Definition: MdtSegmentFitter.cxx:55
MuonR4::MdtSegmentFitter::updateHitSummary
bool updateHitSummary(SegmentFitResult &fitResult) const
Brief updates the hit summary from the contributing hits.
Definition: MdtSegmentFitter.cxx:89
MuonR4::MdtSegmentFitter::calculateStripResiduals
void calculateStripResiduals(const LineWithPartials &line, const CalibratedSpacePoint &spacePoint, ResidualWithPartials &residual) const
Calculates the residual together with hte correspdonding derivatives for strip measurements.
Definition: MdtSegmentFitter.cxx:297
MuonR4::MdtSegmentFitter::LineWithPartials::nPars
static constexpr unsigned nPars
Free parameters of the line (x0,y0,theta,phi)
Definition: MdtSegmentFitter.h:56
MuonR4::MdtSegmentFitter::UpdateStatus
UpdateStatus
Definition: MdtSegmentFitter.h:172
ClusterSeg::residual
@ residual
Definition: ClusterNtuple.h:20
MuonR4::multiply
Amg::Vector3D multiply(const CalibratedSpacePoint::Covariance_t &mat, const Amg::Vector3D &v)
Multiplies a 3D vector with the covariance matrix which can be either 2x2 or 3x3 matrix.
Definition: MuonSpectrometer/MuonPhaseII/Event/MuonSpacePoint/src/UtilFunctions.cxx:40
MuonR4::MdtSegmentFitter::Parameters
SegmentFit::Parameters Parameters
Definition: MdtSegmentFitter.h:22
mergePhysValFiles.start
start
Definition: DataQuality/DataQualityUtils/scripts/mergePhysValFiles.py:13
MuonR4::MdtSegmentFitter::UpdateStatus::outOfBounds
@ outOfBounds
Amg::y
@ y
Definition: GeoPrimitives.h:35
xAOD::UncalibMeasType::MMClusterType
@ MMClusterType
deg
#define deg
Definition: SbPolyhedron.cxx:17
UploadAMITag.l
list l
Definition: UploadAMITag.larcaf.py:157
MuonR4::ISpacePointCalibrator::calibrate
virtual CalibSpacePointPtr calibrate(const EventContext &ctx, const SpacePoint *spacePoint, const Amg::Vector3D &seedPosInChamb, const Amg::Vector3D &seedDirInChamb, const double timeDelay) const =0
Calibrates a single space point.
ATH_MSG_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition: AthMsgStreamMacros.h:28
MuonR4::ISpacePointCalibrator::driftVelocity
virtual double driftVelocity(const EventContext &ctx, const CalibratedSpacePoint &spacePoint) const =0
Returns the drift velocity for a given drift-circle space point.
xAOD::UncalibMeasType::sTgcStripType
@ sTgcStripType
MuonR4::MdtSegmentFitter::Config::noMoveIter
unsigned int noMoveIter
How many iterations with changes below tolerance.
Definition: MdtSegmentFitter.h:42
MuonR4::SegmentFitHelpers::chiSqTerm
double chiSqTerm(const Amg::Vector3D &posInChamber, const Amg::Vector3D &dirInChamber, const SpacePoint &measurement, MsgStream &msg)
Calculates the chi2 contribuation to a linear segment line from an uncalibrated measurement.
Definition: SegmentFitHelperFunctions.cxx:28
MuonR4::SegmentFit::makeLine
std::pair< Amg::Vector3D, Amg::Vector3D > makeLine(const Parameters &pars)
Returns the parsed parameters into an Eigen line parametrization.
Definition: SegmentFitterEventData.cxx:26
AthMessaging::m_imsg
std::atomic< IMessageSvc * > m_imsg
MessageSvc pointer.
Definition: AthMessaging.h:135
MuonR4::MdtSegmentFitter::Config::tolerance
double tolerance
Gradient cut off.
Definition: MdtSegmentFitter.h:30
MuonR4::toString
std::string toString(const CalibratedSpacePoint::Covariance_t &mat)
Returns the matrix in string.
Definition: MuonSpectrometer/MuonPhaseII/Event/MuonSpacePoint/src/UtilFunctions.cxx:75
dq_defect_bulk_create_defects.line
line
Definition: dq_defect_bulk_create_defects.py:27
Athena::getMessageSvc
IMessageSvc * getMessageSvc(bool quiet=false)
Definition: getMessageSvc.cxx:20
config
Definition: PhysicsAnalysis/AnalysisCommon/AssociationUtils/python/config.py:1
MuonR4::SegmentFit::ParamDefs::phi
@ phi
xAOD::UncalibMeasType::TgcStripType
@ TgcStripType
MuonR4::MdtSegmentFitter::Config::calibrator
const ISpacePointCalibrator * calibrator
Pointer to the calibrator tool.
Definition: MdtSegmentFitter.h:40
xAOD::phi
setEt phi
Definition: TrigEMCluster_v1.cxx:29
MuonR4::MdtSegmentFitter::Config::ranges
RangeArray ranges
Definition: MdtSegmentFitter.h:48
MuonR4::SegmentFit::AxisDefs::phi
@ phi
MuonR4::SegmentFitHelpers::driftSign
int driftSign(const Amg::Vector3D &posInChamber, const Amg::Vector3D &dirInChamber, const SpacePoint &uncalibHit, MsgStream &msg)
Calculates whether a segement line travereses the tube measurement on the left (-1) or right (1) side...
Definition: SegmentFitHelperFunctions.cxx:239
MuonR4::SegmentFit::AxisDefs::t0
@ t0
python.utils.AtlRunQueryDQUtils.p
p
Definition: AtlRunQueryDQUtils.py:209
Amg::toString
std::string toString(const Translation3D &translation, int precision=4)
GeoPrimitvesToStringConverter.
Definition: GeoPrimitivesToStringConverter.h:40
TrigConf::MSGTC::Level
Level
Definition: Trigger/TrigConfiguration/TrigConfBase/TrigConfBase/MsgStream.h:21
MuonR4::ISpacePointCalibrator::driftAcceleration
virtual double driftAcceleration(const EventContext &ctx, const CalibratedSpacePoint &spacePoint) const =0
Returns the drift acceleration for a given drift-circle space point.
MuonR4::MdtSegmentFitter::UpdateStatus::noChange
@ noChange
AthMessaging::msgLvl
bool msgLvl(const MSG::Level lvl) const
Test the output level.
Definition: AthMessaging.h:151
TrigVtx::gamma
@ gamma
Definition: TrigParticleTable.h:27
MuonR4::inverse
CalibratedSpacePoint::Covariance_t inverse(const CalibratedSpacePoint::Covariance_t &mat)
Inverts the parsed matrix.
Definition: MuonSpectrometer/MuonPhaseII/Event/MuonSpacePoint/src/UtilFunctions.cxx:65
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
MuonR4::MdtSegmentFitter::MeasCov_t
CalibratedSpacePoint::Covariance_t MeasCov_t
Updates the chi2, its Gradient & Hessian from the measurement residual.
Definition: MdtSegmentFitter.h:163
AmgVector
AmgVector(4) T2BSTrackFilterTool
Definition: T2BSTrackFilterTool.cxx:114
Amg::x
@ x
Definition: GeoPrimitives.h:34
MuonR4::State
CalibratedSpacePoint::State State
Definition: SpacePointCalibrator.cxx:33
chi2
double chi2(TH1 *h0, TH1 *h1)
Definition: comparitor.cxx:525
MuonR4::CalibratedSpacePoint::State
State
State flag to distinguish different space point states.
Definition: CalibratedSpacePoint.h:24
MuonR4::MdtSegmentFitter::Config::useSecOrderDeriv
bool useSecOrderDeriv
Switch toggling whether the second order derivative shall be included.
Definition: MdtSegmentFitter.h:36
AthMessaging
Class to provide easy MsgStream access and capabilities.
Definition: AthMessaging.h:55
MuonR4::MdtSegmentFitter::Config::nParsOutOfBounds
unsigned int nParsOutOfBounds
Abort the fit as soon as more than n parameters leave the fit range.
Definition: MdtSegmentFitter.h:38
xAOD::Other
@ Other
MuonR4::SegmentFit::ParamDefs::x0
@ x0
jobOption.theta
theta
Definition: jobOption.ParticleGun_fwd_sequence.py:13
MuonR4::MdtSegmentFitter::calculateWireResiduals
void calculateWireResiduals(const LineWithPartials &line, const CalibratedSpacePoint &spacePoint, ResidualWithPartials &residual) const
Calculates the residuals together with the corresponding derivatives for a drift-circle measurement.
Definition: MdtSegmentFitter.cxx:167
WriteCalibToCool.swap
swap
Definition: WriteCalibToCool.py:94
HitType
HitType
Definition: FPGATrackSimHit.h:38
MuonR4::SegmentFit::AxisDefs::eta
@ eta
dot.dot
def dot(G, fn, nodesToHighlight=[])
Definition: dot.py:5
AthMessaging::msg
MsgStream & msg() const
The standard message stream.
Definition: AthMessaging.h:164
compute_lumi.denom
denom
Definition: compute_lumi.py:76
MuonR4::SegmentFit::ParamDefs::time
@ time
MuonR4::CalibratedSpacePoint::State::Valid
@ Valid
MuonR4::MdtSegmentFitter::updateDerivatives
void updateDerivatives(const ResidualWithPartials &fitMeas, const MeasCov_t &measCovariance, AmgVector(5)&gradient, AmgSymMatrix(5)&hessian, double &chi2, int startPar) const
Definition: MdtSegmentFitter.cxx:584
MuonR4::MdtSegmentFitter::m_cfg
Config m_cfg
Definition: MdtSegmentFitter.h:111
MuonR4::SegmentFit::ParamDefs::y0
@ y0
MuonR4::MdtSegmentFitter::Config::nMaxCalls
unsigned int nMaxCalls
How many calls shall be executed.
Definition: MdtSegmentFitter.h:28
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:240
createCoolChannelIdFile.par
par
Definition: createCoolChannelIdFile.py:28
plotBeamSpotMon.b
b
Definition: plotBeamSpotMon.py:76
MuonR4::SegmentFit::toInt
constexpr int toInt(const ParamDefs p)
Definition: MuonHoughDefs.h:42
MuonR4::MdtSegmentFitter::Config::reCalibrate
bool reCalibrate
Switch toggling whether the calibrator shall be called at each iteration.
Definition: MdtSegmentFitter.h:34
MuonR4::MdtSegmentFitter::UpdateStatus::allOkay
@ allOkay
Amg::Vector3D
Eigen::Matrix< double, 3, 1 > Vector3D
Definition: GeoPrimitives.h:47
Amg::dirFromAngles
Amg::Vector3D dirFromAngles(const double phi, const double theta)
Constructs a direction vector from the azimuthal & polar angles.
Definition: GeoPrimitivesHelpers.h:299
a
TList * a
Definition: liststreamerinfos.cxx:10
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition: AthMsgStreamMacros.h:32
python.Constants.INFO
int INFO
Definition: Control/AthenaCommon/python/Constants.py:15
MuonR4::MdtSegmentFitter::Config::doTimeFit
bool doTimeFit
Switch toggling whether the T0 shall be fitted or not.
Definition: MdtSegmentFitter.h:32
AthMessaging::m_nm
std::string m_nm
Message source name.
Definition: AthMessaging.h:129
CxxUtils::sincos
Helper to simultaneously calculate sin and cos of the same angle.
Definition: sincos.h:39
Muon::IMuonIdHelperSvc
Interface for Helper service that creates muon Identifiers and can be used to print Identifiers.
Definition: IMuonIdHelperSvc.h:27
Amg::signedDistance
double signedDistance(const Amg::Vector3D &posA, const Amg::Vector3D &dirA, const Amg::Vector3D &posB, const Amg::Vector3D &dirB)
Calculates the signed distance between two lines in 3D space.
Definition: GeoPrimitivesHelpers.h:329
python.Constants.VERBOSE
int VERBOSE
Definition: Control/AthenaCommon/python/Constants.py:13
AthMessaging::initMessaging
void initMessaging() const
Initialize our message level and MessageSvc.
Definition: AthMessaging.cxx:39
PowhegControl_ttFCNC_NLO.params
params
Definition: PowhegControl_ttFCNC_NLO.py:226
AthMessaging::m_msg_tls
boost::thread_specific_ptr< MsgStream > m_msg_tls
MsgStream instance (a std::cout like with print-out levels)
Definition: AthMessaging.h:132
MuonR4::AmgSymMatrix
const AmgSymMatrix(2) &SpacePoint
Definition: MuonSpectrometer/MuonPhaseII/Event/MuonSpacePoint/src/SpacePoint.cxx:90
pow
constexpr int pow(int base, int exp) noexcept
Definition: ap_fixedTest.cxx:15
MuonR4::MdtSegmentFitter::recalibrate
bool recalibrate(const EventContext &ctx, SegmentFitResult &fitResult) const
Recalibrate the measurements participating in the fit & shift them into the centre-of gravity frame.
Definition: MdtSegmentFitter.cxx:62
MuonR4::SegmentFit::ParamDefs::nPars
@ nPars
xAOD::UncalibMeasType::RpcStripType
@ RpcStripType
xAOD::UncalibMeasType::MdtDriftCircleType
@ MdtDriftCircleType
mag
Scalar mag() const
mag method
Definition: AmgMatrixBasePlugin.h:26
MuonR4::MdtSegmentFitter::ParamDefs
SegmentFit::ParamDefs ParamDefs
Definition: MdtSegmentFitter.h:21
MuonR4::SegmentFit::ParamDefs::theta
@ theta
fitman.k
k
Definition: fitman.py:528
generate::Zero
void Zero(TH1D *hin)
Definition: generate.cxx:32
python.SystemOfUnits.ms
float ms
Definition: SystemOfUnits.py:148
MuonR4::MdtSegmentFitter::updateLinePartials
static void updateLinePartials(const Parameters &fitPars, LineWithPartials &line)
Updates the line parameters together with its first & second order partial derivatives.
Definition: MdtSegmentFitter.cxx:125
MuonR4::contract
double contract(const CalibratedSpacePoint::Covariance_t &mat, const Amg::Vector3D &a, const Amg::Vector3D &b)
Definition: MuonSpectrometer/MuonPhaseII/Event/MuonSpacePoint/src/UtilFunctions.cxx:13