diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MSVVariablesFactory.h b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MSVVariablesFactory.h index 72f34061cef2cc704aa1f8c36335b99e8c0ffa1c..419986f708a8e850b4444fc440b161ff89597107 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MSVVariablesFactory.h +++ b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MSVVariablesFactory.h @@ -1,7 +1,7 @@ // -*- c++ -*- /* - Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ #ifndef BTAGTOOL_MSVVARIABLESFACTORY_C @@ -20,38 +20,33 @@ class StoreGateSvc; namespace Analysis { - - //static const InterfaceID IID_JetFitterVariablesFactory("Analysis::JetFitterVariablesFactory", 1, 0); - - class MSVVariablesFactory : public AthAlgTool , virtual public IMSVVariablesFactory { public: - MSVVariablesFactory(const std::string& name, const std::string& n, const IInterface* p); - virtual ~MSVVariablesFactory(); - + virtual ~MSVVariablesFactory() = default; + virtual StatusCode initialize() override; virtual StatusCode finalize() override; - virtual StatusCode fillMSVVariables(const xAOD::Jet &, xAOD::BTagging* BTag, const Trk::VxSecVKalVertexInfo* myInfoVKal, xAOD::VertexContainer* btagVertex, const xAOD::Vertex& PV, std::string basename) const override ; - virtual StatusCode createMSVContainer(const xAOD::Jet &, const Trk::VxSecVKalVertexInfo* myInfoVKal, xAOD::VertexContainer* btagVertex, const xAOD::Vertex& PV) const override; - - void setOrigin(const xAOD::Vertex* priVtx); - // static const InterfaceID& interfaceID() { return IID_JetFitterVariablesFactory; }; - + virtual StatusCode fillMSVVariables + (const xAOD::Jet &, xAOD::BTagging* BTag, + const Trk::VxSecVKalVertexInfo* myInfoVKal, + xAOD::VertexContainer* btagVertex, const xAOD::Vertex& PV, + std::string basename) const override ; + virtual StatusCode createMSVContainer + (const xAOD::Jet &, const Trk::VxSecVKalVertexInfo* myInfoVKal, + xAOD::VertexContainer* btagVertex, const xAOD::Vertex& PV) const override; private: double get3DSignificance(const xAOD::Vertex* priVertex, std::vector<const xAOD::Vertex*>& secVertex, const Amg::Vector3D jetDirection) const; - //const xAOD::Vertex* m_priVtx; }; - //inline void MSVVariablesFactory::setOrigin(const xAOD::Vertex* priVtx) { m_priVtx = priVtx; } }//end Analysis namespace diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MV2Tag.h b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MV2Tag.h index d9f904e498b6d7781b0d4e40df0017c041564632..7932c8c8d85ac0e0ff8b9668f0d6b74e9d0ef6e5 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MV2Tag.h +++ b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MV2Tag.h @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ #ifndef JETTAGTOOLS_MV2TAG_H @@ -19,6 +19,7 @@ #include <list> #include "MVAUtils/BDT.h" #include "JetTagCalibration/JetTagCalibCondData.h" +#include "TLorentzVector.h" namespace Analysis { @@ -32,9 +33,8 @@ namespace Analysis { /** Implementations of the methods defined in the abstract base class */ - virtual ~MV2Tag(); + virtual ~MV2Tag() = default; virtual StatusCode initialize()override ; - virtual StatusCode finalize() override; virtual void assignProbability(xAOD::BTagging* BTag, @@ -44,17 +44,13 @@ namespace Analysis { private: - std::string m_taggerName; std::string m_taggerNameBase; // unique name for regular and flip versions - std::string m_treeName; std::string m_varStrName; /** Key of calibration data: */ SG::ReadCondHandleKey<JetTagCalibCondData> m_readKey{this, "HistosKey", "JetTagCalibHistosKey", "Key of input (derived) JetTag calibration data"}; bool m_forceMV2CalibrationAlias; - std::string m_decTagName; std::string m_MV2CalibAlias; - std::string m_MV2cXX; std::string m_xAODBaseName; std::map<std::string, double > m_defaultvals; @@ -74,29 +70,9 @@ namespace Analysis { (m_runModus=1) where already made reference histograms are read.*/ std::string m_runModus; //!< 0=Do not read histos, 1=Read referece histos (analysis mode) - /** Storage for the primary vertex. Can be removed when JetTag provides origin(). */ - // this pointer does not need to be deleted in the destructor (because it - // points to something in storegate) - //const xAOD::Vertex* m_priVtx; - - /** reader to define the MVA algorithms */ - std::list<std::string> m_undefinedReaders; // keep track of undefined readers to prevent too many warnings. - - - std::string m_ip2d_infosource; - std::string m_ip3d_infosource; - std::string m_sv1_infosource; - std::string m_sv0_infosource; - std::string m_jftNN_infosource; - std::string m_jfprob_infosource; - std::string m_softmuon_infosource; - std::string m_trainingConfig; - float d0sgn_wrtJet(const TLorentzVector& jet, const TLorentzVector& trk, float d0sig) const; float z0sgn_wrtJet(float trackTheta, float trackZ0, float jetEta) const; - //void setInputVariables(xAOD::Jet& jetToTag, xAOD::BTagging* BTag);//for future - //void ClearInputs(); - //void PrintInputs(); + std::vector<float> CreateVariables (const std::map<std::string, double> &inputs, const std::vector<std::string>& inputVars) const; diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MultiSVTag.h b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MultiSVTag.h index c37c61d50fc2d5453952f7d66c2da3b1e73a7555..7b5faa879d53fe6de453e62897225bde6e826d00 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MultiSVTag.h +++ b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MultiSVTag.h @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ /*************************************************************************** @@ -13,9 +13,6 @@ #include "GaudiKernel/ToolHandle.h" #include "JetTagTools/ITagTool.h" #include "JetTagCalibration/JetTagCalibCondData.h" -#include "GeoPrimitives/GeoPrimitives.h" -#include "xAODTracking/TrackParticle.h" -#include "xAODTracking/TrackParticleContainer.h" #include <string> #include <vector> #include <map> @@ -29,9 +26,8 @@ namespace Analysis { public: MultiSVTag(const std::string&,const std::string&,const IInterface*); - virtual ~MultiSVTag(); + virtual ~MultiSVTag() = default; virtual StatusCode initialize() override; - virtual StatusCode finalize() override; virtual StatusCode tagJet(const xAOD::Vertex& priVtx, const xAOD::Jet& jetToTag, @@ -41,30 +37,19 @@ namespace Analysis virtual void finalizeHistos() override; private: - std::string m_taggerName; std::string m_taggerNameBase; // - std::string m_treeName; std::string m_varStrName; /** Key of calibration data: */ SG::ReadCondHandleKey<JetTagCalibCondData> m_readKey{this, "HistosKey", "JetTagCalibHistosKey", "Key of input (derived) JetTag calibration data"}; - std::string m_MultiSV; // std::string m_runModus; - std::string m_refType; - - int m_warnCounter; - std::vector<std::string> m_jetCollectionList; - std::vector<std::string> m_hypotheses; bool m_doForcedCalib; std::string m_ForcedCalibName; std::string m_secVxFinderName; - std::string m_xAODBaseName; - //... - //... - std::string m_sv0_infosource; + std::string m_sv1_infosource; struct Vars diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MultivariateTagManager.h b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MultivariateTagManager.h index 22bdbfc75b207f1847e2fa457f4e7c5d4b679474..3f129d6d6220f326fd36190eadb11f6fdc94e2fe 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MultivariateTagManager.h +++ b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/MultivariateTagManager.h @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ #ifndef BTAGTOOL_MULTIVARIATETAGMANAGER_C @@ -32,12 +32,10 @@ namespace Analysis { const std::string&, const IInterface*); - virtual ~MultivariateTagManager(){}; + virtual ~MultivariateTagManager() = default; virtual StatusCode initialize() override; - virtual StatusCode finalize() override; - virtual void finalizeHistos() override {}; - + virtual void finalizeHistos() override {}; virtual StatusCode tagJet(const xAOD::Vertex& priVtx, const xAOD::Jet& jetToTag, diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/NewLikelihoodTool.h b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/NewLikelihoodTool.h index 4fae388c94ec1b57acd5761bc27cecfc169cb06a..06686364a6bc17a3826ddc0ef3e3a5b8dc4a2a85 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/NewLikelihoodTool.h +++ b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/NewLikelihoodTool.h @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ #ifndef JETTAGTOOLS_NEWLIKELIHOODMULTIDTOOL_H @@ -11,7 +11,6 @@ ********************************************************/ #include "AthenaBaseComps/AthAlgTool.h" -#include "GaudiKernel/ToolHandle.h" #include "JetTagTools/LikelihoodComponents.h" #include "JetTagCalibration/JetTagCalibCondData.h" #include <string> @@ -28,10 +27,9 @@ class NewLikelihoodTool : public AthAlgTool { public: NewLikelihoodTool(const std::string&,const std::string&,const IInterface*); - virtual ~NewLikelihoodTool(); + virtual ~NewLikelihoodTool() = default; virtual StatusCode initialize() override; - virtual StatusCode finalize() override; static const InterfaceID& interfaceID() { return IID_NewLikelihoodTool; }; diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/SVForIPTool.h b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/SVForIPTool.h index d6598a8f147fcc3bbe5cc90edffea7add4539210..d5b7e6254c6a89c00ee91a67c0b76b7cde4f540e 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/SVForIPTool.h +++ b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/SVForIPTool.h @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ #ifndef JetTagTools_SVForIPTool_H @@ -15,10 +15,8 @@ @author Giacinto Piacquadio (giacinto.piacquadio AT physik.uni-freiburg.de) ********************************************************/ -#include "GeoPrimitives/GeoPrimitives.h" #include "AthenaBaseComps/AthAlgTool.h" #include "xAODTracking/Vertex.h" -#include "xAODTracking/VertexContainer.h" #include "xAODBTagging/BTagging.h" #include <vector> @@ -37,18 +35,14 @@ namespace Analysis { SVForIPTool(const std::string& name, const std::string& n, const IInterface* p); - virtual ~SVForIPTool(); - - virtual StatusCode initialize() override; - virtual StatusCode finalize() override; + virtual ~SVForIPTool() = default; - /**Method to get the B flight direction from the secondary vertex info */ + /**Method to get the B flight direction from the secondary vertex info */ void getDirectionFromSecondaryVertexInfo(Amg::Vector3D & SvxDirection, bool & canUseSvxDirection, xAOD::BTagging * BTag, const std::string & secVxFinderName, const xAOD::Vertex & priVtx) const; - // const Trk::RecVertex & priVtx); /**Method to get the tracks from V0 from the secondary vertex info */ diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/SVTag.h b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/SVTag.h index 8b8a34c3e4cb3b1c7b20a9e0425c35a572509cf0..31c128b08f8db1fa3971686f441d365d49c05f02 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/SVTag.h +++ b/PhysicsAnalysis/JetTagging/JetTagTools/JetTagTools/SVTag.h @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ /*************************************************************************** @@ -9,22 +9,15 @@ #ifndef JETTAGTOOLS_SVTAG_H #define JETTAGTOOLS_SVTAG_H -#include "GeoPrimitives/GeoPrimitives.h" #include "AthenaBaseComps/AthAlgTool.h" #include "GaudiKernel/ToolHandle.h" #include "JetTagTools/ITagTool.h" #include "xAODJet/Jet.h" -#include "xAODTracking/TrackParticle.h" -#include "xAODTracking/TrackParticleContainer.h" +#include "xAODTracking/Vertex.h" + #include <vector> #include <map> - -//namespace xAOD { class TrackParticle; class TrackParticleContainer; } -//namespace Trk { class RecVertex;} -//class Jet; -class StoreGateSvc; - namespace Analysis { class SVInfo; @@ -43,7 +36,7 @@ namespace Analysis const xAOD::Jet& jetToTag, xAOD::BTagging& BTag, const std::string &jetName) const override; - virtual void finalizeHistos() override; + virtual void finalizeHistos() override {}; private: @@ -51,9 +44,7 @@ namespace Analysis double get3DSignificance(const xAOD::Vertex& priVertex, std::vector<const xAOD::Vertex*>& secVertex, const Amg::Vector3D jetDirection) const; - // double get3DSignificance(const Trk::RecVertex & priVertex, - // std::vector<const Trk::RecVertex*> & secVertex, - // const Amg::Vector3D jetDirection); + double get3DSignificanceCorr(const xAOD::Vertex& priVertex, std::vector<const xAOD::Vertex*>& secVertex, const Amg::Vector3D jetDirection) const; @@ -71,11 +62,6 @@ namespace Analysis float m_pTjetmin; bool m_checkOverflows; double m_purificationDeltaR; - bool m_UseBinInterpol; - - /** information to persistify: */ - // std::string m_originalTPCollectionName; - // const xAOD::TrackParticleContainer* m_originalTPCollection; /** just print some info at the beginning */ void printParameterSettings(); diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/python/MV2TagConfig.py b/PhysicsAnalysis/JetTagging/JetTagTools/python/MV2TagConfig.py index 7df05c9137e788f5a553efe95b32c32a21e304cf..b2078c56076874e22fcc16e92142711acc4bcd11 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/python/MV2TagConfig.py +++ b/PhysicsAnalysis/JetTagging/JetTagTools/python/MV2TagConfig.py @@ -1,4 +1,4 @@ -# Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration +# Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration from AthenaConfiguration.ComponentAccumulator import ComponentAccumulator from AthenaConfiguration.ComponentFactory import CompFactory @@ -11,7 +11,6 @@ def MV2TagCfg( flags, name = 'MV2c10', scheme = '', useBTagFlagsDefaults = True, The following options have BTaggingFlags defaults: Runmodus default: BTagging.RunModus - taggerName default: "MV2c10" taggerNameBase default: "MV2c10" forceMV2CalibrationAlias default: BTaggingFlags.ForceMV2CalibrationAlias MV2CalibAlias default: BTaggingFlags.MV2CalibAlias @@ -31,7 +30,6 @@ def MV2TagCfg( flags, name = 'MV2c10', scheme = '', useBTagFlagsDefaults = True, if useBTagFlagsDefaults: defaults = { 'Runmodus' : flags.BTagging.RunModus, - 'taggerName' : basename, 'taggerNameBase' : basename, 'forceMV2CalibrationAlias' : ForceMV2CalibrationAlias, 'MV2CalibAlias' : MV2CalibAlias, diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/src/MSVVariablesFactory.cxx b/PhysicsAnalysis/JetTagging/JetTagTools/src/MSVVariablesFactory.cxx index 1280ef51383bc14b1fc230b6aadabe9608c81dde..74c4e36d1d59691da402fad095003405c94d9773 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/src/MSVVariablesFactory.cxx +++ b/PhysicsAnalysis/JetTagging/JetTagTools/src/MSVVariablesFactory.cxx @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ ///////////////////////////////////////////////////////////////////////////////////////////////////// @@ -12,18 +12,12 @@ /// the variables per vertex to MSV. /// /////////////////////////////////////////////////////////////////////////////////////////////////////// -//#include "VxJetVertex/VxVertexOnJetAxis.h" #include "VxVertex/VxTrackAtVertex.h" -//#include "VxJetVertex/VxJetCandidate.h" -//#include "VxJetVertex/SelectedTracksInJet.h" #include "TrkParameters/TrackParameters.h" -//#include "TrkNeutralParameters/MeasuredNeutralPerigee.h" #include "VxSecVertex/VxSecVertexInfo.h" #include "VxSecVertex/VxSecVKalVertexInfo.h" -//#include "TrkTrackLink/ITrackLink.h" - #include <TMath.h> #include "CLHEP/Vector/LorentzVector.h" @@ -33,137 +27,130 @@ #include "GeoPrimitives/GeoPrimitivesHelpers.h" -#include <vector> -// #include "xAODTracking/Vertex.h" #include "xAODTracking/TrackParticle.h" #include "xAODBTagging/SecVtxHelper.h" + +#include <vector> #include <string> namespace Analysis { MSVVariablesFactory::MSVVariablesFactory(const std::string& name, - const std::string& n, const IInterface* p): + const std::string& n, + const IInterface* p): AthAlgTool(name, n,p) -// m_secVxFinderName("InDetVKalVxInJetTool"), - { -// declareProperty("secVxFinderName",m_secVxFinderName); declareInterface<IMSVVariablesFactory>(this); } -///////////////////////////////////////////////////////////////////////////////////// -/// Destructor - check up memory allocation -/// delete any memory allocation on the heap - - MSVVariablesFactory::~MSVVariablesFactory() {} - -StatusCode MSVVariablesFactory::initialize() { - ATH_MSG_DEBUG(" Initialization of MSVVariablesFactory succesfull"); - return StatusCode::SUCCESS; -} + StatusCode MSVVariablesFactory::initialize() { + ATH_MSG_DEBUG(" Initialization of MSVVariablesFactory succesfull"); + return StatusCode::SUCCESS; + } -StatusCode MSVVariablesFactory::finalize() { - ATH_MSG_DEBUG(" Finalization of MSVVariablesFactory succesfull"); - return StatusCode::SUCCESS; -} + StatusCode MSVVariablesFactory::finalize() { + ATH_MSG_DEBUG(" Finalization of MSVVariablesFactory succesfull"); + return StatusCode::SUCCESS; + } - StatusCode MSVVariablesFactory::createMSVContainer(const xAOD::Jet &myJet, const Trk::VxSecVKalVertexInfo* myVertexInfoVKal, xAOD::VertexContainer* VertexContainer,const xAOD::Vertex& PrimaryVtx) const { + StatusCode MSVVariablesFactory::createMSVContainer + (const xAOD::Jet &myJet, const Trk::VxSecVKalVertexInfo* myVertexInfoVKal, + xAOD::VertexContainer* VertexContainer,const xAOD::Vertex& PrimaryVtx) const { + Amg::Vector3D jet_V3(myJet.p4().Px(), myJet.p4().Py(), myJet.p4().Pz()); float jetenergy=0.; const xAOD::Vertex* priVtx = &PrimaryVtx; std::vector< ElementLink< xAOD::VertexContainer > > MSVVertexLinks; const std::vector<xAOD::Vertex*> myVertices = myVertexInfoVKal->vertices(); - if(myVertices.size() == 0){ + if(myVertices.empty()){ ATH_MSG_DEBUG("#BTAG# no MSV vertices...fill default values only... "); - xAOD::Vertex* Vertex = new xAOD::Vertex(); - VertexContainer->push_back(Vertex); - xAOD::SecVtxHelper::setVertexMass(Vertex, -9.); - xAOD::SecVtxHelper::setEnergyFraction(Vertex, -9.); - xAOD::SecVtxHelper::setVtxNtrk(Vertex, -9); - xAOD::SecVtxHelper::setVtxpt(Vertex, -9.); - xAOD::SecVtxHelper::setVtxeta(Vertex, -9.); - xAOD::SecVtxHelper::setVtxphi(Vertex, -9.); - xAOD::SecVtxHelper::setVtxnormDist(Vertex, -9.); + xAOD::Vertex* vertex = new xAOD::Vertex(); + VertexContainer->push_back(vertex); + xAOD::SecVtxHelper::setVertexMass(vertex, -9.); + xAOD::SecVtxHelper::setEnergyFraction(vertex, -9.); + xAOD::SecVtxHelper::setVtxNtrk(vertex, -9); + xAOD::SecVtxHelper::setVtxpt(vertex, -9.); + xAOD::SecVtxHelper::setVtxeta(vertex, -9.); + xAOD::SecVtxHelper::setVtxphi(vertex, -9.); + xAOD::SecVtxHelper::setVtxnormDist(vertex, -9.); return StatusCode::SUCCESS; - } - std::vector<xAOD::Vertex*>::const_iterator verticesBegin = myVertexInfoVKal->vertices().begin(); - std::vector<xAOD::Vertex*>::const_iterator verticesEnd = myVertexInfoVKal->vertices().end(); jetenergy = myVertexInfoVKal->energyTrkInJet(); - for (std::vector<xAOD::Vertex*>::const_iterator verticesIter=verticesBegin; verticesIter!=verticesEnd;++verticesIter) { - - xAOD::Vertex* Vertex = *verticesIter; - VertexContainer->push_back(Vertex); + for (const auto& vertex : myVertexInfoVKal->vertices()){ + VertexContainer->push_back(vertex); //additional info per vertex double sumpx = 0.0; double sumpy = 0.0; double sumpz = 0.0; double sume = 0.0; - const std::vector<ElementLink<xAOD::TrackParticleContainer> > myTrackLinks = (*verticesIter)->trackParticleLinks(); - if (myTrackLinks.size()==0) { + const std::vector<ElementLink<xAOD::TrackParticleContainer> > myTrackLinks = vertex->trackParticleLinks(); + if (myTrackLinks.empty()) { ATH_MSG_WARNING("#BTAG# No Track Links attached to the track at the sec vertex... "); } int npsec = 0; - const std::vector<Trk::VxTrackAtVertex> myTracks=(*verticesIter)->vxTrackAtVertex(); - if (myTracks.size()!=0) { + const std::vector<Trk::VxTrackAtVertex> myTracks=vertex->vxTrackAtVertex(); + if (!myTracks.empty()) { npsec=myTracks.size(); - const std::vector<Trk::VxTrackAtVertex>::const_iterator tracksBegin=myTracks.begin(); - const std::vector<Trk::VxTrackAtVertex>::const_iterator tracksEnd=myTracks.end(); - for (std::vector<Trk::VxTrackAtVertex>::const_iterator tracksIter=tracksBegin; tracksIter!=tracksEnd;++tracksIter) { - const Trk::Perigee* perigee = dynamic_cast<const Trk::Perigee*>((*tracksIter).perigeeAtVertex()); + for (const auto& track : myTracks) { + const Trk::Perigee* perigee = dynamic_cast<const Trk::Perigee*>(track.perigeeAtVertex()); if(perigee){ sumpx += perigee->momentum().x(); sumpy += perigee->momentum().y(); sumpz += perigee->momentum().z(); - sume +=sqrt(perigee->momentum().mag()*perigee->momentum().mag() + 139.5702*139.5702 ); + sume += std::hypot(perigee->momentum().mag(), 139.5702); }else{ - ATH_MSG_WARNING("#BTAG# perigee for VxTrackAtVertex not found"); + ATH_MSG_WARNING("#BTAG# perigee for VxTrackAtVertex not found"); } } } + CLHEP::HepLorentzVector vtxp4(sumpx,sumpy,sumpz,sume); - float efrac = (jetenergy>0)?vtxp4.e()/jetenergy:0; - xAOD::SecVtxHelper::setVertexMass(Vertex, vtxp4.m()); - xAOD::SecVtxHelper::setEnergyFraction(Vertex, efrac); - xAOD::SecVtxHelper::setVtxNtrk(Vertex, npsec); - xAOD::SecVtxHelper::setVtxpt(Vertex, vtxp4.perp()); - xAOD::SecVtxHelper::setVtxeta(Vertex, vtxp4.eta()); - xAOD::SecVtxHelper::setVtxphi(Vertex, vtxp4.phi()); + float efrac = (jetenergy>0) ? vtxp4.e()/jetenergy : 0; + xAOD::SecVtxHelper::setVertexMass(vertex, vtxp4.m()); + xAOD::SecVtxHelper::setEnergyFraction(vertex, efrac); + xAOD::SecVtxHelper::setVtxNtrk(vertex, npsec); + xAOD::SecVtxHelper::setVtxpt(vertex, vtxp4.perp()); + xAOD::SecVtxHelper::setVtxeta(vertex, vtxp4.eta()); + xAOD::SecVtxHelper::setVtxphi(vertex, vtxp4.phi()); ATH_MSG_DEBUG("#BTAG# mass per vertex = "<<vtxp4.m()); - double localdistnrm=0; + double localdistnrm = 0; std::vector<const xAOD::Vertex*> vecVtxHolder; - vecVtxHolder.push_back(*verticesIter); + vecVtxHolder.push_back(vertex); - ATH_MSG_DEBUG("Factory PVX x = " << priVtx->x() << " y = " << priVtx->y() << " z = " << priVtx->z()); + ATH_MSG_DEBUG("Factory PVX x = " << priVtx->x() << " y = " << priVtx->y() << " z = " << priVtx->z()); if (priVtx) { - localdistnrm=get3DSignificance(priVtx, vecVtxHolder, Amg::Vector3D(myJet.p4().Px(),myJet.p4().Py(),myJet.p4().Pz())); + localdistnrm = get3DSignificance(priVtx, vecVtxHolder, jet_V3); } else { ATH_MSG_WARNING("#BTAG# Tagging requested, but no primary vertex supplied."); - localdistnrm=0.; } - xAOD::SecVtxHelper::setVtxnormDist(Vertex, localdistnrm); + xAOD::SecVtxHelper::setVtxnormDist(vertex, localdistnrm); //track links, - Vertex->setTrackParticleLinks(myTrackLinks); + vertex->setTrackParticleLinks(myTrackLinks); } //end loop vertexcontainer return StatusCode::SUCCESS; } - StatusCode MSVVariablesFactory::fillMSVVariables(const xAOD::Jet &myJet, xAOD::BTagging* BTag, const Trk::VxSecVKalVertexInfo* myVertexInfoVKal, xAOD::VertexContainer* VertexContainer,const xAOD::Vertex& PrimaryVtx, std::string basename) const { - //... + StatusCode MSVVariablesFactory::fillMSVVariables + (const xAOD::Jet &myJet, xAOD::BTagging* BTag, + const Trk::VxSecVKalVertexInfo* myVertexInfoVKal, + xAOD::VertexContainer* VertexContainer, const xAOD::Vertex& PrimaryVtx, + std::string basename) const { + + Amg::Vector3D jet_V3(myJet.p4().Px(), myJet.p4().Py(), myJet.p4().Pz()); int nvsec = 0; - float jetenergy=0.; + float jetenergy = 0.; int n2t = 0; - float distnrm=0.; + float distnrm = 0.; const xAOD::Vertex* priVtx = &PrimaryVtx; std::vector< ElementLink< xAOD::VertexContainer > > MSVVertexLinks; const std::vector<xAOD::Vertex*> myVertices = myVertexInfoVKal->vertices(); - if(myVertices.size() == 0){ + if(myVertices.empty()){ ATH_MSG_DEBUG("#BTAG# no MSV vertices...fill default values only... "); BTag->setVariable<int>(basename, "N2Tpair", n2t); BTag->setVariable<float>(basename, "energyTrkInJet", jetenergy); @@ -171,20 +158,17 @@ StatusCode MSVVariablesFactory::finalize() { BTag->setVariable<float>(basename, "normdist", distnrm); BTag->setVariable<std::vector<ElementLink<xAOD::VertexContainer> > >(basename, "vertices", MSVVertexLinks); BTag->setDynVxELName(basename, "vertices"); - xAOD::Vertex* Vertex = new xAOD::Vertex(); - VertexContainer->push_back(Vertex); - xAOD::SecVtxHelper::setVertexMass(Vertex, -9.); - xAOD::SecVtxHelper::setEnergyFraction(Vertex, -9.); - xAOD::SecVtxHelper::setVtxNtrk(Vertex, -9); - xAOD::SecVtxHelper::setVtxpt(Vertex, -9.); - xAOD::SecVtxHelper::setVtxeta(Vertex, -9.); - xAOD::SecVtxHelper::setVtxphi(Vertex, -9.); - xAOD::SecVtxHelper::setVtxnormDist(Vertex, -9.); + xAOD::Vertex* vertex = new xAOD::Vertex(); + VertexContainer->push_back(vertex); + xAOD::SecVtxHelper::setVertexMass(vertex, -9.); + xAOD::SecVtxHelper::setEnergyFraction(vertex, -9.); + xAOD::SecVtxHelper::setVtxNtrk(vertex, -9); + xAOD::SecVtxHelper::setVtxpt(vertex, -9.); + xAOD::SecVtxHelper::setVtxeta(vertex, -9.); + xAOD::SecVtxHelper::setVtxphi(vertex, -9.); + xAOD::SecVtxHelper::setVtxnormDist(vertex, -9.); return StatusCode::SUCCESS; - } - std::vector<xAOD::Vertex*>::const_iterator verticesBegin = myVertexInfoVKal->vertices().begin(); - std::vector<xAOD::Vertex*>::const_iterator verticesEnd = myVertexInfoVKal->vertices().end(); jetenergy = myVertexInfoVKal->energyTrkInJet(); n2t = myVertexInfoVKal->n2trackvertices(); @@ -192,65 +176,61 @@ StatusCode MSVVariablesFactory::finalize() { BTag->setVariable<float>(basename, "energyTrkInJet", jetenergy); std::vector<const xAOD::Vertex*> vecVertices; - for (std::vector<xAOD::Vertex*>::const_iterator verticesIter=verticesBegin; verticesIter!=verticesEnd;++verticesIter) { - - xAOD::Vertex* Vertex = *verticesIter; - VertexContainer->push_back(Vertex); + for (const auto& vertex : myVertexInfoVKal->vertices()) { + VertexContainer->push_back(vertex); //additional info per vertex - vecVertices.push_back(*verticesIter); + vecVertices.push_back(vertex); double sumpx = 0.0; double sumpy = 0.0; double sumpz = 0.0; double sume = 0.0; - const std::vector<ElementLink<xAOD::TrackParticleContainer> > myTrackLinks = (*verticesIter)->trackParticleLinks(); - if (myTrackLinks.size()==0) { + const std::vector<ElementLink<xAOD::TrackParticleContainer> > myTrackLinks = vertex->trackParticleLinks(); + if (myTrackLinks.empty()) { ATH_MSG_WARNING("#BTAG# No Track Links attached to the track at the sec vertex... "); } int npsec = 0; - const std::vector<Trk::VxTrackAtVertex> myTracks=(*verticesIter)->vxTrackAtVertex(); - if (myTracks.size()!=0) { + const std::vector<Trk::VxTrackAtVertex> myTracks=vertex->vxTrackAtVertex(); + if (!myTracks.empty()) { npsec=myTracks.size(); - const std::vector<Trk::VxTrackAtVertex>::const_iterator tracksBegin=myTracks.begin(); - const std::vector<Trk::VxTrackAtVertex>::const_iterator tracksEnd=myTracks.end(); - for (std::vector<Trk::VxTrackAtVertex>::const_iterator tracksIter=tracksBegin; tracksIter!=tracksEnd;++tracksIter) { - const Trk::Perigee* perigee = dynamic_cast<const Trk::Perigee*>((*tracksIter).perigeeAtVertex()); + for (const auto& track : myTracks) { + const Trk::Perigee* perigee = dynamic_cast<const Trk::Perigee*>(track.perigeeAtVertex()); if(perigee){ sumpx += perigee->momentum().x(); sumpy += perigee->momentum().y(); sumpz += perigee->momentum().z(); - sume +=sqrt(perigee->momentum().mag()*perigee->momentum().mag() + 139.5702*139.5702 ); + sume += std::hypot(perigee->momentum().mag(), 139.5702); }else{ - ATH_MSG_WARNING("#BTAG# perigee for VxTrackAtVertex not found"); + ATH_MSG_WARNING("#BTAG# perigee for VxTrackAtVertex not found"); } } } + CLHEP::HepLorentzVector vtxp4(sumpx,sumpy,sumpz,sume); - float efrac = (jetenergy>0)?vtxp4.e()/jetenergy:0; - xAOD::SecVtxHelper::setVertexMass(Vertex, vtxp4.m()); - xAOD::SecVtxHelper::setEnergyFraction(Vertex, efrac); - xAOD::SecVtxHelper::setVtxNtrk(Vertex, npsec); - xAOD::SecVtxHelper::setVtxpt(Vertex, vtxp4.perp()); - xAOD::SecVtxHelper::setVtxeta(Vertex, vtxp4.eta()); - xAOD::SecVtxHelper::setVtxphi(Vertex, vtxp4.phi()); + float efrac = (jetenergy>0) ? vtxp4.e()/jetenergy : 0; + xAOD::SecVtxHelper::setVertexMass(vertex, vtxp4.m()); + xAOD::SecVtxHelper::setEnergyFraction(vertex, efrac); + xAOD::SecVtxHelper::setVtxNtrk(vertex, npsec); + xAOD::SecVtxHelper::setVtxpt(vertex, vtxp4.perp()); + xAOD::SecVtxHelper::setVtxeta(vertex, vtxp4.eta()); + xAOD::SecVtxHelper::setVtxphi(vertex, vtxp4.phi()); ATH_MSG_DEBUG("#BTAG# mass per vertex = "<<vtxp4.m()); - double localdistnrm=0; + double localdistnrm = 0; std::vector<const xAOD::Vertex*> vecVtxHolder; - vecVtxHolder.push_back(*verticesIter); + vecVtxHolder.push_back(vertex); - ATH_MSG_DEBUG("Factory PVX x = " << priVtx->x() << " y = " << priVtx->y() << " z = " << priVtx->z()); + ATH_MSG_DEBUG("Factory PVX x = " << priVtx->x() << " y = " << priVtx->y() << " z = " << priVtx->z()); if (priVtx) { - localdistnrm=get3DSignificance(priVtx, vecVtxHolder, Amg::Vector3D(myJet.p4().Px(),myJet.p4().Py(),myJet.p4().Pz())); + localdistnrm = get3DSignificance(priVtx, vecVtxHolder, jet_V3); } else { ATH_MSG_WARNING("#BTAG# Tagging requested, but no primary vertex supplied."); - localdistnrm=0.; } - xAOD::SecVtxHelper::setVtxnormDist(Vertex, localdistnrm); + xAOD::SecVtxHelper::setVtxnormDist(vertex, localdistnrm); //track links, - Vertex->setTrackParticleLinks(myTrackLinks); + vertex->setTrackParticleLinks(myTrackLinks); ElementLink< xAOD::VertexContainer> linkBTagVertex; - linkBTagVertex.toContainedElement(*VertexContainer, Vertex); + linkBTagVertex.toContainedElement(*VertexContainer, vertex); MSVVertexLinks.push_back(linkBTagVertex); } //end loop vertexcontainer @@ -258,7 +238,7 @@ StatusCode MSVVariablesFactory::finalize() { BTag->setDynVxELName(basename, "vertices"); if (priVtx) { - distnrm=get3DSignificance(priVtx, vecVertices, Amg::Vector3D(myJet.p4().Px(),myJet.p4().Py(),myJet.p4().Pz())); + distnrm = get3DSignificance(priVtx, vecVertices, jet_V3); } else { ATH_MSG_WARNING("#BTAG# Tagging requested, but no primary vertex supplied."); distnrm=0.; @@ -267,42 +247,38 @@ StatusCode MSVVariablesFactory::finalize() { BTag->setVariable<int>(basename, "nvsec", nvsec); BTag->setVariable<float>(basename, "normdist", distnrm); - //... return StatusCode::SUCCESS; - } - double MSVVariablesFactory::get3DSignificance(const xAOD::Vertex* priVertex, - std::vector<const xAOD::Vertex*>& secVertex, - const Amg::Vector3D jetDirection) const { + double MSVVariablesFactory::get3DSignificance + (const xAOD::Vertex* priVertex, + std::vector<const xAOD::Vertex*>& secVertex, + const Amg::Vector3D jetDirection) const { + if(!secVertex.size()) return 0; std::vector<Amg::Vector3D> positions; std::vector<AmgSymMatrix(3)> weightMatrices; - std::vector<const xAOD::Vertex*>::const_iterator secEnd = secVertex.end(); - for (std::vector<const xAOD::Vertex*>::const_iterator secIter = secVertex.begin(); secIter != secEnd; ++secIter){ - positions.push_back((*secIter)->position()); - weightMatrices.push_back((*secIter)->covariancePosition().inverse()); - } Amg::Vector3D weightTimesPosition(0.,0.,0.); AmgSymMatrix(3) sumWeights; sumWeights.setZero(); - int count=0; - for (std::vector<const xAOD::Vertex*>::const_iterator secIter = secVertex.begin(); secIter != secEnd; ++secIter) { - weightTimesPosition+=(weightMatrices[count])*positions[count]; - sumWeights+=(weightMatrices[count]); - ++count; + for (const auto& vertex : secVertex) { + positions.push_back(vertex->position()); + weightMatrices.push_back(vertex->covariancePosition().inverse()); + weightTimesPosition += weightMatrices.back() * positions.back(); + sumWeights += weightMatrices.back(); } + bool invertible; AmgSymMatrix(3) meanCovariance; meanCovariance.setZero(); sumWeights.computeInverseWithCheck(meanCovariance, invertible); - if (! invertible) { + if (!invertible) { ATH_MSG_WARNING("#BTAG# Could not invert sum of sec vtx matrices"); return 0.; } - Amg::Vector3D meanPosition=meanCovariance*weightTimesPosition; + Amg::Vector3D meanPosition = meanCovariance * weightTimesPosition; AmgSymMatrix(3) covariance = meanCovariance + priVertex->covariancePosition(); double Lx = meanPosition[0]-priVertex->position().x(); @@ -315,11 +291,11 @@ StatusCode MSVVariablesFactory::finalize() { double dLdLy = Ly * inv_decaylength; double dLdLz = Lz * inv_decaylength; double decaylength_err = sqrt(dLdLx*dLdLx*covariance(0,0) + - dLdLy*dLdLy*covariance(1,1) + - dLdLz*dLdLz*covariance(2,2) + - 2.*dLdLx*dLdLy*covariance(0,1) + - 2.*dLdLx*dLdLz*covariance(0,2) + - 2.*dLdLy*dLdLz*covariance(1,2)); + dLdLy*dLdLy*covariance(1,1) + + dLdLz*dLdLz*covariance(2,2) + + 2.*dLdLx*dLdLy*covariance(0,1) + + 2.*dLdLx*dLdLz*covariance(0,2) + + 2.*dLdLy*dLdLz*covariance(1,2)); double decaylength_significance = 0.; if (decaylength_err != 0.) decaylength_significance = decaylength/decaylength_err; @@ -328,7 +304,6 @@ StatusCode MSVVariablesFactory::finalize() { return decaylength_significance; - } }//end Analysis namespace diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/src/MV2Tag.cxx b/PhysicsAnalysis/JetTagging/JetTagTools/src/MV2Tag.cxx index 9f5f50ebd32723a504c680d5897a13dbd1b6a4f8..c535eec838557b4a3dafc4e09cb631ee9f63b078 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/src/MV2Tag.cxx +++ b/PhysicsAnalysis/JetTagging/JetTagTools/src/MV2Tag.cxx @@ -1,14 +1,10 @@ /* - Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ -#include "GaudiKernel/IToolSvc.h" -#include "xAODTracking/TrackParticle.h" - #include "JetTagTools/MV2Tag.h" #include "xAODBTagging/BTagging.h" -#include "xAODJet/Jet.h" #include <fstream> #include <sstream> @@ -16,19 +12,13 @@ #include <iostream> #include <limits> #include <cmath> - -#include "xAODTracking/Vertex.h" -#include "xAODTracking/VertexContainer.h" - -#include "AthenaKernel/Units.h" -#include <fstream> -#include <algorithm> #include <utility> #include <vector> #include <map> #include <list> #include <math.h> /* hypot */ +#include "AthenaKernel/Units.h" using Athena::Units::GeV; @@ -47,18 +37,15 @@ namespace Analysis { // force MV2 to always use a calibration derived from MV2CalibAlias jet collection declareProperty("forceMV2CalibrationAlias", m_forceMV2CalibrationAlias = true); - declareProperty("MV2CalibAlias", m_MV2CalibAlias = "AntiKt4TopoEM"); + declareProperty("MV2CalibAlias", m_MV2CalibAlias = "AntiKt4EMTopo"); // global configuration: declareProperty("Runmodus", m_runModus); - //declareProperty("DecorateMvaInputs", m_decorateBTaggingObj=false); declareProperty("xAODBaseName", m_xAODBaseName);//"MV2c20" or etc. // which calibration folder to use declareProperty("taggerNameBase", m_taggerNameBase = "MV2"); - declareProperty("taggerName", m_taggerName = "MV2"); - declareProperty("decTagName", m_decTagName = "MV2_inputs"); declareProperty("defaultvals", m_defaultvals ); declareProperty("MVTMvariableNames", m_MVTM_name_translations ); @@ -67,17 +54,11 @@ namespace Analysis { } - MV2Tag::~MV2Tag() { - - } - - StatusCode MV2Tag::initialize() { m_disableAlgo=false; m_warnCounter=0; - m_treeName = "BDT"; m_varStrName = "variables"; // prepare readKey for calibration data: @@ -87,16 +68,10 @@ namespace Analysis { m_MVTM_name_backtrans[p.second] = p.first; } - //m_egammaBDTs.clear(); return StatusCode::SUCCESS; } - StatusCode MV2Tag::finalize() { - ATH_MSG_DEBUG("#BTAG# Finalizing MV2."); - return StatusCode::SUCCESS; - } - void MV2Tag::assignProbability(xAOD::BTagging *BTag, const std::map<std::string, double> &inputs, const std::string& assigned_jet_author) const diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/src/MultiSVTag.cxx b/PhysicsAnalysis/JetTagging/JetTagTools/src/MultiSVTag.cxx index 4d41247ace14ed278b52fd8905c106ec56bb5ab9..92d549a87f4c8af8cf21488958d16ed3bacb5112 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/src/MultiSVTag.cxx +++ b/PhysicsAnalysis/JetTagging/JetTagTools/src/MultiSVTag.cxx @@ -1,29 +1,18 @@ /* - Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ /*************************************************************************** MultiSVTag.cxx ***************************************************************************/ #include "JetTagTools/MultiSVTag.h" -#include "GaudiKernel/IToolSvc.h" -#include "Navigation/NavigationToken.h" -#include "GaudiKernel/ITHistSvc.h" -#include "JetTagTools/HistoHelperRoot.h" - -#include "VxSecVertex/VxSecVertexInfo.h" -#include "VxSecVertex/VxSecVKalVertexInfo.h" -#include "TrkLinks/LinkToXAODTrackParticle.h" + #include "xAODJet/Jet.h" -#include "xAODTracking/TrackParticle.h" +#include "xAODTracking/TrackParticleContainer.h" #include "xAODTracking/Vertex.h" #include "xAODTracking/VertexContainer.h" #include "xAODBTagging/SecVtxHelper.h" -#include "GeoPrimitives/GeoPrimitivesHelpers.h" -#include "CLHEP/Vector/LorentzVector.h" -#include "VxVertex/RecVertex.h" -#include "VxVertex/VxTrackAtVertex.h" #include <fstream> #include <algorithm> #include <utility> @@ -45,26 +34,14 @@ namespace Analysis m_runModus("analysis") { declareProperty("Runmodus", m_runModus= "analysis"); - declareProperty("jetCollectionList", m_jetCollectionList); declareProperty("useForcedCalibration", m_doForcedCalib = false); declareProperty("ForcedCalibrationName", m_ForcedCalibName = "AntiKt4TopoEM");//Cone4H1Tower declareProperty("SecVxFinderName",m_secVxFinderName); declareProperty("taggerNameBase",m_taggerNameBase = "MultiSVbb1"); - declareProperty("taggerName", m_taggerName = "MultiSVbb1"); - declareProperty("xAODBaseName",m_xAODBaseName); - declareProperty("inputSV0SourceName", m_sv0_infosource = "SV0"); declareProperty("inputSV1SourceName", m_sv1_infosource = "SV1"); } - MultiSVTag::~MultiSVTag() { - } - StatusCode MultiSVTag::initialize() { - // define tagger name: - - m_warnCounter=0; - - m_treeName = "BDT"; m_varStrName = "variables"; // prepare readKey for calibration data: @@ -72,10 +49,6 @@ namespace Analysis return StatusCode::SUCCESS; } - StatusCode MultiSVTag::finalize(){ - return StatusCode::SUCCESS; - } - StatusCode MultiSVTag::tagJet(const xAOD::Vertex& priVtx, const xAOD::Jet& jetToTag, xAOD::BTagging& BTag, @@ -150,24 +123,23 @@ namespace Analysis std::vector<float> v_vtxy = std::vector<float>(10,0); std::vector<float> v_vtxz = std::vector<float>(10,0); // loop in msv vertices - if(msvVertices.size()>0){ - const std::vector<ElementLink<xAOD::VertexContainer> >::const_iterator verticesEnd = msvVertices.end(); - for(std::vector<ElementLink<xAOD::VertexContainer> >::const_iterator vtxIter=msvVertices.begin(); vtxIter!=verticesEnd; ++vtxIter){ + if(!msvVertices.empty()){ + for(const auto& vtx : msvVertices){ if(msvVertices.size()>=10) continue; - float mass = xAOD::SecVtxHelper::VertexMass(**vtxIter); - float efrc = xAOD::SecVtxHelper::EnergyFraction(**vtxIter); - int ntrk = xAOD::SecVtxHelper::VtxNtrk(**vtxIter); - float pt = xAOD::SecVtxHelper::Vtxpt(**vtxIter); - float eta = xAOD::SecVtxHelper::Vtxeta(**vtxIter); - float phi = xAOD::SecVtxHelper::Vtxphi(**vtxIter); - float dls = xAOD::SecVtxHelper::VtxnormDist(**vtxIter); - float x = (**vtxIter)->x(); - float y = (**vtxIter)->y(); - float z = (**vtxIter)->z(); + float mass = xAOD::SecVtxHelper::VertexMass(*vtx); + float efrc = xAOD::SecVtxHelper::EnergyFraction(*vtx); + int ntrk = xAOD::SecVtxHelper::VtxNtrk(*vtx); + float pt = xAOD::SecVtxHelper::Vtxpt(*vtx); + float eta = xAOD::SecVtxHelper::Vtxeta(*vtx); + float phi = xAOD::SecVtxHelper::Vtxphi(*vtx); + float dls = xAOD::SecVtxHelper::VtxnormDist(*vtx); + float x = (*vtx)->x(); + float y = (*vtx)->y(); + float z = (*vtx)->z(); TLorentzVector svp4; svp4.SetPtEtaPhiM(pt,eta,phi,mass); //if(jp4.DeltaR(svp4)>0.4) continue; vars.m_summass += mass; - const std::vector<ElementLink<xAOD::TrackParticleContainer> > svTrackLinks = (**vtxIter)->trackParticleLinks(); + const std::vector<ElementLink<xAOD::TrackParticleContainer> > svTrackLinks = (*vtx)->trackParticleLinks(); if(svTrackLinks.size()>1){ nvtx2trk++; } @@ -196,10 +168,11 @@ namespace Analysis int SV1ntrk = 0; std::vector< ElementLink< xAOD::VertexContainer > > SV1Vertice; status &= BTag.variable<std::vector<ElementLink<xAOD::VertexContainer> > >(m_sv1_infosource, "vertices", SV1Vertice); - if (SV1Vertice.size()>0 && SV1Vertice[0].isValid()){ - status &= BTag.taggerInfo(SV1ntrk, xAOD::BTagInfo::SV1_NGTinSvx); - vars.m_diffntrkSV1 = all_trks - SV1ntrk; - }else{ vars.m_diffntrkSV1 = all_trks; + if (!SV1Vertice.empty() && SV1Vertice[0].isValid()){ + status &= BTag.taggerInfo(SV1ntrk, xAOD::BTagInfo::SV1_NGTinSvx); + vars.m_diffntrkSV1 = all_trks - SV1ntrk; + }else{ + vars.m_diffntrkSV1 = all_trks; } vars.m_diffntrkSV0 = vars.m_diffntrkSV1; @@ -261,8 +234,7 @@ namespace Analysis // distances: max mass vertex to PV, and mx2 to max vertex: if(ivm1>=0&&ivm2>=0) { - vars.m_mx12_2d12 = TMath::Sqrt( (v_vtxx[ivm2] - v_vtxx[ivm1]) * (v_vtxx[ivm2] - v_vtxx[ivm1]) - + (v_vtxy[ivm2] - v_vtxy[ivm1]) * (v_vtxy[ivm2] - v_vtxy[ivm1]) ); + vars.m_mx12_2d12 = std::hypot( v_vtxx[ivm2] - v_vtxx[ivm1], v_vtxy[ivm2] - v_vtxy[ivm1] ); vars.m_mx12_DR = sv1p3.DeltaR(sv2p3); vars.m_mx12_Angle = sv1p3.Angle(sv2p3); @@ -313,28 +285,28 @@ namespace Analysis bool &badVariableFound, std::vector<float*> &inputPointers) { - for (unsigned ivar=0; ivar<inputVars.size(); ivar++) { - if (inputVars.at(ivar)=="pt" ) { inputPointers.push_back(&m_jetpt ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="Nvtx" ) { inputPointers.push_back(&m_nvtx ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="MaxEfrc" ) { inputPointers.push_back(&m_maxefrc ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="sumMass" ) { inputPointers.push_back(&m_summass ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="totalntrk" ) { inputPointers.push_back(&m_totalntrk ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="diffntrkSV0" ) { inputPointers.push_back(&m_diffntrkSV0 ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="diffntrkSV1" ) { inputPointers.push_back(&m_diffntrkSV1 ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="normDist" ) { inputPointers.push_back(&m_normDist ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="maxVtxMass" ) { inputPointers.push_back(&m_mmax_mass ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="maxSecVtxMass" ) { inputPointers.push_back(&m_mmx2_mass ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="EfrcmaxVtxMass" ) { inputPointers.push_back(&m_mmax_efrc ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="EfrcmaxSecVtxMass") { inputPointers.push_back(&m_mmx2_efrc ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="dlsmaxVtxMass" ) { inputPointers.push_back(&m_mmax_dist ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="dlsmaxSecVtxMass" ) { inputPointers.push_back(&m_mmx2_dist ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="dRmaxVtxMassj" ) { inputPointers.push_back(&m_mmax_DRjet ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="dRmaxSecVtxMassj" ) { inputPointers.push_back(&m_mmx2_DRjet ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="d2Mass12" ) { inputPointers.push_back(&m_mx12_2d12 ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="DRMass12" ) { inputPointers.push_back(&m_mx12_DR ) ; nConfgVar++; } - else if (inputVars.at(ivar)=="AngleMass12" ) { inputPointers.push_back(&m_mx12_Angle ) ; nConfgVar++; } + for (const auto& var : inputVars) { + if (var=="pt" ) { inputPointers.push_back(&m_jetpt ) ; nConfgVar++; } + else if (var=="Nvtx" ) { inputPointers.push_back(&m_nvtx ) ; nConfgVar++; } + else if (var=="MaxEfrc" ) { inputPointers.push_back(&m_maxefrc ) ; nConfgVar++; } + else if (var=="sumMass" ) { inputPointers.push_back(&m_summass ) ; nConfgVar++; } + else if (var=="totalntrk" ) { inputPointers.push_back(&m_totalntrk ) ; nConfgVar++; } + else if (var=="diffntrkSV0" ) { inputPointers.push_back(&m_diffntrkSV0 ) ; nConfgVar++; } + else if (var=="diffntrkSV1" ) { inputPointers.push_back(&m_diffntrkSV1 ) ; nConfgVar++; } + else if (var=="normDist" ) { inputPointers.push_back(&m_normDist ) ; nConfgVar++; } + else if (var=="maxVtxMass" ) { inputPointers.push_back(&m_mmax_mass ) ; nConfgVar++; } + else if (var=="maxSecVtxMass" ) { inputPointers.push_back(&m_mmx2_mass ) ; nConfgVar++; } + else if (var=="EfrcmaxVtxMass" ) { inputPointers.push_back(&m_mmax_efrc ) ; nConfgVar++; } + else if (var=="EfrcmaxSecVtxMass") { inputPointers.push_back(&m_mmx2_efrc ) ; nConfgVar++; } + else if (var=="dlsmaxVtxMass" ) { inputPointers.push_back(&m_mmax_dist ) ; nConfgVar++; } + else if (var=="dlsmaxSecVtxMass" ) { inputPointers.push_back(&m_mmx2_dist ) ; nConfgVar++; } + else if (var=="dRmaxVtxMassj" ) { inputPointers.push_back(&m_mmax_DRjet ) ; nConfgVar++; } + else if (var=="dRmaxSecVtxMassj" ) { inputPointers.push_back(&m_mmx2_DRjet ) ; nConfgVar++; } + else if (var=="d2Mass12" ) { inputPointers.push_back(&m_mx12_2d12 ) ; nConfgVar++; } + else if (var=="DRMass12" ) { inputPointers.push_back(&m_mx12_DR ) ; nConfgVar++; } + else if (var=="AngleMass12" ) { inputPointers.push_back(&m_mx12_Angle ) ; nConfgVar++; } else { - msg << MSG::WARNING << "#BTAG# \""<<inputVars.at(ivar)<<"\" <- This variable found in xml/calib-file does not match to any variable declared in MultiSV... the algorithm will be 'disabled'." << endmsg; + msg << MSG::WARNING << "#BTAG# \""<<var<<"\" <- This variable found in xml/calib-file does not match to any variable declared in MultiSV... the algorithm will be 'disabled'." << endmsg; badVariableFound=true; } } diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/src/MultivariateTagManager.cxx b/PhysicsAnalysis/JetTagging/JetTagTools/src/MultivariateTagManager.cxx index fb58bfef5784e964cd4d6549955f8bf04cb414d2..3b88ccb258a14fca3bc6f0e7acfe7488f5e88397 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/src/MultivariateTagManager.cxx +++ b/PhysicsAnalysis/JetTagging/JetTagTools/src/MultivariateTagManager.cxx @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ ////////////////////////////////////////////////////////////////////////////// @@ -14,14 +14,9 @@ /// ///////////////////////////////////////////////////////////////////////////// -#include "GaudiKernel/ToolHandle.h" -#include "GaudiKernel/IToolSvc.h" - #include "JetTagTools/MultivariateTagManager.h" #include "JetTagTools/BTagVariables.h" -#include "TObjString.h" - #include <fstream> #include <sstream> #include <algorithm> @@ -61,8 +56,8 @@ namespace Analysis { { if (m_MultivariateTaggerHandleArray.empty()) { - ATH_MSG_ERROR("MVTagToolList is empty"); - return StatusCode::FAILURE; + ATH_MSG_ERROR("MVTagToolList is empty"); + return StatusCode::FAILURE; } if ( m_MultivariateTaggerHandleArray.retrieve().isFailure() ) { @@ -86,12 +81,6 @@ namespace Analysis { return StatusCode::SUCCESS; } - StatusCode MultivariateTagManager::finalize() { // all taken care of in destructor - ATH_MSG_INFO(" #BTAG# Finalization of MultivariateTagManager successfull" ); - return StatusCode::SUCCESS; - } - - // _______________________________________________________________________ // MultivariateTagManager functions @@ -153,16 +142,16 @@ namespace Analysis { float sm_mu_z0 = NAN; float sm_ID_qOverP = NAN; - BTag.variable<float>(m_softmuon_infosource, "mu_pt" , sm_mu_pt ); + BTag.variable<float>(m_softmuon_infosource, "mu_pt" , sm_mu_pt ); if(!std::isnan(sm_mu_pt) && sm_mu_pt>0){ - BTag.variable<float>(m_softmuon_infosource, "dR" , sm_dR ); - BTag.variable<float>(m_softmuon_infosource, "qOverPratio" , sm_qOverPratio ); - BTag.variable<float>(m_softmuon_infosource, "mombalsignif" , sm_mombalsignif ); - BTag.variable<float>(m_softmuon_infosource, "scatneighsignif" , sm_scatneighsignif); - BTag.variable<float>(m_softmuon_infosource, "pTrel" , sm_pTrel ); - BTag.variable<float>(m_softmuon_infosource, "mu_d0" , sm_mu_d0 ); - BTag.variable<float>(m_softmuon_infosource, "mu_z0" , sm_mu_z0 ); - BTag.variable<float>(m_softmuon_infosource, "ID_qOverP" , sm_ID_qOverP ); + BTag.variable<float>(m_softmuon_infosource, "dR" , sm_dR ); + BTag.variable<float>(m_softmuon_infosource, "qOverPratio" , sm_qOverPratio ); + BTag.variable<float>(m_softmuon_infosource, "mombalsignif" , sm_mombalsignif ); + BTag.variable<float>(m_softmuon_infosource, "scatneighsignif" , sm_scatneighsignif); + BTag.variable<float>(m_softmuon_infosource, "pTrel" , sm_pTrel ); + BTag.variable<float>(m_softmuon_infosource, "mu_d0" , sm_mu_d0 ); + BTag.variable<float>(m_softmuon_infosource, "mu_z0" , sm_mu_z0 ); + BTag.variable<float>(m_softmuon_infosource, "ID_qOverP" , sm_ID_qOverP ); }else{ sm_mu_pt= NAN; } @@ -177,21 +166,20 @@ namespace Analysis { inputs[btagvar::SM_MU_Z0] = sm_mu_z0; inputs[btagvar::SM_ID_QOVERP] = sm_ID_qOverP; - } + void MultivariateTagManager::fill_trkSum(var_map& inputs, xAOD::BTagging& BTag) const { float trkSum_ntrk = NAN; float trkSum_sPt = NAN; float trkSum_vPt = NAN; float trkSum_vAbsEta =NAN; - trkSum_ntrk = BTag.isAvailable<unsigned>("trkSum_ntrk") ? BTag.auxdata<unsigned>("trkSum_ntrk") : NAN; trkSum_sPt = BTag.isAvailable<float >("trkSum_SPt" ) ? BTag.auxdata<float >("trkSum_SPt" ) : NAN; if (!std::isnan(trkSum_ntrk)){ trkSum_vPt = BTag.isAvailable<float>("trkSum_VPt" ) ? BTag.auxdata<float>("trkSum_VPt" ) : NAN; - trkSum_vAbsEta= BTag.isAvailable<float>("trkSum_VEta") ? fabs(BTag.auxdata<float>("trkSum_VEta")) : NAN; + trkSum_vAbsEta= BTag.isAvailable<float>("trkSum_VEta") ? std::abs(BTag.auxdata<float>("trkSum_VEta")) : NAN; } inputs[btagvar::TRKSUM_NTRK] = trkSum_ntrk; @@ -199,7 +187,6 @@ namespace Analysis { inputs[btagvar::TRKSUM_VPT] = trkSum_vPt; inputs[btagvar::TRKSUM_ABSETA] = trkSum_vAbsEta; - } void MultivariateTagManager::fill_jetfitter(var_map& inputs, xAOD::BTagging& BTag) const { @@ -218,7 +205,8 @@ namespace Analysis { float jf_sig3d = NAN; // check if we have vertices - int jf_nvtx_tmp(INT_MISSING), jf_nvtx1t_tmp(INT_MISSING); bool jfitter_ok(false); + int jf_nvtx_tmp(INT_MISSING), jf_nvtx1t_tmp(INT_MISSING); + bool jfitter_ok(false); std::vector< ElementLink< xAOD::BTagVertexContainer > > jf_vertices; BTag.variable<std::vector<ElementLink<xAOD::BTagVertexContainer> > >(m_jftNN_infosource, "JFvertices", jf_vertices); if("JetFitter" == m_jftNN_infosource) { @@ -229,32 +217,35 @@ namespace Analysis { BTag.variable<int>(m_jftNN_infosource, "nVTX", jf_nvtx_tmp ); BTag.variable<int>(m_jftNN_infosource, "nSingleTracks", jf_nvtx1t_tmp); } - if(jf_vertices.size()>0 && jf_vertices[0].isValid() && (jf_nvtx_tmp > 0 || jf_nvtx1t_tmp > 0)) jfitter_ok = true; + + if(!jf_vertices.empty() && jf_vertices[0].isValid() && + (jf_nvtx_tmp > 0 || jf_nvtx1t_tmp > 0)) jfitter_ok = true; if(jfitter_ok) { // Get values from the xAOD if("JetFitter" == m_jftNN_infosource) { // check if JetFitter is known by the xAOD? - BTag.taggerInfo(jf_nvtx, xAOD::BTagInfo::JetFitter_nVTX); - BTag.taggerInfo(jf_nvtx1t, xAOD::BTagInfo::JetFitter_nSingleTracks); - BTag.taggerInfo(jf_ntrkAtVx, xAOD::BTagInfo::JetFitter_nTracksAtVtx); - BTag.taggerInfo(jf_n2tv, xAOD::BTagInfo::JetFitter_N2Tpair); - BTag.taggerInfo(jf_efrc, xAOD::BTagInfo::JetFitter_energyFraction); - BTag.taggerInfo(jf_mass, xAOD::BTagInfo::JetFitter_mass); - BTag.taggerInfo(jf_sig3d, xAOD::BTagInfo::JetFitter_significance3d); - BTag.taggerInfo(jf_dphi, xAOD::BTagInfo::JetFitter_deltaphi); - BTag.taggerInfo(jf_deta, xAOD::BTagInfo::JetFitter_deltaeta); + BTag.taggerInfo(jf_nvtx, xAOD::BTagInfo::JetFitter_nVTX); + BTag.taggerInfo(jf_nvtx1t, xAOD::BTagInfo::JetFitter_nSingleTracks); + BTag.taggerInfo(jf_ntrkAtVx, xAOD::BTagInfo::JetFitter_nTracksAtVtx); + BTag.taggerInfo(jf_n2tv, xAOD::BTagInfo::JetFitter_N2Tpair); + BTag.taggerInfo(jf_efrc, xAOD::BTagInfo::JetFitter_energyFraction); + BTag.taggerInfo(jf_mass, xAOD::BTagInfo::JetFitter_mass); + BTag.taggerInfo(jf_sig3d, xAOD::BTagInfo::JetFitter_significance3d); + BTag.taggerInfo(jf_dphi, xAOD::BTagInfo::JetFitter_deltaphi); + BTag.taggerInfo(jf_deta, xAOD::BTagInfo::JetFitter_deltaeta); } else { // get variables explicitely - BTag.variable<int>(m_jftNN_infosource, "nVTX", jf_nvtx); - BTag.variable<int>(m_jftNN_infosource, "nSingleTracks", jf_nvtx1t); - BTag.variable<int>(m_jftNN_infosource, "nTracksAtVtx", jf_ntrkAtVx); - BTag.variable<int>(m_jftNN_infosource, "N2Tpair", jf_n2tv); - BTag.variable<float>(m_jftNN_infosource, "energyFraction", jf_efrc); - BTag.variable<float>(m_jftNN_infosource, "mass", jf_mass); - BTag.variable<float>(m_jftNN_infosource, "significance3d", jf_sig3d); - BTag.variable<float>(m_jftNN_infosource, "deltaphi", jf_dphi); - BTag.variable<float>(m_jftNN_infosource, "deltaeta", jf_deta); + BTag.variable<int>(m_jftNN_infosource, "nVTX", jf_nvtx); + BTag.variable<int>(m_jftNN_infosource, "nSingleTracks", jf_nvtx1t); + BTag.variable<int>(m_jftNN_infosource, "nTracksAtVtx", jf_ntrkAtVx); + BTag.variable<int>(m_jftNN_infosource, "N2Tpair", jf_n2tv); + BTag.variable<float>(m_jftNN_infosource, "energyFraction", jf_efrc); + BTag.variable<float>(m_jftNN_infosource, "mass", jf_mass); + BTag.variable<float>(m_jftNN_infosource, "significance3d", jf_sig3d); + BTag.variable<float>(m_jftNN_infosource, "deltaphi", jf_dphi); + BTag.variable<float>(m_jftNN_infosource, "deltaeta", jf_deta); } + // NOTE: no need to check for NAN here, it should do the right thing // http://en.cppreference.com/w/cpp/numeric/math/hypot#Error_handling jf_dR = std::hypot(jf_dphi,jf_deta); @@ -295,9 +286,7 @@ namespace Analysis { std::vector<float> weightBofTracksIP2D; BTag.variable<std::vector<float> >(m_ip2d_infosource, "weightBofTracks", weightBofTracksIP2D); - int ntrk_ip2 = weightBofTracksIP2D.size(); - - if(ntrk_ip2>0) { + if(!weightBofTracksIP2D.empty()) { if( m_ip2d_infosource == "IP2D" ) { ip2d_pb = BTag.IP2D_pb(); @@ -362,8 +351,7 @@ namespace Analysis { std::vector<float> weightBofTracksIP3D; BTag.variable<std::vector<float> >(m_ip3d_infosource, "weightBofTracks", weightBofTracksIP3D); - int ntrk_ip3= weightBofTracksIP3D.size(); - if(ntrk_ip3>0) { + if(!weightBofTracksIP3D.empty()) { if( m_ip3d_infosource == "IP3D" ) { ip3d_pb = BTag.IP3D_pb(); ip3d_pc = BTag.IP3D_pc(); @@ -429,25 +417,23 @@ namespace Analysis { std::vector< ElementLink< xAOD::VertexContainer > > myVertices_SV0; BTag.variable<std::vector<ElementLink<xAOD::VertexContainer> > >(m_sv1_infosource, "vertices", myVertices_SV0); - if ( myVertices_SV0.size() > 0 && myVertices_SV0[0].isValid() ) { + if ( !myVertices_SV0.empty() && myVertices_SV0[0].isValid() ) { // if we found a vertex, then sv0 is okay to use sv0_ok = true; } if (sv0_ok) { if (m_sv0_infosource == "SV0") { - BTag.taggerInfo(sv0_mass, xAOD::BTagInfo::SV0_masssvx); - BTag.taggerInfo(sv0_efrc, xAOD::BTagInfo::SV0_efracsvx); - BTag.taggerInfo(sv0_n2t, xAOD::BTagInfo::SV0_N2Tpair); - BTag.taggerInfo(sv0_ntrkv, xAOD::BTagInfo::SV0_NGTinSvx); - //BTag.taggerInfo(sv0_sig3d, xAOD::BTagInfo::SV0_normdist); + BTag.taggerInfo(sv0_mass, xAOD::BTagInfo::SV0_masssvx); + BTag.taggerInfo(sv0_efrc, xAOD::BTagInfo::SV0_efracsvx); + BTag.taggerInfo(sv0_n2t, xAOD::BTagInfo::SV0_N2Tpair); + BTag.taggerInfo(sv0_ntrkv, xAOD::BTagInfo::SV0_NGTinSvx); } else { - BTag.variable<float>(m_sv0_infosource, "masssvx", sv0_mass); - BTag.variable<float>(m_sv0_infosource, "efracsvx", sv0_efrc); - BTag.variable<int>(m_sv0_infosource, "N2Tpair", sv0_n2t); - BTag.variable<int>(m_sv0_infosource, "NGTinSvx", sv0_ntrkv); - //BTag.variable<float>(m_sv0_infosource, "significance3D", sv0_sig3d); + BTag.variable<float>(m_sv0_infosource, "masssvx", sv0_mass); + BTag.variable<float>(m_sv0_infosource, "efracsvx", sv0_efrc); + BTag.variable<int>(m_sv0_infosource, "N2Tpair", sv0_n2t); + BTag.variable<int>(m_sv0_infosource, "NGTinSvx", sv0_ntrkv); } BTag.variable<float>(m_sv0_infosource, "significance3D", sv0_sig3d); @@ -486,7 +472,7 @@ namespace Analysis { bool sv1_ok(false); std::vector< ElementLink< xAOD::VertexContainer > > myVertices_SV1; BTag.variable<std::vector<ElementLink<xAOD::VertexContainer> > >(m_sv1_infosource, "vertices", myVertices_SV1); - if ( myVertices_SV1.size() > 0 && myVertices_SV1[0].isValid() ) { + if ( !myVertices_SV1.empty() && myVertices_SV1[0].isValid() ) { // if we found a vertex, then sv1 is okay to use sv1_ok = true; } @@ -561,7 +547,7 @@ namespace Analysis { inputs[btagvar::SV1_L3D] = sv1_L3d; inputs[btagvar::SV1_SIG3D] = sv1_sig3d; inputs[btagvar::SV1_DR] = sv1_dR; - inputs[btagvar::SV1_DISTMATLAY] = sv1_distmatlay; + inputs[btagvar::SV1_DISTMATLAY] = sv1_distmatlay; } void MultivariateTagManager::fill_arbitrary_aux_data( diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/src/NewLikelihoodTool.cxx b/PhysicsAnalysis/JetTagging/JetTagTools/src/NewLikelihoodTool.cxx index 5baf6d9e333e6a78a0015cd38b27094967fa065b..0bbcf09695e80bfffae019ee27c1124ca832413b 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/src/NewLikelihoodTool.cxx +++ b/PhysicsAnalysis/JetTagging/JetTagTools/src/NewLikelihoodTool.cxx @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ #include "JetTagTools/NewLikelihoodTool.h" @@ -42,16 +42,8 @@ namespace Analysis { m_vetoSmoothingOf.push_back("/Sip3D"); } - NewLikelihoodTool::~NewLikelihoodTool() { - } - StatusCode NewLikelihoodTool::initialize() { ATH_CHECK(m_readKey.initialize()); - - return StatusCode::SUCCESS; - } - - StatusCode NewLikelihoodTool::finalize() { return StatusCode::SUCCESS; } @@ -61,10 +53,10 @@ namespace Analysis { void NewLikelihoodTool::printStatus() const { msg(MSG::INFO) << "#BTAG# - hypotheses : "; - for(unsigned int ih=0;ih<m_hypotheses.size();ih++) msg(MSG::INFO) << m_hypotheses[ih] << ", "; + for(const auto& hypo : m_hypotheses) msg(MSG::INFO) << hypo << ", "; msg(MSG::INFO) << endmsg; msg(MSG::INFO) << "#BTAG# - histograms : " << endmsg; - for(unsigned int ih=0;ih<m_histograms.size();ih++) msg(MSG::INFO) << m_histograms[ih] << endmsg; + for(const auto& histo : m_histograms) msg(MSG::INFO) << histo << endmsg; } std::vector<std::string> NewLikelihoodTool::gradeList(const std::string& histoName) const { @@ -169,8 +161,8 @@ namespace Analysis { if(norm) { // check if smoothing of histogram is not vetoed: bool veto = false; - for(unsigned int iv=0; iv<m_vetoSmoothingOf.size(); iv++) { - if(hname.find(m_vetoSmoothingOf[iv])!=std::string::npos) { + for(const auto& v : m_vetoSmoothingOf) { + if(hname.find(v)!=std::string::npos) { veto = true; ATH_MSG_VERBOSE("#BTAG# Smoothing of " << hname << " is vetoed !"); break; @@ -208,10 +200,12 @@ namespace Analysis { for(int iz=1; iz<=Nz; iz++){ double content=dc_tmp->Integral(1,Nx,1,Ny,iz,iz,""); if(content==0.)content=Nz; double dnorm=total/content/Nz; - for(int ix=1; ix<=Nx; ix++){for(int iy=1; iy<=Ny; iy++){ - double cbin=dc_tmp->GetBinContent(ix,iy,iz)*dnorm; cbin= cbin>0. ? cbin : 0.1; //Protection against empty bins - dc_tmp->SetBinContent(ix,iy,iz, cbin); - }} + for(int ix=1; ix<=Nx; ix++){ + for(int iy=1; iy<=Ny; iy++){ + double cbin=dc_tmp->GetBinContent(ix,iy,iz)*dnorm; cbin= cbin>0. ? cbin : 0.1; //Protection against empty bins + dc_tmp->SetBinContent(ix,iy,iz, cbin); + } + } } HistoHelperRoot::smoothASH3D(dc_tmp, m3d1, m3d1, m3d3, msgLvl(MSG::DEBUG)); } @@ -237,32 +231,30 @@ namespace Analysis { } ATH_MSG_VERBOSE("#BTAG# -- lhVarVal size= " << lhVariableValues.size()); // loop on Tracks in the Jet (IP) / Vertices in the Jet (SV) - for (unsigned int iel = 0; iel<lhVariableValues.size(); iel++) { - ATH_MSG_VERBOSE( "#BTAG# -- element " << iel << " " - << lhVariableValues[iel].name ); - int ncompo = lhVariableValues[iel].composites.size(); - ATH_MSG_VERBOSE( "#BTAG# -- element " << iel << " " - << lhVariableValues[iel].name + for (const auto& value : lhVariableValues) { + ATH_MSG_VERBOSE( "#BTAG# -- element " << value.name ); + int ncompo = value.composites.size(); + ATH_MSG_VERBOSE( "#BTAG# -- element " << value.name << " has " << ncompo << " composites." ); // loop on variables that make up the Tag, e.g. // one 1D for IP2D, one 2D for IP3D, one 1D and one 2D for SV1, one 3D for SV2 - for (int icompo = 0;icompo<ncompo;icompo++) { + for (const auto& compo : value.composites) { double sum(0.); std::vector<double> tmpVector; - std::string histName = lhVariableValues[iel].composites[icompo].name; - int idim = lhVariableValues[iel].composites[icompo].atoms.size(); - ATH_MSG_VERBOSE( "#BTAG# -- composite " << icompo << " histo= " + std::string histName = compo.name; + int idim = compo.atoms.size(); + ATH_MSG_VERBOSE( "#BTAG# -- composite histo= " << histName << " dim= " << idim ); - for (unsigned int ihyp = 0 ; ihyp < m_hypotheses.size(); ++ihyp) { - TH1* tmpHisto = this->prepareHistogram(m_hypotheses[ihyp],histName); + for (const auto& hypo : m_hypotheses) { + TH1* tmpHisto = this->prepareHistogram(hypo,histName); if(tmpHisto) { if(1==idim) { - double valuex = lhVariableValues[iel].composites[icompo].atoms[0].value; + double valuex = compo.atoms[0].value; int binx = (tmpHisto->GetXaxis())->FindBin(valuex); if(valuex >= tmpHisto->GetXaxis()->GetXmax()) binx = tmpHisto->GetXaxis()->GetNbins(); if(valuex <= tmpHisto->GetXaxis()->GetXmin()) binx = 1; double tmp = tmpHisto->GetBinContent(binx); - if( msgLvl(MSG::VERBOSE) ) msg(MSG::VERBOSE) << "#BTAG# For hypothesis= " << m_hypotheses[ihyp] + if( msgLvl(MSG::VERBOSE) ) msg(MSG::VERBOSE) << "#BTAG# For hypothesis= " << hypo << " (1D) actual value= " << valuex << " --> bin= " << binx << " f = " << tmp; if(m_interpolate) { @@ -288,8 +280,8 @@ namespace Analysis { sum += tmp; } if(2==idim) { - double valuex = lhVariableValues[iel].composites[icompo].atoms[0].value; - double valuey = lhVariableValues[iel].composites[icompo].atoms[1].value; + double valuex = compo.atoms[0].value; + double valuey = compo.atoms[1].value; int binx = (tmpHisto->GetXaxis())->FindBin(valuex); int biny = (tmpHisto->GetYaxis())->FindBin(valuey); if(valuex >= tmpHisto->GetXaxis()->GetXmax()) binx = tmpHisto->GetXaxis()->GetNbins(); @@ -297,7 +289,7 @@ namespace Analysis { if(valuey >= tmpHisto->GetYaxis()->GetXmax()) biny = tmpHisto->GetYaxis()->GetNbins(); if(valuey <= tmpHisto->GetYaxis()->GetXmin()) biny = 1; double tmp = tmpHisto->GetBinContent(binx, biny); - if( msgLvl(MSG::VERBOSE) ) msg(MSG::VERBOSE) << "#BTAG# For hypothesis= " << m_hypotheses[ihyp] + if( msgLvl(MSG::VERBOSE) ) msg(MSG::VERBOSE) << "#BTAG# For hypothesis= " << hypo << " (2D) actual value= " << valuex << " " << valuey << " --> bin= " << binx << " " << biny << " f = " << tmp; if(m_interpolate) { @@ -324,9 +316,9 @@ namespace Analysis { sum += tmp; } if(3==idim) { - double valuex = lhVariableValues[iel].composites[icompo].atoms[0].value; - double valuey = lhVariableValues[iel].composites[icompo].atoms[1].value; - double valuez = lhVariableValues[iel].composites[icompo].atoms[2].value; + double valuex = compo.atoms[0].value; + double valuey = compo.atoms[1].value; + double valuez = compo.atoms[2].value; int binx = (tmpHisto->GetXaxis())->FindBin(valuex); int biny = (tmpHisto->GetYaxis())->FindBin(valuey); int binz = (tmpHisto->GetZaxis())->FindBin(valuez); @@ -337,7 +329,7 @@ namespace Analysis { if(valuez >= tmpHisto->GetZaxis()->GetXmax()) binz = tmpHisto->GetZaxis()->GetNbins(); if(valuez <= tmpHisto->GetZaxis()->GetXmin()) binz = 1; double tmp = tmpHisto->GetBinContent(binx, biny, binz); - if( msgLvl(MSG::VERBOSE) ) msg(MSG::VERBOSE) << "#BTAG# For hypothesis= " << m_hypotheses[ihyp] + if( msgLvl(MSG::VERBOSE) ) msg(MSG::VERBOSE) << "#BTAG# For hypothesis= " << hypo << " (3D) actual value= " << valuex << " " << valuey << " " << valuez << " --> bin= " << binx << " " << biny @@ -372,11 +364,10 @@ namespace Analysis { } } // endloop on hypotheses (B,U,C..) unsigned int classCount(0); - for( std::vector<double>::iterator itr3 = tmpVector.begin(); - itr3 != tmpVector.end(); ++itr3 ) { + for (const auto& f : tmpVector) { if(sum != 0.) { if( msgLvl(MSG::DEBUG) ) msg(MSG::DEBUG) << "#BTAG# sum of pX = " << sum << endmsg; - double p = (*itr3); + double p = f; if(m_normalizedProb) p /= sum; probDensityPerEventClassAllVariables[classCount] *= p; } else { diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/src/SVForIPTool.cxx b/PhysicsAnalysis/JetTagging/JetTagTools/src/SVForIPTool.cxx index 3f456bd2cd7ab61e413a978f97273feca9283c78..8cf2a2521ae2593fb8c9ce90cab4e59be8f488d8 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/src/SVForIPTool.cxx +++ b/PhysicsAnalysis/JetTagging/JetTagTools/src/SVForIPTool.cxx @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ #include "JetTagTools/SVForIPTool.h" @@ -16,20 +16,6 @@ namespace Analysis { declareInterface<SVForIPTool>(this); } - SVForIPTool::~SVForIPTool() {} - - StatusCode SVForIPTool::initialize() { - - ATH_MSG_INFO(" Initialization of SVForIPTool succesful"); - return StatusCode::SUCCESS; - } - - StatusCode SVForIPTool::finalize() { - - ATH_MSG_INFO(" Finalization of SVForIPTool succesful"); - return StatusCode::SUCCESS; - } - void SVForIPTool::getDirectionFromSecondaryVertexInfo(Amg::Vector3D & SvxDirection, bool & canUseSvxDirection, xAOD::BTagging* BTag, @@ -39,40 +25,30 @@ namespace Analysis { std::vector< ElementLink< xAOD::VertexContainer > > myVertices; BTag->variable<std::vector<ElementLink<xAOD::VertexContainer> > >(secVxFinderName, "vertices", myVertices); - if (myVertices.size() == 0) { + if (myVertices.empty()) { ATH_MSG_DEBUG(" No secondary vertex found for getting the B flight direction (for the IP sign calculation)"); } else { - if (myVertices[0].isValid()) - { - canUseSvxDirection=true; - SvxDirection=(*myVertices[0])->position()-priVtx.position(); - ATH_MSG_VERBOSE(" Get direction from InDetVKalVertex: phi: " << SvxDirection.phi() << - " theta: " << SvxDirection.theta() ); - } - else - { - ATH_MSG_WARNING("SVX info seems usable, but no SVX available !!!"); - } - } + if (myVertices[0].isValid()) { + canUseSvxDirection=true; + SvxDirection=(*myVertices[0])->position()-priVtx.position(); + ATH_MSG_VERBOSE(" Get direction from InDetVKalVertex: phi: " << SvxDirection.phi() << + " theta: " << SvxDirection.theta() ); + } else { + ATH_MSG_WARNING("SVX info seems usable, but no SVX available !!!"); + } + } } - - - + + void SVForIPTool::getTrkFromV0FromSecondaryVertexInfo(std::vector<const xAOD::TrackParticle*> & TrkFromV0, xAOD::BTagging* BTag, const std::string & secVxFinderName) const { std::vector<ElementLink<xAOD::TrackParticleContainer> > TrkFromV0_ELs; - BTag->variable<std::vector<ElementLink<xAOD::TrackParticleContainer> > >(secVxFinderName, "badTracksIP", TrkFromV0_ELs); - - for (std::vector<ElementLink<xAOD::TrackParticleContainer> >::iterator itr=TrkFromV0_ELs.begin();itr!=TrkFromV0_ELs.end();++itr) - { - if (itr->isValid()) - { - TrkFromV0.push_back((**itr)); - } - } + for (const auto& link : TrkFromV0_ELs) { + if (link.isValid()) TrkFromV0.push_back(*link); + } } }//end namespace diff --git a/PhysicsAnalysis/JetTagging/JetTagTools/src/SVTag.cxx b/PhysicsAnalysis/JetTagging/JetTagTools/src/SVTag.cxx index cc7cbde3c31ae404351276cf48e94cf7b2c52640..de328422a46279b29539b79312f097a9f90ae734 100644 --- a/PhysicsAnalysis/JetTagging/JetTagTools/src/SVTag.cxx +++ b/PhysicsAnalysis/JetTagging/JetTagTools/src/SVTag.cxx @@ -1,5 +1,5 @@ /* - Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration + Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration */ /*************************************************************************** @@ -7,16 +7,12 @@ ***************************************************************************/ #include "JetTagTools/SVTag.h" -#include "GaudiKernel/IToolSvc.h" -#include "Navigation/NavigationToken.h" #include "JetTagTools/NewLikelihoodTool.h" #include "GaudiKernel/ITHistSvc.h" #include "JetTagTools/HistoHelperRoot.h" #include "JetTagTools/LikelihoodComponents.h" #include "ParticleJetTools/JetFlavourInfo.h" -#include "xAODTracking/TrackParticle.h" -#include "xAODTracking/Vertex.h" #include "xAODTracking/VertexContainer.h" #include <string> @@ -39,9 +35,7 @@ namespace Analysis declareProperty("LikelihoodTool", m_likelihoodTool); declareProperty("checkOverflows", m_checkOverflows = false); declareProperty("purificationDeltaR", m_purificationDeltaR = 0.8); - declareProperty("UseBinInterpol", m_UseBinInterpol = false); - // declareProperty("originalTPCollectionName", m_originalTPCollectionName = "InDetTrackParticles"); declareProperty("jetCollectionList", m_jetCollectionList); declareProperty("useForcedCalibration", m_doForcedCalib = false); declareProperty("ForcedCalibrationName", m_ForcedCalibName = "Cone4H1Tower"); @@ -102,19 +96,18 @@ namespace Analysis ATH_MSG_VERBOSE("#BTAG# Booking histos..."); std::vector<double> xb; double xbi[10] = {1., 2., 3., 4., 5., 6., 8., 10., 20., 50.}; - for(uint ijc=0;ijc<m_jetCollectionList.size();ijc++) { - + for(const auto& jetColl : m_jetCollectionList) { for(uint ih=0;ih<m_hypotheses.size();ih++) { // SV1 if (m_SVmode == "SV1") { - std::string hDir = "/RefFile/SV1/"+m_jetCollectionList[ijc]+"/"+m_hypotheses[ih]+"/"; + std::string hDir = "/RefFile/SV1/"+jetColl+"/"+m_hypotheses[ih]+"/"; m_histoHelper->bookHisto(hDir+"N2T", "Number of Good Two Track Vertices",9,xbi); m_histoHelper->bookHisto(hDir+"N2TEffSV1", "Number of Good Two Track Vertices",9,xbi); m_histoHelper->bookHisto(hDir+"N2TNormSV1", "Number of Good Two Track Vertices",30,0.,30.); m_histoHelper->bookHisto(hDir+"BidimME", "(E fraction)**0.7 vs Mass/(Mass+1)" ,50,0.218261406,1.,50,0.,1.); m_histoHelper->bookHisto(hDir+"DRJPVSV", "DeltaR between jet axis and (PV,SV) axis",100,0.,0.5); } else if (m_SVmode == "SV2") { - std::string hDir = "/RefFile/SV2/"+m_jetCollectionList[ijc]+"/"+m_hypotheses[ih]+"/"; + std::string hDir = "/RefFile/SV2/"+jetColl+"/"+m_hypotheses[ih]+"/"; // SV2 m_histoHelper->bookHisto(hDir+"N2TEffSV2", "Number of Good Two Track Vertices",9,xbi); m_histoHelper->bookHisto(hDir+"N2TNormSV2", "Number of Good Two Track Vertices",30,0.,30.); @@ -122,7 +115,7 @@ namespace Analysis else m_histoHelper->bookHisto(hDir+"TridimMEN2T"," ln(N) vs (E fraction)**0.7 vs Mass/(Mass+1)",20,0.,1.,20,0.,1.,7,0.,3.8); if(ih==0) { // Control with SV2 - hDir = "/RefFile/SV2/"+m_jetCollectionList[ijc]+"/controlSV/"; + hDir = "/RefFile/SV2/"+jetColl+"/controlSV/"; m_histoHelper->bookHisto(hDir+"eta","eta",60,-3.,3.); m_histoHelper->bookHisto(hDir+"phi","phi",64,-3.2,3.2); m_histoHelper->bookHisto(hDir+"pt","pt",50,0.,300.); @@ -175,13 +168,14 @@ namespace Analysis ATH_MSG_VERBOSE("#BTAG# Using jet type " << author << " for calibrations."); /* The jet */ + Amg::Vector3D jetDir(jetToTag.p4().Px(),jetToTag.p4().Py(),jetToTag.p4().Pz()); double jeteta = jetToTag.eta(), jetphi = jetToTag.phi(), jetpt = jetToTag.pt(); ATH_MSG_VERBOSE("#BTAG# Jet properties : eta = " << jeteta << " phi = " << jetphi << " pT = " <<jetpt/m_c_mom); // Fill control histograms if (m_runModus=="reference" && m_SVmode == "SV2") { - if (fabs(jeteta) <= 2.5) { + if (std::abs(jeteta) <= 2.5) { m_histoHelper->fillHisto("/RefFile/SV2/"+author+"/controlSV/eta",(double)jeteta); m_histoHelper->fillHisto("/RefFile/SV2/"+author+"/controlSV/phi",(double)jetphi); m_histoHelper->fillHisto("/RefFile/SV2/"+author+"/controlSV/pt",(double)jetpt/m_c_mom); @@ -190,7 +184,7 @@ namespace Analysis // // Get the SV info // - float ambtot = -1., xratio = -1., distnrm = 0., drJPVSV = 0., Lxy = -100., L3d = -100.; + float ambtot = -1., xratio = -1., distnrm = 0., drJPVSV = 0., Lxy = -100., L3d = -100.; int NSVPair = -1; float distnrmCorr=0.; @@ -199,9 +193,8 @@ namespace Analysis std::vector< ElementLink< xAOD::VertexContainer > > myVertices; // don't check the following status BTag.variable<std::vector<ElementLink<xAOD::VertexContainer> > >(m_secVxFinderName, "vertices", myVertices); - // BTag.auxdata<std::vector<ElementLink<xAOD::VertexContainer> > >(m_secVxFinderName+"_vertices"); - if (myVertices.size()>0) { + if (!myVertices.empty()) { status &= BTag.variable<float>(m_secVxFinderName, "masssvx", ambtot);// mass in MeV ambtot/=m_c_mom; @@ -215,11 +208,9 @@ namespace Analysis // DR between Jet axis and PV-SV axis // For the time being computed only for Single Vertex... if (myVertices[0].isValid()) { - // if (myVertices[0]!=0) { const xAOD::Vertex* firstVertex = *(myVertices[0]); //FIXME ugly hack to get a Amg::Vector3D out of a CLHEP::HepLorentzVector - Amg::Vector3D jetDir(jetToTag.p4().Px(),jetToTag.p4().Py(),jetToTag.p4().Pz()); const Amg::Vector3D PVposition = priVtx.position(); const Amg::Vector3D position = firstVertex->position(); Amg::Vector3D PvSvDir( position.x() - PVposition.x(), @@ -232,30 +223,26 @@ namespace Analysis if ( m_isFlipped ) drJPVSV = drJPVSV_2; // for negative tags ATH_MSG_VERBOSE("#BTAG# DRJPVSV regular="<<drJPVSV_1<<" flipped="<<drJPVSV_2<<" chosen="<<drJPVSV); - Lxy=sqrt(pow(PvSvDir(0,0),2)+pow(PvSvDir(1,0),2)); - L3d=sqrt(pow(PvSvDir(0,0),2)+pow(PvSvDir(1,0),2)+pow(PvSvDir(2,0),2)); + Lxy = std::hypot(PvSvDir(0,0), PvSvDir(1,0)); + L3d = std::hypot(PvSvDir(0,0), PvSvDir(1,0), PvSvDir(2,0)); }else{ ATH_MSG_VERBOSE("#BTAG# No secondary vertex."); } std::vector<const xAOD::Vertex*> vecVertices; - const std::vector<ElementLink<xAOD::VertexContainer> >::const_iterator verticesEnd = myVertices.end(); - for (std::vector<ElementLink<xAOD::VertexContainer> >::const_iterator verticesIter = myVertices.begin(); - verticesIter != verticesEnd; ++verticesIter) { - if (!verticesIter->isValid()) { - ATH_MSG_WARNING("#BTAG# Secondary vertex from InDetVKalVxInJetFinder has zero pointer. Skipping... "); + for (const auto& link : myVertices) { + if (!link.isValid()) { + ATH_MSG_WARNING("#BTAG# Secondary vertex from InDetVKalVxInJetFinder has zero pointer. Skipping... "); continue; - } - vecVertices.push_back(**verticesIter); + } + vecVertices.push_back(*link); } if (myVertices[0].isValid()) { const xAOD::Vertex* myVert = *myVertices[0]; - distnrm=get3DSignificance(priVtx, vecVertices, - Amg::Vector3D(jetToTag.p4().Px(),jetToTag.p4().Py(),jetToTag.p4().Pz())); + distnrm = get3DSignificance(priVtx, vecVertices, jetDir); ATH_MSG_VERBOSE("#BTAG# SVX x = " << myVert->position().x() << " y = " << myVert->position().y() << " z = " << myVert->position().z()); - distnrmCorr=get3DSignificanceCorr(priVtx, vecVertices, - Amg::Vector3D(jetToTag.p4().Px(),jetToTag.p4().Py(),jetToTag.p4().Pz())); + distnrmCorr = get3DSignificanceCorr(priVtx, vecVertices, jetDir); } else { ATH_MSG_VERBOSE("#BTAG# No vertex. Cannot calculate normalized distance."); distnrm=0.; @@ -281,7 +268,8 @@ namespace Analysis instanceName = prefix; instanceName += posfix; } -//AA: Move to filling xAOD::BTagging + + //AA: Move to filling xAOD::BTagging // note that this block is filling things other than likelihood, // it should not be conditional on m_save_probabilities if (m_runModus=="analysis") { @@ -319,12 +307,14 @@ namespace Analysis ATH_MSG_VERBOSE("#BTAG# SV mode = " << m_SVmode); if (m_SVmode != "SV0" ) { - float ambtotp = ambtot > 0. ? ambtot/(1.+ambtot): 0.; + float ambtotp = ambtot > 0. ? ambtot/(1.+ambtot) : 0.; float xratiop = xratio > 0. ? (float)pow(xratio,m_expos) : 0.; - float trfJetPt=log(jetToTag.pt()/20000.); if(trfJetPt<0.)trfJetPt=0.01; if(trfJetPt>4.8)trfJetPt=4.79; + float trfJetPt = log(jetToTag.pt()/20000.); + if(trfJetPt<0.) trfJetPt=0.01; + if(trfJetPt>4.8) trfJetPt=4.79; std::string pref = ""; if (m_runModus=="reference") { - if (jetpt >= m_pTjetmin && fabs(jeteta) <= 2.5) { + if (jetpt >= m_pTjetmin && std::abs(jeteta) <= 2.5) { int label = xAOD::jetFlavourLabel(&jetToTag); double deltaRtoClosestB = 999.;//, deltaRtoClosestC = 999.; if (jetToTag.getAttribute("TruthLabelDeltaR_B",deltaRtoClosestB)) { @@ -350,10 +340,9 @@ namespace Analysis m_ncjet++; } } + if (pref == "B" || pref == "C" || pref == "U") { std::string hDir = "/RefFile/"+m_SVmode+"/"+author+"/"+pref+"/"; - //std::cout<<"SVTAG tag histohelper: " << m_histoHelper << std::endl; - //m_histoHelper->print(); if (m_SVmode == "SV1") m_histoHelper->fillHisto(hDir+"N2TNormSV1",(float)NSVPair); if (m_SVmode == "SV2") m_histoHelper->fillHisto(hDir+"N2TNormSV2",(float)NSVPair); if (NSVPair > 0 && ambtot > 0.) { @@ -366,8 +355,8 @@ namespace Analysis } if (m_SVmode == "SV2") { m_histoHelper->fillHisto(hDir+"N2TEffSV2",(float)NSVPair); - if(m_usePtSV2)m_histoHelper->fillHisto(hDir+"TridimMENPt",ambtotp,xratiop,trfJetPt); - else m_histoHelper->fillHisto(hDir+"TridimMEN2T",ambtotp,xratiop,log((float)NSVPair)); + if(m_usePtSV2) m_histoHelper->fillHisto(hDir+"TridimMENPt",ambtotp,xratiop,trfJetPt); + else m_histoHelper->fillHisto(hDir+"TridimMEN2T",ambtotp,xratiop,log((float)NSVPair)); } } } @@ -474,9 +463,6 @@ namespace Analysis return StatusCode::SUCCESS; } - void SVTag::finalizeHistos() { - } - void SVTag::printParameterSettings() { ATH_MSG_INFO("#BTAG# " << name() << "Parameter settings "); ATH_MSG_INFO("#BTAG# I am in " << m_runModus << " modus."); @@ -487,65 +473,34 @@ namespace Analysis double SVTag::get3DSignificance(const xAOD::Vertex& priVertex, std::vector<const xAOD::Vertex*>& secVertex, const Amg::Vector3D jetDirection) const { - // double SVTag::get3DSignificance(const Trk::RecVertex & priVertex, - // std::vector<const Trk::RecVertex* > & secVertex, - // const Amg::Vector3D jetDirection) { std::vector<Amg::Vector3D> positions; std::vector<AmgSymMatrix(3)> weightMatrices; - - std::vector<const xAOD::Vertex*>::const_iterator secEnd = secVertex.end(); - - for (std::vector<const xAOD::Vertex*>::const_iterator secIter = secVertex.begin(); - secIter != secEnd; ++secIter) - { - // Amg::Vector3D position = (*secIter)->position(); - // Amg::Vector3D position; - // position[0]=(*secIter)->position().x(); - // position[1]=(*secIter)->position().y(); - // position[2]=(*secIter)->position().z(); - // positions.push_back(position); - positions.push_back((*secIter)->position()); - - weightMatrices.push_back((*secIter)->covariancePosition().inverse()); - } - // If multiple secondary vertices were reconstructed, then a common (weighted) position will be used // in the signed decay length significance calculation - Amg::Vector3D weightTimesPosition(0.,0.,0.); AmgSymMatrix(3) sumWeights; sumWeights.setZero(); - int count=0; - for (std::vector<const xAOD::Vertex*>::const_iterator secIter = secVertex.begin(); - secIter != secEnd; ++secIter) { - - weightTimesPosition+=(weightMatrices[count])*positions[count]; - sumWeights+=(weightMatrices[count]); - ++count; + for (const auto& vtx : secVertex) { + positions.push_back(vtx->position()); + weightMatrices.push_back(vtx->covariancePosition().inverse()); + weightTimesPosition += weightMatrices.back()*positions.back(); + sumWeights += weightMatrices.back(); } // now we have the sum of the weights, let's invert this matrix to get the mean covariance matrix - -// int failed(0); - // return value of m.inverse() on a non invertible matrix is undefined -// AmgSymMatrix(3) meanCovariance=sumWeights.inverse(); -// if (failed!=0) { -// ATH_MSG_ERROR("#BTAG# Could not invert sum of sec vtx matrices"); -// return 0.; -// } bool invertible; AmgSymMatrix(3) meanCovariance; meanCovariance.setZero(); sumWeights.computeInverseWithCheck(meanCovariance, invertible); - if (! invertible) { + if (!invertible) { ATH_MSG_WARNING("#BTAG# Could not invert sum of sec vtx matrices"); return 0.; } // calculate the weighted mean secondary vertex position - Amg::Vector3D meanPosition=meanCovariance*weightTimesPosition; + Amg::Vector3D meanPosition = meanCovariance*weightTimesPosition; // add the mean covariance matrix of the secondary vertices to that of the primary vertex // this is the covariance matrix for the decay length @@ -560,19 +515,19 @@ namespace Analysis double Lz = meanPosition[2]-priVertex.position().z(); const double decaylength = sqrt(Lx*Lx + Ly*Ly + Lz*Lz); - if(decaylength==0.)return 0.; //Safety + if(decaylength==0.) return 0.; //Safety const double inv_decaylength = 1. / decaylength; double dLdLx = Lx * inv_decaylength; double dLdLy = Ly * inv_decaylength; double dLdLz = Lz * inv_decaylength; - double decaylength_err2 = dLdLx*dLdLx*covariance(0,0) + - dLdLy*dLdLy*covariance(1,1) + - dLdLz*dLdLz*covariance(2,2) + - 2.*dLdLx*dLdLy*covariance(0,1) + - 2.*dLdLx*dLdLz*covariance(0,2) + - 2.*dLdLy*dLdLz*covariance(1,2); - if(decaylength_err2<=0.)return 0.; //Something is wrong + double decaylength_err2 = (dLdLx*dLdLx*covariance(0,0) + + dLdLy*dLdLy*covariance(1,1) + + dLdLz*dLdLz*covariance(2,2) + + 2.*dLdLx*dLdLy*covariance(0,1) + + 2.*dLdLx*dLdLz*covariance(0,2) + + 2.*dLdLy*dLdLz*covariance(1,2)); + if(decaylength_err2<=0.) return 0.; //Something is wrong double decaylength_err = sqrt(decaylength_err2); double decaylength_significance = 0.; @@ -593,26 +548,23 @@ namespace Analysis bool success=true; AmgSymMatrix(3) Wgt; - for (auto & svrt : secVertex) + for (const auto & svrt : secVertex) { Amg::Vector3D SVmPV = svrt->position()-priVertex.position(); - AmgSymMatrix(3) SVmPVCov=svrt->covariancePosition()+priVertex.covariancePosition(); + AmgSymMatrix(3) SVmPVCov = svrt->covariancePosition()+priVertex.covariancePosition(); SVmPVCov.computeInverseWithCheck(Wgt, success); if( !success || Wgt(0,0)<=0. || Wgt(1,1)<=0. || Wgt(2,2)<=0. )continue; //Inversion failure - double significance=SVmPV.transpose()*Wgt*SVmPV; + double significance = SVmPV.transpose()*Wgt*SVmPV; if(significance <= 0.) continue; //Something is still wrong! - significance=std::sqrt(significance); - if(SVmPV.dot(jetDirection)<0.)significance *= -1.; + significance = std::sqrt(significance); + if(SVmPV.dot(jetDirection)<0.) significance *= -1.; Sig3D.push_back(significance); } - if(Sig3D.size()==0)return 0.; + if(Sig3D.size()==0) return 0.; return *std::max_element(Sig3D.begin(),Sig3D.end()); - //return std::accumulate(Sig3D.begin(),Sig3D.end(),0.); } - - }