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) {
412 std::ranges::transform(crot->
containedROTs(), std::back_inserter(clusters),
413 [](const auto& rot){
414 return rot.get();
415 });
416 } else {
417 const Muon::RpcClusterOnTrack* rpc = dynamic_cast<const Muon::RpcClusterOnTrack*>(meas);
419 }
422 auto pos = rpcPrdsPerChamber.find(chamberId);
423 if (pos == rpcPrdsPerChamber.end()) {
424 if (measuresPhi)
425 rpcPrdsPerChamber[chamberId] = std::make_tuple(pars, clusters, RpcClVec());
426 else
427 rpcPrdsPerChamber[chamberId] = std::make_tuple(pars, RpcClVec(), clusters);
428 } else {
429 RpcClVec& clVec = measuresPhi ? std::get<1>(
pos->second) : std::
get<2>(
pos->
second);
431 }
433 const Muon::CscClusterOnTrack* csc = dynamic_cast<const Muon::CscClusterOnTrack*>(meas);
434
438
439 float ix =
pars->position().x();
440 float iy =
pars->position().y();
441 float iz =
pars->position().z();
443 float er = -1;
446 float propTime = 0;
448 candidate.stauHits.push_back(MuGirlNS::StauHit(tech, time + tof, ix, iy, iz, id, ie, er, sh, isEta, propTime));
449 }
450 }
451
455
456 std::vector<const Muon::MuonClusterOnTrack*> calibratedClusters;
457 for (const auto* cluster : clusters) {
459 if (cl) calibratedClusters.push_back(cl);
460 }
461 if (calibratedClusters.empty()) return;
462
463 Muon::IMuonHitTimingTool::TimingResult
result =
m_hitTimingTool->calculateTimingResult(calibratedClusters);
464 for (
const auto* cl : calibratedClusters)
delete cl;
465 if (!
result.valid)
return;
466
467 Identifier
id =
clusters.front()->identify();
468
472 float ix =
pars.position().x();
473 float iy =
pars.position().y();
474 float iz =
pars.position().z();
481 float propTime = 0;
485 << " diff " << std::abs(beta - betaSeed));
486
488
489 hits.push_back(Muon::TimePointBetaFitter::Hit(distance, time, er));
490 candidate.stauHits.push_back(MuGirlNS::StauHit(tech, time + tof, ix, iy, iz, id, ie, er, sh, isEta, propTime));
491 };
492
493
494 RpcClPerChMap::const_iterator chit = rpcPrdsPerChamber.begin();
495 RpcClPerChMap::const_iterator chit_end = rpcPrdsPerChamber.end();
497
498 for (; chit != chit_end; ++chit) {
500 const RpcClVec& phiClusters = std::get<1>(chit->second);
501 const RpcClVec& etaClusters = std::get<2>(chit->second);
502 insertRpcs(*pars, phiClusters, candidate, hits);
503 insertRpcs(*pars, etaClusters, candidate, hits);
504 }
505
506
508 Muon::MuonDriftCircleErrorStrategy calibrationStrategy(bits);
511
512 TrkDriftCircleMath::DCSLFitter mdtFitter;
513 TrkDriftCircleMath::SegmentFinder segmentFinder;
517
518 MdtChamberLayerData::const_iterator mit = mdtDataPerChamberLayer.begin();
519 MdtChamberLayerData::const_iterator mit_end = mdtDataPerChamberLayer.end();
520 for (; mit != mit_end; ++mit) {
522 << " hits " << mit->second.size());
523 if (mit->second.size() < 4) continue;
524
525
526 const MuonGM::MdtReadoutElement* detEl = mit->second.front().second->detectorElement();
527 const Trk::PlaneSurface* surf =
dynamic_cast<const Trk::PlaneSurface*
>(&detEl->
surface());
528 if (!surf) {
529 ATH_MSG_WARNING(
"MdtReadoutElement should always have a PlaneSurface as reference surface");
530 continue;
531 }
533
534
537
538
539 Trk::LocalDirection seedLocDir;
542 TrkDriftCircleMath::LocVec2D seedPos(seedLocPos.y(), seedLocPos.z());
543 TrkDriftCircleMath::Line seedLine(seedPos, seedLocDir.
angleYZ());
545
546 std::vector<std::pair<std::shared_ptr<const Muon::MdtDriftCircleOnTrack>,
const Trk::TrackParameters*>> indexLookUp;
548 for (const auto& entry : mit->second) {
550 const Muon::MdtDriftCircleOnTrack& mdt = *
entry.second;
552
553 std::unique_ptr<const Muon::MdtDriftCircleOnTrack> calibratedMdt(
555 if (!calibratedMdt) {
557 continue;
558 }
562 <<
" r_track " <<
pars.parameters()[
Trk::locR] <<
" residual "
564
565
569
570
572 TrkDriftCircleMath::LocVec2D lpos(locPos.y(), locPos.z());
573
574 double r = std::abs(calibratedMdt->driftRadius());
576
577
580
581
583 TrkDriftCircleMath::DCOnTrack dcOnTrack(dc, 1., 1.);
584
585 dcs.push_back(dcOnTrack);
586 indexLookUp.emplace_back(std::move(calibratedMdt), &pars);
588 }
589
590
591 for (
unsigned int i = 0;
i < dcs.size(); ++
i) {
597 continue;
598 }
599 TrkDriftCircleMath::Segment segment =
result;
601 unsigned int ndofFit = segment.
ndof();
602 if (ndofFit < 1) continue;
603 double chi2NdofSegmentFit = segment.
chi2() / ndofFit;
604 bool hasDropHit = false;
605 unsigned int dropDepth = 0;
606 if (!segmentFinder.
dropHits(segment, hasDropHit, dropDepth)) {
607 ATH_MSG_DEBUG(
"DropHits failed, fit chi2/ndof " << chi2NdofSegmentFit);
608 if (
msgLvl(MSG::VERBOSE)) {
611 segmentFinder.
dropHits(segment, hasDropHit, dropDepth);
613 }
614 continue;
615 }
616 if (i >= segment.
dcs().size())
continue;
617 TrkDriftCircleMath::TransformToLine toLine(segment.
line());
618 const TrkDriftCircleMath::DCOnTrack& dc = segment.
dcs()[
i];
622 double rline = toLine.toLineY(dc.
position());
624 if (index < 0 || index >= (int)indexLookUp.size()) {
625 ATH_MSG_WARNING(
" lookup of TrackParameters and MdtDriftCircleOnTrack failed " << index <<
" range: 0 - "
626 << indexLookUp.size() - 1);
627 continue;
628 }
630 const Muon::MdtDriftCircleOnTrack* mdt = indexLookUp[dc.
index()].first.get();
632
633
634 std::shared_ptr<const Muon::MdtDriftCircleOnTrack> calibratedMdt(
636 if (!calibratedMdt.get()) {
638 continue;
639 }
642
643 float ix =
pars->position().x();
644 float iy =
pars->position().y();
645 float iz =
pars->position().z();
647 float er = -1;
650 float propTime = 0;
652
655
656
657 float driftTime = calibratedMdt->driftTime();
660 auto data = mdtCalibConstants->getCalibData(
id, msgStream());
661 const auto& rtRelation =
data->rtRelation;
662 float drdt = rtRelation->rt()->driftVelocity(driftTime);
663 float rres = rtRelation->rtRes()->resolution(driftTime);
664 float tres = rres / drdt;
665 float TlocR = rtRelation->tr()->driftTime(std::abs(locR)).value_or(0.);
666 float trackTimeRes = errR / drdt;
667 float tofShiftFromBeta = 0.;
668 er = std::sqrt(tres * tres + trackTimeRes * trackTimeRes);
673
676
678 msg(MSG::DEBUG) <<
m_idHelperSvc->toString(
id) << std::setprecision(2) <<
" segment: after fit " << std::setw(5)
679 << chi2NdofSegmentFit << " ndof " << std::setw(2) << ndofFit;
680 if (segment.
ndof() != ndofFit)
681 msg(MSG::DEBUG) <<
" after outlier " << std::setw(5) << chi2NdofSegmentFit <<
" ndof " << std::setw(2) << ndofFit;
682 msg(MSG::DEBUG) <<
" driftR " << std::setw(4) << dc.
r() <<
" rline " << std::setw(5) << rline <<
" residual "
683 << std::setw(5) <<
res <<
" pull " << std::setw(4) <<
pull <<
" time " << std::setw(3) <<
time
684 <<
" beta" << std::setw(2) <<
beta <<
" err " << std::setw(3) << er <<
" intrinsic " << std::setw(3)
685 <<
tres <<
" track " << std::setw(3) << trackTimeRes;
686 if (!isSelected)
msg(MSG::DEBUG) <<
" outlier";
687 msg(MSG::DEBUG) << std::setprecision(5) <<
endmsg;
688 }
689
690 if (!isSelected) {
691
692 candidate.stauHits.emplace_back(MuGirlNS::StauHit(
MuGirlNS::MDTT_STAU_HIT, time + tof, ix, iy, iz,
id, ie, er, sh, isEta, propTime,
false));
696 candidate.stauMDTHitExtras.emplace_back(MuGirlNS::StauMDTHitExtra(iadc, irdrift));
697 }
698 continue;
699 }
700
701
702 hits.emplace_back(distance, time, er);
703 candidate.stauHits.emplace_back(
MuGirlNS::MDTT_STAU_HIT, time + tof, ix, iy, iz,
id, ie, er, sh, isEta, propTime,
true);
705 candidate.stauMDTHitExtras.emplace_back(MuGirlNS::StauMDTHitExtra(iadc, irdrift));
706 }
707 }
708 }
709
710 Muon::TimePointBetaFitter::FitResult betaFitResult =
fitter.fitWithOutlierLogic(hits);
711 ATH_MSG_DEBUG(
" extractTimeMeasurementsFromTrack: extracted " << candidate.stauHits.size() <<
" time measurements "
712 <<
" status fit " << betaFitResult.
status <<
" beta "
713 << betaFitResult.
beta <<
" chi2/ndof " << betaFitResult.
chi2PerDOF());
714
715 candidate.finalBetaFitResult = betaFitResult;
716 }
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< std::unique_ptr< 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
void combinedTrack(long int ICH, double *pv0, double *covi, double effectiveBMAG, double *par, double *covo)
ParametersBase< TrackParametersDim, Charged > TrackParameters
float chi2PerDOF() const
chi2/ndof, return 0 if ndof == 0 or status == 0
float beta
status flag (0 = failed, 1 = ok)
1.16.1