ATLAS Offline Software
MSVVariablesFactory.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
3 */
4 
15 #include "VxVertex/VxTrackAtVertex.h"
16 
17 #include "TrkParameters/TrackParameters.h"
18 #include "VxSecVertex/VxSecVertexInfo.h"
19 #include "VxSecVertex/VxSecVKalVertexInfo.h"
20 
21 #include <TMath.h>
22 
23 #include "CLHEP/Vector/LorentzVector.h"
24 #include "GeoPrimitives/GeoPrimitives.h"
25 
26 #include "JetTagTools/MSVVariablesFactory.h"
27 
28 #include "GeoPrimitives/GeoPrimitivesHelpers.h"
29 
30 #include "xAODTracking/Vertex.h"
31 #include "xAODTracking/TrackParticle.h"
32 #include "xAODBTagging/SecVtxHelper.h"
33 
34 #include <vector>
35 #include <string>
36 
37 namespace Analysis {
38 
39  MSVVariablesFactory::MSVVariablesFactory(const std::string& name,
40  const std::string& n,
41  const IInterface* p):
42  AthAlgTool(name, n,p)
43  {
44  declareInterface<IMSVVariablesFactory>(this);
45  }
46 
47  StatusCode MSVVariablesFactory::initialize() {
48  ATH_MSG_DEBUG(" Initialization of MSVVariablesFactory succesfull");
49  return StatusCode::SUCCESS;
50  }
51 
52  StatusCode MSVVariablesFactory::finalize() {
53  ATH_MSG_DEBUG(" Finalization of MSVVariablesFactory succesfull");
54  return StatusCode::SUCCESS;
55  }
56 
57  StatusCode MSVVariablesFactory::createMSVContainer
58  (const xAOD::Jet &myJet, const Trk::VxSecVKalVertexInfo* myVertexInfoVKal,
59  xAOD::VertexContainer* VertexContainer,const xAOD::Vertex& PrimaryVtx) const {
60 
61  Amg::Vector3D jet_V3(myJet.p4().Px(), myJet.p4().Py(), myJet.p4().Pz());
62  float jetenergy=0.;
63  const xAOD::Vertex* priVtx = &PrimaryVtx;
64  std::vector< ElementLink< xAOD::VertexContainer > > MSVVertexLinks;
65  const std::vector<xAOD::Vertex*> myVertices = myVertexInfoVKal->vertices();
66  if(myVertices.empty()){
67  ATH_MSG_DEBUG("#BTAG# no MSV vertices...fill default values only... ");
68  xAOD::Vertex* vertex = new xAOD::Vertex();
69  VertexContainer->push_back(vertex);
70  xAOD::SecVtxHelper::setVertexMass(vertex, -9.);
71  xAOD::SecVtxHelper::setEnergyFraction(vertex, -9.);
72  xAOD::SecVtxHelper::setVtxNtrk(vertex, -9);
73  xAOD::SecVtxHelper::setVtxpt(vertex, -9.);
74  xAOD::SecVtxHelper::setVtxeta(vertex, -9.);
75  xAOD::SecVtxHelper::setVtxphi(vertex, -9.);
76  xAOD::SecVtxHelper::setVtxnormDist(vertex, -9.);
77  return StatusCode::SUCCESS;
78  }
79 
80  jetenergy = myVertexInfoVKal->energyTrkInJet();
81 
82  for (const auto& vertex : myVertexInfoVKal->vertices()){
83  VertexContainer->push_back(vertex);
84  //additional info per vertex
85  double sumpx = 0.0;
86  double sumpy = 0.0;
87  double sumpz = 0.0;
88  double sume = 0.0;
89  const std::vector<ElementLink<xAOD::TrackParticleContainer> > myTrackLinks = vertex->trackParticleLinks();
90  if (myTrackLinks.empty()) {
91  ATH_MSG_WARNING("#BTAG# No Track Links attached to the track at the sec vertex... ");
92  }
93  int npsec = 0;
94  const std::vector<Trk::VxTrackAtVertex> myTracks=vertex->vxTrackAtVertex();
95  if (!myTracks.empty()) {
96  npsec=myTracks.size();
97  for (const auto& track : myTracks) {
98  const Trk::Perigee* perigee = dynamic_cast<const Trk::Perigee*>(track.perigeeAtVertex());
99  if(perigee){
100  sumpx += perigee->momentum().x();
101  sumpy += perigee->momentum().y();
102  sumpz += perigee->momentum().z();
103  sume += std::hypot(perigee->momentum().mag(), 139.5702);
104  }else{
105  ATH_MSG_WARNING("#BTAG# perigee for VxTrackAtVertex not found");
106  }
107  }
108  }
109 
110  CLHEP::HepLorentzVector vtxp4(sumpx,sumpy,sumpz,sume);
111  float efrac = (jetenergy>0) ? vtxp4.e()/jetenergy : 0;
112  xAOD::SecVtxHelper::setVertexMass(vertex, vtxp4.m());
113  xAOD::SecVtxHelper::setEnergyFraction(vertex, efrac);
114  xAOD::SecVtxHelper::setVtxNtrk(vertex, npsec);
115  xAOD::SecVtxHelper::setVtxpt(vertex, vtxp4.perp());
116  xAOD::SecVtxHelper::setVtxeta(vertex, vtxp4.eta());
117  xAOD::SecVtxHelper::setVtxphi(vertex, vtxp4.phi());
118 
119  ATH_MSG_DEBUG("#BTAG# mass per vertex = "<<vtxp4.m());
120  double localdistnrm = 0;
121  std::vector<const xAOD::Vertex*> vecVtxHolder;
122  vecVtxHolder.push_back(vertex);
123 
124  ATH_MSG_DEBUG("Factory PVX x = " << priVtx->x() << " y = " << priVtx->y() << " z = " << priVtx->z());
125  if (priVtx) {
126  localdistnrm = get3DSignificance(priVtx, vecVtxHolder, jet_V3);
127  } else {
128  ATH_MSG_WARNING("#BTAG# Tagging requested, but no primary vertex supplied.");
129  }
130  xAOD::SecVtxHelper::setVtxnormDist(vertex, localdistnrm);
131  //track links,
132  vertex->setTrackParticleLinks(myTrackLinks);
133 
134  } //end loop vertexcontainer
135 
136  return StatusCode::SUCCESS;
137  }
138 
139  StatusCode MSVVariablesFactory::fillMSVVariables
140  (const xAOD::Jet &myJet, xAOD::BTagging* BTag,
141  const Trk::VxSecVKalVertexInfo* myVertexInfoVKal,
142  xAOD::VertexContainer* VertexContainer, const xAOD::Vertex& PrimaryVtx,
143  std::string basename) const {
144 
145  Amg::Vector3D jet_V3(myJet.p4().Px(), myJet.p4().Py(), myJet.p4().Pz());
146  int nvsec = 0;
147  float jetenergy = 0.;
148  int n2t = 0;
149  float distnrm = 0.;
150  const xAOD::Vertex* priVtx = &PrimaryVtx;
151  std::vector< ElementLink< xAOD::VertexContainer > > MSVVertexLinks;
152  const std::vector<xAOD::Vertex*> myVertices = myVertexInfoVKal->vertices();
153  if(myVertices.empty()){
154  ATH_MSG_DEBUG("#BTAG# no MSV vertices...fill default values only... ");
155  BTag->setVariable<int>(basename, "N2Tpair", n2t);
156  BTag->setVariable<float>(basename, "energyTrkInJet", jetenergy);
157  BTag->setVariable<int>(basename, "nvsec", nvsec);
158  BTag->setVariable<float>(basename, "normdist", distnrm);
159  BTag->setVariable<std::vector<ElementLink<xAOD::VertexContainer> > >(basename, "vertices", MSVVertexLinks);
160  BTag->setDynVxELName(basename, "vertices");
161  xAOD::Vertex* vertex = new xAOD::Vertex();
162  VertexContainer->push_back(vertex);
163  xAOD::SecVtxHelper::setVertexMass(vertex, -9.);
164  xAOD::SecVtxHelper::setEnergyFraction(vertex, -9.);
165  xAOD::SecVtxHelper::setVtxNtrk(vertex, -9);
166  xAOD::SecVtxHelper::setVtxpt(vertex, -9.);
167  xAOD::SecVtxHelper::setVtxeta(vertex, -9.);
168  xAOD::SecVtxHelper::setVtxphi(vertex, -9.);
169  xAOD::SecVtxHelper::setVtxnormDist(vertex, -9.);
170  return StatusCode::SUCCESS;
171  }
172 
173  jetenergy = myVertexInfoVKal->energyTrkInJet();
174  n2t = myVertexInfoVKal->n2trackvertices();
175  BTag->setVariable<int>(basename, "N2Tpair", n2t);
176  BTag->setVariable<float>(basename, "energyTrkInJet", jetenergy);
177 
178  std::vector<const xAOD::Vertex*> vecVertices;
179  for (const auto& vertex : myVertexInfoVKal->vertices()) {
180  VertexContainer->push_back(vertex);
181  //additional info per vertex
182  vecVertices.push_back(vertex);
183  double sumpx = 0.0;
184  double sumpy = 0.0;
185  double sumpz = 0.0;
186  double sume = 0.0;
187  const std::vector<ElementLink<xAOD::TrackParticleContainer> > myTrackLinks = vertex->trackParticleLinks();
188  if (myTrackLinks.empty()) {
189  ATH_MSG_WARNING("#BTAG# No Track Links attached to the track at the sec vertex... ");
190  }
191  int npsec = 0;
192  const std::vector<Trk::VxTrackAtVertex> myTracks=vertex->vxTrackAtVertex();
193  if (!myTracks.empty()) {
194  npsec=myTracks.size();
195  for (const auto& track : myTracks) {
196  const Trk::Perigee* perigee = dynamic_cast<const Trk::Perigee*>(track.perigeeAtVertex());
197  if(perigee){
198  sumpx += perigee->momentum().x();
199  sumpy += perigee->momentum().y();
200  sumpz += perigee->momentum().z();
201  sume += std::hypot(perigee->momentum().mag(), 139.5702);
202  }else{
203  ATH_MSG_WARNING("#BTAG# perigee for VxTrackAtVertex not found");
204  }
205  }
206  }
207 
208  CLHEP::HepLorentzVector vtxp4(sumpx,sumpy,sumpz,sume);
209  float efrac = (jetenergy>0) ? vtxp4.e()/jetenergy : 0;
210  xAOD::SecVtxHelper::setVertexMass(vertex, vtxp4.m());
211  xAOD::SecVtxHelper::setEnergyFraction(vertex, efrac);
212  xAOD::SecVtxHelper::setVtxNtrk(vertex, npsec);
213  xAOD::SecVtxHelper::setVtxpt(vertex, vtxp4.perp());
214  xAOD::SecVtxHelper::setVtxeta(vertex, vtxp4.eta());
215  xAOD::SecVtxHelper::setVtxphi(vertex, vtxp4.phi());
216 
217  ATH_MSG_DEBUG("#BTAG# mass per vertex = "<<vtxp4.m());
218  double localdistnrm = 0;
219  std::vector<const xAOD::Vertex*> vecVtxHolder;
220  vecVtxHolder.push_back(vertex);
221 
222  ATH_MSG_DEBUG("Factory PVX x = " << priVtx->x() << " y = " << priVtx->y() << " z = " << priVtx->z());
223  if (priVtx) {
224  localdistnrm = get3DSignificance(priVtx, vecVtxHolder, jet_V3);
225  } else {
226  ATH_MSG_WARNING("#BTAG# Tagging requested, but no primary vertex supplied.");
227  }
228  xAOD::SecVtxHelper::setVtxnormDist(vertex, localdistnrm);
229  //track links,
230  vertex->setTrackParticleLinks(myTrackLinks);
231 
232  ElementLink< xAOD::VertexContainer> linkBTagVertex;
233  linkBTagVertex.toContainedElement(*VertexContainer, vertex);
234  MSVVertexLinks.push_back(linkBTagVertex);
235  } //end loop vertexcontainer
236 
237  BTag->setVariable<std::vector<ElementLink<xAOD::VertexContainer> > >(basename, "vertices", MSVVertexLinks);
238  BTag->setDynVxELName(basename, "vertices");
239 
240  if (priVtx) {
241  distnrm = get3DSignificance(priVtx, vecVertices, jet_V3);
242  } else {
243  ATH_MSG_WARNING("#BTAG# Tagging requested, but no primary vertex supplied.");
244  distnrm=0.;
245  }
246  nvsec = vecVertices.size();
247  BTag->setVariable<int>(basename, "nvsec", nvsec);
248  BTag->setVariable<float>(basename, "normdist", distnrm);
249 
250  return StatusCode::SUCCESS;
251 
252  }
253 
254  double MSVVariablesFactory::get3DSignificance
255  (const xAOD::Vertex* priVertex,
256  std::vector<const xAOD::Vertex*>& secVertex,
257  const Amg::Vector3D jetDirection) const {
258 
259  if(!secVertex.size()) return 0;
260  std::vector<Amg::Vector3D> positions;
261  std::vector<AmgSymMatrix(3)> weightMatrices;
262  Amg::Vector3D weightTimesPosition(0.,0.,0.);
263  AmgSymMatrix(3) sumWeights;
264  sumWeights.setZero();
265 
266  for (const auto& vertex : secVertex) {
267  positions.push_back(vertex->position());
268  weightMatrices.push_back(vertex->covariancePosition().inverse());
269  weightTimesPosition += weightMatrices.back() * positions.back();
270  sumWeights += weightMatrices.back();
271  }
272 
273  bool invertible;
274  AmgSymMatrix(3) meanCovariance;
275  meanCovariance.setZero();
276  sumWeights.computeInverseWithCheck(meanCovariance, invertible);
277  if (!invertible) {
278  ATH_MSG_WARNING("#BTAG# Could not invert sum of sec vtx matrices");
279  return 0.;
280  }
281  Amg::Vector3D meanPosition = meanCovariance * weightTimesPosition;
282  AmgSymMatrix(3) covariance = meanCovariance + priVertex->covariancePosition();
283 
284  double Lx = meanPosition[0]-priVertex->position().x();
285  double Ly = meanPosition[1]-priVertex->position().y();
286  double Lz = meanPosition[2]-priVertex->position().z();
287 
288  const double decaylength = sqrt(Lx*Lx + Ly*Ly + Lz*Lz);
289  const double inv_decaylength = 1. / decaylength;
290  double dLdLx = Lx * inv_decaylength;
291  double dLdLy = Ly * inv_decaylength;
292  double dLdLz = Lz * inv_decaylength;
293  double decaylength_err = sqrt(dLdLx*dLdLx*covariance(0,0) +
294  dLdLy*dLdLy*covariance(1,1) +
295  dLdLz*dLdLz*covariance(2,2) +
296  2.*dLdLx*dLdLy*covariance(0,1) +
297  2.*dLdLx*dLdLz*covariance(0,2) +
298  2.*dLdLy*dLdLz*covariance(1,2));
299 
300  double decaylength_significance = 0.;
301  if (decaylength_err != 0.) decaylength_significance = decaylength/decaylength_err;
302  double L_proj_jetDir = jetDirection.x()*Lx + jetDirection.y()*Ly + jetDirection.z()*Lz;
303  if (L_proj_jetDir < 0.) decaylength_significance *= -1.;
304 
305  return decaylength_significance;
306 
307  }
308 
309 }//end Analysis namespace
xAOD::Vertex_v1::x
float x() const
Returns the x position.
Analysis::MSVVariablesFactory::finalize
virtual StatusCode finalize() override
Definition: MSVVariablesFactory.cxx:52
TrackParameters.h
Trk::VxSecVKalVertexInfo::n2trackvertices
int n2trackvertices() const
Definition: VxSecVKalVertexInfo.h:177
xAOD::Vertex
Vertex_v1 Vertex
Define the latest version of the vertex class.
Definition: Event/xAOD/xAODTracking/xAODTracking/Vertex.h:16
Analysis::MSVVariablesFactory::fillMSVVariables
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
Definition: MSVVariablesFactory.cxx:140
Trk::ParametersT
Dummy class used to allow special convertors to be called for surfaces owned by a detector element.
Definition: EMErrorDetail.h:25
Trk::VxSecVertexInfo::vertices
const std::vector< xAOD::Vertex * > & vertices() const
Definition: VxSecVertexInfo.cxx:100
xAOD::SecVtxHelper::setVtxNtrk
void setVtxNtrk(xAOD::Vertex *, int value)
Definition: SecVtxHelper.cxx:48
Analysis::MSVVariablesFactory::MSVVariablesFactory
MSVVariablesFactory(const std::string &name, const std::string &n, const IInterface *p)
Definition: MSVVariablesFactory.cxx:39
xAOD::Vertex_v1::position
const Amg::Vector3D & position() const
Returns the 3-pos.
SecVtxHelper.h
xAOD::SecVtxHelper::setEnergyFraction
void setEnergyFraction(xAOD::Vertex *, float value)
Definition: SecVtxHelper.cxx:33
xAOD::VertexContainer
VertexContainer_v1 VertexContainer
Definition of the current "Vertex container version".
Definition: VertexContainer.h:14
AmgSymMatrix
#define AmgSymMatrix(dim)
Definition: EventPrimitives.h:50
GeoPrimitives.h
python.utils.AtlRunQueryDQUtils.p
p
Definition: AtlRunQueryDQUtils.py:210
beamspotman.n
n
Definition: beamspotman.py:731
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
VxTrackAtVertex.h
VxSecVertexInfo.h
xAOD::Vertex_v1::z
float z() const
Returns the z position.
xAOD::SecVtxHelper::setVtxpt
void setVtxpt(xAOD::Vertex *, float value)
Definition: SecVtxHelper.cxx:62
xAOD::BTagging_v1
Definition: BTagging_v1.h:39
DataVector
Derived DataVector<T>.
Definition: DataVector.h:581
Vertex.h
Trk::VxSecVKalVertexInfo::energyTrkInJet
double energyTrkInJet() const
Definition: VxSecVKalVertexInfo.h:165
xAOD::SecVtxHelper::setVtxeta
void setVtxeta(xAOD::Vertex *, float value)
Definition: SecVtxHelper.cxx:76
MSVVariablesFactory.h
Analysis
The namespace of all packages in PhysicsAnalysis/JetTagging.
Definition: BTaggingCnvAlg.h:20
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:221
xAOD::SecVtxHelper::setVtxphi
void setVtxphi(xAOD::Vertex *, float value)
Definition: SecVtxHelper.cxx:90
xAOD::SecVtxHelper::setVtxnormDist
void setVtxnormDist(xAOD::Vertex *, float value)
Definition: SecVtxHelper.cxx:104
DataVector::push_back
value_type push_back(value_type pElem)
Add an element to the end of the collection.
Amg::Vector3D
Eigen::Matrix< double, 3, 1 > Vector3D
Definition: GeoPrimitives.h:47
Trk::VxSecVKalVertexInfo
Definition: VxSecVKalVertexInfo.h:70
TrackParticle.h
Analysis::MSVVariablesFactory::initialize
virtual StatusCode initialize() override
Definition: MSVVariablesFactory.cxx:47
xAOD::Jet_v1
Class describing a jet.
Definition: Jet_v1.h:57
Trk::vertex
@ vertex
Definition: MeasurementType.h:21
xAOD::Jet_v1::p4
virtual FourMom_t p4() const
The full 4-momentum of the particle.
Definition: Jet_v1.cxx:71
VxSecVKalVertexInfo.h
xAOD::Vertex_v1
Class describing a Vertex.
Definition: Vertex_v1.h:42
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition: AthMsgStreamMacros.h:32
GeoPrimitivesHelpers.h
xAOD::Vertex_v1::y
float y() const
Returns the y position.
xAOD::track
@ track
Definition: TrackingPrimitives.h:512
xAODType::BTag
@ BTag
The object is a b-tagging object.
Definition: ObjectType.h:60
AthAlgTool
Definition: AthAlgTool.h:26
xAOD::SecVtxHelper::setVertexMass
void setVertexMass(xAOD::Vertex *, float value)
Definition: SecVtxHelper.cxx:18
Analysis::MSVVariablesFactory::get3DSignificance
double get3DSignificance(const xAOD::Vertex *priVertex, std::vector< const xAOD::Vertex * > &secVertex, const Amg::Vector3D jetDirection) const
Definition: MSVVariablesFactory.cxx:255
Analysis::MSVVariablesFactory::createMSVContainer
virtual StatusCode createMSVContainer(const xAOD::Jet &, const Trk::VxSecVKalVertexInfo *myInfoVKal, xAOD::VertexContainer *btagVertex, const xAOD::Vertex &PV) const override
Definition: MSVVariablesFactory.cxx:58
beamspotman.basename
basename
Definition: beamspotman.py:640
Generated on Thu Nov 7 2024 21:21:31 for ATLAS Offline Software by   doxygen 1.8.18