588 if(!metJet || !
jets) {
590 <<
"MET (" << metJet <<
") or "
591 <<
"jet collection (" <<
jets <<
").");
592 return StatusCode::FAILURE;
596 ATH_MSG_ERROR(
"MET Association Helper isn't associated with a MissingETAssociationMap!");
597 return StatusCode::FAILURE;
600 ATH_MSG_WARNING(
"Requested soft track element links, but no track selection tool supplied.");
605 ATH_MSG_WARNING(
"Incomplete association map received. Cannot rebuild MET.");
606 ATH_MSG_WARNING(
"Note: METMaker should only be run on events containing at least one PV");
607 return StatusCode::SUCCESS;
610 if(doJetJVT &&
m_JvtWP ==
"None"){
616 if(!metSoftClus && !metSoftTrk) {
617 ATH_MSG_WARNING(
"Neither soft cluster nor soft track term has been supplied!");
618 return StatusCode::SUCCESS;
622 dec_constitObjLinks(*metSoftClus) = std::vector<iplink_t>(0);
624 ATH_MSG_ERROR(
"Soft cluster term provided without a core term!");
625 return StatusCode::FAILURE;
629 <<
", mpy " << coreSoftClus->
mpy()
630 <<
" sumet " << coreSoftClus->
sumet());
631 *metSoftClus += *coreSoftClus;
635 if(softConst && acc_softConst.isAvailable(*coreSoftClus)) {
636 for(
const auto& constit : acc_softConst(*coreSoftClus)) {
637 softConst->push_back(*constit);
639 ATH_MSG_DEBUG(softConst->size() <<
" soft constituents from core term");
643 dec_constitObjLinks(*metSoftTrk) = std::vector<iplink_t>(0);
645 ATH_MSG_ERROR(
"Soft track term provided without a core term!");
646 return StatusCode::FAILURE;
650 <<
", mpy " << coreSoftTrk->
mpy()
651 <<
" sumet " << coreSoftTrk->
sumet());
652 *metSoftTrk += *coreSoftTrk;
654 for(
const auto& constit : acc_softConst(*coreSoftTrk)) {
655 softConst->push_back(*constit);
657 ATH_MSG_DEBUG(softConst->size() <<
" soft constituents from trk core term");
661 dec_constitObjLinks(*metJet) = std::vector<iplink_t>(0);
662 dec_constitObjWeights(*metJet) = std::vector<float>(0);
663 std::vector<iplink_t>& uniqueLinks = dec_constitObjLinks(*metJet);
664 std::vector<float>& uniqueWeights = dec_constitObjWeights(*metJet);
665 uniqueLinks.reserve(
jets->size());
666 uniqueWeights.reserve(
jets->size());
667 std::vector<iplink_t> softJetLinks;
668 std::vector<float> softJetWeights;
669 bool originalInputs =
jets->empty() ? false : !acc_originalObject.isAvailable(*
jets->front());
682 return StatusCode::FAILURE;
694 if(!assoc || assoc->isMisc()){
701 if (acc_nominalObject.isAvailable(*
jet)){
707 ATH_MSG_ERROR(
"No nominal calibrated jet available for jet " <<
jet->index() <<
". Cannot simplify overlap removal!");
712 bool JVT_reject(
false);
713 bool isMuFSRJet(
false);
722 ATH_MSG_VERBOSE(
"Jet " << (JVT_reject ?
"fails" :
"passes") <<
" JVT selection");
729 bool caloverlap =
false;
730 caloverlap = calvec.
ce()>0;
731 ATH_MSG_DEBUG(
"Jet " <<
jet->index() <<
" is " << ( caloverlap ?
"" :
"non-") <<
"overlapping");
734 for(
const auto&
object : assoc->objects()) {
743 constjet = assoc->getAlternateConstVec();
746 double denom = (assoc->hasAlternateConstVec() ? assoc->getAlternateConstVec() :
jet->jetP4(
"JetConstitScaleMomentum")).
E();
751 double jpx = constjet.Px();
752 double jpy = constjet.Py();
753 double jpt = constjet.Pt();
754 double opx = jpx - calvec.
cpx();
755 double opy = jpy - calvec.
cpy();
761 met_muonEloss = (*metCont)[
"MuonEloss"];
763 ATH_MSG_WARNING(
"Attempted to apply muon Eloss correction, but corresponding MET term does not exist!");
764 return StatusCode::FAILURE;
768 float total_eloss(0);
770 std::vector<const xAOD::Muon*> muons_in_jet;
771 std::vector<const xAOD::Electron*> electrons_in_jet;
772 bool passJetForEl=
false;
774 if(!acc_ghostMuons.isAvailable(*
jet)){
776 return StatusCode::FAILURE;
778 for(
const auto&
el : acc_ghostMuons(*
jet)) {
780 ATH_MSG_ERROR(
"Invalid element link to ghost muon! Quitting.");
781 return StatusCode::FAILURE;
783 muons_in_jet.push_back(
static_cast<const xAOD::Muon*
>(*
el));
786 for(
const auto&
obj : assoc->objects()) {
792 if(acc_originalObject.isAvailable(*mu_test)) mu_test =
static_cast<const xAOD::Muon*
>(*acc_originalObject(*mu_test));
794 muons_in_jet.push_back(mu_test);
801 if(acc_originalObject.isAvailable(*el_test)) el_test =
static_cast<const xAOD::Electron*
>(*acc_originalObject(*el_test));
802 if(
helper.objSelected(assoc,el_test)){
803 if(el_test->pt()>90.0e3) {
804 electrons_in_jet.push_back(el_test);
812 float jet_ORtrk_sumpt = assoc->overlapTrkVec(
helper).
sumpt();
813 float jet_all_trk_pt = initialTrkMom.
sumpt();
814 float jet_unique_trk_pt = jet_all_trk_pt - jet_ORtrk_sumpt;
817 for(
const auto& elec : electrons_in_jet) {
818 el_calvec += assoc->calVec(elec);
819 el_trkvec += assoc->trkVec(elec);
821 float el_cal_pt = el_calvec.
cpt();
822 float el_trk_pt = el_trkvec.
cpt();
824 <<
" jetalltrk: " << jet_all_trk_pt
825 <<
" jetORtrk: " << jet_ORtrk_sumpt
826 <<
" electrk-jetORtrk: " << (el_trk_pt-jet_ORtrk_sumpt)
827 <<
" elec cal: " << el_cal_pt
828 <<
" jetalltrk-electrk: " << (jet_all_trk_pt-el_trk_pt)
829 <<
" jetalltrk-jetORtrk: " << (jet_all_trk_pt-jet_ORtrk_sumpt) );
833 if(el_trk_pt>1
e-9 && jet_unique_trk_pt>10.0
e3) passJetForEl=
true;
836 for(
const xAOD::Muon* mu_in_jet : muons_in_jet) {
837 if (!mu_in_jet)
continue;
838 float mu_Eloss = acc_Eloss(*mu_in_jet);
843 float mu_id_pt = mu_in_jet->trackParticle(xAOD::Muon::InnerDetectorTrackParticle) ? mu_in_jet->trackParticle(xAOD::Muon::InnerDetectorTrackParticle)->pt() : 0.;
844 float jet_trk_sumpt = acc_trksumpt.isAvailable(*
jet) && this->
getPV() ? acc_trksumpt(*
jet)[this->
getPV()->
index()] : 0.;
847 if(0.9999*mu_id_pt>jet_trk_sumpt)
848 jet_trk_sumpt+=mu_id_pt;
849 float jet_trk_N = acc_trkN.isAvailable(*
jet) && this->
getPV() ? acc_trkN(*
jet)[this->
getPV()->
index()] : 0.;
851 ATH_MSG_VERBOSE(
"Jet has pt " <<
jet->pt() <<
", trk sumpt " << jet_trk_sumpt <<
", trk N " << jet_trk_N);
861 float mu_id_pt = mu_in_jet->trackParticle(xAOD::Muon::InnerDetectorTrackParticle) ? mu_in_jet->trackParticle(xAOD::Muon::InnerDetectorTrackParticle)->pt() : 0.;
862 float jet_trk_sumpt = acc_trksumpt.isAvailable(*
jet) && this->
getPV() ? acc_trksumpt(*
jet)[this->
getPV()->
index()] : 0.;
864 if(0.9999*mu_id_pt>jet_trk_sumpt)
865 jet_trk_sumpt+=mu_id_pt;
866 float jet_trk_N = acc_trkN.isAvailable(*
jet) && this->
getPV() ? acc_trkN(*
jet)[this->
getPV()->
index()] : 0.;
869 if (acc_psf.isAvailable(*
jet)){
870 jet_psE = acc_psf(*
jet);
871 }
else if (acc_sampleE.isAvailable(*
jet)){
872 jet_psE = acc_sampleE(*
jet)[0] + acc_sampleE(*
jet)[4];
874 ATH_MSG_ERROR(
"Jet PS fraction or sampling energy must be available to calculate MET with doSetMuonJetEMScale");
875 return StatusCode::FAILURE;
880 ATH_MSG_VERBOSE(
"Jet has trk sumpt " << jet_trk_sumpt <<
", trk N " << jet_trk_N <<
", PS E " << jet_psE <<
", width " << acc_width(*
jet) <<
", emfrac " << acc_emf(*
jet));
883 ATH_MSG_VERBOSE(
"Jet is from muon -- set to EM scale and subtract Eloss.");
886 ATH_MSG_VERBOSE(
"Jet e: " << constjet.E() <<
", mu Eloss: " << mu_Eloss);
887 float elosscorr = mu_Eloss >= constjet.e() ? 0. : 1.-mu_Eloss/constjet.e();
892 ATH_MSG_VERBOSE(
" Jet eloss factor " << elosscorr <<
", final pt: " << sqrt(opx*opx+opy*opy));
899 switch(mu_in_jet->energyLossType()) {
900 case xAOD::Muon::Parametrized:
901 case xAOD::Muon::MOP:
902 case xAOD::Muon::Tail:
903 case xAOD::Muon::FSRcandidate:
904 case xAOD::Muon::NotIsolated:
908 total_eloss += mu_Eloss;
909 muons_selflags |= (1<<assoc->findIndex(mu_in_jet));
917 for(
size_t iKey = 0; iKey < assoc->sizeCal(); iKey++) {
918 bool selector = (muons_selflags & assoc->calkey()[iKey]);
919 if(
selector) mu_calovec += assoc->calVec(iKey);
923 if(
m_muEloss) mu_calovec *= std::max<float>(0.,1-(total_eloss/mu_calovec.
ce()));
928 ATH_MSG_VERBOSE(
"Jet " <<
jet->index() <<
" const pT after OR " << sqrt(opx*opx+opy*opy));
929 opx += mu_calovec.
cpx();
930 opy += mu_calovec.
cpy();
931 double opt = sqrt( opx*opx+opy*opy );
932 ATH_MSG_VERBOSE(
"Jet " <<
jet->index() <<
" const pT diff after OR readding muon clusters " <<
opt-jpt);
933 double uniquefrac = 1. - (calvec.
ce() - mu_calovec.
ce()) / constjet.E();
934 ATH_MSG_VERBOSE(
"Jet constscale px, py, pt, E = " << jpx <<
", " << jpy <<
", " << jpt <<
", " << constjet.E() );
936 ATH_MSG_VERBOSE(
"Jet OR px, py, pt, E = " << opx <<
", " << opy <<
", " <<
opt <<
", " << constjet.E() - calvec.
ce() );
940 ATH_MSG_ERROR(
"Attempted to apply muon Eloss correction, but corresponding MET term does not exist!");
941 return StatusCode::FAILURE;
943 met_muonEloss->add(opx,opy,
opt);
947 if(selected && !JVT_reject) {
951 if (!tracksForHardJets) {
953 metJet->
add(jpx,jpy,jpt);
961 if(!tracksForHardJets) {
966 double jesF =
jet->pt() / jpt;
967 metJet->
add(opx*jesF,opy*jesF,
opt*jesF);
971 metJet->
add(uniquefrac*jpx,uniquefrac*jpy,uniquefrac*jpt);
976 metJet->
add(uniquefrac*
jet->px(),uniquefrac*
jet->py(),uniquefrac*
jet->pt());
995 uniqueLinks.push_back( jetLink );
996 uniqueWeights.push_back( uniquefrac );
998 if(metSoftClus && !JVT_reject) {
1002 softJetLinks.push_back( jetLink );
1003 softJetWeights.push_back( uniquefrac );
1004 metSoftClus->
add(opx,opy,
opt);
1013 for(
size_t iConst=0; iConst<
jet->numConstituents(); ++iConst) {
1014 const IParticle* constit =
jet->rawConstituent(iConst);
1015 softConst->push_back(constit);
1021 if(!metSoftTrk || (hardJet && !tracksForHardJets))
continue;
1027 if(jettrkvec.
ce()>1
e-9) {
1028 jpx = jettrkvec.
cpx();
1029 jpy = jettrkvec.
cpy();
1030 jpt = jettrkvec.
sumpt();
1031 jettrkvec -= trkvec;
1032 opx = jettrkvec.
cpx();
1033 opy = jettrkvec.
cpy();
1035 ATH_MSG_VERBOSE(
"Jet track px, py, sumpt = " << jpx <<
", " << jpy <<
", " << jpt );
1038 opx = opy =
opt = 0;
1041 if (hardJet) metJet->
add(opx,opy,
opt);
1043 metSoftTrk->
add(opx,opy,
opt);
1046 softJetLinks.push_back( jetLink );
1047 softJetWeights.push_back( uniquefrac );
1056 std::vector<const IParticle*> jettracks;
1058 for(
size_t iConst=0; iConst<jettracks.size(); ++iConst) {
1066 ATH_MSG_DEBUG(
"Number of selected jets: " << dec_constitObjLinks(*metJet).size());
1069 dec_constitObjLinks(*metSoftTrk) = softJetLinks;
1070 ATH_MSG_DEBUG(
"Number of softtrk jets: " << dec_constitObjLinks(*metSoftTrk).size());
1074 dec_constitObjLinks(*metSoftClus) = softJetLinks;
1075 ATH_MSG_DEBUG(
"Number of softclus jets: " << dec_constitObjLinks(*metSoftClus).size());
1078 if(softConst)
ATH_MSG_DEBUG(softConst->size() <<
" soft constituents from core term + jets");
1081 if(!assoc)
return StatusCode::SUCCESS;
1087 double opx = trkvec.
cpx();
1088 double opy = trkvec.
cpy();
1089 double osumpt = trkvec.
sumpt();
1090 ATH_MSG_VERBOSE(
"Misc track px, py, sumpt = " << opx <<
", " << opy <<
", " << osumpt );
1091 metSoftTrk->
add(opx,opy,osumpt);
1093 <<
", mpy " << metSoftTrk->
mpy()
1094 <<
" sumet " << metSoftTrk->
sumet());
1100 float total_eloss(0.);
1103 double opx = calvec.
cpx();
1104 double opy = calvec.
cpy();
1105 double osumpt = calvec.
sumpt();
1106 for(
const auto&
obj : assoc->objects()) {
1109 if(acc_originalObject.isAvailable(*mu_test)) mu_test =
static_cast<const xAOD::Muon*
>(*acc_originalObject(*mu_test));
1111 float mu_Eloss = acc_Eloss(*mu_test);
1112 switch(mu_test->energyLossType()) {
1113 case xAOD::Muon::Parametrized:
1114 case xAOD::Muon::MOP:
1115 case xAOD::Muon::Tail:
1116 case xAOD::Muon::FSRcandidate:
1117 case xAOD::Muon::NotIsolated:
1121 total_eloss += mu_Eloss;
1122 muons_selflags |= (1<<assoc->findIndex(mu_test));
1132 for(
size_t iKey = 0; iKey < assoc->sizeCal(); iKey++) {
1133 bool selector = (muons_selflags & assoc->calkey()[iKey]);
1135 <<
" this calvec E: " << assoc->calVec(iKey).ce());
1136 if(
selector) mu_calovec += assoc->calVec(iKey);
1139 mu_calovec *= std::max<float>(0.,1-(total_eloss/mu_calovec.
ce()));
1140 opx += mu_calovec.
cpx();
1141 opy += mu_calovec.
cpy();
1142 osumpt += mu_calovec.
sumpt();
1146 ATH_MSG_VERBOSE(
"Misc cluster px, py, sumpt = " << opx <<
", " << opy <<
", " << osumpt );
1147 metSoftClus->
add(opx,opy,osumpt);
1149 <<
", mpy " << metSoftClus->
mpy()
1150 <<
" sumet " << metSoftClus->
sumet());
1153 return StatusCode::SUCCESS;
1.8.18