CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
TemplatedTrackCountingComputer.h
Go to the documentation of this file.
1 #ifndef ImpactParameter_TemplatedTrackCountingComputer_h
2 #define ImpactParameter_TemplatedTrackCountingComputer_h
3 
8 #include "Math/GenVector/VectorUtil.h"
10 
11 template <class Container, class Base>
13 public:
14  using Tokens = void;
15 
17 
19  m_minIP = parameters.existsAs<double>("minimumImpactParameter")
20  ? parameters.getParameter<double>("minimumImpactParameter")
21  : -1;
22  m_useSignedIPSig = parameters.existsAs<bool>("useSignedImpactParameterSig")
23  ? parameters.getParameter<bool>("useSignedImpactParameterSig")
24  : true;
25  m_nthTrack = parameters.getParameter<int>("nthTrack");
26  m_ipType = parameters.getParameter<int>("impactParameterType");
27  m_deltaR = parameters.getParameter<double>("deltaR");
28  m_cutMaxDecayLen = parameters.getParameter<double>("maximumDecayLength"); //used
29  m_cutMaxDistToAxis = parameters.getParameter<double>("maximumDistanceToJetAxis"); //used
30  //
31  // access track quality class; "any" takes everything
32  //
33  std::string trackQualityType = parameters.getParameter<std::string>("trackQualityClass"); //used
35  m_useAllQualities = false;
36  if (trackQualityType == "any" || trackQualityType == "Any" || trackQualityType == "ANY")
37  m_useAllQualities = true;
38 
39  uses("ipTagInfos");
40 
42  parameters.existsAs<bool>("useVariableJTA") ? parameters.getParameter<bool>("useVariableJTA") : false;
43  if (useVariableJTA_) {
44  varJTApars = {parameters.getParameter<double>("a_dR"),
45  parameters.getParameter<double>("b_dR"),
46  parameters.getParameter<double>("a_pT"),
47  parameters.getParameter<double>("b_pT"),
48  parameters.getParameter<double>("min_pT"),
49  parameters.getParameter<double>("max_pT"),
50  parameters.getParameter<double>("min_pT_dRcut"),
51  parameters.getParameter<double>("max_pT_dRcut"),
52  parameters.getParameter<double>("max_pT_trackPTcut")};
53  }
54  }
55 
56  float discriminator(const TagInfoHelper& ti) const override {
57  const TagInfo& tkip = ti.get<TagInfo>();
58  std::multiset<float> significances = orderedSignificances(tkip);
59  std::multiset<float>::reverse_iterator nth = significances.rbegin();
60  for (int i = 0; i < m_nthTrack - 1 && nth != significances.rend(); i++)
61  nth++;
62  if (nth != significances.rend())
63  return *nth;
64  else
65  return -100.;
66  }
67 
68 protected:
69  std::multiset<float> orderedSignificances(const TagInfo& tkip) const {
70  const std::vector<reco::btag::TrackIPData>& impactParameters((tkip.impactParameterData()));
71  const Container& tracks(tkip.selectedTracks());
72  std::multiset<float> significances;
73  int i = 0;
74  if (tkip.primaryVertex().isNull()) {
75  return std::multiset<float>();
76  }
77 
78  GlobalPoint pv(tkip.primaryVertex()->position().x(),
79  tkip.primaryVertex()->position().y(),
80  tkip.primaryVertex()->position().z());
81 
82  for (std::vector<reco::btag::TrackIPData>::const_iterator it = impactParameters.begin();
83  it != impactParameters.end();
84  ++it, i++) {
85  if (fabs(impactParameters[i].distanceToJetAxis.value()) < m_cutMaxDistToAxis && // distance to JetAxis
86  (impactParameters[i].closestToJetAxis - pv).mag() < m_cutMaxDecayLen && // max decay len
87  (m_useAllQualities == true ||
88  reco::btag::toTrack(tracks[i])->quality(m_trackQuality)) && // use selected track qualities
89  (fabs(((m_ipType == 0) ? it->ip3d : it->ip2d).value()) > m_minIP) // minimum impact parameter
90  ) {
91  //calculate the signed or un-signed significance
92  float signed_sig = ((m_ipType == 0) ? it->ip3d : it->ip2d).significance();
93  float unsigned_sig = fabs(signed_sig);
94  float significance = (m_useSignedIPSig) ? signed_sig : unsigned_sig;
95 
96  if (useVariableJTA_) {
97  if (tkip.variableJTA(varJTApars)[i])
98  significances.insert(significance);
99  } else // no using variable JTA, use the default method
100  if (m_deltaR <= 0 ||
101  ROOT::Math::VectorUtil::DeltaR((*tkip.jet()).p4().Vect(), (*tracks[i]).momentum()) < m_deltaR)
102  significances.insert(significance);
103  }
104  }
105 
106  return significances;
107  }
108 
111 
112  double m_minIP;
114 
116  int m_ipType;
117  double m_deltaR;
122 };
123 
124 #endif // ImpactParameter_TemplatedTrackCountingComputer_h
bool quality(const TrackQuality) const
Track quality.
Definition: TrackBase.h:552
TrackQuality
track quality
Definition: TrackBase.h:150
const T & get(unsigned int index=0) const
const Container & selectedTracks() const
Definition: IPTagInfo.h:99
reco::TrackBase::TrackQuality m_trackQuality
const reco::Track * toTrack(const reco::TrackBaseRef &t)
Definition: IPTagInfo.h:24
TemplatedTrackCountingComputer(const edm::ParameterSet &parameters)
TEMPL(T2) struct Divides void
Definition: Factorize.h:24
void uses(unsigned int id, const std::string &label)
def pv(vc)
Definition: MetAnalyzer.py:7
const std::vector< btag::TrackIPData > & impactParameterData() const
Definition: IPTagInfo.h:90
reco::IPTagInfo< Container, Base > TagInfo
bool isNull() const
Checks for null.
Definition: Ref.h:235
static TrackQuality qualityByName(const std::string &name)
Definition: TrackBase.cc:126
auto const & tracks
cannot be loose
reco::btag::variableJTAParameters varJTApars
float discriminator(const TagInfoHelper &ti) const override
significance
Definition: met_cff.py:18
std::vector< bool > variableJTA(const btag::variableJTAParameters &params) const
Definition: IPTagInfo.h:210
const edm::Ref< VertexCollection > & primaryVertex() const
Definition: IPTagInfo.h:133
std::multiset< float > orderedSignificances(const TagInfo &tkip) const
edm::AssociationVector< reco::JetRefBaseProd, Values > Container