diff --git a/DataQuality/DataQualityTools/src/DQTGlobalWZFinderTool.cxx b/DataQuality/DataQualityTools/src/DQTGlobalWZFinderTool.cxx
index 2822353b4169dfe290394aeaef6fece92f0c6a56..5c562b3f3244c93cb62933bc134c606d295e8aff 100644
--- a/DataQuality/DataQualityTools/src/DQTGlobalWZFinderTool.cxx
+++ b/DataQuality/DataQualityTools/src/DQTGlobalWZFinderTool.cxx
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 // ********************************************************************
@@ -32,7 +32,7 @@
 #include "xAODJet/Jet.h"
 #include "xAODJet/JetContainer.h"
 
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 #include "IsolationSelection/IIsolationSelectionTool.h"
 
 #include <vector>
diff --git a/MuonSpectrometer/MuonValidation/MuonDQA/MuonPhysValMonitoring/src/MuonPhysValMonitoringTool.h b/MuonSpectrometer/MuonValidation/MuonDQA/MuonPhysValMonitoring/src/MuonPhysValMonitoringTool.h
index 3cb505815b2c325cf483d51099d98e55f515eda6..7af96ecbccbf0f0266dc905ed3080ea9214161c6 100644
--- a/MuonSpectrometer/MuonValidation/MuonDQA/MuonPhysValMonitoring/src/MuonPhysValMonitoringTool.h
+++ b/MuonSpectrometer/MuonValidation/MuonDQA/MuonPhysValMonitoring/src/MuonPhysValMonitoringTool.h
@@ -1,7 +1,5 @@
-///////////////////////// -*- C++ -*- /////////////////////////////
-
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 // MuonPhysValMonitoringTool.cxx 
@@ -30,7 +28,7 @@
 #include "xAODEventInfo/EventInfo.h"
 
 // Tools
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 #include "TrigDecisionTool/TrigDecisionTool.h"
 #include "TrkToolInterfaces/ITrackSelectorTool.h"
 #include "IsolationSelection/IIsolationSelectionTool.h"
diff --git a/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataMonAlg.h b/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataMonAlg.h
index ad308cfb65fb1fd4580bbb2a5254273bc3c7adbf..e8862c82192328dc3f376683d9fa2686afdb2d9d 100755
--- a/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataMonAlg.h
+++ b/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataMonAlg.h
@@ -21,7 +21,7 @@
 #include "GaudiKernel/ToolHandle.h" 
 
 //Helper Includes
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 #include "MdtRawDataMonitoring/MuonChamberIDSelector.h"
 #include "MdtRawDataMonitoring/MDTNoisyTubes.h"
 #include "MdtRawDataMonitoring/MDTChamber.h"
diff --git a/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataValAlg.h b/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataValAlg.h
index 41245dd5ad3362690cc4eea1cf16ef245372f191..3393a01813d1e749dc1436af2b80a2b4036fb31b 100755
--- a/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataValAlg.h
+++ b/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataMonitoring/MdtRawDataValAlg.h
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 ///////////////////////////////////////////////////////////////////////////////////////////
@@ -19,7 +19,7 @@
 #include "GaudiKernel/ToolHandle.h" 
 
 //Helper Includes
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 #include "MdtRawDataMonitoring/MuonChamberIDSelector.h"
 #include "MdtRawDataMonitoring/MDTMonGroupStruct.h"
 #include "MdtRawDataMonitoring/MDTNoisyTubes.h"
diff --git a/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/RpcRawDataMonitoring/CMakeLists.txt b/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/RpcRawDataMonitoring/CMakeLists.txt
index b69196e295efd99acf4861b019d36591ce7bafcc..7027a39ea495e82fa54a3c87a1e23de9b510bbb3 100644
--- a/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/RpcRawDataMonitoring/CMakeLists.txt
+++ b/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/RpcRawDataMonitoring/CMakeLists.txt
@@ -39,7 +39,7 @@ atlas_depends_on_subdirs( PUBLIC
                           Tracking/TrkEvent/TrkTrack
                           Trigger/TrigConfiguration/TrigConfL1Data
                           Trigger/TrigT1/TrigT1Result
-                          PhysicsAnalysis/MuonID/MuonSelectorTools )
+                          PhysicsAnalysis/Interfaces/MuonAnalysisInterfaces )
 
 # External dependencies:
 find_package( Eigen )
