File indexing completed on 2025-01-31 02:19:52
0001 #ifndef ImpactParameter_TemplatedTrackCountingComputer_h
0002 #define ImpactParameter_TemplatedTrackCountingComputer_h
0003
0004 #include "DataFormats/TrackReco/interface/Track.h"
0005 #include "DataFormats/VertexReco/interface/Vertex.h"
0006 #include "DataFormats/BTauReco/interface/TrackCountingTagInfo.h"
0007 #include "DataFormats/BTauReco/interface/IPTagInfo.h"
0008 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0009 #include "Math/GenVector/VectorUtil.h"
0010 #include "RecoBTau/JetTagComputer/interface/JetTagComputer.h"
0011
0012 template <class Container, class Base>
0013 class TemplatedTrackCountingComputer : public JetTagComputer {
0014 public:
0015 using Tokens = void;
0016
0017 typedef reco::IPTagInfo<Container, Base> TagInfo;
0018
0019 TemplatedTrackCountingComputer(const edm::ParameterSet& parameters) {
0020 m_minIP = parameters.getParameter<double>("minimumImpactParameter");
0021 m_useSignedIPSig = parameters.getParameter<bool>("useSignedImpactParameterSig");
0022 m_nthTrack = parameters.getParameter<int>("nthTrack");
0023 m_ipType = parameters.getParameter<int>("impactParameterType");
0024 m_deltaR = parameters.getParameter<double>("deltaR");
0025 m_cutMaxDecayLen = parameters.getParameter<double>("maximumDecayLength");
0026 m_cutMaxDistToAxis = parameters.getParameter<double>("maximumDistanceToJetAxis");
0027
0028
0029
0030 std::string trackQualityType = parameters.getParameter<std::string>("trackQualityClass");
0031 m_trackQuality = reco::TrackBase::qualityByName(trackQualityType);
0032 m_useAllQualities = false;
0033 if (trackQualityType == "any" || trackQualityType == "Any" || trackQualityType == "ANY")
0034 m_useAllQualities = true;
0035
0036 uses("ipTagInfos");
0037
0038 useVariableJTA_ = parameters.getParameter<bool>("useVariableJTA");
0039 if (useVariableJTA_) {
0040 varJTApars = {parameters.getParameter<double>("a_dR"),
0041 parameters.getParameter<double>("b_dR"),
0042 parameters.getParameter<double>("a_pT"),
0043 parameters.getParameter<double>("b_pT"),
0044 parameters.getParameter<double>("min_pT"),
0045 parameters.getParameter<double>("max_pT"),
0046 parameters.getParameter<double>("min_pT_dRcut"),
0047 parameters.getParameter<double>("max_pT_dRcut"),
0048 parameters.getParameter<double>("max_pT_trackPTcut")};
0049 }
0050 }
0051
0052 float discriminator(const TagInfoHelper& ti) const override {
0053 const TagInfo& tkip = ti.get<TagInfo>();
0054 std::multiset<float> significances = orderedSignificances(tkip);
0055 std::multiset<float>::reverse_iterator nth = significances.rbegin();
0056 for (int i = 0; i < m_nthTrack - 1 && nth != significances.rend(); i++)
0057 nth++;
0058 if (nth != significances.rend())
0059 return *nth;
0060 else
0061 return -100.;
0062 }
0063
0064 static void fillPSetDescription(edm::ParameterSetDescription& desc) {
0065 desc.add<double>("minimumImpactParameter", -1.);
0066 desc.add<bool>("useSignedImpactParameterSig", true);
0067 desc.add<int>("nthTrack", -1);
0068 desc.add<int>("impactParameterType", 0)->setComment("0 = 3D, 1 = 2D");
0069 desc.add<double>("deltaR", -1.0)->setComment("maximum deltaR of track to jet. If -ve just use cut from JTA");
0070 desc.add<double>("maximumDecayLength", 999999.0);
0071 desc.add<double>("maximumDistanceToJetAxis", 999999.0);
0072 desc.add<std::string>("trackQualityClass", "any");
0073 desc.add<bool>("useVariableJTA", false);
0074 desc.add<double>("a_dR", -0.001053);
0075 desc.add<double>("b_dR", 0.6263);
0076 desc.add<double>("a_pT", 0.005263);
0077 desc.add<double>("b_pT", 0.3684);
0078 desc.add<double>("min_pT", 120);
0079 desc.add<double>("max_pT", 500);
0080 desc.add<double>("min_pT_dRcut", 0.5);
0081 desc.add<double>("max_pT_dRcut", 0.1);
0082 desc.add<double>("max_pT_trackPTcut", 3.);
0083 }
0084
0085 protected:
0086 std::multiset<float> orderedSignificances(const TagInfo& tkip) const {
0087 const std::vector<reco::btag::TrackIPData>& impactParameters((tkip.impactParameterData()));
0088 const Container& tracks(tkip.selectedTracks());
0089 std::multiset<float> significances;
0090 int i = 0;
0091 if (tkip.primaryVertex().isNull()) {
0092 return std::multiset<float>();
0093 }
0094
0095 GlobalPoint pv(tkip.primaryVertex()->position().x(),
0096 tkip.primaryVertex()->position().y(),
0097 tkip.primaryVertex()->position().z());
0098
0099 for (std::vector<reco::btag::TrackIPData>::const_iterator it = impactParameters.begin();
0100 it != impactParameters.end();
0101 ++it, i++) {
0102 if (fabs(impactParameters[i].distanceToJetAxis.value()) < m_cutMaxDistToAxis &&
0103 (impactParameters[i].closestToJetAxis - pv).mag() < m_cutMaxDecayLen &&
0104 (m_useAllQualities == true ||
0105 reco::btag::toTrack(tracks[i])->quality(m_trackQuality)) &&
0106 (fabs(((m_ipType == 0) ? it->ip3d : it->ip2d).value()) > m_minIP)
0107 ) {
0108
0109 float signed_sig = ((m_ipType == 0) ? it->ip3d : it->ip2d).significance();
0110 float unsigned_sig = fabs(signed_sig);
0111 float significance = (m_useSignedIPSig) ? signed_sig : unsigned_sig;
0112
0113 if (useVariableJTA_) {
0114 if (tkip.variableJTA(varJTApars)[i])
0115 significances.insert(significance);
0116 } else
0117 if (m_deltaR <= 0 ||
0118 ROOT::Math::VectorUtil::DeltaR((*tkip.jet()).p4().Vect(), (*tracks[i]).momentum()) < m_deltaR)
0119 significances.insert(significance);
0120 }
0121 }
0122
0123 return significances;
0124 }
0125
0126 bool useVariableJTA_;
0127 reco::btag::variableJTAParameters varJTApars;
0128
0129 double m_minIP;
0130 bool m_useSignedIPSig;
0131
0132 int m_nthTrack;
0133 int m_ipType;
0134 double m_deltaR;
0135 double m_cutMaxDecayLen;
0136 double m_cutMaxDistToAxis;
0137 reco::TrackBase::TrackQuality m_trackQuality;
0138 bool m_useAllQualities;
0139 };
0140
0141 #endif