587 if(!metJet || !
jets) {
589 <<
"MET (" << metJet <<
") or "
590 <<
"jet collection (" <<
jets <<
").");
591 return StatusCode::FAILURE;
595 ATH_MSG_ERROR(
"MET Association Helper isn't associated with a MissingETAssociationMap!");
596 return StatusCode::FAILURE;
599 ATH_MSG_WARNING(
"Requested soft track element links, but no track selection tool supplied.");
604 ATH_MSG_WARNING(
"Incomplete association map received. Cannot rebuild MET.");
605 ATH_MSG_WARNING(
"Note: METMaker should only be run on events containing at least one PV");
606 return StatusCode::SUCCESS;
609 if(doJetJVT &&
m_JvtWP ==
"None"){
615 if(!metSoftClus && !metSoftTrk) {
616 ATH_MSG_WARNING(
"Neither soft cluster nor soft track term has been supplied!");
617 return StatusCode::SUCCESS;
621 dec_constitObjLinks(*metSoftClus) = std::vector<iplink_t>(0);
623 ATH_MSG_ERROR(
"Soft cluster term provided without a core term!");
624 return StatusCode::FAILURE;
628 <<
", mpy " << coreSoftClus->
mpy()
629 <<
" sumet " << coreSoftClus->
sumet());
630 *metSoftClus += *coreSoftClus;
634 if(softConst && acc_softConst.isAvailable(*coreSoftClus)) {
635 for(
const auto& constit : acc_softConst(*coreSoftClus)) {
636 softConst->push_back(*constit);
638 ATH_MSG_DEBUG(softConst->size() <<
" soft constituents from core term");
642 dec_constitObjLinks(*metSoftTrk) = std::vector<iplink_t>(0);
644 ATH_MSG_ERROR(
"Soft track term provided without a core term!");
645 return StatusCode::FAILURE;
649 <<
", mpy " << coreSoftTrk->
mpy()
650 <<
" sumet " << coreSoftTrk->
sumet());
651 *metSoftTrk += *coreSoftTrk;
653 for(
const auto& constit : acc_softConst(*coreSoftTrk)) {
654 softConst->push_back(*constit);
656 ATH_MSG_DEBUG(softConst->size() <<
" soft constituents from trk core term");
660 dec_constitObjLinks(*metJet) = std::vector<iplink_t>(0);
661 dec_constitObjWeights(*metJet) = std::vector<float>(0);
662 std::vector<iplink_t>& uniqueLinks = dec_constitObjLinks(*metJet);
663 std::vector<float>& uniqueWeights = dec_constitObjWeights(*metJet);
664 uniqueLinks.reserve(
jets->size());
665 uniqueWeights.reserve(
jets->size());
666 std::vector<iplink_t> softJetLinks;
667 std::vector<float> softJetWeights;
668 bool originalInputs =
jets->empty() ? false : !acc_originalObject.isAvailable(*
jets->front());
681 return StatusCode::FAILURE;
693 if(!assoc || assoc->isMisc()){
700 if (acc_nominalObject.isAvailable(*
jet)){
706 ATH_MSG_ERROR(
"No nominal calibrated jet available for jet " <<
jet->index() <<
". Cannot simplify overlap removal!");
711 bool JVT_reject(
false);
712 bool isMuFSRJet(
false);
721 ATH_MSG_VERBOSE(
"Jet " << (JVT_reject ?
"fails" :
"passes") <<
" JVT selection");
728 bool caloverlap =
false;
729 caloverlap = calvec.
ce()>0;
730 ATH_MSG_DEBUG(
"Jet " <<
jet->index() <<
" is " << ( caloverlap ?
"" :
"non-") <<
"overlapping");
733 for(
const auto&
object : assoc->objects()) {
742 constjet = assoc->getAlternateConstVec();
745 double denom = (assoc->hasAlternateConstVec() ? assoc->getAlternateConstVec() :
jet->jetP4(
"JetConstitScaleMomentum")).
E();
750 double jpx = constjet.Px();
751 double jpy = constjet.Py();
752 double jpt = constjet.Pt();
753 double opx = jpx - calvec.
cpx();
754 double opy = jpy - calvec.
cpy();
760 met_muonEloss = (*metCont)[
"MuonEloss"];
762 ATH_MSG_WARNING(
"Attempted to apply muon Eloss correction, but corresponding MET term does not exist!");
763 return StatusCode::FAILURE;
767 float total_eloss(0);
769 std::vector<const xAOD::Muon*> muons_in_jet;
770 std::vector<const xAOD::Electron*> electrons_in_jet;
771 bool passJetForEl=
false;
773 if(!acc_ghostMuons.isAvailable(*
jet)){
775 return StatusCode::FAILURE;
777 for(
const auto&
el : acc_ghostMuons(*
jet)) {
779 ATH_MSG_ERROR(
"Invalid element link to ghost muon! Quitting.");
780 return StatusCode::FAILURE;
782 muons_in_jet.push_back(
static_cast<const xAOD::Muon*
>(*
el));
785 for(
const auto&
obj : assoc->objects()) {
791 if(acc_originalObject.isAvailable(*mu_test)) mu_test =
static_cast<const xAOD::Muon*
>(*acc_originalObject(*mu_test));
793 muons_in_jet.push_back(mu_test);
800 if(acc_originalObject.isAvailable(*el_test)) el_test =
static_cast<const xAOD::Electron*
>(*acc_originalObject(*el_test));
801 if(
helper.objSelected(assoc,el_test)){
802 if(el_test->pt()>90.0e3) {
803 electrons_in_jet.push_back(el_test);
811 float jet_ORtrk_sumpt = assoc->overlapTrkVec(
helper).
sumpt();
812 float jet_all_trk_pt = initialTrkMom.
sumpt();
813 float jet_unique_trk_pt = jet_all_trk_pt - jet_ORtrk_sumpt;
816 for(
const auto& elec : electrons_in_jet) {
817 el_calvec += assoc->calVec(elec);
818 el_trkvec += assoc->trkVec(elec);
820 float el_cal_pt = el_calvec.
cpt();
821 float el_trk_pt = el_trkvec.
cpt();
823 <<
" jetalltrk: " << jet_all_trk_pt
824 <<
" jetORtrk: " << jet_ORtrk_sumpt
825 <<
" electrk-jetORtrk: " << (el_trk_pt-jet_ORtrk_sumpt)
826 <<
" elec cal: " << el_cal_pt
827 <<
" jetalltrk-electrk: " << (jet_all_trk_pt-el_trk_pt)
828 <<
" jetalltrk-jetORtrk: " << (jet_all_trk_pt-jet_ORtrk_sumpt) );
832 if(el_trk_pt>1
e-9 && jet_unique_trk_pt>10.0
e3) passJetForEl=
true;
835 for(
const xAOD::Muon* mu_in_jet : muons_in_jet) {
836 if (!mu_in_jet)
continue;
837 float mu_Eloss = acc_Eloss(*mu_in_jet);
842 float mu_id_pt = mu_in_jet->trackParticle(xAOD::Muon::InnerDetectorTrackParticle) ? mu_in_jet->trackParticle(xAOD::Muon::InnerDetectorTrackParticle)->pt() : 0.;
843 float jet_trk_sumpt = acc_trksumpt.isAvailable(*
jet) && this->
getPV() ? acc_trksumpt(*
jet)[this->
getPV()->
index()] : 0.;
846 if(0.9999*mu_id_pt>jet_trk_sumpt)
847 jet_trk_sumpt+=mu_id_pt;
848 float jet_trk_N = acc_trkN.isAvailable(*
jet) && this->
getPV() ? acc_trkN(*
jet)[this->
getPV()->
index()] : 0.;
850 ATH_MSG_VERBOSE(
"Jet has pt " <<
jet->pt() <<
", trk sumpt " << jet_trk_sumpt <<
", trk N " << jet_trk_N);
860 float mu_id_pt = mu_in_jet->trackParticle(xAOD::Muon::InnerDetectorTrackParticle) ? mu_in_jet->trackParticle(xAOD::Muon::InnerDetectorTrackParticle)->pt() : 0.;
861 float jet_trk_sumpt = acc_trksumpt.isAvailable(*
jet) && this->
getPV() ? acc_trksumpt(*
jet)[this->
getPV()->
index()] : 0.;
863 if(0.9999*mu_id_pt>jet_trk_sumpt)
864 jet_trk_sumpt+=mu_id_pt;
865 float jet_trk_N = acc_trkN.isAvailable(*
jet) && this->
getPV() ? acc_trkN(*
jet)[this->
getPV()->
index()] : 0.;
868 if (acc_psf.isAvailable(*
jet)){
869 jet_psE = acc_psf(*
jet);
870 }
else if (acc_sampleE.isAvailable(*
jet)){
871 jet_psE = acc_sampleE(*
jet)[0] + acc_sampleE(*
jet)[4];
873 ATH_MSG_ERROR(
"Jet PS fraction or sampling energy must be available to calculate MET with doSetMuonJetEMScale");
874 return StatusCode::FAILURE;
879 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));
882 ATH_MSG_VERBOSE(
"Jet is from muon -- set to EM scale and subtract Eloss.");
885 ATH_MSG_VERBOSE(
"Jet e: " << constjet.E() <<
", mu Eloss: " << mu_Eloss);
886 float elosscorr = mu_Eloss >= constjet.e() ? 0. : 1.-mu_Eloss/constjet.e();
891 ATH_MSG_VERBOSE(
" Jet eloss factor " << elosscorr <<
", final pt: " << sqrt(opx*opx+opy*opy));
898 switch(mu_in_jet->energyLossType()) {
899 case xAOD::Muon::Parametrized:
900 case xAOD::Muon::MOP:
901 case xAOD::Muon::Tail:
902 case xAOD::Muon::FSRcandidate:
903 case xAOD::Muon::NotIsolated:
907 total_eloss += mu_Eloss;
908 muons_selflags |= (1<<assoc->findIndex(mu_in_jet));
916 for(
size_t iKey = 0; iKey < assoc->sizeCal(); iKey++) {
917 bool selector = (muons_selflags & assoc->calkey()[iKey]);
918 if(
selector) mu_calovec += assoc->calVec(iKey);
922 if(
m_muEloss) mu_calovec *= std::max<float>(0.,1-(total_eloss/mu_calovec.
ce()));
927 ATH_MSG_VERBOSE(
"Jet " <<
jet->index() <<
" const pT after OR " << sqrt(opx*opx+opy*opy));
928 opx += mu_calovec.
cpx();
929 opy += mu_calovec.
cpy();
930 double opt = sqrt( opx*opx+opy*opy );
931 ATH_MSG_VERBOSE(
"Jet " <<
jet->index() <<
" const pT diff after OR readding muon clusters " <<
opt-jpt);
932 double uniquefrac = 1. - (calvec.
ce() - mu_calovec.
ce()) / constjet.E();
933 ATH_MSG_VERBOSE(
"Jet constscale px, py, pt, E = " << jpx <<
", " << jpy <<
", " << jpt <<
", " << constjet.E() );
935 ATH_MSG_VERBOSE(
"Jet OR px, py, pt, E = " << opx <<
", " << opy <<
", " <<
opt <<
", " << constjet.E() - calvec.
ce() );
939 ATH_MSG_ERROR(
"Attempted to apply muon Eloss correction, but corresponding MET term does not exist!");
940 return StatusCode::FAILURE;
942 met_muonEloss->add(opx,opy,
opt);
946 if(selected && !JVT_reject) {
950 if (!tracksForHardJets) {
952 metJet->
add(jpx,jpy,jpt);
960 if(!tracksForHardJets) {
965 double jesF =
jet->pt() / jpt;
966 metJet->
add(opx*jesF,opy*jesF,
opt*jesF);
970 metJet->
add(uniquefrac*jpx,uniquefrac*jpy,uniquefrac*jpt);
975 metJet->
add(uniquefrac*
jet->px(),uniquefrac*
jet->py(),uniquefrac*
jet->pt());
994 uniqueLinks.push_back( jetLink );
995 uniqueWeights.push_back( uniquefrac );
997 if(metSoftClus && !JVT_reject) {
1001 softJetLinks.push_back( jetLink );
1002 softJetWeights.push_back( uniquefrac );
1003 metSoftClus->
add(opx,opy,
opt);
1012 for(
size_t iConst=0; iConst<
jet->numConstituents(); ++iConst) {
1013 const IParticle* constit =
jet->rawConstituent(iConst);
1014 softConst->push_back(constit);
1020 if(!metSoftTrk || (hardJet && !tracksForHardJets))
continue;
1026 if(jettrkvec.
ce()>1
e-9) {
1027 jpx = jettrkvec.
cpx();
1028 jpy = jettrkvec.
cpy();
1029 jpt = jettrkvec.
sumpt();
1030 jettrkvec -= trkvec;
1031 opx = jettrkvec.
cpx();
1032 opy = jettrkvec.
cpy();
1034 ATH_MSG_VERBOSE(
"Jet track px, py, sumpt = " << jpx <<
", " << jpy <<
", " << jpt );
1037 opx = opy =
opt = 0;
1040 if (hardJet) metJet->
add(opx,opy,
opt);
1042 metSoftTrk->
add(opx,opy,
opt);
1045 softJetLinks.push_back( jetLink );
1046 softJetWeights.push_back( uniquefrac );
1055 std::vector<const IParticle*> jettracks;
1057 for(
size_t iConst=0; iConst<jettracks.size(); ++iConst) {
1065 ATH_MSG_DEBUG(
"Number of selected jets: " << dec_constitObjLinks(*metJet).size());
1068 dec_constitObjLinks(*metSoftTrk) = softJetLinks;
1069 ATH_MSG_DEBUG(
"Number of softtrk jets: " << dec_constitObjLinks(*metSoftTrk).size());
1073 dec_constitObjLinks(*metSoftClus) = softJetLinks;
1074 ATH_MSG_DEBUG(
"Number of softclus jets: " << dec_constitObjLinks(*metSoftClus).size());
1077 if(softConst)
ATH_MSG_DEBUG(softConst->size() <<
" soft constituents from core term + jets");
1080 if(!assoc)
return StatusCode::SUCCESS;
1086 double opx = trkvec.
cpx();
1087 double opy = trkvec.
cpy();
1088 double osumpt = trkvec.
sumpt();
1089 ATH_MSG_VERBOSE(
"Misc track px, py, sumpt = " << opx <<
", " << opy <<
", " << osumpt );
1090 metSoftTrk->
add(opx,opy,osumpt);
1092 <<
", mpy " << metSoftTrk->
mpy()
1093 <<
" sumet " << metSoftTrk->
sumet());
1099 float total_eloss(0.);
1102 double opx = calvec.
cpx();
1103 double opy = calvec.
cpy();
1104 double osumpt = calvec.
sumpt();
1105 for(
const auto&
obj : assoc->objects()) {
1108 if(acc_originalObject.isAvailable(*mu_test)) mu_test =
static_cast<const xAOD::Muon*
>(*acc_originalObject(*mu_test));
1110 float mu_Eloss = acc_Eloss(*mu_test);
1111 switch(mu_test->energyLossType()) {
1112 case xAOD::Muon::Parametrized:
1113 case xAOD::Muon::MOP:
1114 case xAOD::Muon::Tail:
1115 case xAOD::Muon::FSRcandidate:
1116 case xAOD::Muon::NotIsolated:
1120 total_eloss += mu_Eloss;
1121 muons_selflags |= (1<<assoc->findIndex(mu_test));
1131 for(
size_t iKey = 0; iKey < assoc->sizeCal(); iKey++) {
1132 bool selector = (muons_selflags & assoc->calkey()[iKey]);
1134 <<
" this calvec E: " << assoc->calVec(iKey).ce());
1135 if(
selector) mu_calovec += assoc->calVec(iKey);
1138 mu_calovec *= std::max<float>(0.,1-(total_eloss/mu_calovec.
ce()));
1139 opx += mu_calovec.
cpx();
1140 opy += mu_calovec.
cpy();
1141 osumpt += mu_calovec.
sumpt();
1145 ATH_MSG_VERBOSE(
"Misc cluster px, py, sumpt = " << opx <<
", " << opy <<
", " << osumpt );
1146 metSoftClus->
add(opx,opy,osumpt);
1148 <<
", mpy " << metSoftClus->
mpy()
1149 <<
" sumet " << metSoftClus->
sumet());
1152 return StatusCode::SUCCESS;
1.8.18