263 {
264
265 SG::ReadCondHandle<MuonCalib::MdtCalibDataContainer> mdtCalibConstants{
m_calibDbKey, ctx};
266 if (!mdtCalibConstants.
isValid()) {
268 throw std::runtime_error("Failed to retrieve calibration constants");
269 }
270 ATH_MSG_VERBOSE(
"extractTimeMeasurementsFromTrack for candidate: beta seed " << candidate.betaSeed.beta);
272 if (!combinedTrack) return;
273
274
275 float betaSeed = candidate.betaFitResult.beta;
276
277
278 Muon::TimePointBetaFitter
fitter;
280
281
283 if (!states) {
284 ATH_MSG_WARNING(
" track without states, cannot extractTimeMeasurementsFromTrack ");
285 return;
286 }
287
289
290
291 typedef std::vector<const Muon::MuonClusterOnTrack*> RpcClVec;
292 using RpcClPerChMap = std::map<Identifier, std::tuple<const Trk::TrackParameters*, RpcClVec, RpcClVec>>;
293 RpcClPerChMap rpcPrdsPerChamber;
294
295 using MdtTubeData = std::pair<const Trk::TrackParameters*, const Muon::MdtDriftCircleOnTrack*>;
296 using MdtTubeDataVec = std::vector<MdtTubeData>;
297 using MdtChamberLayerData = std::map<int, MdtTubeDataVec>;
298 MdtChamberLayerData mdtDataPerChamberLayer;
299
300
303 for (; tsit != tsit_end; ++tsit) {
305 if (!pars) continue;
306
307
308 const Trk::MeasurementBase* meas = (*tsit)->measurementOnTrack();
310
311
314
315
317
318 const Muon::MdtDriftCircleOnTrack* mdt = dynamic_cast<const Muon::MdtDriftCircleOnTrack*>(meas);
319 if (!mdt) continue;
320
323 if (chIndexWithBIR ==
toInt(ChIndex::BIL)) {
325 if (
stName[2] ==
'R') { chIndexWithBIR += 1000; }
326 }
327 mdtDataPerChamberLayer[chIndexWithBIR].push_back(std::make_pair(pars, mdt));
328 } else {
332
333 float ix =
pars->position().x();
334 float iy =
pars->position().y();
335 float iz =
pars->position().z();
337 float er = -1;
340 float propTime = 0;
342
343
347 auto data = mdtCalibConstants->getCalibData(
id, msgStream());
348 const auto& rtRelation =
data->rtRelation;
349 float drdt = rtRelation->rt()->driftVelocity(driftTime);
350 float rres = rtRelation->rtRes()->resolution(driftTime);
351 float tres = rres / drdt;
352 float TlocR = rtRelation->tr()->driftTime(std::abs(locR)).value_or(0.);
353 float trackTimeRes = errR / drdt;
354 float tofShiftFromBeta =
calculateTof(betaSeed, distance) - tof;
355 er = std::sqrt(tres * tres + trackTimeRes * trackTimeRes);
360
362 std::unique_ptr<const Trk::TrackParameters> unbiasedPars(
364 if (unbiasedPars) {
365 float locRu = unbiasedPars->parameters()[
Trk::locR];
366 float TlocRu = rtRelation->tr()->driftTime(std::abs(locRu)).value_or(0.);
367 float errRu = unbiasedPars->covariance() ?
Amg::error(*unbiasedPars->covariance(),
Trk::locR) : 0.3;
368 float trackTimeResu = errRu / drdt;
371 er = std::sqrt(tres * tres + trackTimeResu * trackTimeResu);
373 ATH_MSG_VERBOSE(
" Got unbiased parameters: r " << locR <<
" ur " << locRu <<
" err " << errR <<
" uerr "
374 << errRu << " terr " << trackTimeRes << " terru "
375 << trackTimeResu);
376 }
377 }
378
380 << " TlocR " << TlocR << " diff " << driftTime - TlocR << " tofShift " << tofShiftFromBeta
381 << " time " << time << " err " << er << " intrinsic " << tres << " track " << trackTimeRes);
382
385 << beta << " diff " << std::abs(beta - betaSeed));
387
388 candidate.stauHits.emplace_back(MuGirlNS::StauHit(tech, time + tof, ix, iy, iz, id, ie, er, sh, isEta, propTime, false));
392 candidate.stauMDTHitExtras.emplace_back(MuGirlNS::StauMDTHitExtra(iadc, irdrift));
393 }
394 continue;
395 }
396
397
398 hits.emplace_back(Muon::TimePointBetaFitter::Hit(distance, time, er));
399 candidate.stauHits.emplace_back(MuGirlNS::StauHit(tech, time + tof, ix, iy, iz, id, ie, er, sh, isEta, propTime, true));
400
404 candidate.stauMDTHitExtras.emplace_back(MuGirlNS::StauMDTHitExtra(iadc, irdrift));
405 }
406 }
408
409 std::vector<const Muon::MuonClusterOnTrack*>
clusters;
410 const Muon::CompetingMuonClustersOnTrack* crot = dynamic_cast<const Muon::CompetingMuonClustersOnTrack*>(meas);
411 if (crot) {
413 } else {
414 const Muon::RpcClusterOnTrack* rpc = dynamic_cast<const Muon::RpcClusterOnTrack*>(meas);
416 }
419 auto pos = rpcPrdsPerChamber.find(chamberId);
420 if (pos == rpcPrdsPerChamber.end()) {
421 if (measuresPhi)
422 rpcPrdsPerChamber[chamberId] = std::make_tuple(pars, clusters, RpcClVec());
423 else
424 rpcPrdsPerChamber[chamberId] = std::make_tuple(pars, RpcClVec(), clusters);
425 } else {
426 RpcClVec& clVec = measuresPhi ? std::get<1>(
pos->second) : std::
get<2>(
pos->
second);
428 }
430 const Muon::CscClusterOnTrack* csc = dynamic_cast<const Muon::CscClusterOnTrack*>(meas);
431
435
436 float ix =
pars->position().x();
437 float iy =
pars->position().y();
438 float iz =
pars->position().z();
440 float er = -1;
443 float propTime = 0;
445 candidate.stauHits.push_back(MuGirlNS::StauHit(tech, time + tof, ix, iy, iz, id, ie, er, sh, isEta, propTime));
446 }
447 }
448
452
453 std::vector<const Muon::MuonClusterOnTrack*> calibratedClusters;
454 for (const auto* cluster : clusters) {
456 if (cl) calibratedClusters.push_back(cl);
457 }
458 if (calibratedClusters.empty()) return;
459
460 Muon::IMuonHitTimingTool::TimingResult
result =
m_hitTimingTool->calculateTimingResult(calibratedClusters);
461 for (
const auto* cl : calibratedClusters)
delete cl;
462 if (!
result.valid)
return;
463
464 Identifier
id =
clusters.front()->identify();
465
469 float ix =
pars.position().x();
470 float iy =
pars.position().y();
471 float iz =
pars.position().z();
478 float propTime = 0;
482 << " diff " << std::abs(beta - betaSeed));
483
485
486 hits.push_back(Muon::TimePointBetaFitter::Hit(distance, time, er));
487 candidate.stauHits.push_back(MuGirlNS::StauHit(tech, time + tof, ix, iy, iz, id, ie, er, sh, isEta, propTime));
488 };
489
490
491 RpcClPerChMap::const_iterator chit = rpcPrdsPerChamber.begin();
492 RpcClPerChMap::const_iterator chit_end = rpcPrdsPerChamber.end();
494
495 for (; chit != chit_end; ++chit) {
497 const RpcClVec& phiClusters = std::get<1>(chit->second);
498 const RpcClVec& etaClusters = std::get<2>(chit->second);
499 insertRpcs(*pars, phiClusters, candidate, hits);
500 insertRpcs(*pars, etaClusters, candidate, hits);
501 }
502
503
505 Muon::MuonDriftCircleErrorStrategy calibrationStrategy(bits);
508
509 TrkDriftCircleMath::DCSLFitter mdtFitter;
510 TrkDriftCircleMath::SegmentFinder segmentFinder;
514
515 MdtChamberLayerData::const_iterator mit = mdtDataPerChamberLayer.begin();
516 MdtChamberLayerData::const_iterator mit_end = mdtDataPerChamberLayer.end();
517 for (; mit != mit_end; ++mit) {
519 << " hits " << mit->second.size());
520 if (mit->second.size() < 4) continue;
521
522
523 const MuonGM::MdtReadoutElement* detEl = mit->second.front().second->detectorElement();
524 const Trk::PlaneSurface* surf =
dynamic_cast<const Trk::PlaneSurface*
>(&detEl->
surface());
525 if (!surf) {
526 ATH_MSG_WARNING(
"MdtReadoutElement should always have a PlaneSurface as reference surface");
527 continue;
528 }
530
531
534
535
536 Trk::LocalDirection seedLocDir;
539 TrkDriftCircleMath::LocVec2D seedPos(seedLocPos.y(), seedLocPos.z());
540 TrkDriftCircleMath::Line seedLine(seedPos, seedLocDir.
angleYZ());
542
543 std::vector<std::pair<std::shared_ptr<const Muon::MdtDriftCircleOnTrack>,
const Trk::TrackParameters*>> indexLookUp;
545 for (const auto& entry : mit->second) {
547 const Muon::MdtDriftCircleOnTrack& mdt = *
entry.second;
549
550 std::unique_ptr<const Muon::MdtDriftCircleOnTrack> calibratedMdt(
552 if (!calibratedMdt) {
554 continue;
555 }
559 <<
" r_track " <<
pars.parameters()[
Trk::locR] <<
" residual "
561
562
566
567
569 TrkDriftCircleMath::LocVec2D lpos(locPos.y(), locPos.z());
570
571 double r = std::abs(calibratedMdt->driftRadius());
573
574
577
578
580 TrkDriftCircleMath::DCOnTrack dcOnTrack(dc, 1., 1.);
581
582 dcs.push_back(dcOnTrack);
583 indexLookUp.emplace_back(std::move(calibratedMdt), &pars);
585 }
586
587
588 for (
unsigned int i = 0;
i < dcs.size(); ++
i) {
594 continue;
595 }
596 TrkDriftCircleMath::Segment segment =
result;
598 unsigned int ndofFit = segment.
ndof();
599 if (ndofFit < 1) continue;
600 double chi2NdofSegmentFit = segment.
chi2() / ndofFit;
601 bool hasDropHit = false;
602 unsigned int dropDepth = 0;
603 if (!segmentFinder.
dropHits(segment, hasDropHit, dropDepth)) {
604 ATH_MSG_DEBUG(
"DropHits failed, fit chi2/ndof " << chi2NdofSegmentFit);
605 if (
msgLvl(MSG::VERBOSE)) {
608 segmentFinder.
dropHits(segment, hasDropHit, dropDepth);
610 }
611 continue;
612 }
613 if (i >= segment.
dcs().size())
continue;
614 TrkDriftCircleMath::TransformToLine toLine(segment.
line());
615 const TrkDriftCircleMath::DCOnTrack& dc = segment.
dcs()[
i];
619 double rline = toLine.toLineY(dc.
position());
621 if (index < 0 || index >= (int)indexLookUp.size()) {
622 ATH_MSG_WARNING(
" lookup of TrackParameters and MdtDriftCircleOnTrack failed " << index <<
" range: 0 - "
623 << indexLookUp.size() - 1);
624 continue;
625 }
627 const Muon::MdtDriftCircleOnTrack* mdt = indexLookUp[dc.
index()].first.get();
629
630
631 std::shared_ptr<const Muon::MdtDriftCircleOnTrack> calibratedMdt(
633 if (!calibratedMdt.get()) {
635 continue;
636 }
639
640 float ix =
pars->position().x();
641 float iy =
pars->position().y();
642 float iz =
pars->position().z();
644 float er = -1;
647 float propTime = 0;
649
652
653
654 float driftTime = calibratedMdt->driftTime();
657 auto data = mdtCalibConstants->getCalibData(
id, msgStream());
658 const auto& rtRelation =
data->rtRelation;
659 float drdt = rtRelation->rt()->driftVelocity(driftTime);
660 float rres = rtRelation->rtRes()->resolution(driftTime);
661 float tres = rres / drdt;
662 float TlocR = rtRelation->tr()->driftTime(std::abs(locR)).value_or(0.);
663 float trackTimeRes = errR / drdt;
664 float tofShiftFromBeta = 0.;
665 er = std::sqrt(tres * tres + trackTimeRes * trackTimeRes);
670
673
675 msg(MSG::DEBUG) <<
m_idHelperSvc->toString(
id) << std::setprecision(2) <<
" segment: after fit " << std::setw(5)
676 << chi2NdofSegmentFit << " ndof " << std::setw(2) << ndofFit;
677 if (segment.
ndof() != ndofFit)
678 msg(MSG::DEBUG) <<
" after outlier " << std::setw(5) << chi2NdofSegmentFit <<
" ndof " << std::setw(2) << ndofFit;
679 msg(MSG::DEBUG) <<
" driftR " << std::setw(4) << dc.
r() <<
" rline " << std::setw(5) << rline <<
" residual "
680 << std::setw(5) <<
res <<
" pull " << std::setw(4) <<
pull <<
" time " << std::setw(3) <<
time
681 <<
" beta" << std::setw(2) <<
beta <<
" err " << std::setw(3) << er <<
" intrinsic " << std::setw(3)
682 <<
tres <<
" track " << std::setw(3) << trackTimeRes;
683 if (!isSelected)
msg(MSG::DEBUG) <<
" outlier";
684 msg(MSG::DEBUG) << std::setprecision(5) <<
endmsg;
685 }
686
687 if (!isSelected) {
688
689 candidate.stauHits.emplace_back(MuGirlNS::StauHit(
MuGirlNS::MDTT_STAU_HIT, time + tof, ix, iy, iz,
id, ie, er, sh, isEta, propTime,
false));
693 candidate.stauMDTHitExtras.emplace_back(MuGirlNS::StauMDTHitExtra(iadc, irdrift));
694 }
695 continue;
696 }
697
698
699 hits.emplace_back(distance, time, er);
700 candidate.stauHits.emplace_back(
MuGirlNS::MDTT_STAU_HIT, time + tof, ix, iy, iz,
id, ie, er, sh, isEta, propTime,
true);
702 candidate.stauMDTHitExtras.emplace_back(MuGirlNS::StauMDTHitExtra(iadc, irdrift));
703 }
704 }
705 }
706
707 Muon::TimePointBetaFitter::FitResult betaFitResult =
fitter.fitWithOutlierLogic(hits);
708 ATH_MSG_DEBUG(
" extractTimeMeasurementsFromTrack: extracted " << candidate.stauHits.size() <<
" time measurements "
709 <<
" status fit " << betaFitResult.
status <<
" beta "
710 << betaFitResult.
beta <<
" chi2/ndof " << betaFitResult.
chi2PerDOF());
711
712 candidate.finalBetaFitResult = betaFitResult;
713 }
char data[hepevt_bytes_allocation_ATLAS]
std::pair< std::vector< unsigned int >, bool > res
DataModel_detail::const_iterator< DataVector > const_iterator
virtual const Trk::Surface & surface() const override final
Return surface associated with this detector element.
virtual Amg::Transform3D GlobalToAmdbLRSTransform() const
const std::vector< const MuonClusterOnTrack * > & containedROTs() const
returns the vector of SCT_ClusterOnTrack objects .
virtual const CscPrepData * prepRawData() const override final
Returns the CscPrepData - is a CscPrepData in this scope.
double time() const
Returns the time.
double driftRadius() const
Returns the value of the drift radius.
double driftTime() const
Returns the value of the drift time used to obtain the drift radius.
virtual const Trk::StraightLineSurface & associatedSurface() const override final
Returns the surface on which this measurement was taken.
virtual const MdtPrepData * prepRawData() const override final
Returns the PrepRawData used to create this corrected measurement.
int adc() const
Returns the ADC (typically range is 0 to 250)
@ Segment
Treating a segment or a track.
void residual(double res)
set residual
void errorTrack(double error)
set track error
virtual bool fit(Segment &result, const Line &line, const DCOnTrackVec &dcs, double t0Seed=-99999.) const
unsigned int index() const
double r() const
access to drift radius
const LocVec2D & position() const
access to local position
double dr() const
access to error drift radius
@ InTime
drift time too small to be compatible with drift spectrum
void setMaxDropDepth(int max)
bool dropHits(Segment &segment, bool &hasDroppedHit, unsigned int &dropDepth) const
void setDeltaCut(double cut)
void debugLevel(int debugLevel)
void setChi2DropCut(double chi2)
void hitsOnTrack(unsigned int hitsOnTrack)
const Line & line() const
const DCOnTrackVec & dcs() const
unsigned int ndof() const
double angleYZ() const
access method for angle of local YZ projection
const LocalParameters & localParameters() const
Interface method to get the LocalParameters.
const Amg::MatrixX & localCovariance() const
Interface method to get the localError.
const Amg::Vector3D & momentum() const
Access method for the momentum.
const Amg::Vector3D & position() const
Access method for the position.
virtual Intersection straightLineIntersection(const Amg::Vector3D &pos, const Amg::Vector3D &dir, bool forceDir, Trk::BoundaryCheck bchk) const override final
fast straight line intersection schema - standard: provides closest intersection and (signed) path le...
void globalToLocalDirection(const Amg::Vector3D &glodir, Trk::LocalDirection &locdir) const
This method transforms the global direction to a local direction wrt the plane.
Identifier identify() const
return the identifier -extends MeasurementBase
virtual void localToGlobal(const Amg::Vector2D &locp, const Amg::Vector3D &mom, Amg::Vector3D &glob) const override final
Specified for StraightLineSurface: LocalToGlobal method without dynamic memory allocation.
@ Outlier
This TSoS contains an outlier, that is, it contains a MeasurementBase/RIO_OnTrack which was not used ...
const std::string selection
T * get(TKey *tobj)
get a TObject* from a TKey* (why can't a TObject be a TKey?)
Eigen::Affine3d Transform3D
Eigen::Matrix< double, 2, 1 > Vector2D
Eigen::Matrix< double, 3, 1 > Vector3D
const std::string & stName(StIndex index)
convert StIndex into a string
const std::string & chName(ChIndex index)
convert ChIndex into a string
ChIndex
enum to classify the different chamber layers in the muon spectrometer
std::bitset< 23 > MuonDriftCircleErrorStrategyInput
std::vector< bool > HitSelection
std::vector< DCOnTrack > DCOnTrackVec
DataVector< const Trk::TrackStateOnSurface > TrackStates
ParametersBase< TrackParametersDim, Charged > TrackParameters
void combinedTrack(long int ICH, double *pv0, double *covi, double BMAG, double *par, double *covo)
float chi2PerDOF() const
chi2/ndof, return 0 if ndof == 0 or status == 0
float beta
status flag (0 = failed, 1 = ok)