execute execute the combined muon FEX.
Execute() method.
457{
458 using namespace xAOD;
459
460
461
462 auto timer = Monitored::Timer(
"TIME_execute");
463
464 auto ptMS = Monitored::Scalar("PtMS", -9999.);
465 auto etaMS = Monitored::Scalar(
"EtaMS", -9999.);
466 auto phiMS = Monitored::Scalar("PhiMS", -9999.);
467 auto zetaMS = Monitored::Scalar("ZetaMS", -9999.);
468
469 auto ptID = Monitored::Scalar("PtID", -9999.);
470 auto etaID = Monitored::Scalar("EtaID", -9999.);
471 auto phiID = Monitored::Scalar(
"PhiID", -9999.);
472 auto zetaID = Monitored::Scalar("ZetaID", -9999.);
473
474 auto ptMC = Monitored::Scalar("PtMC", -9999.);
475 auto dEta = Monitored::Scalar(
"DEta", -9999.);
476 auto dPhi = Monitored::Scalar(
"DPhi", -9999.);
477 auto dZeta = Monitored::Scalar("DZeta", -9999.);
478 auto dR = Monitored::Scalar("DeltaR", -9999.);
479
480 auto ptFL = Monitored::Scalar("PtFL", -9999.);
481 auto etaFL = Monitored::Scalar("EtaFL", -9999.);
482 auto phiFL = Monitored::Scalar("PhiFL", -9999.);
483
484 auto efficiency = Monitored::Scalar<int>(
"Efficiency", -1);
485 auto StrategyMC = Monitored::Scalar<int>("StrategyMC", -1);
486 auto ErrorFlagMC = Monitored::Scalar<int>("ErrorFlagMC", 0);
487 auto MatchFlagMC = Monitored::Scalar<int>("MatchFlagMC", 0);
488
489
490 auto mon = Monitored::Group(
m_monTool, timer, ptMS, etaMS, phiMS, zetaMS,
491 ptID, etaID, phiID, zetaID,
492 ptMC, dEta, dPhi, dZeta, dR,
493 ptFL, etaFL, phiFL,
494 efficiency, StrategyMC, ErrorFlagMC, MatchFlagMC);
495
496
498
499
502
504 const AtlasFieldCacheCondObj* fieldCondObj{*readHandle};
505 if (fieldCondObj == nullptr) {
507 return StatusCode::FAILURE;
508 }
509 MagField::AtlasFieldCache fieldCache;
511
514 ATH_MSG_DEBUG(
"=========== Magnetic Field Status ========== " );
517 ATH_MSG_DEBUG(
" ---> Solenoid : " << ((solenoidOn) ?
"ON" :
"OFF") );
518 ATH_MSG_DEBUG(
" ---> Toroid : " << ((toroidOn) ?
"ON" :
"OFF") );
519
520
521
522 int usealgo = 0;
523 if (solenoidOn) {
524 if (toroidOn) {
526 else usealgo = 1;
527 } else {
528 usealgo = 2;
529 }
530 } else {
531 if (toroidOn) usealgo = 3;
532 else usealgo = 4;
533 }
534
536
537 ATH_MSG_DEBUG(
"=========== Matching windows g4 ========== " );
546
547
549 ATH_CHECK( muonCBColl.record (std::make_unique<xAOD::L2CombinedMuonContainer>(),
550 std::make_unique<xAOD::L2CombinedMuonAuxContainer>()) );
551
552
554
555
556
557 if (muonColl->size() == 0) {
558
559
560
561
562
563 ATH_MSG_DEBUG(
" L2StandAloneMuonContainer empty -> stop processing RoI, no match" );
564 return StatusCode::SUCCESS;
565
566 }
567
568
569
570
571
572 for(
uint i_muonSA = 0; i_muonSA < muonColl->size(); i_muonSA++){
573
575
580
581
582 ElementLink<xAOD::L2StandAloneMuonContainer> muonSAEL(*muonColl, i_muonSA);
584
585
586 if (muonSA->
pt() == 0.) {
587 ErrorFlagMC = 1;
588 if (usealgo == 2 || usealgo == 4) {
589 ATH_MSG_DEBUG(
" L2StandAloneMuon pt = 0 --> using angular match" );
590
591 } else {
592 ATH_MSG_DEBUG(
" L2StandAloneMuon pt = 0 --> stop processing RoI, no match" );
594 muonCBColl->push_back(muonCB);
595 return StatusCode::SUCCESS;
596 }
597 }
598
599
600
601
602
606 double eta_ms = -999.;
607 double phi_ms = -999.;
608 double zeta_ms = -999.;
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
633 eta_ms = muonSA->
etaMS();
634 phi_ms = muonSA->
phiMS();
635 zeta_ms = muonSA->
zMS();
636
638 <<
" / eta = " <<
eta
639 <<
" / phi = " <<
phi
640 << " / etaMS = " << eta_ms
641 << " / phiMS = " << phi_ms
642 << " / zMS = " << zeta_ms);
643
644
645
646
647
649 if (!idTrackParticles.isValid()){
650 ATH_MSG_DEBUG(
" Failed to get xAOD::TrackParticleContainer --> no match" );
651 ErrorFlagMC = 2;
653 muonCBColl->push_back(muonCB);
654 return StatusCode::SUCCESS;
655 }
656 if (idTrackParticles->size() < 1){
657 ATH_MSG_DEBUG(
" xAOD::TrackParticleContainer has 0 tracks --> no match" );
658 ErrorFlagMC = 2;
660 muonCBColl->push_back(muonCB);
661 return StatusCode::SUCCESS;
662 }
663
664 ATH_MSG_DEBUG(
" Found xAOD::TrackParticleContainer with size: " << idTrackParticles->size() );
665
666
667 double ptinv_comb = 0.;
668 double ptres_comb = 0.001;
669 double chi2_comb = 1.0e33;
670 int ndof_comb = -1;
671 double best_dr = 1.0e33;
672 double best_deta = 1.0e33;
673 double best_dphi = 1.0e33;
674 bool has_match = false;
675 int imatch = -1;
676 int matched_trk_idx = -1;
677 int matched_trk_idx_tmp = -1;
678
679
680
681
682
683 ExtrapolationResult extr{};
686 }
687
688 for(const auto trkit:(*idTrackParticles)) {
689
690 matched_trk_idx_tmp++;
691 double ptinv_tmp = 0.;
692 double ptres_tmp = 0.001;
693 double deta_tmp = 0.;
694 double dphi_tmp = 0.;
695 double dr_tmp = 0.;
696 double chi2_tmp = 0.;
697 int ndof_tmp = 0;
698 bool has_match_tmp = false;
699 int imatch_tmp = -1;
700
701
702 double phi_id = (trkit)->
phi();
703 double eta_id = (trkit)->
eta();
704 double qoverp_id = (trkit)->
qOverP();
705 double q_id = ((trkit)->qOverP() > 0 ? 1.0 : -1.0);
706 double pt_id = (trkit)->
pt() * q_id;
707
711 double tanthetaov2 = std::fabs(std::exp(-
eta));
712 double e_eta_id = std::fabs(0.5 * (1./tanthetaov2 + tanthetaov2) * e_theta_id);
713
714 double e_qoverpt_id = e_qoverp_id;
715 double theta_id = (trkit)->
theta();
716 if (std::sin(theta_id) != 0) e_qoverpt_id /= std::fabs(std::sin(theta_id));
717
719 << " with pt (GeV) = " << pt_id / Gaudi::Units::GeV
720 << ", q = " << q_id
721 << ", eta = " << eta_id
722 << ", phi = " << phi_id
723 << ", th = " << theta_id
724 << ", ephi = " << e_phi_id
725 << ", eth = " << e_theta_id
726 << ", eeta = " << e_eta_id
727 << ", ip = " << qoverp_id
728 << ", eip = " << e_qoverp_id
729 << ", eipt = " << e_qoverpt_id );
730
731 if (usealgo != 3) {
732 if ((std::fabs(pt_id) / Gaudi::Units::GeV) <
m_PtMinTrk)
continue;
733 }
735
737
738 if (usealgo > 0) {
739
740
741
742
743 imatch_tmp =
drptMatch(muonSA, pt_id, eta_id, phi_id, usealgo, ptinv_tmp, ptres_tmp, deta_tmp, dphi_tmp, chi2_tmp);
744
745 if (imatch_tmp == 0) has_match_tmp = true;
746 } else {
748 imatch_tmp =
mfMatch(muonSA, eta_id, phi_id, pt_id, q_id, ptinv_tmp, ptres_tmp, deta_tmp, dphi_tmp, chi2_tmp, ndof_tmp);
749 if (imatch_tmp == 0) has_match_tmp = true;
750 } else {
751 imatch_tmp =
g4Match(extr, eta_id, phi_id, pt_id, q_id, e_eta_id, e_phi_id, e_qoverpt_id, ptinv_tmp, ptres_tmp, deta_tmp, dphi_tmp, chi2_tmp, ndof_tmp);
752 if (imatch_tmp == 0) has_match_tmp = true;
753 }
754 }
755
756 imatch = imatch_tmp;
757 if (!has_match_tmp) continue;
758
759
760 dr_tmp = std::sqrt(deta_tmp * deta_tmp + dphi_tmp * dphi_tmp);
761
762 if (chi2_tmp <= chi2_comb) {
763 has_match = true;
764 chi2_comb = chi2_tmp;
765 ndof_comb = ndof_tmp;
766 ptinv_comb = ptinv_tmp;
767 ptres_comb = ptres_tmp;
768 best_deta = deta_tmp;
769 best_dphi = dphi_tmp;
770 best_dr = dr_tmp;
771 imatch = imatch_tmp;
772 matched_trk_idx = matched_trk_idx_tmp;
773 }
774
775 }
776
777
778 if (usealgo == 0 && std::fabs(pt) >= 6.) {
782 zetaMS = zeta_ms;
786 if (ptMS > 100.) ptMS = 101.5;
787 if (ptMS < -100.) ptMS = -101.5;
788 if (ptFL > 100.) ptFL = 101.5;
789 if (ptFL < -100.) ptFL = -101.5;
790 }
791
792 if (!has_match) {
793 ErrorFlagMC = 3;
794 MatchFlagMC = imatch;
795 if (usealgo == 0 && std::fabs(pt) >= 6.)
efficiency = 0;
799 muonCBColl->push_back(muonCB);
800 return StatusCode::SUCCESS;
801 }
802
803
804 ElementLink<xAOD::TrackParticleContainer> idtrkEL(*idTrackParticles, matched_trk_idx);
806
808 double q_id = ((muonCB->
idTrack()->
qOverP()) > 0 ? 1.0 : -1.0);
811
812 double zPos_id = muonCB->
idTrack()->
z0();
813
815
816
817 MatchFlagMC = imatch;
818 ptFL = -9999.;
819 etaFL = -9999.;
820 phiFL = -9999.;
821 if (usealgo == 0 && std::fabs(pt) >= 6.) {
823 ptID = pt_id / Gaudi::Units::GeV;
824 etaID = eta_id;
826 zetaID = zPos_id;
827 ptMC = 1. / (ptinv_comb * Gaudi::Units::GeV);
828 dZeta = zeta_ms - zPos_id;
831 dR = best_dr;
832
833 if (ptMC > 100.) ptMC = 101.5;
834 if (ptMC < -100.) ptMC = -101.5;
835 if (ptID > 100.) ptID = 101.5;
836 if (ptID < -100.) ptID = -101.5;
837 }
838
840
841
843 << " usealgo / IdPt (GeV) / muonPt (GeV) / CombPt (GeV) / chi2 / ndof: " << " / " << usealgo << " / " << pt_id*q_id / Gaudi::Units::GeV << " / " << prt_pt << " / " << 1. / ptinv_comb / Gaudi::Units::GeV << " / " << chi2_comb << " / " << ndof_comb );
844
845 muonCB->
setPt(std::fabs(1. / ptinv_comb));
848
849 float mcq = -1.0;
850 if (ptinv_comb > 0) mcq = 1.0;
851
853
854 float mcresu = std::fabs(ptres_comb / (ptinv_comb * ptinv_comb));
855 ATH_MSG_DEBUG(
" SigmaPt (GeV) is: " << mcresu / Gaudi::Units::GeV );
857
860
861 muonCBColl->push_back(muonCB);
862 }
863 return StatusCode::SUCCESS;
864}
Scalar theta() const
theta method
#define ATH_CHECK
Evaluate an expression and check for errors.
void getInitializedCache(MagField::AtlasFieldCache &cache) const
get B field cache for evaluation as a function of 2-d or 3-d position.
bool solenoidOn() const
status of the magnets
void makePrivateStore()
Create a new (empty) private store for this object.
Gaudi::Property< double > m_Chi2Weight_g4
Scale factor for the Chi2 1/pt resolutions (MS and ID) (Barrel)
Gaudi::Property< double > m_WinEta_g4
Number of sigmas for the Eta matching window LUT backextrapolator (Barrel)
Gaudi::Property< double > m_WeightPhi_g4
Scale factor for the Phi matching window in LUT backextrapolator.
SG::ReadHandleKey< xAOD::TrackParticleContainer > m_TrackParticleContainerKey
Gaudi::Property< bool > m_assumeSolenoidOff
flag to assume B_Solenoid=0 anyway
int mfMatch(const xAOD::L2StandAloneMuon *feature, double, double, double, double, double &, double &, double &, double &, double &, int &) const
Gaudi::Property< double > m_WeightEta_g4
Scale factor for the Eta matching window in LUT backextrapolator.
SG::WriteHandleKey< xAOD::L2CombinedMuonContainer > m_outputCBmuonCollKey
int g4Match(const ExtrapolationResult &extr, double, double, double, double, double, double, double, double &, double &, double &, double &, double &, int &) const
SG::ReadCondHandleKey< AtlasFieldCacheCondObj > m_fieldCacheCondObjInputKey
Gaudi::Property< double > m_PtMinTrk
Min Pt to select the ID track for matching.
Gaudi::Property< double > m_WinPhi_g4
Number of sigmas for the Phi matching window LUT backextrapolator (Barrel)
ExtrapolationResult getExtrapolatedMuon(const EventContext &ctx, const xAOD::L2StandAloneMuon *feature) const
Gaudi::Property< double > m_EtaMaxTrk
Max abs(eta) to select the ID track for matching.
Gaudi::Property< bool > m_useBackExtrapolatorG4
flag to switch between G4 and LUT based back-extrapolation
Gaudi::Property< double > m_WinPhi_EC_g4
Number of sigmas for the Phi matching window LUT backextrapolator (EndCaps)
Gaudi::Property< int > m_AlgoStrategy
muComb matching strategy: 0: auto select best option 1: simplified R,(Pt) matching
ToolHandle< GenericMonitoringTool > m_monTool
Gaudi::Property< bool > m_assumeToroidOff
flag to assume B_Toroid=0 anyway
Gaudi::Property< double > m_WinEta_EC_g4
Number of sigmas for the Eta matching window LUT backextrapolator (EndCaps)
SG::ReadHandleKey< xAOD::L2StandAloneMuonContainer > m_muonCollKey
void setPt(float pt)
Set the transverse momentum ( ) of the muon.
void setSigmaPt(float value)
set sigma combined Pt
void setStrategy(int value)
set algorithm strategy flag
void setErrorFlag(int value)
set algorithm error flag
void setMatchFlag(int value)
set algorithm match flag
void setIdTrackLink(const ElementLink< xAOD::TrackParticleContainer > &link)
set ID track used to make the CB muon
void setMuSATrackLink(const ElementLink< xAOD::L2StandAloneMuonContainer > &link)
set SA muon used to make the CB muon
void setPhi(float phi)
Set the azimuthal angle ( ) of the muon.
void setEta(float eta)
Set the pseudorapidity ( ) of the muon.
void setCharge(float value)
set seeding muon charge
const xAOD::TrackParticle * idTrack() const
Get the ID track as a bare pointer.
virtual double phi() const
The azimuthal angle ( ) of the particle.
float zMS() const
Get the Z at muon spectrometer.
virtual double eta() const
The pseudorapidity ( ) of the particle.
float z0() const
Returns the parameter.
virtual double phi() const override final
The azimuthal angle ( ) of the particle (has range to .)
virtual double pt() const override final
The transverse momentum ( ) of the particle.
virtual double eta() const override final
The pseudorapidity ( ) of the particle.
void efficiency(std::vector< double > &bins, std::vector< double > &values, const std::vector< std::string > &files, const std::string &histname, const std::string &tplotname, const std::string &label="")
SG::ReadCondHandle< T > makeHandle(const SG::ReadCondHandleKey< T > &key, const EventContext &ctx=Gaudi::Hive::currentContext())
bool dPhi(const xAOD::TauJet &tau, const xAOD::CaloVertexedTopoCluster &cluster, float &out)
bool dEta(const xAOD::TauJet &tau, const xAOD::CaloVertexedTopoCluster &cluster, float &out)
timer(name, disabled=False)
xAOD::ParametersCovMatrix_t definingParametersCovMatrix(std::span< const float > covMatrixDiag, std::span< const float > covMatrixOffDiag, bool &valid)
L2CombinedMuon_v1 L2CombinedMuon
Define the latest version of the muon CB class.
L2StandAloneMuon_v2 L2StandAloneMuon
Define the latest version of the muon SA class.
1.14.0