 |
ATLAS Offline Software
|
#include <MdtSegmentFitter.h>
|
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 ¤tPars, 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...
|
|
|
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...
|
|
Definition at line 19 of file MdtSegmentFitter.h.
◆ HitType
◆ HitVec
◆ 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
-
fitMeas | Residual of the prediction from the measurement together with its derivatives |
measCovariance | Covariance matrix of the measurement |
gradient | Reference to the gradient of the chi2 which is updated by the function |
hessian | Reference to the Hessian matrix of the chi2 which is updated by the function |
chi2 | Reference to the chi2 value itself which is updated by the function |
startPar | Number of the maximum fit parameter to consider |
Definition at line 163 of file MdtSegmentFitter.h.
◆ ParamDefs
◆ Parameters
◆ UpdateStatus
◆ MdtSegmentFitter()
MuonR4::MdtSegmentFitter::MdtSegmentFitter |
( |
const std::string & |
name, |
|
|
Config && |
config |
|
) |
| |
Standard constructor.
- Parameters
-
name | Name to be printed in the messaging |
config | Fit configuration parameters |
Definition at line 51 of file MdtSegmentFitter.cxx.
◆ blockCovariance()
template<unsigned int nDim>
Definition at line 606 of file MdtSegmentFitter.cxx.
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);
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();
◆ calculateStripResiduals()
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
-
line | Current segment line together with the first & second order derivatives |
spacePoint | Reference to the measurement to calculate the residual to |
residual | Output residual object to save the actual residual and the derivatives |
Update the residual accordingly
Definition at line 297 of file MdtSegmentFitter.cxx.
304 ? hit.spacePoint()->planeNormal() : Amg::Vector3D::UnitZ();
305 const double planeOffSet = normal.dot(hitPos);
307 const double normDot = normal.dot(
line.dir);
308 if (std::abs(normDot) < std::numeric_limits<double>::epsilon()){
321 const double travelledDist = (planeOffSet -
line.pos.dot(normal)) / normDot;
325 residual.residual.block<2,1>(0,0) = (planeIsect - hitPos).block<2,1>(0,0);
327 for (
unsigned fitPar = 0; fitPar <
line.gradient.size(); ++fitPar) {
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);
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);
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);
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);
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);
◆ calculateWireResiduals()
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
-
line | Current segment line together with the first & second order derivatives |
spacePoint | Reference to the measurement to calculate the residual to |
residual | Output 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.
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();
189 resObj.invProjLenSq = 1. / resObj.projDirLenSq;
190 resObj.invProjLen = std::sqrt(resObj.invProjLenSq);
192 resObj.projDir = (
line.dir - resObj.projIntoWirePlane*hitDir) * resObj.invProjLen;
196 const double lineDist = resObj.projDir.cross(hitDir).dot(hitMinSeg);
197 const double resVal = (lineDist - hit.driftRadius());
199 ATH_MSG_VERBOSE(
"Mdt drift radius: "<<hit.driftRadius()<<
" distance: "<<lineDist<<
";"
203 if (hit.dimension() == 2) {
204 resObj.residual[
toInt(
AxisDefs::phi)] = (hitMinSeg.dot(
line.dir)*resObj.projIntoWirePlane - hitMinSeg.dot(hitDir)) * resObj.invProjLenSq;
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);
222 if (hit.dimension() == 2) {
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;
231 const double partialDist = - resObj.projDir.cross(hitDir).dot(
line.gradient[param]);
233 if (hit.dimension() == 2) {
235 line.gradient[param].dot(hitDir)) * resObj.invProjLenSq;
252 for (
int param1 = param; param1>=0; --param1) {
256 const int lineIdx = vecIdxFromSymMat<nLinePars>(param, param1);
257 const int resIdx = vecIdxFromSymMat<toInt(ParamDefs::nPars)>(param, param1);
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;
269 const double partialSqDist = projDirPartSq.cross(hitDir).dot(hitMinSeg);
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;
285 const double partialSqDist = - resObj.partProjDir[param].cross(hitDir).dot(
line.gradient[param1]);
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;
◆ fitSegment()
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.
380 std::stringstream hitStream{};
381 const auto [startPos, startDir] =
makeLine(startPars);
382 for (
const HitType& hit : calibHits) {
388 hitStream<<
", channel dir: "<<
Amg::toString(hit->directionInChamber())<<std::endl;
391 <<
", plane location: "<<
Amg::toString(localToGlobal)<<std::endl<<hitStream.str());
395 SegmentFitResult fitResult{};
396 fitResult.segmentPars = startPars;
398 fitResult.calibMeasurements = std::move(calibHits);
404 const auto [segPos, segDir] =
makeLine(startPars);
411 LineWithPartials segmentLine{};
415 unsigned int noChangeIter{0};
417 ATH_MSG_DEBUG(
"Iteration: "<<fitResult.nIter<<
" parameters: "<<
toString(fitResult.segmentPars)<<
", chi2: "<<fitResult.chi2);
431 for (
HitType& hit : fitResult.calibMeasurements) {
437 ATH_MSG_DEBUG(
"Update chi2 from measurement "<<(hit->spacePoint() ? idHelperSvc->toString(hit->spacePoint()->identify())
441 switch (hit->type()) {
445 const double dSign = (hit->driftRadius() > 0 ? 1. : -1.);
449 hit->setDriftRadius(dSign*hit->driftRadius());
453 if (!fitResult.timeFit) {
469 if (!fitResult.timeFit || !hit->measuresTime()) {
475 const double totFlightDist = (localToGlobal * planeIsect).
mag() * c_inv;
480 ATH_MSG_VERBOSE(
"Partial derivative of "<<idHelperSvc->toString(hit->spacePoint()->identify())
482 <<
" "<<std::endl<<
toString(hit->covariance())<<std::endl<<
" w.r.t "
490 const auto &measurement = *hit->spacePoint()->primaryMeasurement();
491 ATH_MSG_WARNING(
"MdtSegmentFitter() - Unsupported measurment type" <<
typeid(measurement).
name());
495 ATH_MSG_VERBOSE(
"Update derivatives for hit "<< (hit->spacePoint() ? idHelperSvc->toString(hit->spacePoint()->identify()) :
"beamspot"));
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);
508 fitResult.converged =
true;
509 ATH_MSG_VERBOSE(
"Fit converged after "<<fitResult.nIter<<
" iterations with "<<fitResult.chi2);
515 if (!fitResult.nPhiMeas && !fitResult.timeFit) {
516 paramUpdate = updateParameters<2>(fitResult.segmentPars, prevPars, gradient, prevGrad, hessian);
517 }
else if (!fitResult.nPhiMeas && fitResult.timeFit) {
525 paramUpdate = updateParameters<3>(fitResult.segmentPars, prevPars, gradient, prevGrad, hessian);
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);
538 fitResult.converged =
true;
549 fitResult.nDoF-=fitResult.timeFit;
552 const auto [segPos, segDir] = fitResult.makeLine();
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();
563 for (
const HitType& hit : fitResult.calibMeasurements) {
564 hit->setDriftRadius(std::abs(hit->driftRadius()));
568 if (!fitResult.nPhiMeas&& !fitResult.timeFit) {
569 blockCovariance<2>(std::move(hessian), fitResult.segmentParErrs);
570 }
else if (!fitResult.nPhiMeas && fitResult.timeFit) {
573 blockCovariance<3>(std::move(hessian), fitResult.segmentParErrs);
576 }
else if (fitResult.nPhiMeas) {
577 blockCovariance<4>(std::move(hessian), fitResult.segmentParErrs);
578 }
else if (fitResult.nPhiMeas && fitResult.timeFit) {
579 blockCovariance<5>(std::move(hessian), fitResult.segmentParErrs);
◆ 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.
◆ 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.
◆ 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
-
lvl | The message level to test against |
- Returns
- boolean Indicating if messages at given level will be printed
- Return values
-
true | Messages at level "lvl" will be printed |
Definition at line 151 of file AthMessaging.h.
◆ 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
-
ctx | EventContext to fetch the calibration constants |
fitResult | Current fit state |
centerOfGravity | Vector to save the center of gravity position later on. |
invCov | Vector 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.
65 const auto [segPos, segDir] =
makeLine(fitResult.segmentPars);
67 fitResult.calibMeasurements =
m_cfg.
calibrator->
calibrate(ctx, std::move(fitResult.calibMeasurements), segPos, segDir,
74 if (fitResult.timeFit && fitResult.nDoF <= 1) {
75 fitResult.timeFit =
false;
76 ATH_MSG_DEBUG(
"Switch of the time fit because nDoF: "<<fitResult.nDoF);
79 fitResult.calibMeasurements =
m_cfg.
calibrator->
calibrate(ctx, std::move(fitResult.calibMeasurements), segPos, segDir,
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;
◆ setLevel()
void AthMessaging::setLevel |
( |
MSG::Level |
lvl | ) |
|
|
inherited |
◆ updateDerivatives()
Definition at line 584 of file MdtSegmentFitter.cxx.
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]);
596 const int symMatIdx = vecIdxFromSymMat<ResidualWithPartials::nPars>(
p,
k);
597 hessian(
p,
k)+=
contract(invCov, covRes, fitMeas.hessian[symMatIdx]);
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));
◆ updateDriftSigns()
Update the signs of the measurement.
Definition at line 55 of file MdtSegmentFitter.cxx.
57 for (std::unique_ptr<CalibratedSpacePoint>& meas : fitResult.calibMeasurements) {
◆ 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
-
fitResult | Reference 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.
92 fitResult.nPhiMeas = fitResult.nDoF = fitResult.nTimeMeas = 0;
93 for (
const HitType& hit : fitResult.calibMeasurements) {
98 fitResult.nPhiMeas+= hit->measuresPhi();
99 fitResult.nDoF+= hit->measuresPhi();
100 fitResult.nDoF+= hit->measuresEta();
104 fitResult.nTimeMeas+=hit->measuresTime();
106 if (!fitResult.nDoF) {
110 if (!fitResult.nPhiMeas) {
113 const double nHits = fitResult.calibMeasurements.size();
115 std::ranges::for_each(fitResult.calibMeasurements,[&avgX, nHits](
const HitType& hit) {
116 return avgX += hit->positionInChamber().x() / nHits;
121 fitResult.nDoF = fitResult.nDoF - 2 - (fitResult.nPhiMeas > 0 ? 2 : 0);
◆ updateLinePartials()
Updates the line parameters together with its first & second order partial derivatives.
- Parameters
-
fitPars | Set of segment parameters in the current iteration |
line | Reference 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.
◆ updateParameters()
template<unsigned int nDim>
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
-
currentPars | Best segment estimator parameters |
previousPars | Segment estimator parameters from the last iteration |
currGrad | Gradient of the chi2 from this iteration |
prevGrad | Gradient of the chi2 from the previous iteration |
hessian | Hessian 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.
626 AmgSymMatrix(nDim) miniHessian = currentHessian.block<nDim, nDim>(0,0);
628 <<
", hessian ("<<miniHessian.determinant()<<
")"<<std::endl<<miniHessian);
630 double changeMag{0.};
631 if (std::abs(miniHessian.determinant()) > detCutOff) {
632 prevPars.block<nDim,1>(0,0) = currPars.block<nDim,1>(0,0);
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 -"
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))
646 ATH_MSG_VERBOSE(
"Hessian determinant invalid. Try deepest descent - \nprev parameters: "
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();
654 unsigned int nOutOfBound{0};
655 for (
unsigned int p = 0;
p< nDim; ++
p) {
656 double& parValue{currPars[
p]};
657 if constexpr(nDim == 3) {
◆ ATLAS_THREAD_SAFE
std::atomic_flag m_initialized AthMessaging::ATLAS_THREAD_SAFE = ATOMIC_FLAG_INIT |
|
mutableprivateinherited |
◆ m_cfg
Config MuonR4::MdtSegmentFitter::m_cfg {} |
|
private |
◆ m_imsg
std::atomic<IMessageSvc*> AthMessaging::m_imsg { nullptr } |
|
mutableprivateinherited |
◆ m_lvl
std::atomic<MSG::Level> AthMessaging::m_lvl { MSG::NIL } |
|
mutableprivateinherited |
◆ 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 |
The documentation for this class was generated from the following files:
std::atomic< MSG::Level > m_lvl
Current logging level.
void updateDriftSigns(const Amg::Vector3D &segPos, const Amg::Vector3D &segDir, SegmentFitResult &fitRes) const
Update the signs of the measurement.
bool updateHitSummary(SegmentFitResult &fitResult) const
Brief updates the hit summary from the contributing hits.
void calculateStripResiduals(const LineWithPartials &line, const CalibratedSpacePoint &spacePoint, ResidualWithPartials &residual) const
Calculates the residual together with hte correspdonding derivatives for strip measurements.
static constexpr unsigned nPars
Free parameters of the line (x0,y0,theta,phi)
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.
SegmentFit::Parameters Parameters
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.
#define ATH_MSG_VERBOSE(x)
virtual double driftVelocity(const EventContext &ctx, const CalibratedSpacePoint &spacePoint) const =0
Returns the drift velocity for a given drift-circle space point.
unsigned int noMoveIter
How many iterations with changes below tolerance.
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.
std::pair< Amg::Vector3D, Amg::Vector3D > makeLine(const Parameters &pars)
Returns the parsed parameters into an Eigen line parametrization.
std::atomic< IMessageSvc * > m_imsg
MessageSvc pointer.
double tolerance
Gradient cut off.
std::string toString(const CalibratedSpacePoint::Covariance_t &mat)
Returns the matrix in string.
IMessageSvc * getMessageSvc(bool quiet=false)
const ISpacePointCalibrator * calibrator
Pointer to the calibrator tool.
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...
std::string toString(const Translation3D &translation, int precision=4)
GeoPrimitvesToStringConverter.
virtual double driftAcceleration(const EventContext &ctx, const CalibratedSpacePoint &spacePoint) const =0
Returns the drift acceleration for a given drift-circle space point.
bool msgLvl(const MSG::Level lvl) const
Test the output level.
CalibratedSpacePoint::Covariance_t inverse(const CalibratedSpacePoint::Covariance_t &mat)
Inverts the parsed matrix.
CalibratedSpacePoint::Covariance_t MeasCov_t
Updates the chi2, its Gradient & Hessian from the measurement residual.
CalibratedSpacePoint::State State
double chi2(TH1 *h0, TH1 *h1)
State
State flag to distinguish different space point states.
bool useSecOrderDeriv
Switch toggling whether the second order derivative shall be included.
Class to provide easy MsgStream access and capabilities.
unsigned int nParsOutOfBounds
Abort the fit as soon as more than n parameters leave the fit range.
void calculateWireResiduals(const LineWithPartials &line, const CalibratedSpacePoint &spacePoint, ResidualWithPartials &residual) const
Calculates the residuals together with the corresponding derivatives for a drift-circle measurement.
def dot(G, fn, nodesToHighlight=[])
MsgStream & msg() const
The standard message stream.
void updateDerivatives(const ResidualWithPartials &fitMeas, const MeasCov_t &measCovariance, AmgVector(5)&gradient, AmgSymMatrix(5)&hessian, double &chi2, int startPar) const
unsigned int nMaxCalls
How many calls shall be executed.
constexpr int toInt(const ParamDefs p)
bool reCalibrate
Switch toggling whether the calibrator shall be called at each iteration.
Eigen::Matrix< double, 3, 1 > Vector3D
Amg::Vector3D dirFromAngles(const double phi, const double theta)
Constructs a direction vector from the azimuthal & polar angles.
#define ATH_MSG_WARNING(x)
bool doTimeFit
Switch toggling whether the T0 shall be fitted or not.
std::string m_nm
Message source name.
Helper to simultaneously calculate sin and cos of the same angle.
Interface for Helper service that creates muon Identifiers and can be used to print Identifiers.
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.
void initMessaging() const
Initialize our message level and MessageSvc.
boost::thread_specific_ptr< MsgStream > m_msg_tls
MsgStream instance (a std::cout like with print-out levels)
const AmgSymMatrix(2) &SpacePoint
constexpr int pow(int base, int exp) noexcept
bool recalibrate(const EventContext &ctx, SegmentFitResult &fitResult) const
Recalibrate the measurements participating in the fit & shift them into the centre-of gravity frame.
Scalar mag() const
mag method
SegmentFit::ParamDefs ParamDefs
static void updateLinePartials(const Parameters &fitPars, LineWithPartials &line)
Updates the line parameters together with its first & second order partial derivatives.
double contract(const CalibratedSpacePoint::Covariance_t &mat, const Amg::Vector3D &a, const Amg::Vector3D &b)