diff --git a/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/RpcRawDataMonitoring/RpcRawDataMonitoring/RPCMonitorAlgorithm.h b/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/RpcRawDataMonitoring/RpcRawDataMonitoring/RPCMonitorAlgorithm.h
index 31b15f6fa84b3d979daede68b345c0443fc213d9..4761a5c831f944c5f9a5dcc79d87f996ac95641b 100644
--- a/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/RpcRawDataMonitoring/RpcRawDataMonitoring/RPCMonitorAlgorithm.h
+++ b/MuonSpectrometer/MuonValidation/MuonDQA/MuonRawDataMonitoring/RpcRawDataMonitoring/RpcRawDataMonitoring/RPCMonitorAlgorithm.h
@@ -13,7 +13,7 @@
 #include "MuonIdHelpers/RpcIdHelper.h"
 #include "MuonRDO/RpcPadContainer.h"
 #include "xAODTrigger/MuonRoIContainer.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
 class RPCMonitorAlgorithm : public AthMonitorAlgorithm
 {
diff --git a/MuonSpectrometer/MuonValidation/MuonDQA/MuonTrackMonitoring/MuonTrackMonitoring/MuonGenericTracksMon.h b/MuonSpectrometer/MuonValidation/MuonDQA/MuonTrackMonitoring/MuonTrackMonitoring/MuonGenericTracksMon.h
index 3611eade8ad6ec3cf699f32ab7e348e401bf33f2..17d4edbf92196c5b1a4de9e942b5c368a1d9a59e 100644
--- a/MuonSpectrometer/MuonValidation/MuonDQA/MuonTrackMonitoring/MuonTrackMonitoring/MuonGenericTracksMon.h
+++ b/MuonSpectrometer/MuonValidation/MuonDQA/MuonTrackMonitoring/MuonTrackMonitoring/MuonGenericTracksMon.h
@@ -38,7 +38,7 @@
 #include "FourMomUtils/P4Helpers.h"
 
 #include "MuonRecHelperTools/IMuonEDMHelperSvc.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 #include "MuonResonanceTools/IMuonResonanceSelectionTool.h"
 #include "MuonResonanceTools/IMuonResonancePairingTool.h"
 
diff --git a/PhysicsAnalysis/JetMissingEtID/MissingEtDQA/src/PhysValMET.h b/PhysicsAnalysis/JetMissingEtID/MissingEtDQA/src/PhysValMET.h
index 50165308691f5ffa4772349acffa13b6e2630331..a42ffd4fca707251c519fed52575456d1a0c29eb 100644
--- a/PhysicsAnalysis/JetMissingEtID/MissingEtDQA/src/PhysValMET.h
+++ b/PhysicsAnalysis/JetMissingEtID/MissingEtDQA/src/PhysValMET.h
@@ -1,7 +1,5 @@
-///////////////////////// -*- C++ -*- /////////////////////////////
-
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 // PhysValMET.h 
@@ -20,7 +18,7 @@
 // Local includes
 #include "AthenaMonitoring/ManagedMonitorToolBase.h"
 
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 #include "EgammaAnalysisInterfaces/IAsgElectronLikelihoodTool.h"
 #include "EgammaAnalysisInterfaces/IAsgPhotonIsEMSelector.h"
 #include "JetInterface/IJetUpdateJvt.h"
@@ -36,9 +34,6 @@
 class IMETMaker;
 class IAsgElectronLikelihoodTool;
 class IAsgPhotonIsEMSelector;
-namespace CP {
-  class IMuonSelectionTool;
-}
 namespace TauAnalysisTools {
   class ITauSelectionTool;
 }
diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/src/JetVertexCharge.cxx b/PhysicsAnalysis/JetTagging/JetTagTools/src/JetVertexCharge.cxx
index 0ac88cbf34a03d3daa1b2f620a5423646157dca9..70a3911b36a4984e4e874eea62efe538c85a9423 100644
--- a/PhysicsAnalysis/JetTagging/JetTagTools/src/JetVertexCharge.cxx
+++ b/PhysicsAnalysis/JetTagging/JetTagTools/src/JetVertexCharge.cxx
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 /***************************************************************************
@@ -25,7 +25,7 @@
 
 #include "xAODMuon/MuonContainer.h"
 
-#include "MuonSelectorTools/IMuonSelectionTool.h" 
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h" 
 #include "MuonMomentumCorrections/IMuonCalibrationAndSmearingTool.h"
 #include "TMVA/Reader.h"
 #include "TList.h"
diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/src/SoftMuonTag.cxx b/PhysicsAnalysis/JetTagging/JetTagTools/src/SoftMuonTag.cxx
index 7b0dc066c70fc0f15c14db2be43a54de0a474205..e5843e9e1f0982c203e0e962f24d6a9c05b4d49c 100644
--- a/PhysicsAnalysis/JetTagging/JetTagTools/src/SoftMuonTag.cxx
+++ b/PhysicsAnalysis/JetTagging/JetTagTools/src/SoftMuonTag.cxx
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 /********************************************************************
@@ -32,7 +32,7 @@ PURPOSE:  b-tagging based on soft muon identification
 #include "GaudiKernel/IToolSvc.h"
 #include "ITrackToVertex/ITrackToVertex.h"
 #include "TrkVertexFitterInterfaces/ITrackToVertexIPEstimator.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h" 
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h" 
 
 #include "JetTagInfo/TruthInfo.h"
 #include "JetTagInfo/SoftMuonInfo.h"
diff --git a/PhysicsAnalysis/MuonID/MuonIDAnalysis/MuonMomentumCorrections/src/TestMCASTTool.h b/PhysicsAnalysis/MuonID/MuonIDAnalysis/MuonMomentumCorrections/src/TestMCASTTool.h
index 445609b4d094cd07eddffe352a5a8e1cf9cfb614..6d5c03d16d53fd67321e59ff9cd5936355b56e5b 100644
--- a/PhysicsAnalysis/MuonID/MuonIDAnalysis/MuonMomentumCorrections/src/TestMCASTTool.h
+++ b/PhysicsAnalysis/MuonID/MuonIDAnalysis/MuonMomentumCorrections/src/TestMCASTTool.h
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 #ifndef MCAST_TOOLTESTER_H
@@ -11,7 +11,7 @@
 
 // Local include(s):
 #include "MuonMomentumCorrections/IMuonCalibrationAndSmearingTool.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
 // Root include(s)
 #include "TTree.h"
diff --git a/PhysicsAnalysis/MuonID/MuonPerformanceAnalysis/MuonResonanceTools/MuonResonanceTools/IMuonResonanceSelectionTool.h b/PhysicsAnalysis/MuonID/MuonPerformanceAnalysis/MuonResonanceTools/MuonResonanceTools/IMuonResonanceSelectionTool.h
index 07fd22937c93bb39f126b01ce42eadebd63999b6..2274803c62bd348532afe6d567f722c30600828c 100644
--- a/PhysicsAnalysis/MuonID/MuonPerformanceAnalysis/MuonResonanceTools/MuonResonanceTools/IMuonResonanceSelectionTool.h
+++ b/PhysicsAnalysis/MuonID/MuonPerformanceAnalysis/MuonResonanceTools/MuonResonanceTools/IMuonResonanceSelectionTool.h
@@ -1,10 +1,7 @@
 /*
-  Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
-
-// IMuonResonanceSelectionTool.h
-
 #ifndef IMuonResonanceSelectionTool_H
 #define IMuonResonanceSelectionTool_H
 
@@ -14,7 +11,7 @@
 #include "xAODBase/IParticleContainer.h"
 #include "xAODTracking/Vertex.h"
 #include "EventPrimitives/EventPrimitivesHelpers.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 #include "MuonAnalysisInterfaces/IMuonEfficiencyScaleFactors.h"
 #include "MuonMomentumCorrections/IMuonCalibrationAndSmearingTool.h"
 #include "TrigMuonMatching/ITrigMuonMatching.h"
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/CMakeLists.txt b/PhysicsAnalysis/MuonID/MuonSelectorTools/CMakeLists.txt
index 2d2eeb51182fd383a8e3aba1b974abac8d9f9c79..b4e6f41c7dabebb4333f63cbb9b20d4b8affa7bc 100644
--- a/PhysicsAnalysis/MuonID/MuonSelectorTools/CMakeLists.txt
+++ b/PhysicsAnalysis/MuonID/MuonSelectorTools/CMakeLists.txt
@@ -1,4 +1,3 @@
-# $Id: CMakeLists.txt 777718 2016-10-11 16:52:13Z krasznaa $
 ################################################################################
 # Package: MuonSelectorTools
 ################################################################################
@@ -20,6 +19,7 @@ atlas_depends_on_subdirs(
    Control/AthToolSupport/AsgTools
    Event/xAOD/xAODMuon
    PhysicsAnalysis/AnalysisCommon/PATCore
+   PhysicsAnalysis/Interfaces/MuonAnalysisInterfaces
    PRIVATE
    Event/xAOD/xAODCore
    Event/xAOD/xAODEventInfo
@@ -28,20 +28,20 @@ atlas_depends_on_subdirs(
    ${extra_deps} )
 
 # External dependencies:
-find_package( ROOT COMPONENTS Core Hist RIO )
+find_package( ROOT COMPONENTS Core Hist RIO TMVA )
 
 # Libraries in the package:
 atlas_add_library( MuonSelectorToolsLib
    MuonSelectorTools/*.h Root/*.cxx
    PUBLIC_HEADERS MuonSelectorTools
    INCLUDE_DIRS ${ROOT_INCLUDE_DIRS}
-   LINK_LIBRARIES ${ROOT_LIBRARIES} AsgTools xAODMuon PATCoreLib
+   LINK_LIBRARIES ${ROOT_LIBRARIES} AsgTools xAODEventInfo xAODMuon PATCoreLib MuonAnalysisInterfacesLib
    PRIVATE_LINK_LIBRARIES xAODTracking PathResolver )
 
 if( NOT XAOD_STANDALONE )
    atlas_add_component( MuonSelectorTools
       src/*.h src/*.cxx src/components/*.cxx
-      LINK_LIBRARIES AthenaBaseComps GaudiKernel xAODCore xAODMuon
+      LINK_LIBRARIES AthenaBaseComps GaudiKernel xAODCore xAODEventInfo xAODMuon MuonAnalysisInterfacesLib
       MuonSelectorToolsLib )
 endif()
 
@@ -59,6 +59,15 @@ if( XAOD_STANDALONE )
       xAODTracking xAODCore MuonSelectorToolsLib )
 endif()
 
+
+# Test(s) in the package: 
+if( XAOD_STANDALONE )
+   atlas_add_test( ut_MuonSelectorToolsTester_data
+      SOURCES test/ut_MuonSelectorToolsTester_data.cxx
+      INCLUDE_DIRS ${ROOT_INCLUDE_DIRS}
+      LINK_LIBRARIES ${ROOT_LIBRARIES} AsgTools )
+endif()
+
 # Install files from the package:
 atlas_install_joboptions( share/*.py )
 atlas_install_python_modules( python/*.py )
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/IMuonSelectionTool.h b/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/IMuonSelectionTool.h
deleted file mode 100644
index 744b5b1e3c68142b8799f3c3de3cd5873fda40f4..0000000000000000000000000000000000000000
--- a/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/IMuonSelectionTool.h
+++ /dev/null
@@ -1,87 +0,0 @@
-// Dear emacs, this is -*- c++ -*-
-
-/*
-  Copyright (C) 2002-2017, 2019 CERN for the benefit of the ATLAS collaboration
-*/
-
-// $Id: IMuonSelectionTool.h 299883 2014-03-28 17:34:16Z krasznaa $
-#ifndef MUONSELECTORTOOLS_IMUONSELECTIONTOOL_H
-#define MUONSELECTORTOOLS_IMUONSELECTIONTOOL_H
-
-// Framework include(s):
-#include "AsgTools/IAsgTool.h"
-#include "PATCore/AcceptInfo.h"
-#include "PATCore/AcceptData.h"
-
-// EDM include(s):
-#include "xAODMuon/Muon.h"
-
-namespace CP {
-
-   /// Interface for (a) muon selector tool(s)
-   ///
-   /// This is an example of how to define object selection in
-   /// a tool.
-   ///
-   /// @author Attila Krasznahorkay <Attila.Krasznahorkay@cern.ch>
-   ///
-   /// $Revision: 299883 $
-   /// $Date: 2014-03-28 18:34:16 +0100 (Fri, 28 Mar 2014) $
-   ///
-   class IMuonSelectionTool : public virtual asg::IAsgTool {
-
-      /// Declare the interface that the class provides
-      ASG_TOOL_INTERFACE( CP::IMuonSelectionTool )
-
-   public:
-      /// Decide whether the muon in question is a "good muon" or not
-      virtual const asg::AcceptInfo& getAcceptInfo() const = 0;
-
-      /// Decide whether the muon in question is a "good muon" or not
-      virtual asg::AcceptData accept( const xAOD::Muon& mu ) const = 0;
-
-      /// set the passes ID cuts variable of the muon 
-      virtual void setPassesIDCuts( xAOD::Muon& mu ) const = 0;
-
-      /// set the passes high pT cuts variable of the muon 
-      virtual void setPassesHighPtCuts( xAOD::Muon& mu ) const = 0;
-
-      /// set the passes low pT cuts variable of the muon
-      //virtual void setPassesLowPtEfficiencyCuts( xAOD::Muon& mu ) const = 0;
-
-      /// set the passes quality variable of the muon 
-      virtual void setQuality( xAOD::Muon& mu ) const = 0;
-
-      /// Returns true if the muon passes the standard MCP ID cuts. To set the value on the muon, instead call setPassesIDCuts(xAOD::Muon&) const
-      virtual bool passedIDCuts(const xAOD::Muon&) const =0;
-
-      /// Returns true if the muon passes a standardized loose preselection.
-      virtual bool passedMuonCuts(const xAOD::Muon&) const =0;
-      
-      /// Returns true if the track particle passes the standard MCP ID cuts.
-      virtual bool passedIDCuts(const xAOD::TrackParticle&) const=0;
-
-      /// Returns true if the muon passes the standard MCP high pt cuts. To set the value on the muon, instead call setPassesHighPtCuts(xAOD::Muon&) const
-      virtual bool passedHighPtCuts(const xAOD::Muon&) const =0;
-     
-      /// Returns true if the muon passes the standard MCP low pt cuts. To set the value on the muon, instead call setPassesHighPtCuts(xAOD::Muon&) const
-      virtual bool passedLowPtEfficiencyCuts(const xAOD::Muon&) const =0;
-      virtual bool passedLowPtEfficiencyCuts(const xAOD::Muon&, xAOD::Muon::Quality thisMu_quality) const =0;
-     
-      /// Returns true if a CB muon fails a pt- and eta-dependent cut on the relative CB q/p error   
-      virtual bool passedErrorCutCB(const xAOD::Muon&) const=0;
-
-      /// Returns true if a CB muon fails some loose quaility requirements designed to remove pathological tracks 
-      virtual bool isBadMuon(const xAOD::Muon&) const=0;
-
-      /// Returns the quality of the muon. To set the value on the muon, instead call setQuality(xAOD::Muon&) const
-      virtual xAOD::Muon::Quality getQuality( const xAOD::Muon& mu ) const =0;
-
-     /// Returns true if the muon passes additional calo-tag quality cuts
-     virtual bool passedCaloTagQuality (const xAOD::Muon& mu) const = 0;
-
-   }; // class IMuonSelectionTool
-
-} // namespace CP
-
-#endif // CPTOOLTESTS_IMUONSELECTIONTOOL_H
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/MuonSelectionTool.h b/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/MuonSelectionTool.h
index 454632425f04692d39a1ed139286ca14feab238b..e6dca75c59f1854ad8e66e1f2794b5e3351eba6c 100644
--- a/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/MuonSelectionTool.h
+++ b/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/MuonSelectionTool.h
@@ -1,22 +1,22 @@
-// Dear emacs, this is -*- c++ -*-
-
 /*
-  Copyright (C) 2002-2017, 2019 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
-// $Id: MuonSelectionTool.h 299883 2014-03-28 17:34:16Z krasznaa $
 #ifndef MUONSELECTORTOOLS_MUONSELECTIONTOOL_H
 #define MUONSELECTORTOOLS_MUONSELECTIONTOOL_H
 
-// Framework include(s):
 #include "AsgTools/AsgTool.h"
 #include "PATCore/IAsgSelectionTool.h"
+#include "StoreGate/ReadHandleKey.h"
+#include "xAODEventInfo/EventInfo.h"
+
 #include "TFile.h"
 #include "TH2D.h"
-#include "TSystem.h" // Replace with PathResolver   
+#include "TSystem.h" // Replace with PathResolver
+#include "TMVA/Reader.h"
+
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
-// Local include(s):
-#include "MuonSelectorTools/IMuonSelectionTool.h"
 
 namespace CP {
 
@@ -96,6 +96,7 @@ namespace CP {
       /// Returns true if the muon passes the standard MCP low pt cuts. To set the value on the muon, instead call setPassesLowPtEfficiencyCuts(xAOD::Muon&) const
       virtual bool passedLowPtEfficiencyCuts(const xAOD::Muon&) const override;
       virtual bool passedLowPtEfficiencyCuts(const xAOD::Muon&, xAOD::Muon::Quality thisMu_quality) const override;
+      bool passedLowPtEfficiencyMVACut(const xAOD::Muon&) const;
 
       /// Returns true if a CB muon fails a pt- and eta-dependent cut on the relative CB q/p error
       virtual bool passedErrorCutCB(const xAOD::Muon&) const override;
@@ -121,7 +122,6 @@ namespace CP {
      const std::string m_name;
       /// Maximum pseudorapidity for the selected muons
      double m_maxEta;
-     /// xAOD::Muon::Quality m_quality;
      int  m_quality;
      bool m_isSimulation;
      
@@ -135,6 +135,16 @@ namespace CP {
      bool m_PixCutOff;
      bool m_SiHolesCutOff;
      bool m_TurnOffMomCorr;
+     bool m_useAllAuthors;
+     bool m_use2stationMuonsHighPt;
+     bool m_useMVALowPt;
+     
+     SG::ReadHandleKey<xAOD::EventInfo> m_eventInfo{this, "EventInfoContName", "EventInfo", "event info key"};
+
+     std::string m_MVAreaderFile_EVEN_MuidCB;
+     std::string m_MVAreaderFile_ODD_MuidCB;
+     std::string m_MVAreaderFile_EVEN_MuGirl;
+     std::string m_MVAreaderFile_ODD_MuGirl;
 
      /// Checks for each histogram  
      StatusCode getHist( TFile* file, const char* histName, TH2D*& hist );
@@ -150,6 +160,31 @@ namespace CP {
      // possible override for the calibration version
      std::string m_custom_dir;
 
+     //Need run number (or random run number) to apply period-dependent selections.
+     //If selection depends only on data taking year, this can be specified by passing
+     //argument needOnlyCorrectYear=true, in which case the random run number decoration
+     //from the pile-up reweighting tool is not needed.
+     unsigned int getRunNumber(bool needOnlyCorrectYear = false) const;
+
+     //TMVA readers for low-pT working point
+     TMVA::Reader* m_readerE_MUID;
+     TMVA::Reader* m_readerO_MUID;
+     TMVA::Reader* m_readerE_MUGIRL;
+     TMVA::Reader* m_readerO_MUGIRL;
+
+     //TMVA initialize function
+     void PrepareReader(TMVA::Reader* reader);
+
+     //variables for the TMVA readers
+     Float_t *m_lowPTmva_middleHoles;
+     Float_t *m_lowPTmva_muonSeg1ChamberIdx;
+     Float_t *m_lowPTmva_muonSeg2ChamberIdx;
+     Float_t *m_lowPTmva_momentumBalanceSig;
+     Float_t *m_lowPTmva_scatteringCurvatureSig;
+     Float_t *m_lowPTmva_scatteringNeighbourSig;
+     Float_t *m_lowPTmva_energyLoss;
+     Float_t *m_lowPTmva_muonSegmentDeltaEta;
+
    }; // class MuonSelectionTool
 
 } // namespace CP
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/MuonSelectorToolsDict.h b/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/MuonSelectorToolsDict.h
index c50773d474841fbffae7081213396578564c9dbd..6833109af821cc1670f665c76ad6fe0fd46bc543 100644
--- a/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/MuonSelectorToolsDict.h
+++ b/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/MuonSelectorToolsDict.h
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 #ifndef MUONSELECTORTOOLS_MUONSELECTORTOOLSDICT_H
@@ -10,7 +10,6 @@
 #endif // __GCCXML__
 
 // Includes for the dictionary generation:
-#include "MuonSelectorTools/IMuonSelectionTool.h"
-#include "MuonSelectorTools/MuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
 #endif // MUONSELECTORTOOLS_MUONSELECTORTOOLSDICT_H
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/errorcheck.h b/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/errorcheck.h
index 6c22677eec197e6018577ec7d6f3f6bd68c563e5..f8b71d0db3011b4e65be21ccfa347d2082b9306d 100644
--- a/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/errorcheck.h
+++ b/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/errorcheck.h
@@ -1,10 +1,7 @@
-// Dear emacs, this is -*- c++ -*-
-
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
-// $Id: errorcheck.h 299732 2014-03-27 17:41:34Z krasznaa $
 #ifndef CPTOOLTESTS_ERRORCHECK_H
 #define CPTOOLTESTS_ERRORCHECK_H
 
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/selection.xml b/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/selection.xml
index 1b576e58c210c55c2d54fedfa693e98f42f4afe3..8ca690ba1e95211f11d9954e52a75fb4bec6aa8c 100644
--- a/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/selection.xml
+++ b/PhysicsAnalysis/MuonID/MuonSelectorTools/MuonSelectorTools/selection.xml
@@ -1,5 +1,4 @@
 <lcgdict>
-   <class name="CP::IMuonSelectionTool" />
    <class name="CP::MuonSelectionTool" />
    <exclusion>
      <class name="SG::IConstAuxStore"/>
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/Root/MuonSelectionTool.cxx b/PhysicsAnalysis/MuonID/MuonSelectorTools/Root/MuonSelectionTool.cxx
index 5d5cce1aa2a5cc138794b4401d72ed5d44f82cd5..9864a06726815b0efece98af8253157a327ec582 100644
--- a/PhysicsAnalysis/MuonID/MuonSelectorTools/Root/MuonSelectionTool.cxx
+++ b/PhysicsAnalysis/MuonID/MuonSelectorTools/Root/MuonSelectionTool.cxx
@@ -1,10 +1,7 @@
 /*
-  Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
-// $Id: MuonSelectionTool.cxx 299883 2014-03-28 17:34:16Z krasznaa $
-
-// Local include(s):
 #include "MuonSelectorTools/MuonSelectionTool.h"
 #include "xAODTracking/TrackingPrimitives.h"
 #include "PathResolver/PathResolver.h"
@@ -14,7 +11,11 @@ namespace CP {
   MuonSelectionTool::MuonSelectionTool( const std::string& name )
     : asg::AsgTool( name ),
       m_name(name),
-      m_acceptInfo( "MuonSelection" ){
+      m_acceptInfo( "MuonSelection" ),
+      m_readerE_MUID(nullptr),
+      m_readerO_MUID(nullptr),
+      m_readerE_MUGIRL(nullptr),
+      m_readerO_MUGIRL(nullptr) {
     
     declareProperty( "MaxEta", m_maxEta = 2.7 );
     //xAOD::MuonQuality enum {Tight, Medium, Loose, VeryLoose}
@@ -24,20 +25,44 @@ namespace CP {
     declareProperty( "TurnOffMomCorr", m_TurnOffMomCorr = false );
     declareProperty( "CalibrationRelease", m_calibration_version = "PreRec2016_2016-04-13" );
 
+    //for users of high-pT working point to choose whether to include "safe" 2-station muons
+    declareProperty( "Use2stationMuonsHighPt", m_use2stationMuonsHighPt = true );
+
+    //for users of low-pT working point to choose whether to use MVA
+    declareProperty( "UseMVALowPt", m_useMVALowPt = false );
+
+    //MVA configs for low-pT working point
+    declareProperty( "MVAreaderFile_EVEN_MuidCB", m_MVAreaderFile_EVEN_MuidCB = 
+		     "MuonSelectorTools/190118_PrelimLowPtMVA/LowPtMVA_Weights/BDTG_9JAN2019_MuidCB_EVEN.weights.xml");
+    declareProperty( "MVAreaderFile_ODD_MuidCB", m_MVAreaderFile_ODD_MuidCB =
+		     "MuonSelectorTools/190118_PrelimLowPtMVA/LowPtMVA_Weights/BDTG_9JAN2019_MuidCB_ODD.weights.xml");
+    declareProperty( "MVAreaderFile_EVEN_MuGirl", m_MVAreaderFile_EVEN_MuGirl =
+		     "MuonSelectorTools/190118_PrelimLowPtMVA/LowPtMVA_Weights/BDTG_9JAN2019_MuGirl_EVEN.weights.xml");
+    declareProperty( "MVAreaderFile_ODD_MuGirl", m_MVAreaderFile_ODD_MuGirl =
+		     "MuonSelectorTools/190118_PrelimLowPtMVA/LowPtMVA_Weights/BDTG_9JAN2019_MuGirl_ODD.weights.xml");
+
+
     // DEVELOPEMENT MODE: EXPERTS ONLY!!! 
     declareProperty( "ExpertDevelopMode", m_developMode = false );
     // these are for debugging / testing, *not* for general use!
     declareProperty( "CustomInputFolder", m_custom_dir = "");
-    declareProperty( "TrtCutOff", m_TrtCutOff = false );
+    declareProperty( "TrtCutOff", m_TrtCutOff = true );
     declareProperty( "SctCutOff", m_SctCutOff = false );
     declareProperty( "PixCutOff", m_PixCutOff = false );
     declareProperty( "SiHolesCutOff", m_SiHolesCutOff = false );
+    declareProperty( "UseAllAuthors", m_useAllAuthors = true );
     //
     m_tightWP_lowPt_rhoCuts = 0;
     m_tightWP_lowPt_qOverPCuts = 0;
     m_tightWP_mediumPt_rhoCuts = 0;
     m_tightWP_highPt_rhoCuts = 0;
     //
+    m_lowPTmva_middleHoles = new Float_t; m_lowPTmva_muonSeg1ChamberIdx = new Float_t;
+    m_lowPTmva_muonSeg2ChamberIdx = new Float_t; m_lowPTmva_momentumBalanceSig = new Float_t;
+    m_lowPTmva_scatteringCurvatureSig = new Float_t; m_lowPTmva_scatteringNeighbourSig = new Float_t;
+    m_lowPTmva_energyLoss = new Float_t; m_lowPTmva_muonSegmentDeltaEta = new Float_t;
+
+
     ATH_MSG_DEBUG("Creating MuonSelectionTool named "<<m_name);
   }
     
@@ -55,7 +80,11 @@ namespace CP {
       m_SiHolesCutOff( toCopy.m_SiHolesCutOff ),
       m_TurnOffMomCorr(  toCopy.m_TurnOffMomCorr ),
       m_calibration_version( toCopy.m_calibration_version ),
-      m_custom_dir( toCopy.m_custom_dir )
+      m_custom_dir( toCopy.m_custom_dir ),
+      m_readerE_MUID(nullptr),
+      m_readerO_MUID(nullptr),
+      m_readerE_MUGIRL(nullptr),
+      m_readerO_MUGIRL(nullptr)
   {
     //
     m_tightWP_lowPt_rhoCuts = 0;
@@ -63,6 +92,11 @@ namespace CP {
     m_tightWP_mediumPt_rhoCuts = 0;
     m_tightWP_highPt_rhoCuts = 0;
     //
+    m_lowPTmva_middleHoles = new Float_t; m_lowPTmva_muonSeg1ChamberIdx = new Float_t;
+    m_lowPTmva_muonSeg2ChamberIdx = new Float_t; m_lowPTmva_momentumBalanceSig = new Float_t;
+    m_lowPTmva_scatteringCurvatureSig = new Float_t; m_lowPTmva_scatteringNeighbourSig = new Float_t;
+    m_lowPTmva_energyLoss = new Float_t; m_lowPTmva_muonSegmentDeltaEta = new Float_t;
+
     ATH_MSG_DEBUG("Creating copy of MuonSelectionTool named "<<m_name);
   }
   
@@ -85,6 +119,26 @@ namespace CP {
       delete m_tightWP_highPt_rhoCuts;
       m_tightWP_highPt_rhoCuts = 0;
     }
+    //
+    if( m_readerE_MUID ){
+      delete m_readerE_MUID;
+      m_readerE_MUID = nullptr;
+    }
+    if( m_readerO_MUID ){
+      delete m_readerO_MUID;
+      m_readerO_MUID = nullptr;
+    }
+    if( m_readerE_MUGIRL ){
+      delete m_readerE_MUGIRL;
+      m_readerE_MUGIRL = nullptr;
+    }
+    if( m_readerO_MUGIRL ){
+      delete m_readerO_MUGIRL;
+      m_readerO_MUGIRL = nullptr;
+    }
+    //
+    delete m_lowPTmva_middleHoles; delete m_lowPTmva_muonSeg1ChamberIdx; delete m_lowPTmva_muonSeg2ChamberIdx; delete m_lowPTmva_momentumBalanceSig;
+    delete m_lowPTmva_scatteringCurvatureSig; delete m_lowPTmva_scatteringNeighbourSig; delete m_lowPTmva_energyLoss; delete m_lowPTmva_muonSegmentDeltaEta;
   }
   
   StatusCode MuonSelectionTool::initialize() {
@@ -94,26 +148,35 @@ namespace CP {
     ATH_MSG_INFO( "Maximum eta: " << m_maxEta );
     ATH_MSG_INFO( "Muon quality: " << m_quality );
     if ( m_toroidOff ) ATH_MSG_INFO( "!! CONFIGURED FOR TOROID-OFF COLLISIONS !!" );
-    //if ( m_TrtCutOff ) ATH_MSG_WARNING( "!! SWITCHING TRT REQUIREMENTS OFF !! FOR DEVELOPMENT USE ONLY !!" );
     if ( m_SctCutOff ) ATH_MSG_WARNING( "!! SWITCHING SCT REQUIREMENTS OFF !! FOR DEVELOPMENT USE ONLY !!" );
     if ( m_PixCutOff ) ATH_MSG_WARNING( "!! SWITCHING PIXEL REQUIREMENTS OFF !! FOR DEVELOPMENT USE ONLY !!" );
     if ( m_SiHolesCutOff ) ATH_MSG_WARNING( "!! SWITCHING SILICON HOLES REQUIREMENTS OFF !! FOR DEVELOPMENT USE ONLY !!" );
     if (m_custom_dir!="") ATH_MSG_WARNING("!! SETTING UP WITH USER SPECIFIED INPUT LOCATION \""<<m_custom_dir<<"\"!! FOR DEVELOPMENT USE ONLY !! ");
+    if (!m_useAllAuthors) ATH_MSG_WARNING("Not using allAuthors variable as currently missing in many derivations; LowPtEfficiency working point will always return false, but this is expected at the moment. Have a look here: https://twiki.cern.ch/twiki/bin/view/Atlas/MuonSelectionToolR21#New_LowPtEfficiency_working_poin");
+
+    //Print warning to ensure that users including 2-station muons in the high-pT selection are aware of this
+    if (!m_use2stationMuonsHighPt) ATH_MSG_INFO("You have opted select 3-station muons in the high-pT selection! "<<
+        "Please feed 'HighPt3Layers' to the 'WorkingPoint'  property to retrieve the appropiate scale-factors");
+
 
     // Set up the TAccept object:
     m_acceptInfo.addCut( "Eta",
-                         "Selection of muons according to their pseudorapidity" );
+		     "Selection of muons according to their pseudorapidity" );
     m_acceptInfo.addCut( "IDHits",
-                         "Selection of muons according to whether they passed the MCP ID Hit cuts" );
+                     "Selection of muons according to whether they passed the MCP ID Hit cuts" );
     m_acceptInfo.addCut( "Preselection",
-                         "Selection of muons according to their type/author" );
+                     "Selection of muons according to their type/author" );
     m_acceptInfo.addCut( "Quality",
-                         "Selection of muons according to their tightness" );
+		     "Selection of muons according to their tightness" );
     // Sanity check
     if(m_quality>5 ){
       ATH_MSG_ERROR( "Invalid quality (i.e. selection WP) set: " << m_quality << " - it must be an integer between 0 and 5! (0=Tight, 1=Medium, 2=Loose, 3=Veryloose, 4=HighPt, 5=LowPtEfficiency)" );
       return StatusCode::FAILURE;
     }
+    if(m_quality==5 && !m_useAllAuthors){
+      ATH_MSG_ERROR("Cannot use lowPt working point if allAuthors is not available!");
+      return StatusCode::FAILURE;
+    }
 
     // Load Tight WP cut-map
     ATH_MSG_INFO( "Initialising tight working point histograms..." );
@@ -125,15 +188,10 @@ namespace CP {
     	tightWP_rootFile_fullPath = PathResolverFindCalibFile(Form("MuonSelectorTools/%s/muonSelection_tightWPHisto.root",
     			m_calibration_version.c_str()));
     }
-    // HARD-CODED (TESTING ONLY) !!!
-    /*std::string tightWP_rootFile_fullPath = gSystem->ExpandPathName("$ROOTCOREBIN/data/MuonSelectorTools/");
-    m_tightWP_rootFile = "muonSelection_tightWPHisto_2016_03_15.root";
-    tightWP_rootFile_fullPath.append( m_tightWP_rootFile );*/
-    // ! HARD-CODED END
 
     ATH_MSG_INFO( "Reading muon tight working point histograms from " << tightWP_rootFile_fullPath  );
     // 
-    TFile* file = TFile::Open( tightWP_rootFile_fullPath.c_str() ,"READ");
+    std::unique_ptr<TFile> file ( TFile::Open( tightWP_rootFile_fullPath.c_str() ,"READ"));
 
     if( !file->IsOpen() ){
       ATH_MSG_ERROR( "Cannot read tight working point file from " << tightWP_rootFile_fullPath );
@@ -141,13 +199,38 @@ namespace CP {
     }
 
     // Retrieve all the relevant histograms 
-    ATH_CHECK( getHist( file,"tightWP_lowPt_rhoCuts",m_tightWP_lowPt_rhoCuts) ) ;
-    ATH_CHECK( getHist( file,"tightWP_lowPt_qOverPCuts",m_tightWP_lowPt_qOverPCuts) );
-    ATH_CHECK( getHist( file,"tightWP_mediumPt_rhoCuts",m_tightWP_mediumPt_rhoCuts) ) ;
-    ATH_CHECK( getHist( file,"tightWP_highPt_rhoCuts",m_tightWP_highPt_rhoCuts) ) ;
+    ATH_CHECK( getHist( file.get(),"tightWP_lowPt_rhoCuts",m_tightWP_lowPt_rhoCuts) ) ;
+    ATH_CHECK( getHist( file.get(),"tightWP_lowPt_qOverPCuts",m_tightWP_lowPt_qOverPCuts) );
+    ATH_CHECK( getHist( file.get(),"tightWP_mediumPt_rhoCuts",m_tightWP_mediumPt_rhoCuts) ) ;
+    ATH_CHECK( getHist( file.get(),"tightWP_highPt_rhoCuts",m_tightWP_highPt_rhoCuts) ) ;
     // 
     file->Close();
-    delete file;
+  
+
+    //Set up TMVA readers for MVA-based low-pT working point
+    //E and O refer to even and odd event numbers to avoid applying the MVA on events used for training
+    TString weightPath_EVEN_MuidCB = PathResolverFindCalibFile(m_MVAreaderFile_EVEN_MuidCB);
+    TString weightPath_ODD_MuidCB = PathResolverFindCalibFile(m_MVAreaderFile_ODD_MuidCB);
+    TString weightPath_EVEN_MuGirl = PathResolverFindCalibFile(m_MVAreaderFile_EVEN_MuGirl);
+    TString weightPath_ODD_MuGirl = PathResolverFindCalibFile(m_MVAreaderFile_ODD_MuGirl);
+
+    m_readerE_MUID = new TMVA::Reader();
+    PrepareReader( m_readerE_MUID );
+    m_readerE_MUID->BookMVA("BDTG", weightPath_EVEN_MuidCB);
+
+    m_readerO_MUID = new TMVA::Reader();
+    PrepareReader( m_readerO_MUID );
+    m_readerO_MUID->BookMVA("BDTG", weightPath_ODD_MuidCB);
+
+    m_readerE_MUGIRL = new TMVA::Reader();
+    PrepareReader( m_readerE_MUGIRL );
+    m_readerE_MUGIRL->BookMVA("BDTG", weightPath_EVEN_MuGirl);
+
+    m_readerO_MUGIRL = new TMVA::Reader();
+    PrepareReader( m_readerO_MUGIRL );
+    m_readerO_MUGIRL->BookMVA("BDTG", weightPath_ODD_MuGirl);
+
+    ATH_CHECK( m_eventInfo.initialize() );
     
     // Return gracefully:
     return StatusCode::SUCCESS;
@@ -172,6 +255,21 @@ namespace CP {
     return StatusCode::SUCCESS;
   }
 
+
+  void MuonSelectionTool::PrepareReader( TMVA::Reader* reader ) {
+    // ::
+    reader->AddVariable( "momentumBalanceSignificance", m_lowPTmva_momentumBalanceSig );
+    reader->AddVariable( "scatteringCurvatureSignificance", m_lowPTmva_scatteringCurvatureSig );
+    reader->AddVariable( "scatteringNeighbourSignificance", m_lowPTmva_scatteringNeighbourSig );
+    reader->AddVariable( "EnergyLoss", m_lowPTmva_energyLoss );
+    reader->AddVariable( "middleLargeHoles+middleSmallHoles", m_lowPTmva_middleHoles );
+    reader->AddVariable( "muonSegmentDeltaEta", m_lowPTmva_muonSegmentDeltaEta );
+    reader->AddVariable( "muonSeg1ChamberIdx", m_lowPTmva_muonSeg1ChamberIdx );
+    reader->AddVariable( "muonSeg2ChamberIdx", m_lowPTmva_muonSeg2ChamberIdx );
+    // ::
+    return;
+  }
+
   
   const asg::AcceptInfo& MuonSelectionTool::getAcceptInfo() const {
     
@@ -180,7 +278,7 @@ namespace CP {
   
   asg::AcceptData
   MuonSelectionTool::accept( const xAOD::IParticle* p ) const {
-  
+    
     // Check if this is a muon:
     if( p->type() != xAOD::Type::Muon ) {
       ATH_MSG_ERROR( "accept(...) Function received a non-muon" );
@@ -200,11 +298,29 @@ namespace CP {
   
   asg::AcceptData MuonSelectionTool::accept( const xAOD::Muon& mu ) const {
     
+    // Verbose information
+    ATH_MSG_VERBOSE( "-----------------------------------" );
+    ATH_MSG_VERBOSE( "New muon passed to accept function:" );
+    if(mu.muonType() == xAOD::Muon::Combined)
+      ATH_MSG_VERBOSE( "Muon type: combined" );
+    else if(mu.muonType() == xAOD::Muon::MuonStandAlone)
+      ATH_MSG_VERBOSE( "Muon type: stand-alone" );
+    else if(mu.muonType() == xAOD::Muon::SegmentTagged)
+      ATH_MSG_VERBOSE( "Muon type: segment-tagged" );
+    else if(mu.muonType() == xAOD::Muon::CaloTagged)
+      ATH_MSG_VERBOSE( "Muon type: calorimeter-tagged" );
+    else if(mu.muonType() == xAOD::Muon::SiliconAssociatedForwardMuon)
+      ATH_MSG_VERBOSE( "Muon type: silicon-associated forward" );
+    ATH_MSG_VERBOSE( "Muon pT [GeV]: " << mu.pt()/1000. );
+    ATH_MSG_VERBOSE( "Muon eta: " << mu.eta() );
+    ATH_MSG_VERBOSE( "Muon phi: " << mu.phi() );
+
+
     asg::AcceptData acceptData (&m_acceptInfo);
 
     // Do the eta cut:
-    ATH_MSG_VERBOSE( "Muon eta: " << mu.eta() );
     if( std::abs( mu.eta() ) > m_maxEta ) {
+      ATH_MSG_VERBOSE( "Failed eta cut" );
       return acceptData;
     }
     acceptData.setCutResult( "Eta", true );
@@ -225,8 +341,11 @@ namespace CP {
 
     // Passes quality requirements 
     xAOD::Muon::Quality thisMu_quality = getQuality(mu);
-    bool thisMu_highpt = passedHighPtCuts(mu);
-    bool thisMu_lowptE = passedLowPtEfficiencyCuts(mu,thisMu_quality);
+    bool thisMu_highpt=false;
+    thisMu_highpt = passedHighPtCuts(mu);
+    bool thisMu_lowptE=false;
+    thisMu_lowptE = passedLowPtEfficiencyCuts(mu,thisMu_quality);
+    ATH_MSG_VERBOSE( "Summary of quality information for this muon: ");
     ATH_MSG_VERBOSE( "Muon quality: " << thisMu_quality << " passes HighPt: "<< thisMu_highpt << " passes LowPtEfficiency: "<< thisMu_lowptE );
     if(m_quality<4 && thisMu_quality > m_quality){
       return acceptData;
@@ -248,13 +367,17 @@ namespace CP {
   }
   
   xAOD::Muon::Quality MuonSelectionTool::getQuality( const xAOD::Muon& mu ) const {
-    using namespace xAOD;
     
+    ATH_MSG_VERBOSE( "Evaluating muon quality..." );
+
     // Combined muons
     if (mu.muonType()==xAOD::Muon::Combined){
 
+      ATH_MSG_VERBOSE( "Muon is combined" );
+
       // formally switching off STACO for Rel 20.7 (preliminary) 
       if (mu.author()==xAOD::Muon::STACO){
+	ATH_MSG_VERBOSE( "Muon is STACO - return VeryLoose" );
 	return xAOD::Muon::VeryLoose;
       }
 
@@ -265,6 +388,7 @@ namespace CP {
         return xAOD::Muon::VeryLoose;
       }
       if (combinedTrackOutBoundsPrecisionHits>0){
+	ATH_MSG_VERBOSE( "Muon has out-of-bounds precision hits - return VeryLoose" );
         return xAOD::Muon::VeryLoose;
       }
       
@@ -299,16 +423,17 @@ namespace CP {
 	}
 
         float cbPt = mu.pt();
-        float meP  = 1.0 / ( sin(metrack->theta()) / mePt);
-        float idP  = 1.0 / ( sin(idtrack->theta()) / idPt);
+        float meP  = 1.0 / ( std::sin(metrack->theta()) / mePt);
+        float idP  = 1.0 / ( std::sin(idtrack->theta()) / idPt);
 
-        float rho           = fabs( idPt - mePt ) / cbPt;
-        float qOverPsigma  = sqrt( idtrack->definingParametersCovMatrix()(4,4) + metrack->definingParametersCovMatrix()(4,4) );
-        float qOverPsignif  = fabs( (metrack->charge() / meP) - (idtrack->charge() / idP) ) / qOverPsigma;        
+        float rho           = std::fabs( idPt - mePt ) / cbPt;
+        float qOverPsigma  = std::sqrt( idtrack->definingParametersCovMatrix()(4,4) + metrack->definingParametersCovMatrix()(4,4) );
+        float qOverPsignif  = std::fabs( (metrack->charge() / meP) - (idtrack->charge() / idP) ) / qOverPsigma;        
         float reducedChi2 = mu.primaryTrackParticle()->chiSquared()/mu.primaryTrackParticle()->numberDoF();
 
 	//---- FIX FOR CSC ----
-	if( fabs(mu.eta()) > 2.0 ) {
+	if( std::fabs(mu.eta()) > 2.0 ) {
+	  ATH_MSG_VERBOSE( "Recalculating number of precision layers for combined muon" );
 	  nprecisionLayers = 0;
 	  uint8_t innerSmallHits, innerLargeHits, middleSmallHits, middleLargeHits, outerSmallHits, outerLargeHits;
 	  if ( !mu.summaryValue(innerSmallHits, xAOD::MuonSummaryType::innerSmallHits) ||
@@ -326,42 +451,67 @@ namespace CP {
 	  if( outerSmallHits>2  || outerLargeHits>2  ) nprecisionLayers += 1;
 	}
 
+	ATH_MSG_VERBOSE( "Relevant cut variables:" );
+	ATH_MSG_VERBOSE( "number of precision layers = " << (int)nprecisionLayers );
+	ATH_MSG_VERBOSE( "reduced Chi2 = " << reducedChi2 );
+	ATH_MSG_VERBOSE( "qOverP significance = " << qOverPsignif );
+
+
 	// NEW TIGHT WP  
-        if( nprecisionLayers>1 && reducedChi2<8 && fabs(qOverPsignif)<7 ) {
+        if( nprecisionLayers>1 && reducedChi2<8 && std::fabs(qOverPsignif)<7 ) {
 	  if( passTight(mu,rho,qOverPsignif) ) {
+	    ATH_MSG_VERBOSE( "Muon is tight" );
 	    return xAOD::Muon::Tight;
 	  }
         }
 
+	ATH_MSG_VERBOSE( "Muon did not pass requirements for tight combined muon" );
+
         // MEDIUM WP
-        if( (fabs(qOverPsignif) < 7 || m_toroidOff)
-            && (nprecisionLayers > 1 || ( nprecisionLayers == 1 && nprecisionHoleLayers < 2 && fabs(mu.eta())<0.1) ) ) {
+        if( (std::fabs(qOverPsignif) < 7 || m_toroidOff)
+            && (nprecisionLayers > 1 || ( nprecisionLayers == 1 && nprecisionHoleLayers < 2 && std::fabs(mu.eta())<0.1) ) ) {
+	  ATH_MSG_VERBOSE( "Muon is medium" );
           return xAOD::Muon::Medium;
         }
 
+	ATH_MSG_VERBOSE( "Muon did not pass requirements for medium combined muon" );
+
       } else {
+
+	ATH_MSG_VERBOSE( "Muon is missing the ID and/or ME tracks..." );
+
         // CB muons with missing ID or ME track
-        if( (nprecisionLayers > 1 || ( nprecisionLayers == 1 && nprecisionHoleLayers < 2 && fabs(mu.eta())<0.1) ) )
+        if( (nprecisionLayers > 1 || ( nprecisionLayers == 1 && nprecisionHoleLayers < 2 && std::fabs(mu.eta())<0.1) ) ) {
           // In toroid-off data ME/MS tracks often missing - need special treatment  => flagging as "Medium"
           // In toroid-on data ME/MS tracks missing only for <1% of CB muons, mostly MuGirl (to be fixed) => flagging as "Loose"
-          return (m_toroidOff ? xAOD::Muon::Medium : xAOD::Muon::Loose);
+	  if (m_toroidOff) {
+	    ATH_MSG_VERBOSE( "...this is toroid-off data - returning medium" );
+	    return xAOD::Muon::Medium;
+	  }
+	  else {
+	    ATH_MSG_VERBOSE( "...this is not toroid-off data - returning loose" );
+            return xAOD::Muon::Loose;
+	  }
+	}
       }
 
-      // Improvement for Loose targeting low-pT muons
-      /*if (fabs(mu.eta())>0.1 && fabs(mu.eta())<1.05 && nprecisionLayers==1 && 
-	  (mu.author()==xAOD::Muon::MuidCo || 
-	   (mu.author()==xAOD::Muon::MuGirl && mu.isAuthor(xAOD::Muon::MuTagIMO)) || 
-	   mu.author()==xAOD::Muon::MuTagIMO ) ) 
-	   return xAOD::Muon::Loose;*/
+      // Improvement for Loose targeting low-pT muons (pt<7 GeV)
+      if ( mu.pt()/1000.<7. && std::fabs(mu.eta())<1.3 && nprecisionLayers>0 && (mu.author()==xAOD::Muon::MuGirl && mu.isAuthor(xAOD::Muon::MuTagIMO)) ) {
+	ATH_MSG_VERBOSE( "Muon passed selection for loose working point at low pT" );
+	return xAOD::Muon::Loose;
+      }
 
       // didn't pass the set of requirements for a medium or tight combined muon
+      ATH_MSG_VERBOSE( "Did not pass selections for combined muon - returning VeryLoose" );
       return xAOD::Muon::VeryLoose;
     }
     
     // SA muons
     if ( mu.author()==xAOD::Muon::MuidSA ) {
 
-      if(fabs(mu.eta())>2.5){
+      ATH_MSG_VERBOSE( "Muon is stand-alone" );
+
+      if(std::fabs(mu.eta())>2.5){
 
 	uint8_t nprecisionLayers = 0;
 	uint8_t innerSmallHits(0), innerLargeHits(0), middleSmallHits(0), middleLargeHits(0), outerSmallHits(0), outerLargeHits(0); //, nGoodPrecLayers(0);
@@ -377,35 +527,34 @@ namespace CP {
 	  return xAOD::Muon::VeryLoose;
 	}
 
-	// requiring at least 3 good precision layers (for future improvements)                                                                                                    
-	// if( nGoodPrecLayers <3 )                                                                                                                                                
-	//return xAOD::Muon::VeryLoose;                                                                                                                                            
 	// ---- require 3 MS stations ----                                                                                                                                         
 	if( innerSmallHits>1  || innerLargeHits>1  ) nprecisionLayers += 1;
-	if( middleSmallHits>2 && outerSmallHits>2  ) nprecisionLayers += 2;
-	if( middleLargeHits>2 && outerLargeHits>2  ) nprecisionLayers += 2;
-	
-	/*if (!mu.summaryValue(nprecisionLayers, xAOD::SummaryType::numberOfPrecisionLayers) ) {
-	  ATH_MSG_VERBOSE("getQuality - #precision layers missing in Standalone muon! Aborting.");
-	  return xAOD::Muon::VeryLoose;
-	  }*/
+	if( (middleSmallHits>2 && outerSmallHits>2) || (middleLargeHits>2 && outerLargeHits>2) ) nprecisionLayers += 2;
+
+	ATH_MSG_VERBOSE( "number of precision layers = " << (int)nprecisionLayers );
 
 	// 3 station requirement for medium
-	if( nprecisionLayers>2 && !m_toroidOff )  return xAOD::Muon::Medium;
+	if( nprecisionLayers>2 && !m_toroidOff )  {
+	  ATH_MSG_VERBOSE( "Muon is medium" );
+	  return xAOD::Muon::Medium;
+	}
       }
 
       // didn't pass the set of requirements for a medium SA muon
+      ATH_MSG_VERBOSE( "Muon did not pass selection for medium stand-alone muon - return VeryLoose" );
       return xAOD::Muon::VeryLoose;
     }
 
     // SiliconAssociatedForward (SAF) muons
     if ( mu.muonType()==xAOD::Muon::SiliconAssociatedForwardMuon ){
 
+      ATH_MSG_VERBOSE( "Muon is silicon-associated forward muon" );
+
       const xAOD::TrackParticle* cbtrack = mu.trackParticle( xAOD::Muon::CombinedTrackParticle );
       const xAOD::TrackParticle* metrack = mu.trackParticle( xAOD::Muon::ExtrapolatedMuonSpectrometerTrackParticle );
 
       if( cbtrack && metrack ) {
-	if( fabs(cbtrack->eta()) > 2.5 ) { 
+	if( std::fabs(cbtrack->eta()) > 2.5 ) { 
 	  uint8_t nprecisionLayers = 0;
           uint8_t innerSmallHits(0), innerLargeHits(0), middleSmallHits(0), middleLargeHits(0), outerSmallHits(0), outerLargeHits(0); //, nGoodPrecLayers(0);
           if ( !mu.summaryValue(innerSmallHits, xAOD::MuonSummaryType::innerSmallHits) ||
@@ -419,18 +568,12 @@ namespace CP {
             ATH_MSG_WARNING("getQuality - SAF muon with missing MS hits information!!! Returning VeryLoose...");
             return xAOD::Muon::VeryLoose;
           }
-	  // requiring at least 3 good precision layers (for future improvements)
-	  // if( nGoodPrecLayers <3 ) 
-	  //return xAOD::Muon::VeryLoose;
+
 	  // ---- require 3 MS stations ----
-          if( innerSmallHits>1  || innerLargeHits>1  ) nprecisionLayers += 1;
-          if( middleSmallHits>2 && outerSmallHits>2  ) nprecisionLayers += 2;
-          if( middleLargeHits>2 && outerLargeHits>2  ) nprecisionLayers += 2;
+	  if( innerSmallHits>1  || innerLargeHits>1  ) nprecisionLayers += 1;
+	  if( (middleSmallHits>2 && outerSmallHits>2) || (middleLargeHits>2 && outerLargeHits>2) ) nprecisionLayers += 2;
 
-	  /*if (!mu.summaryValue(nprecisionLayers, xAOD::SummaryType::numberOfPrecisionLayers) ) {
-	    ATH_MSG_VERBOSE("getQuality - #precision layers missing in SiliconAssociatedForward muon! Aborting.");
-	    return xAOD::Muon::VeryLoose;
-	    }*/
+	  ATH_MSG_VERBOSE( "number of precision layers = " << (int)nprecisionLayers );
 	  
 	  if( nprecisionLayers >2 && !m_toroidOff  ) {
 	    if (mu.trackParticle( xAOD::Muon::Primary )  == mu.trackParticle( xAOD::Muon::InnerDetectorTrackParticle ) && !m_developMode ){
@@ -438,29 +581,44 @@ namespace CP {
 			    "This is a bug fixed starting with xAODMuon-00-17-07, which should be present in this release. "<<
 			    "Please report this to the Muon CP group!");
 	    }
+	    ATH_MSG_VERBOSE( "Muon is medium" );
 	    return xAOD::Muon::Medium;
 	  }
 	}
       }
 
       // didn't pass the set of requirements for a medium SAF muon
+      ATH_MSG_VERBOSE( "Muon did not pass selection for medium silicon-associated forward muon - return VeryLoose" );
       return xAOD::Muon::VeryLoose;      
     }
     
     // SegmentTagged muons
     if (mu.muonType()==xAOD::Muon::SegmentTagged){
-      if(fabs(mu.eta())<0.1)
+
+      ATH_MSG_VERBOSE( "Muon is segment-tagged" );
+
+      if(std::fabs(mu.eta())<0.1) {
+	ATH_MSG_VERBOSE( "Muon is loose" );
     	return xAOD::Muon::Loose;
-      else
+      }
+      else {
+	ATH_MSG_VERBOSE( "Do not allow segment-tagged muon at |eta| > 0.1 - return VeryLoose" );
     	return xAOD::Muon::VeryLoose;
+      }
     }
 
     // CaloTagged muons
     if (mu.muonType()==xAOD::Muon::CaloTagged ){
-      if(fabs(mu.eta())<0.1 && passedCaloTagQuality(mu))
+
+      ATH_MSG_VERBOSE( "Muon is calorimeter-tagged" );
+
+      if(std::fabs(mu.eta())<0.1 && passedCaloTagQuality(mu)) {
+	ATH_MSG_VERBOSE( "Muon is loose" );
 	return xAOD::Muon::Loose;
+      }
     }
     
+    ATH_MSG_VERBOSE( "Muon did not pass selection for loose/medium/tight for any muon type - return VeryLoose" );
     return xAOD::Muon::VeryLoose;
   }
   
@@ -471,29 +629,15 @@ namespace CP {
   void MuonSelectionTool::setPassesHighPtCuts( xAOD::Muon& mu ) const {    
     mu.setPassesHighPtCuts(passedHighPtCuts(mu));
   }
-
-  /*void MuonSelectionTool::setPassesLowPtEfficiencyCuts( xAOD::Muon& mu ) const {
-    mu.setPassesLowPtEfficiencyCuts(passedLowPtEfficiencyCuts(mu));
-    }*/
   
   bool MuonSelectionTool::passedIDCuts( const xAOD::Muon& mu ) const {
-    //using namespace xAOD;
     //do not apply the ID hit requirements for SA muons for |eta| > 2.5
-    if ( mu.author()==xAOD::Muon::MuidSA && fabs(mu.eta())>2.5 ) {
+    if ( mu.author()==xAOD::Muon::MuidSA && std::fabs(mu.eta())>2.5 ) {
       return true;
     } else if( mu.muonType()==xAOD::Muon::SiliconAssociatedForwardMuon ) {
       const xAOD::TrackParticle* cbtrack = mu.trackParticle( xAOD::Muon::CombinedTrackParticle );
       if( cbtrack ) {
-	if( fabs(cbtrack->eta()) >2.5 ) {
-	  //uint8_t value=0;
-	  // ---
-	  // Postponing this to post-ICHEP recommendations 
-	  // ---
-	  /*if( cbtrack->summaryValue(value, xAOD::SummaryType::numberOfInnermostPixelLayerHits) ) {
-	    // at least one b-layer hit 
-	    if( value>0 ) return true;
-	    }*/
-	  // TODO -COMMENT THIS RETURN FOR ICHEP!!!
+	if( std::fabs(cbtrack->eta()) >2.5 ) {
 	  return true;
 	}
       }
@@ -514,23 +658,39 @@ namespace CP {
     const xAOD::TrackParticle*       metrack = mu.trackParticle( xAOD::Muon::ExtrapolatedMuonSpectrometerTrackParticle );
     const xAOD::TrackParticle*       cbtrack = mu.trackParticle( xAOD::Muon::CombinedTrackParticle );
     // ::
+    // Some spurious muons are found to have negative ME track fit covariance, and are typically poorly reconstructed
+    if (metrack && metrack->definingParametersCovMatrix()(4,4) < 0.0)
+      return true;
+    // ::
     bool IsBadMuon = false;
     if( idtrack && metrack && cbtrack ) {
       // ::
       double qOverP_ID = idtrack->qOverP();
-      double qOverPerr_ID = sqrt( idtrack->definingParametersCovMatrix()(4,4) );
+      double qOverPerr_ID = std::sqrt( idtrack->definingParametersCovMatrix()(4,4) );
       double qOverP_ME = metrack->qOverP();
-      double qOverPerr_ME = sqrt( metrack->definingParametersCovMatrix()(4,4) );
+      double qOverPerr_ME = std::sqrt( metrack->definingParametersCovMatrix()(4,4) );
       double qOverP_CB = cbtrack->qOverP();
-      double qOverPerr_CB = sqrt( cbtrack->definingParametersCovMatrix()(4,4) );
+      double qOverPerr_CB = std::sqrt( cbtrack->definingParametersCovMatrix()(4,4) );
       // ::
       if( m_quality==4 ) { 
 	// recipe for high-pt selection
 	IsBadMuon = !passedErrorCutCB(mu);
+
+	uint8_t nprecisionLayers(0);
+	if ( !mu.summaryValue( nprecisionLayers, xAOD::SummaryType::numberOfPrecisionLayers ) )
+	  ATH_MSG_WARNING("isBadMuon - unable to retrieve numberOfPrecisionLayers");
+	
+	//temporarily apply same recipe as for other working points in addition to CB error
+	//cut for 2-station muons, pending better treatment of ID/MS misalignments
+	if (m_use2stationMuonsHighPt && nprecisionLayers == 2) {
+	  double IdCbRatio = std::fabs( (qOverPerr_ID/qOverP_ID) / (qOverPerr_CB/qOverP_CB) );
+	  double MeCbRatio = std::fabs( (qOverPerr_ME/qOverP_ME) / (qOverPerr_CB/qOverP_CB) );
+	  IsBadMuon = ( IdCbRatio<0.8 || MeCbRatio<0.8 || IsBadMuon );
+	}
       } else {
 	// recipe for other WP
-	double IdCbRatio = fabs( (qOverPerr_ID/qOverP_ID) / (qOverPerr_CB/qOverP_CB) );
-	double MeCbRatio = fabs( (qOverPerr_ME/qOverP_ME) / (qOverPerr_CB/qOverP_CB) );
+	double IdCbRatio = std::fabs( (qOverPerr_ID/qOverP_ID) / (qOverPerr_CB/qOverP_CB) );
+	double MeCbRatio = std::fabs( (qOverPerr_ME/qOverP_ME) / (qOverPerr_CB/qOverP_CB) );
 	IsBadMuon = ( IdCbRatio<0.8 || MeCbRatio<0.8 );
       }	
     } else {
@@ -546,18 +706,41 @@ namespace CP {
 
   bool MuonSelectionTool::passedLowPtEfficiencyCuts( const xAOD::Muon& mu, xAOD::Muon::Quality thisMu_quality ) const {
 
+    ATH_MSG_VERBOSE( "Checking whether muon passes low-pT selection..." );
+
+    if(!m_useAllAuthors) { //no allAuthors, always fail the WP
+      ATH_MSG_VERBOSE( "Do not have allAuthors variable - fail low-pT" );
+      return false; 
+    }
+
     // requiring combined muons
-    if( mu.muonType() != xAOD::Muon::Combined ) return false;
-    if( mu.author()!=xAOD::Muon::MuGirl && mu.author()!=xAOD::Muon::MuidCo ) return false;
+    if( mu.muonType() != xAOD::Muon::Combined ) {
+      ATH_MSG_VERBOSE( "Muon is not combined - fail low-pT" );
+      return false;
+    }
+    if( mu.author()!=xAOD::Muon::MuGirl && mu.author()!=xAOD::Muon::MuidCo ) {
+      ATH_MSG_VERBOSE( "Muon is neither MuGirl nor MuidCo - fail low-pT" );
+      return false;
+    }
 
     // applying Medium selection above pT = 18 GeV 
     if( mu.pt()/1000.>18. ) {
-      if( thisMu_quality <= xAOD::Muon::Medium ) return true;
-      else return false;
+      ATH_MSG_VERBOSE( "pT > 18 GeV - apply medium selection" );
+      if( thisMu_quality <= xAOD::Muon::Medium ) {
+	ATH_MSG_VERBOSE( "Muon passed low-pT selection" );
+	return true;
+      }
+      else {
+	ATH_MSG_VERBOSE( "Muon failed low-pT selection" );
+	return false;
+      }
     }
 
     // requiring Medium in forward regions
-    if( fabs(mu.eta())>1.55 && thisMu_quality > xAOD::Muon::Medium ) return false;
+    if( !m_useMVALowPt && std::fabs(mu.eta())>1.55 && thisMu_quality > xAOD::Muon::Medium ) {
+      ATH_MSG_VERBOSE( "Not using MVA selection, failing low-pT selection due to medium requirement in forward region" );
+      return false;
+    }
     
     // rejection of muons with out-of-bounds hits 
     uint8_t combinedTrackOutBoundsPrecisionHits;
@@ -566,6 +749,7 @@ namespace CP {
       return false;
     }
     if (combinedTrackOutBoundsPrecisionHits>0) {
+      ATH_MSG_VERBOSE( "Muon has out-of-bounds precision hits - fail low-pT" );
       return false;
     }
 
@@ -575,38 +759,131 @@ namespace CP {
       ATH_MSG_WARNING("passedLowPtEfficiencyCuts - #precision layers missing in combined muon! Failing selection...");
       return false;
     }
-    uint nStationsCut = (fabs(mu.eta())>1.3&&fabs(mu.eta())<1.55) ? 2 : 1;
+    uint nStationsCut = (std::fabs(mu.eta())>1.3&&std::fabs(mu.eta())<1.55) ? 2 : 1;
     if ( nprecisionLayers<nStationsCut ) {
+      ATH_MSG_VERBOSE( "number of precision layers = " << (int)nprecisionLayers << " is lower than cut value " << nStationsCut << " - fail low-pT" );
       return false;
     }
 
-    // reject MuGirl muon if not found also by MuTagIMO 
-    if( mu.author()==xAOD::Muon::MuGirl && !mu.isAuthor(xAOD::Muon::MuTagIMO) ) {
-      return false;
+    // reject MuGirl muon if not found also by MuTagIMO
+    if(m_useAllAuthors){
+      if( mu.author()==xAOD::Muon::MuGirl && !mu.isAuthor(xAOD::Muon::MuTagIMO) ) {
+	ATH_MSG_VERBOSE( "MuGirl muon is not confirmed by MuTagIMO - fail low-pT" );
+	return false;
+      }
+    }
+    else return false;
+
+    if (m_useMVALowPt) {
+      ATH_MSG_VERBOSE( "Applying MVA-based selection" );
+      return passedLowPtEfficiencyMVACut(mu);
     }
 
+    ATH_MSG_VERBOSE( "Applying cut-based selection" );
+
     // apply some loose quality requirements 
     float momentumBalanceSignificance(0.), scatteringCurvatureSignificance(0.), scatteringNeighbourSignificance(0.);
     if( !mu.parameter(momentumBalanceSignificance,xAOD::Muon::momentumBalanceSignificance) 
 	|| !mu.parameter(scatteringCurvatureSignificance,xAOD::Muon::scatteringCurvatureSignificance) 
 	|| !mu.parameter(scatteringNeighbourSignificance,xAOD::Muon::scatteringNeighbourSignificance) ) {
-      ATH_MSG_WARNING("passedLowPtEfficiencyCuts - momentum balance, scatternig curvature or neighbour significances missing in combined muon! Failing selection...");
+      ATH_MSG_WARNING("passedLowPtEfficiencyCuts - momentum balance, scattering curvature or neighbour significances missing in combined muon! Failing selection...");
       return false;
     }
-    if( fabs(momentumBalanceSignificance)>3. || fabs(scatteringCurvatureSignificance)>3. || fabs(scatteringNeighbourSignificance)>3. ) {
+
+    ATH_MSG_VERBOSE( "momentum balance significance: " << momentumBalanceSignificance );
+    ATH_MSG_VERBOSE( "scattering curvature significance: " << scatteringCurvatureSignificance );
+    ATH_MSG_VERBOSE( "scattering neighbour significance: " << scatteringNeighbourSignificance );
+
+    if( std::fabs(momentumBalanceSignificance)>3. || std::fabs(scatteringCurvatureSignificance)>3. || std::fabs(scatteringNeighbourSignificance)>3. ) {
+      ATH_MSG_VERBOSE( "Muon failed cut-based low-pT selection" );
       return false;
     }
 
     // passed low pt selection 
+    ATH_MSG_VERBOSE( "Muon passed cut-based low-pT selection" );
     return true;
   }
 
+
+  bool MuonSelectionTool::passedLowPtEfficiencyMVACut( const xAOD::Muon& mu ) const{
+
+    //set values for all BDT input variables from the muon in question
+    uint8_t middleSmallHoles,middleLargeHoles;
+    if( !mu.parameter(*m_lowPTmva_momentumBalanceSig,xAOD::Muon::momentumBalanceSignificance)
+	|| !mu.parameter(*m_lowPTmva_scatteringCurvatureSig,xAOD::Muon::scatteringCurvatureSignificance)
+	|| !mu.parameter(*m_lowPTmva_scatteringNeighbourSig,xAOD::Muon::scatteringNeighbourSignificance)
+	|| !mu.parameter(*m_lowPTmva_energyLoss,xAOD::Muon::EnergyLoss)
+	|| !mu.parameter(*m_lowPTmva_muonSegmentDeltaEta,xAOD::Muon::segmentDeltaEta)
+	|| !mu.summaryValue(middleSmallHoles, xAOD::MuonSummaryType::middleSmallHoles)
+	|| !mu.summaryValue(middleLargeHoles, xAOD::MuonSummaryType::middleLargeHoles) ) {
+      ATH_MSG_WARNING("passedLowPtEfficiencyMVACut - some low-pT MVA input variables missing in combined muon! Failing selection...");
+      return false;
+    }
+
+    *m_lowPTmva_middleHoles = (Float_t)(middleSmallHoles + middleLargeHoles);
+
+    if (mu.nMuonSegments() > 0)
+      *m_lowPTmva_muonSeg1ChamberIdx = mu.muonSegment(0)->chamberIndex();
+    else
+      *m_lowPTmva_muonSeg1ChamberIdx = -9.0;
+
+    if (mu.nMuonSegments() > 1)
+      *m_lowPTmva_muonSeg2ChamberIdx = mu.muonSegment(1)->chamberIndex();
+    else
+      *m_lowPTmva_muonSeg2ChamberIdx = -9.0;
+
+
+    //get event number from event info
+    SG::ReadHandle<xAOD::EventInfo> eventInfo(m_eventInfo);
+
+    //use different trainings for even/odd numbered events
+    TMVA::Reader *reader_MUID, *reader_MUGIRL;
+    if( eventInfo->eventNumber() % 2 == 1) {
+      reader_MUID = m_readerE_MUID;
+      reader_MUGIRL = m_readerE_MUGIRL;
+    } 
+    else {
+      reader_MUID = m_readerO_MUID;
+      reader_MUGIRL = m_readerO_MUGIRL;
+    }
+
+    // get the BDT discriminant response
+    float BDTdiscriminant;
+    
+    if( mu.author() == 1 )
+      BDTdiscriminant = reader_MUID->EvaluateMVA( "BDTG" );
+    else if( mu.author() == 6 || mu.author() == 4 )
+      BDTdiscriminant = reader_MUGIRL->EvaluateMVA( "BDTG" );
+    else {
+      ATH_MSG_WARNING("Invalid author for low-pT MVA, failing selection...");
+      return false;
+    }
+
+    //cut on dicriminant at -0.6
+    if (BDTdiscriminant > -0.6) {
+      ATH_MSG_VERBOSE( "Passed low-pT MVA cut" );
+      return true;
+    }
+    else {
+      ATH_MSG_VERBOSE( "Failed low-pT MVA cut" );
+      return false;
+    }
+  }
+
+
   bool MuonSelectionTool::passedHighPtCuts( const xAOD::Muon& mu ) const {
-    using namespace xAOD;
+
+    ATH_MSG_VERBOSE( "Checking whether muon passes high-pT selection..." );
     
     // :: Request combined muons
-    if( mu.muonType() != xAOD::Muon::Combined ) return false;
-    if( mu.author()==xAOD::Muon::STACO ) return false;
+    if( mu.muonType() != xAOD::Muon::Combined ) {
+      ATH_MSG_VERBOSE( "Muon is not combined - fail high-pT" );
+      return false;
+    }
+    if( mu.author()==xAOD::Muon::STACO ) {
+      ATH_MSG_VERBOSE( "Muon is STACO - fail high-pT" );
+      return false;
+    }
 
     // :: Reject muons with out-of-bounds hits
     uint8_t combinedTrackOutBoundsPrecisionHits;
@@ -615,12 +892,13 @@ namespace CP {
       return false;
     }
     if (combinedTrackOutBoundsPrecisionHits>0){
+      ATH_MSG_VERBOSE( "Muon has out-of-bounds precision hits - fail high-pT" );
       return false;
     }
 
     // :: Access MS hits information 
     uint8_t nprecisionLayers(0), nGoodPrecLayers(0), innerSmallHits(0), innerLargeHits(0), middleSmallHits(0), middleLargeHits(0), 
-      outerSmallHits(0), outerLargeHits(0), extendedSmallHits(0), extendedLargeHits(0), extendedSmallHoles(0), isSmallGoodSectors(0);
+      outerSmallHits(0), outerLargeHits(0), extendedSmallHits(0), extendedLargeHits(0), extendedSmallHoles(0), isSmallGoodSectors(0), cscUnspoiledEtaHits(0);
     if ( !mu.summaryValue( nprecisionLayers, xAOD::SummaryType::numberOfPrecisionLayers ) || 
 	 !mu.summaryValue( nGoodPrecLayers, xAOD::numberOfGoodPrecisionLayers ) ||
 	 !mu.summaryValue(innerSmallHits, xAOD::MuonSummaryType::innerSmallHits) ||
@@ -632,15 +910,14 @@ namespace CP {
 	 !mu.summaryValue(extendedSmallHits, xAOD::MuonSummaryType::extendedSmallHits) ||
 	 !mu.summaryValue(extendedLargeHits, xAOD::MuonSummaryType::extendedLargeHits) ||
 	 !mu.summaryValue(extendedSmallHoles, xAOD::MuonSummaryType::extendedSmallHoles) ||
-	 !mu.summaryValue(isSmallGoodSectors, xAOD::MuonSummaryType::isSmallGoodSectors)
+	 !mu.summaryValue(isSmallGoodSectors, xAOD::MuonSummaryType::isSmallGoodSectors) ||
+	 !mu.summaryValue(cscUnspoiledEtaHits, xAOD::MuonSummaryType::cscUnspoiledEtaHits)
        ){
       ATH_MSG_WARNING("passedHighPtCuts - MS hits information missing!!! Failing High-pT selection...");
       return false;
     }
 
-    //::: Require 3 (good) station muons
-    if( nprecisionLayers < 3 ) return false;
-    //if( nGoodPrecLayers < 3 ) return false; // postponed (further studies needed)
+    ATH_MSG_VERBOSE( "number of precision layers: " << (int)nprecisionLayers );
 
     //::: Apply MS Chamber Vetoes
     // Given according to their eta-phi locations in the muon spectrometer
@@ -664,57 +941,93 @@ namespace CP {
       float phiMS = MS_track->phi();
       float etaCB = CB_track->eta();
 
+      //::: no unspoiled clusters in CSC
+      if( std::fabs(etaMS)>2.0 || std::fabs(etaCB)>2.0  ) {
+	if( cscUnspoiledEtaHits==0 ) {
+	  ATH_MSG_VERBOSE( "Muon has only spoiled CSC clusters - fail high-pT" );
+	  return false;
+	}
+      }
+      
+      // veto bad CSC giving troubles with scale factors
+      if( mu.eta()<-1.899 && std::fabs(mu.phi())<0.211 ) {
+	ATH_MSG_VERBOSE( "Muon is in eta/phi region vetoed due to disabled chambers in MC - fail high-pT" );
+	return false;
+      }
+
       //::: Barrel/Endcap overlap region
-      if( 1.01 < fabs( etaMS ) && fabs( etaMS ) < 1.1 ) return false;
-      if( 1.01 < fabs( etaCB ) && fabs( etaCB ) < 1.1 ) return false;
+      if( (1.01 < std::fabs( etaMS ) && std::fabs( etaMS ) < 1.1) || (1.01 < std::fabs( etaCB ) && std::fabs( etaCB ) < 1.1) ) {
+	ATH_MSG_VERBOSE( "Muon is in barrel/endcap overlap region - fail high-pT" );
+	return false;
+      }
+      
 
       //::: BIS78
       float BIS78_eta[ 2 ] = { 1.05, 1.3 };
       float BIS78_phi[ 8 ] = { 0.21, 0.57, 1.00, 1.33, 1.78, 2.14, 2.57, 2.93 };
 
-      if ( fabs( etaMS ) >= BIS78_eta[ 0 ] && fabs( etaMS ) <= BIS78_eta[ 1 ] ) {
-        if ( ( fabs( phiMS ) >= BIS78_phi[ 0 ] && fabs( phiMS ) <= BIS78_phi[ 1 ] ) 
-	     || ( fabs( phiMS ) >= BIS78_phi[ 2 ] && fabs( phiMS ) <= BIS78_phi[ 3 ] ) 
-	     || ( fabs( phiMS ) >= BIS78_phi[ 4 ] && fabs( phiMS ) <= BIS78_phi[ 5 ] ) 
-	     || ( fabs( phiMS ) >= BIS78_phi[ 6 ] && fabs( phiMS ) <= BIS78_phi[ 7 ] ) 
-	   ) return false;
+      if ( std::fabs( etaMS ) >= BIS78_eta[ 0 ] && std::fabs( etaMS ) <= BIS78_eta[ 1 ] ) {
+        if ( ( std::fabs( phiMS ) >= BIS78_phi[ 0 ] && std::fabs( phiMS ) <= BIS78_phi[ 1 ] ) 
+	     || ( std::fabs( phiMS ) >= BIS78_phi[ 2 ] && std::fabs( phiMS ) <= BIS78_phi[ 3 ] ) 
+	     || ( std::fabs( phiMS ) >= BIS78_phi[ 4 ] && std::fabs( phiMS ) <= BIS78_phi[ 5 ] ) 
+	     || ( std::fabs( phiMS ) >= BIS78_phi[ 6 ] && std::fabs( phiMS ) <= BIS78_phi[ 7 ] ) 
+	     ) {
+	  ATH_MSG_VERBOSE( "Muon is in BIS7/8 eta/phi region - fail high-pT" );
+	  return false;
+	}
+      }
+
+      //::: BMG - only veto in 2017+2018 data and corresponding MC
+      float BMG_eta[ 6 ] = { 0.35, 0.47, 0.68, 0.80, 0.925, 1.04 };
+      float BMG_phi[ 4 ] = { -1.93, -1.765, -1.38, -1.21 };
+      if ( getRunNumber(true) >= 324320 ) {
+	if ( ( std::fabs( etaMS ) >= BMG_eta[ 0 ] && std::fabs( etaMS ) <= BMG_eta[ 1 ] )
+	     || ( std::fabs( etaMS ) >= BMG_eta[ 2 ] && std::fabs( etaMS ) <= BMG_eta[ 3 ] )
+	     || ( std::fabs( etaMS ) >= BMG_eta[ 4 ] && std::fabs( etaMS ) <= BMG_eta[ 5 ] ) ) {
+	  if ( ( phiMS >= BMG_phi[ 0 ] && phiMS <= BMG_phi[ 1 ] ) 
+	       || ( phiMS >= BMG_phi[ 2 ] && phiMS <= BMG_phi[ 3 ] )
+	       ) {
+	    ATH_MSG_VERBOSE( "Muon is in BMG eta/phi region - fail high-pT" );
+	    return false;
+	  }
+	}
       }
 
       //::: BEE
       float BEE_eta[ 2 ] = { 1.440, 1.692 };
       float BEE_phi[ 8 ] = { 0.301, 0.478, 1.086, 1.263, 1.872, 2.049, 2.657, 2.834 };     
-      if ( fabs( etaMS ) >= BEE_eta[ 0 ] && fabs( etaMS ) <= BEE_eta[ 1 ] ) {
-	if ( ( fabs( phiMS ) >= BEE_phi[ 0 ] && fabs( phiMS ) <= BEE_phi[ 1 ] ) 
-	     || ( fabs( phiMS ) >= BEE_phi[ 2 ] && fabs( phiMS ) <= BEE_phi[ 3 ] ) 
-	     || ( fabs( phiMS ) >= BEE_phi[ 4 ] && fabs( phiMS ) <= BEE_phi[ 5 ] ) 
-	     || ( fabs( phiMS ) >= BEE_phi[ 6 ] && fabs( phiMS ) <= BEE_phi[ 7 ] ) 
+      if ( std::fabs( etaMS ) >= BEE_eta[ 0 ] && std::fabs( etaMS ) <= BEE_eta[ 1 ] ) {
+	if ( ( std::fabs( phiMS ) >= BEE_phi[ 0 ] && std::fabs( phiMS ) <= BEE_phi[ 1 ] ) 
+	     || ( std::fabs( phiMS ) >= BEE_phi[ 2 ] && std::fabs( phiMS ) <= BEE_phi[ 3 ] ) 
+	     || ( std::fabs( phiMS ) >= BEE_phi[ 4 ] && std::fabs( phiMS ) <= BEE_phi[ 5 ] ) 
+	     || ( std::fabs( phiMS ) >= BEE_phi[ 6 ] && std::fabs( phiMS ) <= BEE_phi[ 7 ] ) 
 	     ) {
 	  // Muon falls in the BEE eta-phi region: asking for 4 good precision layers
-	  //if( nGoodPrecLayers < 4 ) return false; // postponed (further studies needed) 
-	  if( nprecisionLayers < 4 ) return false;
+	  if( nprecisionLayers < 4 ) {
+	    ATH_MSG_VERBOSE( "Muon is in BEE eta/phi region and does not have 4 precision layers - fail high-pT" );
+	    return false;
+	  }
 	}  
       }
-      if( fabs(etaCB)>1.4 ) {
+      if( std::fabs(etaCB)>1.4 ) {
 	// Veto residual 3-station muons in BEE region due to MS eta/phi resolution effects
 	//if( nGoodPrecLayers<4 && (extendedSmallHits>0||extendedSmallHoles>0) ) return false; // postponed (further studies needed)
-	if( nprecisionLayers<4 && (extendedSmallHits>0||extendedSmallHoles>0) ) return false;
+	if( nprecisionLayers<4 && (extendedSmallHits>0||extendedSmallHoles>0) ) {
+	  ATH_MSG_VERBOSE( "Muon is in BEE eta/phi region and does not have 4 precision layers - fail high-pT" );
+	  return false;
+	}
       }
     } else {
       ATH_MSG_WARNING( "passedHighPtCuts - MS or CB track missing in muon! Failing High-pT selection..." );
       return false;
     }
 
-    // Remove 3-station muons with small-large sectors overlap
-    if( isSmallGoodSectors ) {
-      if( !(innerSmallHits > 2 && middleSmallHits > 2 && (outerSmallHits>2||extendedSmallHits>2)) ) return false;
-    } else {
-      if( !(innerLargeHits > 2 && middleLargeHits > 2 && (outerLargeHits>2||extendedLargeHits>2)) ) return false;
-    }
+
     
     //::: Apply 1/p significance cut
     const xAOD::TrackParticle* idtrack = mu.trackParticle( xAOD::Muon::InnerDetectorTrackParticle );
     const xAOD::TrackParticle* metrack = mu.trackParticle( xAOD::Muon::ExtrapolatedMuonSpectrometerTrackParticle );
-    if( idtrack && metrack ){
+    if( idtrack && metrack && metrack->definingParametersCovMatrix()(4,4)>0 ){
       float mePt = -999999., idPt = -999999.;
       
       if( !m_TurnOffMomCorr ) { // now using corrected ID/MS momenta
@@ -734,16 +1047,69 @@ namespace CP {
 	idPt = idtrack->pt();
       }
       
-      float meP  = 1.0 / ( sin(metrack->theta()) / mePt);
-      float idP  = 1.0 / ( sin(idtrack->theta()) / idPt);
-      
-      float qOverPsigma  = sqrt( idtrack->definingParametersCovMatrix()(4,4) + metrack->definingParametersCovMatrix()(4,4) );
-      float qOverPsignif  = fabs( (metrack->charge() / meP) - (idtrack->charge() / idP) ) / qOverPsigma;        
-      if( fabs(qOverPsignif) > 7 ) return false;
+      float meP  = 1.0 / ( std::sin(metrack->theta()) / mePt);
+      float idP  = 1.0 / ( std::sin(idtrack->theta()) / idPt);
+      float qOverPsigma  = std::sqrt( idtrack->definingParametersCovMatrix()(4,4) + metrack->definingParametersCovMatrix()(4,4) );
+      float qOverPsignif  = std::fabs( (metrack->charge() / meP) - (idtrack->charge() / idP) ) / qOverPsigma;        
+
+      ATH_MSG_VERBOSE( "qOverP significance: " << qOverPsignif );
 
+      if( std::fabs(qOverPsignif) > 7 ) {
+	ATH_MSG_VERBOSE( "Muon failed qOverP significance cut" );
+	return false;
+      }
+    }
+    else {
+      ATH_MSG_VERBOSE( "Muon missing ID or ME tracks - fail high-pT" );
+      return false;
+    }
+
+
+    //Accept good 2-station muons if the user has opted to include these
+    if (m_use2stationMuonsHighPt && nprecisionLayers == 2) {
+
+      //should not accept EM+EO muons due to ID/MS alignment issues
+      if (std::fabs(mu.eta()) > 1.2 && extendedSmallHits < 3 && extendedLargeHits < 3) {
+	ATH_MSG_VERBOSE( "2-station muon with EM+EO - fail high-pT" );
+	return false;
+      }
+
+      //only select muons missing the inner precision layer
+      //apply strict veto on overlap between small and large sectors
+
+      if (innerLargeHits == 0 && middleLargeHits == 0 && outerLargeHits == 0 && extendedLargeHits == 0 
+	  && middleSmallHits > 2 && (outerSmallHits>2||extendedSmallHits>2) ) {
+	ATH_MSG_VERBOSE( "Accepted 2-station muon in small sector" );
+	return true;
+      }
+
+      if (innerSmallHits == 0 && middleSmallHits == 0 && outerSmallHits == 0 && extendedSmallHits == 0 
+	  && middleLargeHits > 2 && (outerLargeHits>2||extendedLargeHits>2) ) {
+	ATH_MSG_VERBOSE( "Accepted 2-station muon in large sector" );
+	return true;
+      }
     }
-    else return false;
 
+    //::: Require 3 (good) station muons
+    if( nprecisionLayers < 3 ) {
+      ATH_MSG_VERBOSE( "Muon has less than 3 precision layers - fail high-pT" );
+      return false;
+    }
+
+    // Remove 3-station muons with small-large sectors overlap
+    if( isSmallGoodSectors ) {
+      if( !(innerSmallHits > 2 && middleSmallHits > 2 && (outerSmallHits>2||extendedSmallHits>2)) ) {
+	ATH_MSG_VERBOSE( "Muon has small/large sectors overlap - fail high-pT" );
+	return false;
+      }
+    } else {
+      if( !(innerLargeHits > 2 && middleLargeHits > 2 && (outerLargeHits>2||extendedLargeHits>2)) ) {
+	ATH_MSG_VERBOSE( "Muon has small/large sectors overlap - fail high-pT" );
+	return false;
+      }
+    }
+    
+    ATH_MSG_VERBOSE( "Muon passed high-pT selection" );
     return true;
   }
 
@@ -751,35 +1117,68 @@ namespace CP {
     // ::
     if( mu.muonType() != xAOD::Muon::Combined ) return false;
     // :: 
-    float fabs_eta = fabs(mu.eta());
-    float p0(8.0), p1(0.034), p2(0.00011);
+    double start_cut = 2.5;
+    double fabs_eta = std::fabs(mu.eta());
+
+    //parametrization of expected q/p error as function of pT
+    double p0(8.0), p1(0.034), p2(0.00011);
     if( fabs_eta>1.05 && fabs_eta<1.3 ) {
       p1=0.036;
       p2=0.00012;
     } else if( fabs_eta>1.3 && fabs_eta<1.7 ) {
       p1=0.051;
       p2=0.00014;
+      start_cut = 2.0;
     } else if( fabs_eta>1.7 && fabs_eta<2.0 ) {
       p1=0.042;
       p2=0.00010;
     } else if( fabs_eta>2.0) {
       p1=0.034;
       p2=0.00013;
+      start_cut = 2.0;
     }
     // :: 
+    uint8_t nprecisionLayers(0);
+    if ( !mu.summaryValue( nprecisionLayers, xAOD::SummaryType::numberOfPrecisionLayers ) )
+      ATH_MSG_WARNING("passedErrorCutCB - unable to retrieve numberOfPrecisionLayers");
+
+    //independent parametrization for 2-station muons
+    if (m_use2stationMuonsHighPt && nprecisionLayers == 2) {
+      p1 = 0.0739568;
+      p2 = 0.00012443;
+      if( fabs_eta>1.05 && fabs_eta<1.3 ) {
+	p1=0.0674484;
+	p2=0.000119879;
+      } else if( fabs_eta>=1.3 && fabs_eta<1.7 ) {
+	p1=0.041669;
+	p2=0.000178349;
+      } else if( fabs_eta>=1.7 && fabs_eta<2.0 ) {
+	p1=0.0488664;
+	p2=0.000137648;
+      } else if( fabs_eta>=2.0) {
+	p1=0.028077;
+	p2=0.000152707;
+      }
+    }
+    // ::
     bool passErrorCutCB = false;
     const xAOD::TrackParticle* cbtrack = mu.trackParticle( xAOD::Muon::CombinedTrackParticle );
     if( cbtrack ) {
       // ::
       double pt_CB = (cbtrack->pt() / 1000. < 5000.) ? cbtrack->pt() / 1000. : 5000.; // GeV
       double qOverP_CB = cbtrack->qOverP();
-      double qOverPerr_CB = sqrt( cbtrack->definingParametersCovMatrix()(4,4) );
+      double qOverPerr_CB = std::sqrt( cbtrack->definingParametersCovMatrix()(4,4) );
       // sigma represents the average expected error at the muon's pt/eta 
-      double sigma = sqrt( pow(p0/pt_CB,2) + pow(p1,2) + pow(p2*pt_CB,2) );
-      // cuttting at 1.8*sigma for pt <=1 TeV, then linearly tightening untill 1*sigma is reached at pt >= 5TeV. 
-      double coefficient = (pt_CB > 1000.) ? (2.0-0.0002*pt_CB) : 1.8;
+      double sigma = std::sqrt( std::pow(p0/pt_CB,2) + std::pow(p1,2) + std::pow(p2*pt_CB,2) );
+      // cutting at 2.5*sigma or 2.0*sigma for pt <=1 TeV depending on eta region, 
+      // then linearly tightening untill 1*sigma is reached at pt >= 5TeV. 
+      double a = (1.0-start_cut)/4000.0;
+      double b = 1.0 - a*5000.0;
+      double coefficient = (pt_CB > 1000.) ? (a*pt_CB + b) : start_cut;
+      if (m_use2stationMuonsHighPt && nprecisionLayers == 2)
+	coefficient = (pt_CB > 1000.) ? (1.2-0.0001*pt_CB) : 1.1; //for 2-station muons 1.1*sigma -> 0.7*sigma @ 5 TeV
       // ::
-      if( fabs(qOverPerr_CB/qOverP_CB) < coefficient*sigma ) {
+      if( std::fabs(qOverPerr_CB/qOverP_CB) < coefficient*sigma ) {
         passErrorCutCB = true;
       }
     }
@@ -794,16 +1193,16 @@ namespace CP {
       else return true;
     }
     // ::
-    if( mu.muonType() == xAOD::Muon::CaloTagged && fabs(mu.eta())<0.105 ) return true;
+    if( mu.muonType() == xAOD::Muon::CaloTagged && std::fabs(mu.eta())<0.105 && passedCaloTagQuality(mu)) return true;
     // ::
-    if( mu.muonType() == xAOD::Muon::SegmentTagged && fabs(mu.eta())<0.105 ) return true;
+    if( mu.muonType() == xAOD::Muon::SegmentTagged && std::fabs(mu.eta())<0.105 ) return true;
     // ::
-    if( mu.author()==xAOD::Muon::MuidSA && fabs(mu.eta())>2.4 ) return true;
+    if( mu.author()==xAOD::Muon::MuidSA && std::fabs(mu.eta())>2.4 ) return true;
     // ::
     if( mu.muonType() == xAOD::Muon::SiliconAssociatedForwardMuon ) {
       const xAOD::TrackParticle* cbtrack = mu.trackParticle( xAOD::Muon::CombinedTrackParticle );
       if( cbtrack ) {
-        if( fabs(cbtrack->eta()) >2.4 ) return true; 
+        if( std::fabs(cbtrack->eta()) >2.4 ) return true; 
       }
     }
     // ::
@@ -811,7 +1210,6 @@ namespace CP {
   }
   
   bool MuonSelectionTool::passedIDCuts( const xAOD::TrackParticle & track ) const {
-    using namespace xAOD;    
     uint8_t value1=0;
     uint8_t value2=0;
 
@@ -869,23 +1267,6 @@ namespace CP {
     //float CaloLRLikelihood = 1.0;
     int CaloMuonIDTag = -20;
 
-    // Rel. 20.1
-    /*try{
-      bool readLR = mu.parameter(CaloLRLikelihood, xAOD::Muon::CaloLRLikelihood);
-      bool readID = mu.parameter(CaloMuonIDTag, xAOD::Muon::CaloMuonIDTag);
-      
-      if (!readLR || !readID){
-	ATH_MSG_VERBOSE("Unable to read CaloTag Quality information! Don't do anything for now, accept the calotag muons as they are!");
-	//just a hack for now -- this information is not available in the DC14 xAODs
-	//return false;
-      }
-      return ( CaloLRLikelihood > 0.9 || CaloMuonIDTag > 10 );
-    }
-
-    catch (SG::ExcBadAuxVar b) {
-      return false;
-    }*/
-
     // Rel 20.7
     try{                                                                                                                                                                                                                                                                                                                    
       bool readID = mu.parameter(CaloMuonIDTag, xAOD::Muon::CaloMuonIDTag);       
@@ -902,11 +1283,11 @@ namespace CP {
 
   bool MuonSelectionTool::passTight( const xAOD::Muon& mu, float rho, float oneOverPSig ) const
   {
-    float symmetric_eta = fabs( mu.eta() );
+    float symmetric_eta = std::fabs( mu.eta() );
     float pt = mu.pt() / 1000.0; // GeV                                                                                                                                                                                                                                                                                                                                     
     // Impose pT and eta cuts; the bounds of the cut maps  
     if( pt < 4.0 || symmetric_eta>2.5 ) return false;
-    ATH_MSG_DEBUG( "Muon is passing tight WP kinemartic cuts with pT,eta " << mu.pt() << "  ,  " << mu.eta()  );
+    ATH_MSG_VERBOSE( "Muon is passing tight WP kinematic cuts with pT,eta " << mu.pt() << "  ,  " << mu.eta()  );
 
     // ** Low pT specific cuts ** //  
     if( pt > 4.0 && pt <= 20.0 ){
@@ -914,15 +1295,15 @@ namespace CP {
       double rhoCut    = m_tightWP_lowPt_rhoCuts->Interpolate( pt, symmetric_eta );
       double qOverPCut = m_tightWP_lowPt_qOverPCuts->Interpolate( pt , symmetric_eta );
 
-      ATH_MSG_DEBUG( "Applying tight WP cuts to a low pt muon with (pt,eta) ( " << pt << " , " << mu.eta() << " ) " );
-      ATH_MSG_DEBUG( "Rho value " << rho << ", required to be less than " << rhoCut  );
-      ATH_MSG_DEBUG( "Momentum significance value " << oneOverPSig << ", required to be less than " << qOverPCut );
+      ATH_MSG_VERBOSE( "Applying tight WP cuts to a low pt muon with (pt,eta) ( " << pt << " , " << mu.eta() << " ) " );
+      ATH_MSG_VERBOSE( "Rho value " << rho << ", required to be less than " << rhoCut  );
+      ATH_MSG_VERBOSE( "Momentum significance value " << oneOverPSig << ", required to be less than " << qOverPCut );
 
       if( rho > rhoCut ) return false;
-      ATH_MSG_DEBUG( "Muon passed tight WP, low pT rho cut!" );
+      ATH_MSG_VERBOSE( "Muon passed tight WP, low pT rho cut!" );
 
       if( oneOverPSig > qOverPCut ) return false;
-      ATH_MSG_DEBUG( "Muon passed tight WP, low pT momentum significance cut" );
+      ATH_MSG_VERBOSE( "Muon passed tight WP, low pT momentum significance cut" );
 
       // Tight muon!
       return true;
@@ -933,12 +1314,12 @@ namespace CP {
     else if ( pt > 20.0 && pt <= 100.0 ) {
       double rhoCut = m_tightWP_mediumPt_rhoCuts->Interpolate( pt , symmetric_eta );
       // 
-      ATH_MSG_DEBUG( "Applying tight WP cuts to a medium pt muon with (pt,eta) (" << pt << "," << mu.eta() << ")" );
-      ATH_MSG_DEBUG( "Rho value " << rho << " required to be less than " << rhoCut );
+      ATH_MSG_VERBOSE( "Applying tight WP cuts to a medium pt muon with (pt,eta) (" << pt << "," << mu.eta() << ")" );
+      ATH_MSG_VERBOSE( "Rho value " << rho << " required to be less than " << rhoCut );
 
       // Apply cut 
       if( rho > rhoCut ) return false;
-      ATH_MSG_DEBUG( "Muon passed tight WP, medium pT rho cut!" );
+      ATH_MSG_VERBOSE( "Muon passed tight WP, medium pT rho cut!" );
 
       // Tight muon! 
       return true;
@@ -947,21 +1328,21 @@ namespace CP {
     // ** High pT specific cuts  
     else if ( pt > 100.0 && pt <= 500.0 ){
       // 
-      ATH_MSG_DEBUG( "Applying tight WP cuts to a high pt muon with (pt,eta) (" << pt << "," << mu.eta() << ")" );
+      ATH_MSG_VERBOSE( "Applying tight WP cuts to a high pt muon with (pt,eta) (" << pt << "," << mu.eta() << ")" );
       // No interpolation, since bins with -1 mean we should cut really loose       
       double rhoCut = m_tightWP_highPt_rhoCuts->GetBinContent( m_tightWP_highPt_rhoCuts->FindBin(pt, symmetric_eta) );
-      ATH_MSG_DEBUG( "Rho value " << rho << ", required to be less than " << rhoCut << " unless -1, in which no cut is applied" );
+      ATH_MSG_VERBOSE( "Rho value " << rho << ", required to be less than " << rhoCut << " unless -1, in which no cut is applied" );
       //
       if( rhoCut <  0.0 ) return true;
       if( rho > rhoCut  ) return false;
-      ATH_MSG_DEBUG( "Muon passd tight WP, high pT rho cut!" );
+      ATH_MSG_VERBOSE( "Muon passd tight WP, high pT rho cut!" );
 
 
       return true;
     }
     // For muons with pT > 500 GeV, no extra cuts                    
     else{
-      ATH_MSG_DEBUG( "Not applying any tight WP cuts to a very high pt muon with (pt,eta) (" << pt << "," << mu.eta() << ")" );
+      ATH_MSG_VERBOSE( "Not applying any tight WP cuts to a very high pt muon with (pt,eta) (" << pt << "," << mu.eta() << ")" );
       return true;
     }
     
@@ -969,6 +1350,47 @@ namespace CP {
     return false;
   }
 
+
+  //need run number (or random run number) to apply period-dependent selections
+  unsigned int MuonSelectionTool::getRunNumber(bool needOnlyCorrectYear) const {
+
+    static const SG::AuxElement::ConstAccessor<unsigned int> acc_rnd("RandomRunNumber");
+
+    SG::ReadHandle<xAOD::EventInfo> eventInfo(m_eventInfo);
+
+    if (!eventInfo->eventType(xAOD::EventInfo::IS_SIMULATION)) {
+      ATH_MSG_DEBUG("The current event is a data event. Return runNumber.");
+      return eventInfo->runNumber();
+    }
+
+    if (!acc_rnd.isAvailable(*eventInfo)) {
+
+      if (needOnlyCorrectYear) {
+	if (eventInfo->runNumber() < 300000) {
+	  ATH_MSG_DEBUG("Random run number not available and this is mc16a, returning dummy 2016 run number.");
+	  return 311071;
+	}
+	else if (eventInfo->runNumber() < 310000) {
+	  ATH_MSG_DEBUG("Random run number not available and this is mc16c/d, returning dummy 2017 run number.");
+	  return 340072;
+	}
+	else {
+	  ATH_MSG_DEBUG("Random run number not available and this is mc16e, returning dummy 2018 run number.");
+	  return 351359;
+	}
+      }
+      else {
+	ATH_MSG_WARNING("Failed to find the RandomRunNumber decoration. Please call the apply() method from the PileupReweightingTool before applying the selection tool. Returning dummy 2017 run number.");
+	return 340072;
+      }
+    } else if (acc_rnd(*eventInfo) == 0) {
+
+      ATH_MSG_DEBUG("Pile up tool has given runNumber 0. Returning dummy 2017 run number.");
+      return 340072;
+    }
+
+    return acc_rnd(*eventInfo);
+  }
   
 } // namespace CP
 
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/python/MuonSelectorCutDefs.py b/PhysicsAnalysis/MuonID/MuonSelectorTools/python/MuonSelectorCutDefs.py
index 15e47261c6c2353d0256d1c572fbf789857b9247..8f140b97342860bde2ce71e22e77886b1ea0dd28 100644
--- a/PhysicsAnalysis/MuonID/MuonSelectorTools/python/MuonSelectorCutDefs.py
+++ b/PhysicsAnalysis/MuonID/MuonSelectorTools/python/MuonSelectorCutDefs.py
@@ -1,4 +1,4 @@
-# Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+# Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 
 ##=============================================================================
 ## Name:        MuonSelectorCutDefs
@@ -10,11 +10,6 @@
 ##
 ##=============================================================================
 
-# import the needed Reflex and ROOT stuff
-import PyCintex
-PyCintex.Cintex.Enable()
-import ROOT
-
 # Import a needed helper
 from PATCore.HelperUtils import *
 
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/src/MuonQualityUpdaterAlg.h b/PhysicsAnalysis/MuonID/MuonSelectorTools/src/MuonQualityUpdaterAlg.h
index a89bed2db858a8bda4e6cc6d6323323cf47c0005..35ec4e58ea694c17cc29ed076e7098622e675bcb 100644
--- a/PhysicsAnalysis/MuonID/MuonSelectorTools/src/MuonQualityUpdaterAlg.h
+++ b/PhysicsAnalysis/MuonID/MuonSelectorTools/src/MuonQualityUpdaterAlg.h
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 #ifndef MUONSELECTORTOOLS_MUONQUALITYUPDATERALG
@@ -7,7 +7,7 @@
 
 #include "AthenaBaseComps/AthAlgorithm.h"
 #include "GaudiKernel/ToolHandle.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
 namespace CP {
 
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/src/MuonSelectionAlg.h b/PhysicsAnalysis/MuonID/MuonSelectorTools/src/MuonSelectionAlg.h
index 010bab0e300c17dfc80947b5bc07d88b8ae4baae..518db20c409a4bdd1c8d2a5883a1219ceb7cf46d 100644
--- a/PhysicsAnalysis/MuonID/MuonSelectorTools/src/MuonSelectionAlg.h
+++ b/PhysicsAnalysis/MuonID/MuonSelectorTools/src/MuonSelectionAlg.h
@@ -1,19 +1,16 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 //Simple alg to wrap the selectiontool
 //outputs a view container of the selected muons, if m_outputMuons different to m_inputMuons
 
-//author: will buttinger <will@cern.ch>
-
-
 #ifndef MUONSELECTORTOOLS_MUONSELECTIONALG
 #define MUONSELECTORTOOLS_MUONSELECTIONALG
 
 #include "AthenaBaseComps/AthAlgorithm.h"
 #include "GaudiKernel/ToolHandle.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
 namespace CP {
 
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/src/components/MuonSelectorTools_entries.cxx b/PhysicsAnalysis/MuonID/MuonSelectorTools/src/components/MuonSelectorTools_entries.cxx
index edb246225158bb249e69a2880610f0b49882b7f7..4e6fe8f3d1bc3c0c832286a87850afa4665b01c7 100644
--- a/PhysicsAnalysis/MuonID/MuonSelectorTools/src/components/MuonSelectorTools_entries.cxx
+++ b/PhysicsAnalysis/MuonID/MuonSelectorTools/src/components/MuonSelectorTools_entries.cxx
@@ -1,13 +1,10 @@
+/*
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
+*/
 #include "MuonSelectorTools/MuonSelectionTool.h"
 #include "../MuonSelectionAlg.h"
 #include "../MuonQualityUpdaterAlg.h"
 
-
-using namespace CP;
-
-DECLARE_COMPONENT( CP::MuonSelectionTool )//DECLARE_COMPONENT( MuonSelectionTool ) is there a difference!?
+DECLARE_COMPONENT( CP::MuonSelectionTool )
 DECLARE_COMPONENT( CP::MuonSelectionAlg )
 DECLARE_COMPONENT( CP::MuonQualityUpdaterAlg )
-
-    
-
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/test/ut_MuonSelectorToolsTester_data.cxx b/PhysicsAnalysis/MuonID/MuonSelectorTools/test/ut_MuonSelectorToolsTester_data.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..512e4958664ca192730fa26a1378b8c1f2a8f131
--- /dev/null
+++ b/PhysicsAnalysis/MuonID/MuonSelectorTools/test/ut_MuonSelectorToolsTester_data.cxx
@@ -0,0 +1,31 @@
+/*
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include <iostream>
+#include <cstdlib>
+
+#include <AsgTools/MessageCheck.h>
+
+int main()
+{
+  using namespace asg::msgUserCode;
+
+  ATH_MSG_INFO ("Unit test for MuonSelectorTools on data");
+
+  // Full `env` makes log file diffs useless.  Check if the input file changed, though - points us quickly to an issue.
+  ATH_MSG_INFO ("Test files");
+  system("env | grep ASG_TEST_FILE_ | sort");
+
+  std::string cmd("MuonSelectorToolsTester $ASG_TEST_FILE_DATA");
+  ATH_MSG_INFO ("Will now run this command: " << cmd);
+  int ret = system(cmd.c_str());
+
+  if (ret != 0) {
+    ATH_MSG_ERROR ("Test failed (return code was " << ret << ")");
+    return 1;
+  }
+
+  ATH_MSG_INFO ("Finished (return code was " << ret << ")");
+  return 0;
+}
diff --git a/PhysicsAnalysis/MuonID/MuonSelectorTools/util/MuonSelectorToolsTester.cxx b/PhysicsAnalysis/MuonID/MuonSelectorTools/util/MuonSelectorToolsTester.cxx
index 5423445d184f741134e04aed70bfb54105af9019..3e025ec8ccc5da7226f7047ac114e5aff208e687 100644
--- a/PhysicsAnalysis/MuonID/MuonSelectorTools/util/MuonSelectorToolsTester.cxx
+++ b/PhysicsAnalysis/MuonID/MuonSelectorTools/util/MuonSelectorToolsTester.cxx
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 /// a simple testing macro for the MuonSelectorTools_xAOD package
@@ -10,6 +10,7 @@
 #include <cstdlib>
 #include <string>
 #include <map>
+#include <iomanip>
 
 // ROOT include(s):
 #include <TFile.h>
@@ -42,13 +43,13 @@
 int main( int argc, char* argv[] ) {
 
 
-	// The application's name:
+   // The application's name:
    const char* APP_NAME = argv[ 0 ];
 
    // Check if we received a file name:
    if( argc < 2 ) {
       Error( APP_NAME, "No file name received!" );
-      Error( APP_NAME, "  Usage: %s [xAOD file name]", APP_NAME );
+      Error( APP_NAME, "  Usage: %s [xAOD file name] [Nevts to process]", APP_NAME );
       return 1;
    }
 
@@ -58,7 +59,7 @@ int main( int argc, char* argv[] ) {
    // Open the input file:
    const TString fileName = argv[ 1 ];
    Info( APP_NAME, "Opening file: %s", fileName.Data() );
-   std::auto_ptr< TFile > ifile( TFile::Open( fileName, "READ" ) );
+   std::unique_ptr< TFile > ifile( TFile::Open( fileName, "READ" ) );
    CHECK( ifile.get() );
 
    // Create a TEvent object:
@@ -76,32 +77,134 @@ int main( int argc, char* argv[] ) {
       }
    }
 
-   CP::MuonSelectionTool m_muonSelection( "MuonSelection" );
+   // Create a TStore object
+   xAOD::TStore store;
 
-   m_muonSelection.msg().setLevel( MSG::INFO );
-   m_muonSelection.setProperty( "MaxEta", 2.7 );
-   m_muonSelection.setProperty( "MuQuality", 3);
-   //m_muonSelection.setProperty( "Author", 12 );
-   //m_muonSelection.initialize();
-   m_muonSelection.setProperty( "TurnOffMomCorr", true );
-   CHECK (m_muonSelection.initialize());
+   // Create a set of selector tools configured for each of the available working points
 
-   int allMuons = 0;
-   int badMuons = 0;
-   int verylooseMuons = 0;
-   int looseMuons = 0;
-   int mediumMuons = 0;
-   int tightMuons = 0;
-   int allgoodMuons = 0;
-   int nCaloTagged = 0;
+   std::vector<std::unique_ptr<CP::MuonSelectionTool> > selectorTools;
+   selectorTools.clear();
+
+   const int Nwp = 6; // number of working points (tool instances)
+
+   std::vector<std::string> WPnames = {"Tight", "Medium", "Loose", "VeryLoose", "HighPt", "LowPt"};
+
+   for (int wp = 0; wp < Nwp; wp++) {
+
+     Info( APP_NAME, "Creating selector tool for working point: %s ...", WPnames[wp].c_str() );
+
+     CP::MuonSelectionTool* muonSelection = new CP::MuonSelectionTool( "MuonSelection_"+WPnames[wp] );
+
+     muonSelection->msg().setLevel( MSG::INFO );
+     CHECK (muonSelection->setProperty( "MaxEta", 2.7 ));
+     CHECK (muonSelection->setProperty( "MuQuality", wp));
+     CHECK (muonSelection->setProperty( "TurnOffMomCorr", true ));
+     CHECK (muonSelection->initialize());
+
+     selectorTools.emplace_back(muonSelection);
+   }
+
+   //done setting up selector tools
+
+   //Create "padding" strings to facilitate easy table view of results
+   std::vector<std::string> padding;
+   padding.clear();
 
-   int nSAMuons = 0;
+   unsigned int maxNameLength = 0;
+   for (int wp = 0; wp < Nwp; wp++)
+     if (WPnames[wp].size() > maxNameLength)
+       maxNameLength = WPnames[wp].size();
 
+   for (int wp = 0; wp < Nwp; wp++) {
+
+     std::string pad = "";
+     for (unsigned int i = 0; i < maxNameLength - WPnames[wp].size(); i++)
+       pad += " ";
+
+     padding.push_back(pad);
+   }
+
+   //string with names of working points for selection results display
+   std::string namesString = "                           ";
+   for (int wp = 0; wp < Nwp; wp++)
+     namesString += WPnames[wp] + "   ";
+
+   //Muon counters
+   int allMuons = 0;
    int nPositive = 0;
    int nNegative = 0;
 
-   int passesCutsTool = 0;
-   int passesCutsxAOD = 0;
+   //Summary information - how many muons passed each working point (with and without vetoing bad muons)
+   int selectedMuons[Nwp];
+   for (int wp = 0; wp < Nwp; wp++)
+     selectedMuons[wp] = 0;
+
+   int selectedMuonsNotBad[Nwp];
+   for (int wp = 0; wp < Nwp; wp++)
+     selectedMuonsNotBad[wp] = 0;
+
+   //Obtain summary information also split by muon type
+   const int Ntype = 5;
+
+   std::string typeNames[Ntype];
+   for (int type = 0; type < Ntype; type++) {
+
+     if(type == xAOD::Muon::Combined)
+       typeNames[type] = "combined";
+     else if(type == xAOD::Muon::MuonStandAlone)
+       typeNames[type] = "stand-alone";
+     else if(type == xAOD::Muon::SegmentTagged)
+       typeNames[type] = "segment-tagged";
+     else if(type == xAOD::Muon::CaloTagged)
+       typeNames[type] = "calo-tagged";
+     else if(type == xAOD::Muon::SiliconAssociatedForwardMuon)
+       typeNames[type] = "forward";
+     else
+       typeNames[type] = "unknown";
+   }
+
+
+   //Muon counters for each type
+   int allMuonsType[Ntype];
+   for (int type = 0; type < Ntype; type++)
+     allMuonsType[type] = 0;
+
+   //Muon counters for muons of each type passing each working point
+   int selectedMuonsType[Ntype][Nwp];
+   for (int type = 0; type < Ntype; type++)
+     for (int wp = 0; wp < Nwp; wp++)
+       selectedMuonsType[type][wp] = 0;
+
+   int selectedMuonsTypeNotBad[Ntype][Nwp];
+   for (int type = 0; type < Ntype; type++)
+     for (int wp = 0; wp < Nwp; wp++)
+       selectedMuonsTypeNotBad[type][wp] = 0;
+
+
+
+   //Obtain summary information also split by muon |eta|
+   const int Neta = 4;
+   double etaCuts[Neta-1] = {1.0, 2.0, 2.5};
+
+   std::string etaRegions = "|eta| < 1.0      1.0 < |eta| < 2.0   2.0 < |eta| < 2.5     |eta| > 2.5";
+
+   //Muon counters for each eta region
+   int allMuonsEta[Neta];
+   for (int eta = 0; eta < Neta; eta++)
+     allMuonsEta[eta] = 0;
+
+   //Muon counters for muons in each eta region passing each working point
+   int selectedMuonsEta[Neta][Nwp];
+   for (int eta = 0; eta < Neta; eta++)
+     for (int wp = 0; wp < Nwp; wp++)
+       selectedMuonsEta[eta][wp] = 0;
+
+   int selectedMuonsEtaNotBad[Neta][Nwp];
+   for (int eta = 0; eta < Neta; eta++)
+     for (int wp = 0; wp < Nwp; wp++)
+       selectedMuonsEtaNotBad[eta][wp] = 0;
+
+
 
    // Loop over the events:
    for( Long64_t entry = 0; entry < entries; ++entry ) {
@@ -109,7 +212,15 @@ int main( int argc, char* argv[] ) {
       // Tell the object which entry to look at:
       event.getEntry( entry );
 
-      int goodMuons = 0;
+      //Counters for selected muons within each event
+      int selectedMuonsEvent[Nwp];
+      for (int wp = 0; wp < Nwp; wp++)
+	selectedMuonsEvent[wp] = 0;
+
+      int selectedMuonsEventNotBad[Nwp];
+      for (int wp = 0; wp < Nwp; wp++)
+	selectedMuonsEventNotBad[wp] = 0;
+
       // Print some event information for fun:
       const xAOD::EventInfo* ei = 0;
       CHECK( event.retrieve( ei, "EventInfo" ) );
@@ -126,228 +237,141 @@ int main( int argc, char* argv[] ) {
       Info( APP_NAME, "Number of muons: %i",
             static_cast< int >( muons->size() ) );
 
+
       xAOD::Muon::Quality my_quality;
-      bool passesIDRequirements = false;
-      
+      bool passesIDRequirements;
+      bool passesPreselectionCuts;
+
+      int muCounter = 0; //muon counter for each event
       
       // Print their properties, using the tools:
       xAOD::MuonContainer::const_iterator mu_itr = muons->begin();
       xAOD::MuonContainer::const_iterator mu_end = muons->end();
       for( ; mu_itr != mu_end; ++mu_itr ) {
 
-	allMuons++;
+	int etaIndex = Neta-1;
+	for (int eta = 0; eta < Neta-1; eta++)
+	  if (std::abs((*mu_itr)->eta()) < etaCuts[eta]) {
+	    etaIndex = eta;
+	    break;
+	  }
 
-	uint8_t PixHits = 0, PixDead = 0, SCTHits = 0, SCTDead = 0, PixHoles = 0, SCTHoles = 0, TRTHits = 0, TRTOut = 0;
-	uint8_t nprecisionLayers = 0, nprecisionHoleLayers = 0;
-	//float momSigBal = 0;
-	//float CaloLRLikelihood = 0;
-	//int CaloMuonIDTag = 0;
+	allMuons++;
+	allMuonsType[(*mu_itr)->muonType()]++;
+	allMuonsEta[etaIndex]++;
+	muCounter++;
 
-	float abseta = std::abs((*mu_itr)->eta());
-	//bool passesTool = false;
-	//bool passesxAOD = false;
+	Info( APP_NAME, "===== Muon number: %i",
+	      static_cast< int >( muCounter ) );
 
 
-	Info( APP_NAME, "Muon pT:             %g ", std::abs((*mu_itr)->pt())/1000.);
-	Info( APP_NAME, "Muon eta:            %g ", std::abs((*mu_itr)->eta()));
-	Info( APP_NAME, "Muon muonType:       %d ", std::abs((*mu_itr)->muonType()));
-	// Info( APP_NAME, "Muon phi:            %g ", std::abs((*mu_itr)->phi()));
-	// Info( APP_NAME, "Muon charge:         %g ", (*mu_itr)->charge());
 	if((*mu_itr)->charge() > 0)
 	  nPositive++;
 	else
 	  nNegative++;
-	// //Info( APP_NAME, "Muon allAuthors:     %d ", std::abs((*mu_itr)->allAuthors()));
-	// Info( APP_NAME, "Muon author:         %d ", std::abs((*mu_itr)->author()));
-	//Info( APP_NAME, "Muon muonType:       %d ", std::abs((*mu_itr)->muonType()));
-	// Info( APP_NAME, "Muon quality:        %d ", std::abs((*mu_itr)->quality()));
-	// Info( APP_NAME, "----------------------------------------------------------");
-	//	if (!(*mu_itr)->parameter(momSigBal, xAOD::Muon_v1::momentumBalanceSignificance))
-	//	  std::cout << "No momSigBal " << std::endl;
-	// Info( APP_NAME, "Muon: momSigBal %f ", momSigBal);
-	//	  if (!(*mu_itr)->parameter(CaloLRLikelihood, xAOD::Muon_v1::CaloLRLikelihood))
-	//	  std::cout << "No caloLRLikelihood " << std::endl;
-	//Info( APP_NAME, "Muon: caloLRLH  %f ", CaloLRLikelihood);
-	//	    if (!(*mu_itr)->parameter(CaloMuonIDTag, xAOD::Muon_v1::CaloMuonIDTag))
-	//	  std::cout << "No caloMuonIDTag " << std::endl;
-	//	Info( APP_NAME, "Muon: caloIDTag %u ", CaloMuonIDTag);
-	
-	
-	// if( (*mu_itr)->muonType() == xAOD::Muon_v1::CaloTagged)	{
-	//   std::cout << "Found a calo-tagged muon! Yay! With abseta " << abseta << std::endl;
-	//   nCaloTagged++;
-	// }
-        passesIDRequirements = m_muonSelection.passedIDCuts(**mu_itr);
-	my_quality = m_muonSelection.getQuality(**mu_itr);	
-	if( (*mu_itr)->muonType() == xAOD::Muon_v1::MuonStandAlone)
-	  std::cout << "SA MUON " << passesIDRequirements << " " << my_quality << std::endl;
-	
-	//std::cout << "quality " << my_quality << " IDHits " << passesIDRequirements << std::endl;
-	//std::cout << "Comparison of the quality: xAOD vs SelectorTool " << (*mu_itr)->quality() << " " <<  my_quality << std::endl;
-	
-	
-	//std::cout << passesxAOD << " "  << passesTool << std::endl;
-	if(abseta < 2.5 && passesIDRequirements && (*mu_itr)->quality() <= xAOD::Muon::Medium){
-	  //std:: cout << "Passes the selection (eta, IDCuts and quality) from the xAOD " << std::endl;
-	  passesCutsxAOD++; //passesxAOD = true;
-	}
-	
-	if(my_quality <= xAOD::Muon::VeryLoose)
-	  verylooseMuons++;
-	if(my_quality <= xAOD::Muon::Loose){
-	  looseMuons++;}
-	if(my_quality <= xAOD::Muon::Medium){
-	  mediumMuons++;}
-	if(my_quality <= xAOD::Muon::Tight)
-	  tightMuons++;
-	if(my_quality <= xAOD::Muon::VeryLoose && !(my_quality <= xAOD::Muon::Loose))
-	  badMuons++;
-			  	
-	if(m_muonSelection.accept(**mu_itr)){
-	  //std::cout << "FAILED ACCEPT FUNCTION " << std::endl;
-	  //std:: cout << "Passes the selection (eta, IDCuts and quality) from the SelectorTool " << std::endl;
-	  passesCutsTool++; //passesTool = true;
-	}
-	
-	//std::cout << passesxAOD << " "  << passesTool << std::endl;
-	
-	// if((passesTool && !passesxAOD) || (!passesTool && passesxAOD)){
-	//   std::cout << "DISCREPANCY!!!! " << std::endl;
-	//   std::cout << "eta " << abseta << " IDCuts " << passesIDRequirements << std::endl;
-	//   std::cout << "Comparison of the quality: xAOD vs SelectorTool " << (*mu_itr)->quality() << " " <<  my_quality << std::endl;
-	if( (*mu_itr)->muonType() == xAOD::Muon_v1::Combined)
-	  std::cout << "combined " << std::endl;
-	else if( (*mu_itr)->muonType() == xAOD::Muon_v1::MuonStandAlone){
-	  nSAMuons++;
-	  std::cout << "stand-alone " << std::endl;}
-	else if( (*mu_itr)->muonType() == xAOD::Muon_v1::SegmentTagged)
-	  std::cout << "segment-tagged " << std::endl;
-	else if( (*mu_itr)->muonType() == xAOD::Muon_v1::CaloTagged)
-	  std::cout << "calo-tagged " << std::endl;
-	else if( (*mu_itr)->muonType() == xAOD::Muon_v1::SiliconAssociatedForwardMuon)
-	  std::cout << "forward " << std::endl;
-	else
-	  std::cout << "no type! " << std::endl;
+
+
+        passesIDRequirements = selectorTools[0]->passedIDCuts(**mu_itr);
+        passesPreselectionCuts = selectorTools[0]->passedMuonCuts(**mu_itr);
+	my_quality = selectorTools[0]->getQuality(**mu_itr);
+
+	//Print some general information about the muon
+	Info( APP_NAME, "Muon pT [GeV]:       %g ", std::abs((*mu_itr)->pt())/1000.);
+	Info( APP_NAME, "Muon eta, phi:       %g, %g ", (*mu_itr)->eta(),(*mu_itr)->phi());
+	Info( APP_NAME, "Muon muonType:       %d (%s)", (*mu_itr)->muonType(), typeNames[(*mu_itr)->muonType()].c_str());
+
+	Info( APP_NAME, "Muon quality (from tool, from xAOD):      %d, %d", my_quality, (*mu_itr)->quality());
+	Info( APP_NAME, "Muon passes cuts (ID hits, preselection): %d, %d", passesIDRequirements, passesPreselectionCuts);
+
+
+	//Check availability of variables and dump a message if they are missing
+	uint8_t PixHits = 0, PixDead = 0, SCTHits = 0, SCTDead = 0, PixHoles = 0, SCTHoles = 0, TRTHits = 0, TRTOut = 0;
+	uint8_t nprecisionLayers = 0, nprecisionHoleLayers = 0, cscUnspoiledEtaHits = 0;
 
 	if (!(*mu_itr)->primaryTrackParticle()->summaryValue(nprecisionLayers, xAOD::numberOfPrecisionLayers))
-	  std::cout << "no nprecisionlayers! " << std::endl;
+	   Info( APP_NAME, "no nprecisionlayers! ");
 	
 	if (!(*mu_itr)->primaryTrackParticle()->summaryValue(nprecisionHoleLayers, xAOD::numberOfPrecisionHoleLayers))
-	  std::cout << "no nprecisionholelayers! " << std::endl;
-	Info( APP_NAME, "Primary track: NPrecisionLayers %i NPrecisionHoleLayers %i ", nprecisionLayers, nprecisionHoleLayers);
-	  // Info( APP_NAME, "Muon: momSigBal %f ", momSigBal);
-
-	// }
-
-	if(!m_muonSelection.accept(**mu_itr)) continue;
+	  Info( APP_NAME, "no nprecisionholelayers! ");
 	
+	if (!(*mu_itr)->summaryValue(cscUnspoiledEtaHits, xAOD::MuonSummaryType::cscUnspoiledEtaHits))
+	  Info( APP_NAME, "no cscUnspoiledEtaHits! ");
+
 	if(!(*mu_itr)->primaryTrackParticle()->summaryValue(PixHits,xAOD::numberOfPixelHits))
 	  Info( APP_NAME, "Missing info!");
-	//Info( APP_NAME, "Primary track: PixHits %i ", PixHits);
 
 	if(!(*mu_itr)->primaryTrackParticle()->summaryValue(PixDead,xAOD::numberOfPixelDeadSensors))
 	  Info( APP_NAME, "Missing info!");
-	//Info( APP_NAME, "Primary track: PixDead %i ", PixDead);
 
 	if(!(*mu_itr)->primaryTrackParticle()->summaryValue(SCTHits,xAOD::numberOfSCTHits))
 	  Info( APP_NAME, "Missing info!");
-	//Info( APP_NAME, "Primary track: SCTHits %i ", SCTHits);
 
 	if(!(*mu_itr)->primaryTrackParticle()->summaryValue(SCTDead,xAOD::numberOfSCTDeadSensors))
 	  Info( APP_NAME, "Missing info!");
-	//Info( APP_NAME, "Primary track: SCTDead %i ", SCTDead);
-
 
 	if(!(*mu_itr)->primaryTrackParticle()->summaryValue(PixHoles,xAOD::numberOfPixelHoles))
 	  Info( APP_NAME, "Missing info!");
-	//Info( APP_NAME, "Primary track: PixHoles %i ", PixHoles);
 
 	if(!(*mu_itr)->primaryTrackParticle()->summaryValue(SCTHoles,xAOD::numberOfSCTHoles))
 	  Info( APP_NAME, "Missing info!");
-	//Info( APP_NAME, "Primary track: SCTHoles %i ", SCTHoles);
 
 	if(!(*mu_itr)->primaryTrackParticle()->summaryValue(TRTHits,xAOD::numberOfTRTHits))
 	  Info( APP_NAME, "Missing info!");
-	//Info( APP_NAME, "Primary track: TRTHits %i ", TRTHits);
 
 	if(!(*mu_itr)->primaryTrackParticle()->summaryValue(TRTOut,xAOD::numberOfTRTOutliers))
 	  Info( APP_NAME, "Missing info!");
-	//Info( APP_NAME, "Primary track: TRTOut %i ", TRTOut);
-
-	//uint8_t totTRT = TRTHits+TRTOut;
-
-	//Info( APP_NAME, "Muon eta:  %g ", std::abs((*mu_itr)->eta()));
-	//Info( APP_NAME, "TotTRT %i > 5; OutCond %i < %g ", totTRT, TRTOut, 0.9*totTRT);
-
-	// if (!((0.1<abseta && abseta<=1.9 && totTRT>5 && TRTOut<(0.9 * totTRT)) || (abseta <= 0.1 || abseta > 1.9)))
-	//   std::cout << "didn't pass the TRT cuts! v1 " << std::endl;
-
-	// if ((0.1<abseta && abseta<=1.9) && !(totTRT>5 && TRTOut<(0.9 * totTRT)))
-	//   std::cout << "didn't pass the TRT cuts! v2 " << std::endl;
-
-
-         // Select "good" muons:
-	//	if( ! m_muonSelection.accept( **mu_itr ) ) std::cout << "didn't pass! " << std::endl;
-
-	//	if( (*mu_itr)->quality() == xAOD::Muon_v1::Tight)
-	//	  std::cout << "tight " << std::endl;
-	//	else if( (*mu_itr)->quality() == xAOD::Muon_v1::Medium)
-	//	  std::cout << "medium " << std::endl;
-	//	else if( (*mu_itr)->quality() == xAOD::Muon_v1::Loose)
-	//	  std::cout << "loose " << std::endl;
-	//	else
-	//	  std::cout << "no quality! " << std::endl;
-
-	//	if( (*mu_itr)->muonType() == xAOD::Muon_v1::Combined)
-	//	  std::cout << "combined " << std::endl;
-	//	else if( (*mu_itr)->muonType() == xAOD::Muon_v1::MuonStandAlone)
-	//	  std::cout << "stand-alone " << std::endl;
-	//	else if( (*mu_itr)->muonType() == xAOD::Muon_v1::SegmentTagged)
-	//	  std::cout << "segment-tagged " << std::endl;
-	//	else if( (*mu_itr)->muonType() == xAOD::Muon_v1::CaloTagged)
-	//	  std::cout << "calo-tagged " << std::endl;
-	//	else if( (*mu_itr)->muonType() == xAOD::Muon_v1::SiliconAssociatedForwardMuon)
-	//	  std::cout << "forward " << std::endl;
-	//	else
-	//	  std::cout << "no type! " << std::endl;
-
-	//	Info(APP_NAME, "MyVerboseType %i ", (*mu_itr)->muonType());
-	//	Info(APP_NAME, "MyVerboseQuality %i ", (*mu_itr)->quality());
-	
-	//	Info(APP_NAME, "MyVerboseGetQuality %i ", m_muonSelection.getQuality( **mu_itr ));
-
-	//	if(m_muonSelection.getQuality( **mu_itr ) != xAOD::Muon_v1::Medium) {
-	  //	  Info(APP_NAME, "DIDNT PASS THE QUALITY CUTS!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!");
-	  //continue;
-	//	}
-	//	if(!(*mu_itr)->primaryTrackParticle()->summaryValue(SCTHoles,xAOD::numberOfSCTHoles))
-	// if (!(*mu_itr)->parameter(momSigBal, xAOD::Muon_v1::momentumBalanceSignificance))
-	//   std::cout << "No momSigBal " << std::endl;
-	
-	// Info( APP_NAME, "Muon: momSigBal %f ", momSigBal);
-
-	goodMuons++;
-	allgoodMuons++;
-
-         // // Print some info about the selected muon:
-         // Info( APP_NAME, "  Selected muon: eta = %g, phi = %g, pt = %g",
-         //       ( *mu_itr )->eta(), ( *mu_itr )->phi(), ( *mu_itr )->pt() );
-         // Info( APP_NAME, "    Primary track: eta = %g, phi = %g, pt = %g",
-         //       ( *mu_itr )->primaryTrackParticle()->eta(),
-         //       ( *mu_itr )->primaryTrackParticle()->phi(),
-         //       ( *mu_itr )->primaryTrackParticle()->pt() );
-
-         // const xAOD::TrackParticle* idtrack =
-         //    ( *mu_itr )->trackParticle( xAOD::Muon::InnerDetectorTrackParticle );
-         // if( idtrack ) {
-         //    Info( APP_NAME, "         ID track: eta = %g, phi = %g, pt = %g",
-         //          idtrack->eta(), idtrack->phi(), idtrack->pt() );
-         // }
 
+
+	//Now, let's check whether the muon passes the different working points and also whether it is flagged as bad
+
+	std::string selectionResults = "Muon selection acceptance:  ";
+	std::string badMuonResults =   "Bad muon flag:              ";
+
+	for (int wp = 0; wp < Nwp; wp++) {
+
+	  if (selectorTools[wp]->accept(*mu_itr)) {
+	    selectedMuons[wp]++;
+	    selectedMuonsEvent[wp]++;
+	    selectedMuonsType[(*mu_itr)->muonType()][wp]++;
+	    selectedMuonsEta[etaIndex][wp]++;
+	    selectionResults += "pass     ";
+
+	    if (!selectorTools[wp]->isBadMuon(**mu_itr)) {
+	      selectedMuonsNotBad[wp]++;
+	      selectedMuonsEventNotBad[wp]++;
+	      selectedMuonsTypeNotBad[(*mu_itr)->muonType()][wp]++;
+	      selectedMuonsEtaNotBad[etaIndex][wp]++;
+	    }
+	  }
+	  else
+	    selectionResults += "fail     ";
+
+	  if (!selectorTools[wp]->isBadMuon(**mu_itr))
+	    badMuonResults += "good     ";
+	  else
+	    badMuonResults += "bad      ";
+	}
+
+	//Print table of selection results for this muon
+	Info( APP_NAME, "%s", namesString.c_str() );
+	Info( APP_NAME, "%s", selectionResults.c_str() );
+	Info( APP_NAME, "%s", badMuonResults.c_str() );
+
+      } //done loop over muons
+
+
+      //Print table of number of selected muons in this event
+      std::string NselectedString = "Number of selected muons:     ";
+      std::string NselectedStringNotBad = "Including bad muon veto:      ";
+      for (int wp = 0; wp < Nwp; wp++) {
+	NselectedString += std::to_string(selectedMuonsEvent[wp]) + "        ";
+	NselectedStringNotBad += std::to_string(selectedMuonsEventNotBad[wp]) + "        ";
       }
 
-      Info( APP_NAME, "Number of good muons: %i",
-	    goodMuons );
+      Info( APP_NAME, "===== Event summary:");
+      Info( APP_NAME, "%s", namesString.c_str() );
+      Info( APP_NAME, "%s", NselectedString.c_str() );
+      Info( APP_NAME, "%s", NselectedStringNotBad.c_str() );
 
       // Close with a message:
       Info( APP_NAME,
@@ -356,23 +380,95 @@ int main( int argc, char* argv[] ) {
             static_cast< int >( ei->eventNumber() ),
             static_cast< int >( ei->runNumber() ),
             static_cast< int >( entry + 1 ) );
-   }
-
-   Info(APP_NAME, "All muons %i and good muons %i " , allMuons, allgoodMuons);
-   
-   Info(APP_NAME, "All muons %i, tight %i, medium %i, loose %i, veryloose %i, bad muons %i " , allMuons, tightMuons, mediumMuons, looseMuons, verylooseMuons, badMuons);
-   
-   Info(APP_NAME, "CaloTagged muons %i Stand-alone muons %i " , nCaloTagged, nSAMuons);
-
-   Info(APP_NAME, "Positive %i and negative muons %i " , nPositive, nNegative);
 
-   Info(APP_NAME, "Good muons xAOD %i and SelectorTool %i " , passesCutsxAOD, passesCutsTool);
+   } //done loop over events
+
+
+   //Now, let's summarize all we have found in the processed events
+
+   Info(APP_NAME, "======================================");
+   Info(APP_NAME, "========= Full run summary ===========");
+   Info(APP_NAME, "======================================");
+
+   Info(APP_NAME, "Processed %i events and %i muons" , static_cast< int >(entries), allMuons);
+
+   Info(APP_NAME, "%i positive and %i negative muons" , nPositive, nNegative);
+
+   Info(APP_NAME, "Selected muons by working point (numbers in parenthesis include bad muon veto):");
+   Info(APP_NAME, "--------------------------");
+   for (int wp = 0; wp < Nwp; wp++)
+     Info(APP_NAME, "%s: %s %i (%i)", WPnames[wp].c_str(), padding[wp].c_str(), selectedMuons[wp], selectedMuonsNotBad[wp]);
+   Info(APP_NAME, "--------------------------");
+
+   //Make table of selected muons by type and working point
+   Info(APP_NAME, "Selected muons by type and working point (numbers in parenthesis include bad muon veto):");
+   Info(APP_NAME, "---------------------------------------------------------------------------------------");
+   for (int l = 0; l < Nwp+2; l++) {
+
+     std::string line = "";
+     if (l == 0) { //line with type names
+       line += "              ";
+       for (int type = 0; type < Ntype; type++)
+	 line += typeNames[type] + "     ";
+     }
+     else if (l == 1) { //line for all muons inclusive
+       line += "All muons:      ";
+       for (int type = 0; type < Ntype; type++) {
+	 std::stringstream ss;
+	 ss << std::left << std::setw(16) << std::to_string(allMuonsType[type]);
+	 line += ss.str();
+       }
+     }
+     else { //lines for each of the working points
+       int wp = l - 2;
+       line += WPnames[wp] + ":" + padding[wp] + "     ";
+       for (int type = 0; type < Ntype; type++) {
+	 std::stringstream ss;
+	 ss << std::left << std::setw(16) << (std::to_string(selectedMuonsType[type][wp]) + " (" + std::to_string(selectedMuonsTypeNotBad[type][wp]) + ")");
+	 line += ss.str();
+       }
+     }
+
+     Info(APP_NAME, "%s", line.c_str());
+   }
+   Info(APP_NAME, "---------------------------------------------------------------------------------------");
+
+
+   //Make table of selected muons by |eta| and working point
+   Info(APP_NAME, "Selected muons by |eta| and working point (numbers in parenthesis include bad muon veto):");
+   Info(APP_NAME, "---------------------------------------------------------------------------------------");
+   for (int l = 0; l < Nwp+2; l++) {
+
+     std::string line = "";
+     if (l == 0) { //line with eta regions
+       line += "              ";
+       line += etaRegions;
+     }
+     else if (l == 1) { //line for all muons inclusive
+       line += "All muons:      ";
+       for (int eta = 0; eta < Neta; eta++) {
+	 std::stringstream ss;
+	 ss << std::left << std::setw(20) << std::to_string(allMuonsEta[eta]);
+	 line += ss.str();
+       }
+     }
+     else { //lines for each of the working points
+       int wp = l - 2;
+       line += WPnames[wp] + ":" + padding[wp] + "     ";
+       for (int eta = 0; eta < Neta; eta++) {
+	 std::stringstream ss;
+	 ss << std::left << std::setw(20) << (std::to_string(selectedMuonsEta[eta][wp]) + " (" + std::to_string(selectedMuonsEtaNotBad[eta][wp]) + ")");
+	 line += ss.str();
+       }
+     }
+
+     Info(APP_NAME, "%s", line.c_str());
+   }
+   Info(APP_NAME, "---------------------------------------------------------------------------------------");
 
    // Needed for Smart Slimming
    xAOD::IOStats::instance().stats().printSmartSlimmingBranchList();
-  
+
    // Return gracefully:
    return 0;
-
-
 }
diff --git a/PhysicsAnalysis/MuonID/MuonTagTools/MuonTagTools/MuonTagTool.h b/PhysicsAnalysis/MuonID/MuonTagTools/MuonTagTools/MuonTagTool.h
index 1817a047eac3f4ac77e8253831fcb5e057c4f7e3..2d5e293b099cc992a6aad25c777df008560d24de 100644
--- a/PhysicsAnalysis/MuonID/MuonTagTools/MuonTagTools/MuonTagTool.h
+++ b/PhysicsAnalysis/MuonID/MuonTagTools/MuonTagTools/MuonTagTool.h
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 #ifndef MUONTAGTOOL_H 
@@ -19,7 +19,7 @@ Purpose : build the Muon Tag objects - MuonTagCollection.h.
 #include <inttypes.h>
 #include "xAODMuon/MuonContainer.h"
 #include "MuonMomentumCorrections/IMuonCalibrationAndSmearingTool.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 #include "IsolationSelection/IIsolationSelectionTool.h"
 #include "xAODEventInfo/EventInfo.h"
 #include <map>
diff --git a/PhysicsAnalysis/MuonID/MuonTagTools/src/MuonTagTool.cxx b/PhysicsAnalysis/MuonID/MuonTagTools/src/MuonTagTool.cxx
index fd9ca1c60a044b35e97bf25630c2cbb014363cef..59d27766db79aefa9e971cc625f942292218b8e8 100644
--- a/PhysicsAnalysis/MuonID/MuonTagTools/src/MuonTagTool.cxx
+++ b/PhysicsAnalysis/MuonID/MuonTagTools/src/MuonTagTool.cxx
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 /*****************************************************************************
@@ -21,7 +21,6 @@ Purpose : create a collection of MuonTag
 #include "TagEvent/MuonAttributeNames.h"
 #include "AnalysisUtils/AnalysisMisc.h"
 #include "AthenaPoolUtilities/AthenaAttributeSpecification.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
 #include "xAODEventInfo/EventInfo.h"
 #include <sstream>
 #include <cmath>
diff --git a/PhysicsAnalysis/PrimaryDPDMaker/PrimaryDPDMaker/DRAW_ZMUMUSkimmingTool.h b/PhysicsAnalysis/PrimaryDPDMaker/PrimaryDPDMaker/DRAW_ZMUMUSkimmingTool.h
index 9092750caa730b84cfef4dff0d064178791ca742..d326797cee142a025a17b0ec18a958ebf3bf38e1 100644
--- a/PhysicsAnalysis/PrimaryDPDMaker/PrimaryDPDMaker/DRAW_ZMUMUSkimmingTool.h
+++ b/PhysicsAnalysis/PrimaryDPDMaker/PrimaryDPDMaker/DRAW_ZMUMUSkimmingTool.h
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 ///////////////////////////////////////////////////////////////////
@@ -16,10 +16,7 @@
 
 // DerivationFramework includes
 #include "DerivationFrameworkInterfaces/ISkimmingTool.h"
-
-namespace CP {
-  class IMuonSelectionTool;
-}
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
 namespace DerivationFramework {
 
diff --git a/PhysicsAnalysis/PrimaryDPDMaker/src/DRAW_ZMUMU_SkimmingTool.cxx b/PhysicsAnalysis/PrimaryDPDMaker/src/DRAW_ZMUMU_SkimmingTool.cxx
index df51466fd310d216e72e6166c808ec04a80671ac..f70549e5dcbd9b0873e73a6d5e265cbc28582edc 100644
--- a/PhysicsAnalysis/PrimaryDPDMaker/src/DRAW_ZMUMU_SkimmingTool.cxx
+++ b/PhysicsAnalysis/PrimaryDPDMaker/src/DRAW_ZMUMU_SkimmingTool.cxx
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 /////////////////////////////////////////////////////////////////
@@ -12,7 +12,6 @@
 
 #include "PrimaryDPDMaker/DRAW_ZMUMUSkimmingTool.h"
 #include "xAODMuon/MuonContainer.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
 #include <vector>
 #include <string>
 
diff --git a/PhysicsAnalysis/SUSYPhys/LongLivedParticleDPDMaker/LongLivedParticleDPDMaker/KinkTrkSingleJetMetFilterTool.h b/PhysicsAnalysis/SUSYPhys/LongLivedParticleDPDMaker/LongLivedParticleDPDMaker/KinkTrkSingleJetMetFilterTool.h
index 8dedf66cd6fc572857455dad113a368d0346f514..738d9bda0b7f8ea8bbd6da24358b7094eb05715a 100644
--- a/PhysicsAnalysis/SUSYPhys/LongLivedParticleDPDMaker/LongLivedParticleDPDMaker/KinkTrkSingleJetMetFilterTool.h
+++ b/PhysicsAnalysis/SUSYPhys/LongLivedParticleDPDMaker/LongLivedParticleDPDMaker/KinkTrkSingleJetMetFilterTool.h
@@ -1,11 +1,7 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
-///////////////////////////////////////////////////////////////////
-// KinkTrkSingleJetFilterTool.h, (c) ATLAS Detector software
-///////////////////////////////////////////////////////////////////
-
 #ifndef DERIVATIONFRAMEWORK_KINKTRKSINGLEJETMETFILTERTOOL_H
 #define DERIVATIONFRAMEWORK_KINKTRKSINGLEJETMETFILTERTOOL_H 1
 
@@ -18,7 +14,7 @@
 // DerivationFramework includes
 #include "DerivationFrameworkInterfaces/ISkimmingTool.h"
 
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
 namespace DerivationFramework {
   class KinkTrkSingleJetMetFilterTool : public AthAlgTool, public ISkimmingTool {
diff --git a/PhysicsAnalysis/SUSYPhys/LongLivedParticleDPDMaker/LongLivedParticleDPDMaker/KinkTrkZmumuTagTool.h b/PhysicsAnalysis/SUSYPhys/LongLivedParticleDPDMaker/LongLivedParticleDPDMaker/KinkTrkZmumuTagTool.h
index c891261094a1a3c98cc8df24cf2bbd0c092d8ecb..c62803d5a68d5cbd64608aad2dc302c2163f7f57 100644
--- a/PhysicsAnalysis/SUSYPhys/LongLivedParticleDPDMaker/LongLivedParticleDPDMaker/KinkTrkZmumuTagTool.h
+++ b/PhysicsAnalysis/SUSYPhys/LongLivedParticleDPDMaker/LongLivedParticleDPDMaker/KinkTrkZmumuTagTool.h
@@ -1,11 +1,7 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
-///////////////////////////////////////////////////////////////////
-// KinkTrkZmumuTagTool.h, (c) ATLAS Detector software
-///////////////////////////////////////////////////////////////////
-
 #ifndef DERIVATIONFRAMEWORK_KINKTRKZMUMUTAGTOOL_H
 #define DERIVATIONFRAMEWORK_KINKTRKZMUMUTAGTOOL_H 1
 
@@ -22,7 +18,7 @@
 
 #include "xAODMuon/MuonContainer.h"
 #include "xAODTracking/TrackParticle.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
 namespace DerivationFramework {
 
diff --git a/PhysicsAnalysis/SUSYPhys/SUSYTools/Root/Muons.cxx b/PhysicsAnalysis/SUSYPhys/SUSYTools/Root/Muons.cxx
index 044287f3a4bd79b96adea2a9d846f781dfdf228c..0447efd2f201aa7d30ce9e9f68f90ec4ed90fad8 100644
--- a/PhysicsAnalysis/SUSYPhys/SUSYTools/Root/Muons.cxx
+++ b/PhysicsAnalysis/SUSYPhys/SUSYTools/Root/Muons.cxx
@@ -1,27 +1,23 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 // This source file implements all of the functions related to <OBJECT>
 // in the SUSYObjDef_xAOD class
 
-// Local include(s):
 #include "SUSYTools/SUSYObjDef_xAOD.h"
 
 #include "xAODBase/IParticleHelpers.h"
 #include "EventPrimitives/EventPrimitivesHelpers.h"
-//#include "xAODPrimitives/IsolationType.h"
 #include "FourMomUtils/xAODP4Helpers.h"
 #include "xAODTracking/TrackParticlexAODHelpers.h"
 #include "AthContainers/ConstDataVector.h"
-//#include "PATInterfaces/SystematicsUtil.h"
 
-//#include "PileupReweighting/IPileupReweightingTool.h"
 #include "AsgAnalysisInterfaces/IPileupReweightingTool.h"
 
 #include "MuonMomentumCorrections/IMuonCalibrationAndSmearingTool.h"
 #include "MuonAnalysisInterfaces/IMuonEfficiencyScaleFactors.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 #include "MuonAnalysisInterfaces/IMuonTriggerScaleFactors.h"
 
 #include "IsolationCorrections/IIsolationCorrectionTool.h"
diff --git a/PhysicsAnalysis/SUSYPhys/SUSYTools/Root/SUSYObjDef_xAOD.cxx b/PhysicsAnalysis/SUSYPhys/SUSYTools/Root/SUSYObjDef_xAOD.cxx
index 3d11b4ae25c3f061d201c3fae630c8b081002ce3..af233c0211df2e33291ddf7fcc97cfce6694c8e2 100644
--- a/PhysicsAnalysis/SUSYPhys/SUSYTools/Root/SUSYObjDef_xAOD.cxx
+++ b/PhysicsAnalysis/SUSYPhys/SUSYTools/Root/SUSYObjDef_xAOD.cxx
@@ -2,7 +2,6 @@
   Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
-// Local include(s):
 #include "SUSYTools/SUSYObjDef_xAOD.h"
 
 #include "xAODBase/IParticleHelpers.h"
@@ -37,7 +36,7 @@
 #include "ElectronPhotonShowerShapeFudgeTool/IElectronPhotonShowerShapeFudgeTool.h"
 #include "ElectronPhotonSelectorTools/IEGammaAmbiguityTool.h"
 
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 #include "MuonMomentumCorrections/IMuonCalibrationAndSmearingTool.h"
 #include "MuonAnalysisInterfaces/IMuonEfficiencyScaleFactors.h"
 #include "MuonAnalysisInterfaces/IMuonTriggerScaleFactors.h"
diff --git a/PhysicsAnalysis/TopPhys/xAOD/TopCPTools/TopCPTools/TopMuonCPTools.h b/PhysicsAnalysis/TopPhys/xAOD/TopCPTools/TopCPTools/TopMuonCPTools.h
index cb3f285d92e93d04d5feabd5ac8499b623857d22..bf62cf5e8d63e1e53c5be178540d0d4e1e4f0027 100644
--- a/PhysicsAnalysis/TopPhys/xAOD/TopCPTools/TopCPTools/TopMuonCPTools.h
+++ b/PhysicsAnalysis/TopPhys/xAOD/TopCPTools/TopCPTools/TopMuonCPTools.h
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 #ifndef TOPCPTOOLS_TOPMUONCPTOOLS_H_
@@ -17,7 +17,7 @@
 
 // Muon include(s):
 #include "MuonMomentumCorrections/IMuonCalibrationAndSmearingTool.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 #include "MuonAnalysisInterfaces/IMuonTriggerScaleFactors.h"
 #include "MuonAnalysisInterfaces/IMuonEfficiencyScaleFactors.h"
 
diff --git a/PhysicsAnalysis/TopPhys/xAOD/TopEventSelectionTools/TopEventSelectionTools/NoBadMuonSelector.h b/PhysicsAnalysis/TopPhys/xAOD/TopEventSelectionTools/TopEventSelectionTools/NoBadMuonSelector.h
index 2d1c9bbc809e24198c37d0aafa13f9068369ac65..c0bef23a8488699ff22ca06ca6a87161b7f6d8e9 100644
--- a/PhysicsAnalysis/TopPhys/xAOD/TopEventSelectionTools/TopEventSelectionTools/NoBadMuonSelector.h
+++ b/PhysicsAnalysis/TopPhys/xAOD/TopEventSelectionTools/TopEventSelectionTools/NoBadMuonSelector.h
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 #ifndef NOBADMUONSELECTOR_H_
@@ -8,7 +8,7 @@
 #include "TopEventSelectionTools/EventSelectorBase.h"
 
 #include "AsgTools/ToolHandle.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
 namespace top {
 
diff --git a/PhysicsAnalysis/TopPhys/xAOD/TopObjectSelectionTools/TopObjectSelectionTools/AntiMuonMC15.h b/PhysicsAnalysis/TopPhys/xAOD/TopObjectSelectionTools/TopObjectSelectionTools/AntiMuonMC15.h
index c675f80eb7e932001e0b875feeef9e8e3bf12d81..359c0bcd48df3c32326c7e1fdb4becbc48943bf4 100644
--- a/PhysicsAnalysis/TopPhys/xAOD/TopObjectSelectionTools/TopObjectSelectionTools/AntiMuonMC15.h
+++ b/PhysicsAnalysis/TopPhys/xAOD/TopObjectSelectionTools/TopObjectSelectionTools/AntiMuonMC15.h
@@ -1,7 +1,7 @@
-/*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
-*/
-
+/*
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
+*/
+
 #ifndef ANTIMUONMC15_H_
 #define ANTIMUONMC15_H_
 
@@ -9,7 +9,7 @@
 #include "TopObjectSelectionTools/IsolationTools.h"
 
 #include "AsgTools/ToolHandle.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
 namespace top {
 
diff --git a/PhysicsAnalysis/TopPhys/xAOD/TopObjectSelectionTools/TopObjectSelectionTools/MuonMC15.h b/PhysicsAnalysis/TopPhys/xAOD/TopObjectSelectionTools/TopObjectSelectionTools/MuonMC15.h
index 49f0040a286df500953ed15bd30538dae0f77947..f0801c397a74b139ae39576185ccb498bf8ab399 100644
--- a/PhysicsAnalysis/TopPhys/xAOD/TopObjectSelectionTools/TopObjectSelectionTools/MuonMC15.h
+++ b/PhysicsAnalysis/TopPhys/xAOD/TopObjectSelectionTools/TopObjectSelectionTools/MuonMC15.h
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 #ifndef MUONMC15_H_
@@ -9,7 +9,7 @@
 #include "TopObjectSelectionTools/IsolationTools.h"
 
 #include "AsgTools/ToolHandle.h"
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
 namespace top {
 
diff --git a/PhysicsAnalysis/TopPhys/xAOD/TopSystematicObjectMaker/TopSystematicObjectMaker/MuonObjectCollectionMaker.h b/PhysicsAnalysis/TopPhys/xAOD/TopSystematicObjectMaker/TopSystematicObjectMaker/MuonObjectCollectionMaker.h
index 5cde1ae054e356a0143d93b37c9efb0371d8ebb3..c3ad4217749624bee4555b687e01e92f3997aba7 100644
--- a/PhysicsAnalysis/TopPhys/xAOD/TopSystematicObjectMaker/TopSystematicObjectMaker/MuonObjectCollectionMaker.h
+++ b/PhysicsAnalysis/TopPhys/xAOD/TopSystematicObjectMaker/TopSystematicObjectMaker/MuonObjectCollectionMaker.h
@@ -1,8 +1,7 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
-// $Id: MuonObjectCollectionMaker.h 799839 2017-03-08 11:07:28Z grancagn $
 #ifndef ANALYSISTOP_TOPSYSTEMATICOBJECTMAKER_MUONOBJECTCOLLECTIONMAKER_H
 #define ANALYSISTOP_TOPSYSTEMATICOBJECTMAKER_MUONOBJECTCOLLECTIONMAKER_H
 
@@ -35,7 +34,7 @@
 #include "MuonMomentumCorrections/IMuonCalibrationAndSmearingTool.h"
 #include "IsolationSelection/IIsolationSelectionTool.h"
 // the following is needed to make sure all muons for which d0sig is calculated are at least Loose
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 
 // Forward declaration(s):
 namespace top{
diff --git a/Reconstruction/DiTauRec/DiTauRec/ElMuFinder.h b/Reconstruction/DiTauRec/DiTauRec/ElMuFinder.h
index e114f613a3cccd99ab0d1b34d8115056326c211e..764dc9a649933b978de7f6568587f2d91cacfd3e 100644
--- a/Reconstruction/DiTauRec/DiTauRec/ElMuFinder.h
+++ b/Reconstruction/DiTauRec/DiTauRec/ElMuFinder.h
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 #ifndef DITAUREC_ELMUFINDER_H
@@ -10,9 +10,6 @@
 #include "StoreGate/ReadHandle.h"
 #include "GaudiKernel/ToolHandle.h"
 
-// #include "MuonSelectorTools/IMuonSelectionTool.h"
-// #include "MuonSelectorTools/errorcheck.h"
-#include "MuonSelectorTools/MuonSelectionTool.h"
 #include "xAODEgamma/ElectronContainer.h"
 #include "xAODMuon/MuonContainer.h"
 
diff --git a/Reconstruction/DiTauRec/src/ElMuFinder.cxx b/Reconstruction/DiTauRec/src/ElMuFinder.cxx
index 1cbc471cb07fa3c8742dbd64bbd94c03e76fc20c..ee7c170aac2fb2ebdc086cf54d34e0a1a085652b 100644
--- a/Reconstruction/DiTauRec/src/ElMuFinder.cxx
+++ b/Reconstruction/DiTauRec/src/ElMuFinder.cxx
@@ -1,8 +1,7 @@
 /*
-  Copyright (C) 2002-2017, 2019 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
-
 #include "DiTauRec/ElMuFinder.h"
 #include "DiTauRec/DiTauToolBase.h"
 
@@ -11,9 +10,6 @@
 #include "xAODEgamma/ElectronContainer.h"
 #include "xAODMuon/MuonContainer.h"
 
-// #include "MuonSelectorTools/errorcheck.h"
-#include "MuonSelectorTools/MuonSelectionTool.h"
-
 #include "tauRecTools/KineUtils.h"
 
 #include "StoreGate/ReadHandle.h"
diff --git a/Reconstruction/MET/METUtilities/src/METMakerAlg.cxx b/Reconstruction/MET/METUtilities/src/METMakerAlg.cxx
index ef9e041e61244028ea4d58810c7269df639fb5e4..51d92a4ee0055257d4f678570991f0a28873392f 100644
--- a/Reconstruction/MET/METUtilities/src/METMakerAlg.cxx
+++ b/Reconstruction/MET/METUtilities/src/METMakerAlg.cxx
@@ -1,9 +1,7 @@
 /*
-  Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
-// METMakerAlg.cxx
-
 #include "METMakerAlg.h"
 #include "METInterface/IMETMaker.h"
 
@@ -12,7 +10,7 @@
 #include "xAODMissingET/MissingETAssociationMap.h"
 #include "xAODMissingET/MissingETAssociationHelper.h"
 
-#include "MuonSelectorTools/IMuonSelectionTool.h"
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h"
 #include "EgammaAnalysisInterfaces/IAsgElectronLikelihoodTool.h"
 #include "EgammaAnalysisInterfaces/IAsgPhotonIsEMSelector.h"
 #include "PATCore/AcceptData.h"
diff --git a/Reconstruction/MuonIdentification/MuonCombinedBaseTools/src/MuonCreatorTool.h b/Reconstruction/MuonIdentification/MuonCombinedBaseTools/src/MuonCreatorTool.h
index fee05160323344efb5d0c14cb5e5632a491cadd6..f4e51d7939b98bc64070550032ab034b0a6b6fc8 100644
--- a/Reconstruction/MuonIdentification/MuonCombinedBaseTools/src/MuonCreatorTool.h
+++ b/Reconstruction/MuonIdentification/MuonCombinedBaseTools/src/MuonCreatorTool.h
@@ -22,7 +22,7 @@
 #include "xAODTracking/TrackParticleContainer.h"
 #include "MuonCombinedToolInterfaces/IMuonMomentumBalanceSignificance.h"
 #include "MuonCombinedToolInterfaces/IMuonScatteringAngleSignificance.h" 
-#include "MuonSelectorTools/IMuonSelectionTool.h" 
+#include "MuonAnalysisInterfaces/IMuonSelectionTool.h" 
 #include "RecoToolInterfaces/IParticleCaloExtensionTool.h"
 #include "TrkParametersIdentificationHelpers/TrackParametersIdHelper.h"
 #include "TrkSegment/SegmentCollection.h"