CMS 3D CMS Logo

TemplatedSecondaryVertexTagInfo.h
Go to the documentation of this file.
1 #ifndef DataFormats_BTauReco_TemplatedSecondaryVertexTagInfo_h
2 #define DataFormats_BTauReco_TemplatedSecondaryVertexTagInfo_h
3 
4 #include <utility>
5 #include <vector>
6 
18 #include <functional>
19 #include <algorithm>
20 
24 
25 
26 
27 
28 namespace reco {
29 namespace btag{
30 inline float weight(const reco::TrackRef & t, const reco::Vertex &v) {return v.trackWeight(t);}
31 inline float weight(const reco::CandidatePtr & c, const reco::VertexCompositePtrCandidate &v) {return std::find(v.daughterPtrVector().begin(),v.daughterPtrVector().end(),c)!=v.daughterPtrVector().end();}
32 
33  struct TrackData {
34  enum Status {
38  };
39 
40  inline bool usedForVertexFit() const
41  { return (int)svStatus >= (int)trackUsedForVertexFit; }
42  inline bool associatedToVertex() const
43  { return (int)svStatus >= (int)trackAssociatedToVertex; }
44  inline bool associatedToVertex(unsigned int index) const
45  { return (int)svStatus == (int)trackAssociatedToVertex + (int)index; }
46 
48  };
49  typedef std::pair<unsigned int, TrackData> IndexedTrackData;
50 
51 }
52 template <class IPTI, class VTX>
54  public:
57 
58  struct VertexData {
59  VTX vertex;
60  Measurement1D dist1d,dist2d, dist3d;
62 
63  // Used by ROOT storage
65  };
66 
67  struct TrackFinder {
68  TrackFinder(const typename IPTI::input_container &tracks,
69  const typename IPTI::input_container::value_type &track) :
70  tracks(tracks), track(track) {}
71 
72  bool operator () (const IndexedTrackData &idt)
73  { return tracks[idt.first] == track; }
74 
75  const typename IPTI::input_container &tracks;
77  };
78 
80  bool operator () (const IndexedTrackData &idt)
81  { return idt.second.associatedToVertex(); }
82  };
83 
86  index(index) {}
87 
88  bool operator () (const IndexedTrackData &idt)
89  { return idt.second.associatedToVertex(index); }
90 
91  unsigned int index;
92  };
93 
94 
95 
96  typedef typename IPTI::input_container input_container;
97 
100 
102  const std::vector<IndexedTrackData> &trackData,
103  const std::vector<VertexData> &svData,
104  unsigned int vertexCandidates,
105  const edm::Ref<std::vector<IPTI> >&);
106 
108  TemplatedSecondaryVertexTagInfo * clone(void) const override {
109  return new TemplatedSecondaryVertexTagInfo(*this);
110  }
111 
113  { return m_trackIPTagInfoRef; }
114 
115  edm::RefToBase<Jet> jet(void) const override
116  { return m_trackIPTagInfoRef->jet(); }
117 
118 // virtual input_container ipTracks(void) const
119 // { return m_trackIPTagInfoRef->tracks(); }
120 
121 //AR TODO is it needed?
122 // const JetTracksAssociationRef &jtaRef(void) const
123 // { return m_trackIPTagInfoRef->jtaRef(); }
124 
125  const VTX &secondaryVertex(unsigned int index) const
126  { return m_svData[index].vertex; }
127 
128  unsigned int nSelectedTracks() const { return m_trackData.size(); }
129  unsigned int nVertexTracks() const;
130  unsigned int nVertexTracks(unsigned int index) const;
131  unsigned int nVertices() const { return m_svData.size(); }
132  unsigned int nVertexCandidates() const { return m_vertexCandidates; }
133 
134  input_container selectedTracks() const;
135  input_container vertexTracks() const;
136  input_container vertexTracks(unsigned int index) const;
137 
138  typename input_container::value_type track(unsigned int index) const;
139  unsigned int findTrack(const typename input_container::value_type &track) const;
140 
141  const TrackData &trackData(unsigned int index) const;
142  const TrackData &trackData(const typename input_container::value_type &track) const;
143 
144  const reco::btag::TrackIPData &trackIPData(unsigned int index) const;
145  const reco::btag::TrackIPData &trackIPData(const typename input_container::value_type &track) const;
146 
147  float trackWeight(unsigned int svIndex, unsigned int trackindex) const;
148  float trackWeight(unsigned int svIndex, const typename input_container::value_type &track) const;
149 
151  flightDistance(unsigned int index, int dim =0) const{
152  if(dim==1) return m_svData[index].dist1d;
153  else if(dim==2) return m_svData[index].dist2d;
154  else return m_svData[index].dist3d;
155  }
156  const GlobalVector &flightDirection(unsigned int index) const
157  { return m_svData[index].direction; }
158  TaggingVariableList taggingVariables() const override;
159 
160  // Used by ROOT storage
162 
163  private:
164  std::vector<IndexedTrackData> m_trackData;
165  std::vector<VertexData> m_svData;
166  unsigned int m_vertexCandidates;
167 
168  edm::Ref<std::vector<IPTI> > m_trackIPTagInfoRef;
169 };
170 
171 
173  const std::vector<IndexedTrackData> &trackData,
174  const std::vector<VertexData> &svData,
175  unsigned int vertexCandidates,
176  const edm::Ref<std::vector<IPTI> > & trackIPTagInfoRef) :
177  m_trackData(trackData),
178  m_svData(svData),
179  m_vertexCandidates(vertexCandidates),
180  m_trackIPTagInfoRef(trackIPTagInfoRef)
181 {
182 }
183 
184 template<class IPTI,class VTX> unsigned int TemplatedSecondaryVertexTagInfo<IPTI,VTX>::nVertexTracks() const
185 {
186  return std::count_if(m_trackData.begin(), m_trackData.end(),
188 }
189 
190 template<class IPTI,class VTX> unsigned int TemplatedSecondaryVertexTagInfo<IPTI,VTX>::nVertexTracks(unsigned int index) const
191 {
192  return std::count_if(m_trackData.begin(), m_trackData.end(),
194 }
195 
196 template<class IPTI,class VTX>
198 {
199  input_container trackRefs;
200  const input_container &trackIPTrackRefs =
201  m_trackIPTagInfoRef->selectedTracks();
202 
203  for(typename std::vector<typename reco::TemplatedSecondaryVertexTagInfo<IPTI,VTX>::IndexedTrackData>::const_iterator iter =
204  m_trackData.begin(); iter != m_trackData.end(); iter++)
205 
206  trackRefs.push_back(trackIPTrackRefs[iter->first]);
207 
208  return trackRefs;
209 }
210 
211 template<class IPTI,class VTX>
213 {
214  input_container trackRefs;
215  const input_container &trackIPTrackRefs =
216  m_trackIPTagInfoRef->selectedTracks();
217 
218  for(typename std::vector<typename reco::TemplatedSecondaryVertexTagInfo<IPTI,VTX>::IndexedTrackData>::const_iterator iter =
219  m_trackData.begin(); iter != m_trackData.end(); iter++)
220 
221  if (iter->second.associatedToVertex())
222  trackRefs.push_back(trackIPTrackRefs[iter->first]);
223 
224  return trackRefs;
225 }
226 
227 template<class IPTI,class VTX>
229 {
230  input_container trackRefs;
231  const input_container &trackIPTrackRefs =
232  m_trackIPTagInfoRef->selectedTracks();
233 
234  for(typename std::vector<typename reco::TemplatedSecondaryVertexTagInfo<IPTI,VTX>::IndexedTrackData>::const_iterator iter =
235  m_trackData.begin(); iter != m_trackData.end(); iter++)
236 
237  if (iter->second.associatedToVertex(index))
238  trackRefs.push_back(trackIPTrackRefs[iter->first]);
239 
240  return trackRefs;
241 }
242 
244 {
245  return m_trackIPTagInfoRef->selectedTracks()[m_trackData[index].first];
246 }
247 
248 template<class IPTI,class VTX> unsigned int TemplatedSecondaryVertexTagInfo<IPTI,VTX>::findTrack(const typename input_container::value_type &track) const
249 {
251  std::find_if(m_trackData.begin(), m_trackData.end(),
252  TrackFinder(m_trackIPTagInfoRef->selectedTracks(),
253  track));
254 
255  if (pos == m_trackData.end())
257  << "Track not found in "
258  " TemplatedSecondaryVertexTagInfo<IPTI,VTX>::findTrack." << std::endl;
259 
260  return pos - m_trackData.begin();
261 }
262 
263 template<class IPTI,class VTX>
265 {
266  return m_trackData[index].second;
267 }
268 
269 template<class IPTI,class VTX>
271 {
272  return m_trackData[findTrack(track)].second;
273 }
274 
275 template<class IPTI,class VTX>
277 {
278  return m_trackIPTagInfoRef->impactParameterData()[
279  m_trackData[index].first];
280 }
281 
282 template<class IPTI,class VTX>
284 {
285  return trackIPData(findTrack(track));
286 }
287 
288 template<class IPTI,class VTX> float TemplatedSecondaryVertexTagInfo<IPTI,VTX>::trackWeight(unsigned int svIndex,
289  const typename input_container::value_type &track) const
290 {
291  return reco::btag::weight(track,m_svData[svIndex].vertex);
292 }
293 
294 template<class IPTI,class VTX> float TemplatedSecondaryVertexTagInfo<IPTI,VTX>::trackWeight(unsigned int svIndex,
295  unsigned int trackIndex) const
296 {
297  return trackWeight(svIndex, track(trackIndex));
298 }
299 
301 {
303 
304  for(typename std::vector<typename TemplatedSecondaryVertexTagInfo<IPTI,VTX>::VertexData>::const_iterator iter = m_svData.begin();
305  iter != m_svData.end(); iter++) {
307  iter->dist1d.value(), true);
309  iter->dist1d.significance(), true);
311  iter->dist2d.value(), true);
313  iter->dist2d.significance(), true);
315  iter->dist3d.value(), true);
317  iter->dist3d.significance(), true);
318 
320  Geom::deltaR(iter->direction, jet()->momentum()), true);
321  }
322 
323  vars.insert(btau::jetNSecondaryVertices, m_vertexCandidates, true);
324  vars.insert(btau::vertexNTracks, nVertexTracks(), true);
325 
326  vars.finalize();
327  return vars;
328 }
329 
330 }
331 #endif // DataFormats_BTauReco_TemplatedSecondaryVertexTagInfo_h
const VTX & secondaryVertex(unsigned int index) const
float weight(const reco::TrackRef &t, const reco::Vertex &v)
edm::RefToBase< Jet > jet(void) const override
returns a polymorphic reference to the tagged jet
#define CMS_CLASS_VERSION(_version_)
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:20
Definition: Matching.h:10
virtual const daughters & daughterPtrVector() const
references to daughtes
const edm::Ref< std::vector< IPTI > > & trackIPTagInfoRef() const
unsigned int findTrack(const typename input_container::value_type &track) const
float trackWeight(unsigned int svIndex, unsigned int trackindex) const
const reco::btag::TrackIPData & trackIPData(unsigned int index) const
const TrackData & trackData(unsigned int index) const
def template(fileName, svg, replaceme="REPLACEME")
Definition: svgfig.py:521
TemplatedSecondaryVertexTagInfo * clone(void) const override
clone
float trackWeight(const TREF &r) const
returns the weight with which a Track has contributed to the vertex-fit.
Definition: Vertex.h:81
input_container::value_type track(unsigned int index) const
const GlobalVector & flightDirection(unsigned int index) const
ROOT::Math::DisplacementVector3D< ROOT::Math::Cartesian3D< float >, ROOT::Math::GlobalCoordinateSystemTag > GlobalVector
vector in glovbal coordinate system
Definition: Vector3D.h:27
TrackFinder(const typename IPTI::input_container &tracks, const typename IPTI::input_container::value_type &track)
TaggingVariableList taggingVariables() const override
returns a description of the extended informations in a TaggingVariableList
bool associatedToVertex(unsigned int index) const
fixed size matrix
HLT enums.
Measurement1D flightDistance(unsigned int index, int dim=0) const
vars
Definition: DeepTauId.cc:77
std::pair< unsigned int, TrackData > IndexedTrackData
float trackWeight(const reco::Vertex &sv, const reco::TransientTrack &track)
void insert(const TaggingVariable &variable, bool delayed=false)