ATLAS Offline Software
Loading...
Searching...
No Matches
Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx
Go to the documentation of this file.
1/*
2 Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
3*/
7// Header include
8#include "NewVrtSecInclusiveTool/NewVrtSecInclusiveTool.h"
9#include "TrkNeutralParameters/NeutralParameters.h"
10#include "TrkTrackSummary/TrackSummary.h"
11#include "TrkVKalVrtFitter/TrkVKalVrtFitter.h"
12#include "AthContainers/ConstAccessor.h"
13#include "CxxUtils/sincos.h"
14//-------------------------------------------------
15#include "TrkGeometry/TrackingGeometry.h"
16#include "TrkGeometry/TrackingVolume.h"
17#include "TrkGeometry/Layer.h"
18#include "TruthUtils/MagicNumbers.h"
19#include "GeoPrimitives/GeoPrimitivesHelpers.h"
20 // Other stuff
21#include <cmath>
22
23
24namespace Rec{
25
26
27 void NewVrtSecInclusiveTool::printWrkSet(const std::vector<WrkVrt> *WrkVrtSet, const std::string &name) const {
28 int nGoodV=0;
29 if(msgLvl(MSG::INFO)){
30 for(int iv=0; iv<(int)WrkVrtSet->size(); iv++) {
31 msg(MSG::INFO)<<name
32 <<"= "<<(*WrkVrtSet)[iv].vertex[0]
33 <<", "<<(*WrkVrtSet)[iv].vertex[1]
34 <<", "<<(*WrkVrtSet)[iv].vertex[2]
35 <<" NTrk="<<(*WrkVrtSet)[iv].selTrk.size()
36 <<" is good="<<std::boolalpha<<(*WrkVrtSet)[iv].Good<<std::noboolalpha
37 <<" Chi2="<<(*WrkVrtSet)[iv].chi2
38 <<" Mass="<<(*WrkVrtSet)[iv].vertexMom.M()
39 <<" detached="<<(*WrkVrtSet)[iv].detachedTrack
40 <<" proj.dist="<<(*WrkVrtSet)[iv].projectedVrt
41 <<" trk=";
42 for(int kk=0; kk<(int)(*WrkVrtSet)[iv].selTrk.size(); kk++) {
43 msg(MSG::INFO)<<", "<<(*WrkVrtSet)[iv].selTrk[kk];}
44 for(int kk=0; kk<(int)(*WrkVrtSet)[iv].selTrk.size(); kk++) {
45 msg(MSG::INFO)<<", "<<momAtVrt((*WrkVrtSet)[iv].trkAtVrt[kk]).Pt();}
46 msg(MSG::INFO)<<endmsg;
47 if((*WrkVrtSet)[iv].Good)nGoodV++;
48 }
49 msg(MSG::INFO)<<name<<" N="<<nGoodV<<endmsg;
50 }
51 }
52
53 /* Technicalities */
54 double NewVrtSecInclusiveTool::projSV_PV(const Amg::Vector3D & SV, const xAOD::Vertex & PV, const TLorentzVector & Direction)
55 {
56 TVector3 SV_PV( SV.x()-PV.x(), SV.y()-PV.y(), SV.z()-PV.z() );
57 return Direction.Vect().Unit()*SV_PV.Unit();
58 }
59
60 double NewVrtSecInclusiveTool::vrtVrtDist(const xAOD::Vertex & primVrt, const Amg::Vector3D & secVrt,
61 const std::vector<double>& secVrtErr, double& signif)
62
63 {
64 double distx = primVrt.x()- secVrt.x();
65 double disty = primVrt.y()- secVrt.y();
66 double distz = primVrt.z()- secVrt.z();
67
68
69 AmgSymMatrix(3) primCovMtx=primVrt.covariancePosition(); //Create
70 primCovMtx(0,0) += secVrtErr[0];
71 primCovMtx(0,1) += secVrtErr[1];
72 primCovMtx(1,0) += secVrtErr[1];
73 primCovMtx(1,1) += secVrtErr[2];
74 primCovMtx(0,2) += secVrtErr[3];
75 primCovMtx(2,0) += secVrtErr[3];
76 primCovMtx(1,2) += secVrtErr[4];
77 primCovMtx(2,1) += secVrtErr[4];
78 primCovMtx(2,2) += secVrtErr[5];
79
80 AmgSymMatrix(3) wgtMtx = primCovMtx.inverse();
81
82 signif = distx*wgtMtx(0,0)*distx
83 +disty*wgtMtx(1,1)*disty
84 +distz*wgtMtx(2,2)*distz
85 +2.*distx*wgtMtx(0,1)*disty
86 +2.*distx*wgtMtx(0,2)*distz
87 +2.*disty*wgtMtx(1,2)*distz;
88 signif=std::sqrt(std::abs(signif));
89 if( signif!=signif ) signif = 0.;
90 return std::sqrt(distx*distx+disty*disty+distz*distz);
91 }
92
93 double NewVrtSecInclusiveTool::vrtVrtDist2D(const xAOD::Vertex & primVrt, const Amg::Vector3D & secVrt,
94 const std::vector<double>& secVrtErr, double& signif)
95
96 {
97 double distx = primVrt.x()- secVrt.x();
98 double disty = primVrt.y()- secVrt.y();
99
100
101 AmgSymMatrix(3) primCovMtx=primVrt.covariancePosition(); //Create
102 AmgSymMatrix(2) covMtx;
103 covMtx(0,0) = primCovMtx(0,0) + secVrtErr[0];
104 covMtx(0,1) = primCovMtx(0,1) + secVrtErr[1];
105 covMtx(1,0) = primCovMtx(1,0) + secVrtErr[1];
106 covMtx(1,1) = primCovMtx(1,1) + secVrtErr[2];
107
108 AmgSymMatrix(2) wgtMtx = covMtx.inverse();
109
110 signif = distx*wgtMtx(0,0)*distx
111 +disty*wgtMtx(1,1)*disty
112 +2.*distx*wgtMtx(0,1)*disty;
113 signif=std::sqrt(std::abs(signif));
114 if( signif!=signif ) signif = 0.;
115 return std::sqrt(distx*distx+disty*disty);
116 }
117
118
119 double NewVrtSecInclusiveTool::vrtVrtDist(const Amg::Vector3D & vrt1, const std::vector<double> & vrtErr1,
120 const Amg::Vector3D & vrt2, const std::vector<double> & vrtErr2)
121
122 {
123 double distx = vrt1.x()- vrt2.x();
124 double disty = vrt1.y()- vrt2.y();
125 double distz = vrt1.z()- vrt2.z();
126
127 AmgSymMatrix(3) primCovMtx; //Create
128 primCovMtx(0,0) = vrtErr1[0]+vrtErr2[0];
129 primCovMtx(0,1) = primCovMtx(1,0) = vrtErr1[1]+vrtErr2[1];
130 primCovMtx(1,1) = vrtErr1[2]+vrtErr2[2];
131 primCovMtx(0,2) = primCovMtx(2,0) = vrtErr1[3]+vrtErr2[3];
132 primCovMtx(1,2) = primCovMtx(2,1) = vrtErr1[4]+vrtErr2[4];
133 primCovMtx(2,2) = vrtErr1[5]+vrtErr2[5];
134
135 AmgSymMatrix(3) wgtMtx = primCovMtx.inverse();
136
137 double signif =
138 distx*wgtMtx(0,0)*distx
139 +disty*wgtMtx(1,1)*disty
140 +distz*wgtMtx(2,2)*distz
141 +2.*distx*wgtMtx(0,1)*disty
142 +2.*distx*wgtMtx(0,2)*distz
143 +2.*disty*wgtMtx(1,2)*distz;
144 signif=std::sqrt(std::abs(signif));
145 if(signif != signif) signif = 0.;
146 return signif;
147 }
148//
149 double NewVrtSecInclusiveTool::PntPntDist(const Amg::Vector3D & Vrt1, const Amg::Vector3D & Vrt2)
150 {
151 double dx = Vrt1.x()- Vrt2.x();
152 double dy = Vrt1.y()- Vrt2.y();
153 double dz = Vrt1.z()- Vrt2.z();
154 return std::sqrt(dx*dx+dy*dy*dz*dz);
155 }
156
157//----------------------------
158// Vertex error along radius
159//----------------------------
160 double NewVrtSecInclusiveTool::vrtRadiusError(const Amg::Vector3D & SecVrt, const std::vector<double> & VrtErr)
161 {
162 double DirX=SecVrt.x(), DirY=SecVrt.y();
163 double Covar = DirX*VrtErr[0]*DirX
164 +2.*DirX*VrtErr[1]*DirY
165 +DirY*VrtErr[2]*DirY;
166 Covar /= DirX*DirX + DirY*DirY;
167 Covar=std::sqrt(std::abs(Covar));
168 if(Covar != Covar) Covar = 0.;
169 return Covar;
170 }
171
172
173 /* Invariant mass calculation for V0 decays*/
174 /* Gives correct mass assignment in case of nonequal masses*/
175
176
177 double NewVrtSecInclusiveTool::massV0(const std::vector< std::vector<double> >& TrkAtVrt,
178 double massP, double massPi )
179
180 {
181 double ap1i=std::abs(TrkAtVrt[0][2]); double ap2i=std::abs(TrkAtVrt[1][2]);
182 CxxUtils::sincos phi1(TrkAtVrt[0][0]);
183 CxxUtils::sincos theta1(TrkAtVrt[0][1]);
184 CxxUtils::sincos phi2(TrkAtVrt[1][0]);
185 CxxUtils::sincos theta2(TrkAtVrt[1][1]);
186 double px = phi1.cs * theta1.sn * ap1i
187 + phi2.cs * theta2.sn * ap2i;
188 double py = phi1.sn * theta1.sn * ap1i
189 + phi2.sn * theta2.sn * ap2i;
190 double pz = theta1.cs * ap1i
191 + theta2.cs * ap2i;
192 double ee= (ap1i > ap2i) ?
193 (std::sqrt(ap1i*ap1i+massP*massP)+std::sqrt(ap2i*ap2i+massPi*massPi)):
194 (std::sqrt(ap2i*ap2i+massP*massP)+std::sqrt(ap1i*ap1i+massPi*massPi));
195 double test=(ee-pz)*(ee+pz)-px*px-py*py;
196 return test>0 ? std::sqrt(test) : 0.;
197 }
198
199
200
201 ROOT::Math::PxPyPzEVector NewVrtSecInclusiveTool::momAtVrt(const std::vector< double >& inpTrk) const
202 {
203 double api=1./std::abs(inpTrk[2]);
204 CxxUtils::sincos phi(inpTrk[0]);
205 CxxUtils::sincos theta(inpTrk[1]);
206 double px = phi.cs * theta.sn * api;
207 double py = phi.sn * theta.sn * api;
208 double pz = theta.cs * api;
209 double ee = std::sqrt( px*px + py*py + pz*pz + m_massPi*m_massPi);
210 return {px,py,pz,ee};
211 }
212
213/*************************************************************************************************************/
214 int NewVrtSecInclusiveTool::getIBLHit(const xAOD::TrackParticle* Part)
215 {
216 uint8_t IBLhit,IBLexp;
217 if(!Part->summaryValue( IBLexp, xAOD::expectInnermostPixelLayerHit) ) IBLexp = 0;
218 if( IBLexp==0 ) return -1;
219 if(!Part->summaryValue( IBLhit, xAOD::numberOfInnermostPixelLayerHits) ) IBLhit = 0;
220 if(IBLhit) return 1;
221 else return 0;
222 }
223 int NewVrtSecInclusiveTool::getBLHit(const xAOD::TrackParticle* Part)
224 {
225 uint8_t BLhit,BLexp;
226 if(!Part->summaryValue( BLexp, xAOD::expectNextToInnermostPixelLayerHit) ) BLexp = 0;
227 if( BLexp==0 ) return -1;
228 if(!Part->summaryValue( BLhit, xAOD::numberOfNextToInnermostPixelLayerHits) ) BLhit = 0;
229 if(BLhit) return 1;
230 else return 0;
231 }
232
233 void NewVrtSecInclusiveTool::getPixelDiscs(const xAOD::TrackParticle* Part, int &d0Hit, int &d1Hit, int &d2Hit)
234 {
235 uint32_t HitPattern=Part->hitPattern();
236 d0Hit=0; if( HitPattern&((1<<Trk::pixelEndCap0)) ) d0Hit=1;
237 d1Hit=0; if( HitPattern&((1<<Trk::pixelEndCap1)) ) d1Hit=1;
238 d2Hit=0; if( HitPattern&((1<<Trk::pixelEndCap2)) ) d2Hit=1;
239 }
240/*************************************************************************************************************/
241
242 int NewVrtSecInclusiveTool::getIdHF(const xAOD::TrackParticle* TP ) const {
243 static const SG::ConstAccessor< ElementLink< xAOD::TruthParticleContainer> >
244 truthParticleLinkAcc ( "truthParticleLink");
245 if( truthParticleLinkAcc.isAvailable(*TP) ) {
246 const ElementLink<xAOD::TruthParticleContainer>& tplink = truthParticleLinkAcc(*TP);
247 if( !tplink.isValid() ) return 0;
248 static const SG::ConstAccessor< float >truthMatchProbabilityAcc ( "truthMatchProbability" );
249 if( truthMatchProbabilityAcc( *TP ) < 0.75 ) return 0;
250 if( HepMC::is_simulation_particle((*tplink))) return 0;
251 //-------------
252 const xAOD::TruthParticle * parent1=getPreviousParent(*tplink);
253 if(!parent1 || !parent1->isHadron()) return 0;
254 if(parent1->isBottomHadron()) return isExcitedHadron(parent1) ? 0 : 1; //Ground state B-hadron
255 if(parent1->isCharmHadron() && !isExcitedHadron(parent1))return 1; //Ground state C-hadron
256 //--------------
257 const xAOD::TruthParticle * parent2=getPreviousParent(parent1);
258 if(!parent2 || !parent2->isHadron()) return 0;
259 if(parent2->isBottomHadron()) return isExcitedHadron(parent2) ? 0 : 1; //Ground state B-hadron
260 if(parent2->isCharmHadron() && !isExcitedHadron(parent2))return 1; //Ground state C-hadron
261 //--------------
262 const xAOD::TruthParticle * parent3=getPreviousParent(parent2);
263 if(!parent3 || !parent3->isHadron()) return 0;
264 if(parent3->isBottomHadron()) return isExcitedHadron(parent3) ? 0 : 1; //Ground state B-hadron
265 if(parent3->isCharmHadron() && !isExcitedHadron(parent3))return 1; //Ground state C-hadron
266 }
267 return 0;
268 }
269 const xAOD::TruthParticle * NewVrtSecInclusiveTool::getPreviousParent(const xAOD::TruthParticle * child) {
270 if( child->hasProdVtx() ){
271 if( child->prodVtx()->nIncomingParticles()==1 ){
272 return *(child->prodVtx()->incomingParticleLinks())[0];
273 }
274 }
275 return nullptr;
276 }
277 bool NewVrtSecInclusiveTool::isExcitedHadron(const xAOD::TruthParticle * tp) {
278 if( !tp->hasProdVtx() )return false;
279 if( !tp->hasDecayVtx() )return false;
280 Amg::Vector3D pvrt(tp->prodVtx()->x(),tp->prodVtx()->y(),tp->prodVtx()->z());
281 Amg::Vector3D dvrt(tp->decayVtx()->x(),tp->decayVtx()->y(),tp->decayVtx()->z());
282 return Amg::distance(pvrt,dvrt) < 0.001*Gaudi::Units::mm ? true : false;
283 }
284
285 bool NewVrtSecInclusiveTool::isDisplaced(const xAOD::TrackParticle * TP) {
286 if(!TP)return false;
287 static const SG::ConstAccessor<ElementLink<xAOD::TruthParticleContainer>> truthParticleLinkAcc( "truthParticleLink");
288 if( truthParticleLinkAcc.isAvailable( *TP ) ) {
289 const ElementLink<xAOD::TruthParticleContainer>& tplink = truthParticleLinkAcc( *TP );
290 if( !tplink.isValid() ) return false;
291 if((*tplink)->hasProdVtx() && (*tplink)->prodVtx()->perp() > 1.*Gaudi::Units::mm) return true;
292 }
293 return false;
294 }
295
296 int NewVrtSecInclusiveTool::getG4Inter(const xAOD::TrackParticle* TP ) {
297 static const SG::ConstAccessor< ElementLink< xAOD::TruthParticleContainer> >
298 truthParticleLinkAcc ( "truthParticleLink");
299 if( truthParticleLinkAcc.isAvailable( *TP ) ) {
300 const ElementLink<xAOD::TruthParticleContainer>& tplink =
301 truthParticleLinkAcc( *TP );
302 if( tplink.isValid() && HepMC::is_simulation_particle((*tplink))) return 1;
303 }
304 return 0;
305 }
306 int NewVrtSecInclusiveTool::getMCPileup(const xAOD::TrackParticle* TP ) {
307 static const SG::ConstAccessor< ElementLink< xAOD::TruthParticleContainer> >
308 truthParticleLinkAcc ( "truthParticleLink");
309 if( truthParticleLinkAcc.isAvailable( *TP ) ) {
310 const ElementLink<xAOD::TruthParticleContainer>& tplink =
311 truthParticleLinkAcc( *TP );
312 if( !tplink.isValid() ) return 1;
313 } else { return 1; }
314 return 0;
315 }
316
317 int NewVrtSecInclusiveTool::getProdVrtBarcode(const xAOD::TrackParticle * TP, float resolLimit) // FIXME barcode-based
318 {
319 int barVrt{HepMC::UNDEFINED_ID};
320 if(!TP)return barVrt;
321 static const SG::ConstAccessor<ElementLink<xAOD::TruthParticleContainer>> truthParticleLinkAcc( "truthParticleLink");
322 if( truthParticleLinkAcc.isAvailable( *TP ) ) {
323 barVrt=-1;
324 const ElementLink<xAOD::TruthParticleContainer>& tplink = truthParticleLinkAcc( *TP );
325 if( !tplink.isValid() ) return barVrt;
326 const xAOD::TruthParticle *tparticle = (*tplink); // truth particle for TrackPartile
327 if(tparticle->hasProdVtx()){ // truth particle has production vertex
328 barVrt = HepMC::uniqueID(tparticle->prodVtx());
329 if(HepMC::is_simulation_vertex(tparticle->prodVtx())) return barVrt; // Geant4 vertex
330 if(tparticle->prodVtx()->nIncomingParticles()!=1) return barVrt; // Safety!
331 const xAOD::TruthParticle *parent = *(tparticle->prodVtx()->incomingParticleLinks())[0];
332 Amg::Vector3D vpos0(tparticle->prodVtx()->x(),tparticle->prodVtx()->y(),tparticle->prodVtx()->z()); //Truth vertex position
333 //
334 //-- Now check parent
335 if(!parent->isHadron()) return barVrt; // Parent is parton not particle! Stop.
336 if(!parent->hasProdVtx()) return barVrt; // Parent particle doesn't have production vertex
337 Amg::Vector3D vpos1(parent->prodVtx()->x(),parent->prodVtx()->y(),parent->prodVtx()->z()); //Truth vertex position
338 if( Amg::distance(vpos0,vpos1) > resolLimit) return barVrt; // Parent vertex is far
339 barVrt = HepMC::uniqueID(parent->prodVtx()); // Else use parent vertex as reference
340 if(HepMC::is_simulation_vertex(parent->prodVtx())) return barVrt; // Geant4 vertex
341 if(parent->prodVtx()->nIncomingParticles()!=1) return barVrt; // Not a decay vertex
342 const xAOD::TruthParticle *grandparent = *(parent->prodVtx()->incomingParticleLinks())[0];
343 //
344 //-- Now check grandparent
345 if(!grandparent->isHadron()) return barVrt; // Grandparent is parton not particle! Stop.
346 if(!grandparent->hasProdVtx()) return barVrt; // Parent particle doesn't have production vertex
347 Amg::Vector3D vpos2(grandparent->prodVtx()->x(),grandparent->prodVtx()->y(),grandparent->prodVtx()->z()); //Truth vertex position
348 if( Amg::distance(vpos0,vpos2) > resolLimit) return barVrt; // Grandparent vertex is far
349 barVrt = HepMC::uniqueID(grandparent->prodVtx()); // Use grandparent vertex as reference
350 if(HepMC::is_simulation_vertex(grandparent->prodVtx())) return barVrt; // Geant4 vertex
351 if(grandparent->prodVtx()->nIncomingParticles()!=1) return barVrt; // Not a decay vertex
352 const xAOD::TruthParticle *biggrandparent = *(grandparent->prodVtx()->incomingParticleLinks())[0];
353 //
354 //-- Now check biggrandparent
355 if(!biggrandparent->isHadron()) return barVrt; // BigGrandparent is parton not particle! Stop.
356 if(!biggrandparent->hasProdVtx()) return barVrt; // Parent particle doesn't have production vertex
357 Amg::Vector3D vpos3(biggrandparent->prodVtx()->x(),biggrandparent->prodVtx()->y(),biggrandparent->prodVtx()->z());
358 if( Amg::distance(vpos0,vpos3) > resolLimit) return barVrt; // Grandparent vertex is far
359 barVrt = HepMC::uniqueID(biggrandparent->prodVtx()); // Use grandparent vertex as reference
360 }
361 }
362 return barVrt;
363 }
364
365 bool NewVrtSecInclusiveTool::checkTrue2TrVrt(const xAOD::TrackParticle * TP1, const xAOD::TrackParticle * TP2, float nearCut)
366 {
367 if( !TP1 || !TP2 )return false;
368 static const SG::ConstAccessor<ElementLink<xAOD::TruthParticleContainer>> truthParticleLinkAcc( "truthParticleLink");
369 if( !truthParticleLinkAcc.isAvailable(*TP1) ) return false;
370 if( !truthParticleLinkAcc.isAvailable(*TP2) ) return false;
371 const ElementLink<xAOD::TruthParticleContainer>& tplink1 = truthParticleLinkAcc( *TP1 );
372 if( !tplink1.isValid() ) return false;
373 const ElementLink<xAOD::TruthParticleContainer>& tplink2 = truthParticleLinkAcc( *TP2 );
374 if( !tplink2.isValid() ) return false;
375 const xAOD::TruthParticle *tpart1 = (*tplink1); // truth particle for TrackPartile1
376 const xAOD::TruthParticle *tpart2 = (*tplink2); // truth particle for TrackPartile2
377 //-----
378 if(!tpart1->hasProdVtx()) return false;
379 if(!tpart2->hasProdVtx()) return false;
380 if(tpart1->prodVtx()->nOutgoingParticles()<2) return false;
381 if(tpart2->prodVtx()->nOutgoingParticles()<2) return false;
382 //-----
383 Amg::Vector3D vpos1(tpart1->prodVtx()->x(),tpart1->prodVtx()->y(),tpart1->prodVtx()->z());
384 Amg::Vector3D vpos2(tpart2->prodVtx()->x(),tpart2->prodVtx()->y(),tpart2->prodVtx()->z());
385 return Amg::distance(vpos1,vpos2)<nearCut ? true : false;
386 }
387
388
389 double NewVrtSecInclusiveTool::distToMatLayerSignificance(Vrt2Tr & Vrt) const
390 {
391 const EventContext& ctx = Gaudi::Hive::currentContext();
392 if(Vrt.fitVertex.perp()<20.) return 1.e9;
393 double normP=1./Vrt.momentum.P();
394 Amg::Vector3D momentumP(Vrt.momentum.Px()*normP,Vrt.momentum.Py()*normP,Vrt.momentum.Pz()*normP);
395 Amg::Vector3D momentumN=-momentumP;
396
397 const Trk::Layer * someLayer = nullptr;
398 const Trk::Layer * nextLayerP = nullptr;
399 const Trk::Layer * nextLayerN = nullptr;
400 const auto *volume = m_extrapolator->trackingGeometry()->lowestTrackingVolume(Vrt.fitVertex);
401 someLayer = volume->associatedLayer(Vrt.fitVertex);
402 const auto *material = someLayer->layerMaterialProperties();
403 if(material){
404 nextLayerP=someLayer;
405 } else {
406 nextLayerP = someLayer->nextLayer(Vrt.fitVertex,momentumP);
407 if(nextLayerP){ if(!nextLayerP->layerMaterialProperties())nextLayerP=nullptr; }
408 nextLayerN = someLayer->nextLayer(Vrt.fitVertex,momentumN);
409 if(nextLayerN){ if(!nextLayerN->layerMaterialProperties())nextLayerN=nullptr; }
410 }
411 momentumP *= 1.e5; //100GeV to have straight trajectory
412 double charge = 1.;
413 const Trk::Perigee pseudoVrtPart(Vrt.fitVertex, momentumP, charge, Vrt.fitVertex);
414
415 const Trk::TrackParameters * extrapParP=nullptr; //along momentum
416 const Trk::TrackParameters * extrapParN=nullptr; //backward
417 if(nextLayerP){ extrapParP = m_extrapolator->extrapolate(ctx, pseudoVrtPart,
418 nextLayerP->surfaceRepresentation(), Trk::anyDirection, false, Trk::nonInteractingMuon).release();}
419 if(nextLayerN){ extrapParN = m_extrapolator->extrapolate(ctx, pseudoVrtPart,
420 nextLayerN->surfaceRepresentation(), Trk::anyDirection, false, Trk::nonInteractingMuon).release();}
421
422 float distanceP=1.e9, distanceN=1.e9;
423 if(extrapParP)distanceP=PntPntDist(extrapParP->position(), Vrt.fitVertex);
424 if(extrapParN)distanceN=PntPntDist(extrapParN->position(), Vrt.fitVertex);
425 if(distanceP==1.e9 && distanceN==1.e9) return 1.e9;
426
427 //std::pair<const Trk::TrackParameters*,const Trk::Layer*> next=
428 // m_extrapolator->extrapolateToNextActiveLayer(pseudoVrtPart,Trk::anyDirection,true,Trk::pion) ;
429
430 double signif=1.e9;
431 std::vector<double> pntCovar={1.e-2,0.,1.e-2,0.,0.,4.e-2};
432 if(distanceP<distanceN)signif=vrtVrtDist(Vrt.fitVertex, Vrt.errorMatrix, extrapParP->position(), pntCovar);
433 else signif=vrtVrtDist(Vrt.fitVertex, Vrt.errorMatrix, extrapParN->position(), pntCovar);
434 delete extrapParP;
435 delete extrapParN;
436 return signif;
437 }
438
439 std::vector<double> NewVrtSecInclusiveTool::estimVrtPos( int nTrk, std::deque<long int> &selTrk, std::map<long int, std::vector<double>> & vrt)
440 {
441 std::vector<double> estimation(3,0.);
442 int ntsel=selTrk.size();
443 for( int i=0; i<ntsel-1; i++){
444 for( int j=i+1; j<ntsel; j++){
445 int k = selTrk[i]<selTrk[j] ? selTrk[i]*nTrk+selTrk[j] : selTrk[j]*nTrk+selTrk[i];
446 estimation[0]+=vrt.at(k)[0];
447 estimation[1]+=vrt[k][1];
448 estimation[2]+=vrt[k][2];
449 } }
450 estimation[0] /= ntsel*(ntsel-1)/2;
451 estimation[1] /= ntsel*(ntsel-1)/2;
452 estimation[2] /= ntsel*(ntsel-1)/2;
453 return estimation;
454 }
455
456} //end namespace
phi
Scalar phi() const
phi method
Definition AmgMatrixBasePlugin.h:67
theta
Scalar theta() const
theta method
Definition AmgMatrixBasePlugin.h:75
endmsg
#define endmsg
Definition AnalysisConfig_Ntuple.cxx:63
charge
double charge(const T &p)
Definition AtlasPID.h:997
ConstAccessor.h
Helper class to provide constant type-safe access to aux data.
AmgSymMatrix
#define AmgSymMatrix(dim)
Definition EventPrimitives.h:50
AthCommonMsg< AlgTool >::msgLvl
bool msgLvl(const MSG::Level lvl) const
Definition AthCommonMsg.h:30
AthCommonMsg< AlgTool >::msg
MsgStream & msg() const
Definition AthCommonMsg.h:24
Rec::NewVrtSecInclusiveTool::m_extrapolator
ToolHandle< Trk::IExtrapolator > m_extrapolator
Definition NewVrtSecInclusiveTool.h:184
Rec::NewVrtSecInclusiveTool::vrtVrtDist2D
static double vrtVrtDist2D(const xAOD::Vertex &primVrt, const Amg::Vector3D &secVrt, const std::vector< double > &vrtErr, double &signif)
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:93
Rec::NewVrtSecInclusiveTool::getProdVrtBarcode
static int getProdVrtBarcode(const xAOD::TrackParticle *TP, float resolLimit=0.1)
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:317
Rec::NewVrtSecInclusiveTool::PntPntDist
static double PntPntDist(const Amg::Vector3D &Vrt1, const Amg::Vector3D &Vrt2)
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:149
Rec::NewVrtSecInclusiveTool::projSV_PV
static double projSV_PV(const Amg::Vector3D &SV, const xAOD::Vertex &PV, const TLorentzVector &Direction)
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:54
Rec::NewVrtSecInclusiveTool::getPixelDiscs
static void getPixelDiscs(const xAOD::TrackParticle *Part, int &d0Hit, int &d1Hit, int &d2Hit)
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:233
Rec::NewVrtSecInclusiveTool::massV0
static double massV0(const std::vector< std::vector< double > > &TrkAtVrt, double massP, double massPi)
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:177
Rec::NewVrtSecInclusiveTool::momAtVrt
ROOT::Math::PxPyPzEVector momAtVrt(const std::vector< double > &inpTrk) const
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:201
Rec::NewVrtSecInclusiveTool::printWrkSet
void printWrkSet(const std::vector< WrkVrt > *WrkSet, const std::string &name) const
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:27
Rec::NewVrtSecInclusiveTool::vrtVrtDist
static double vrtVrtDist(const xAOD::Vertex &primVrt, const Amg::Vector3D &secVrt, const std::vector< double > &vrtErr, double &signif)
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:60
Rec::NewVrtSecInclusiveTool::getPreviousParent
static const xAOD::TruthParticle * getPreviousParent(const xAOD::TruthParticle *child)
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:269
Rec::NewVrtSecInclusiveTool::estimVrtPos
static std::vector< double > estimVrtPos(int nTrk, std::deque< long int > &selTrk, std::map< long int, std::vector< double > > &vrt)
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:439
Rec::NewVrtSecInclusiveTool::vrtRadiusError
static double vrtRadiusError(const Amg::Vector3D &secVrt, const std::vector< double > &vrtErr)
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:160
Rec::NewVrtSecInclusiveTool::checkTrue2TrVrt
static bool checkTrue2TrVrt(const xAOD::TrackParticle *TP1, const xAOD::TrackParticle *TP2, float nearCut=0.1)
Definition Reconstruction/VKalVrt/NewVrtSecInclusiveTool/src/Utilities.cxx:365
SG::ConstAccessor
Helper class to provide constant type-safe access to aux data.
Definition ConstAccessor.h:55
SG::ConstAccessor::isAvailable
bool isAvailable(const ELT &e) const
Test to see if this variable exists in the store.
Trk::Layer
Base Class for a Detector Layer in the Tracking realm.
Definition Layer.h:72
Trk::Layer::surfaceRepresentation
virtual const Surface & surfaceRepresentation() const =0
Transforms the layer into a Surface representation for extrapolation.
Trk::Layer::nextLayer
const Layer * nextLayer(const Amg::Vector3D &gp, const Amg::Vector3D &udir) const
getting the next/previous Layer if registered - unit for direction vector required
Definition Layer.cxx:161
Trk::Layer::layerMaterialProperties
const LayerMaterialProperties * layerMaterialProperties() const
getting the LayerMaterialProperties including full/pre/post update
Trk::ParametersBase::position
const Amg::Vector3D & position() const
Access method for the position.
xAOD::TruthParticle_v1::isBottomHadron
bool isBottomHadron() const
Determine if the PID is that of a b-hadron.
xAOD::TruthParticle_v1::isHadron
bool isHadron() const
Whether the particle is a hadron.
xAOD::TruthParticle_v1::hasProdVtx
bool hasProdVtx() const
Check for a production vertex on this particle.
Definition TruthParticle_v1.cxx:69
xAOD::TruthParticle_v1::prodVtx
const TruthVertex_v1 * prodVtx() const
The production vertex of this particle.
Definition TruthParticle_v1.cxx:75
xAOD::TruthParticle_v1::isCharmHadron
bool isCharmHadron() const
Determine if the PID is that of a c-hadron.
xAOD::TruthVertex_v1::z
float z() const
Vertex longitudinal distance along the beam line form the origin.
xAOD::TruthVertex_v1::y
float y() const
Vertex y displacement.
xAOD::TruthVertex_v1::nOutgoingParticles
size_t nOutgoingParticles() const
Get the number of outgoing particles.
xAOD::TruthVertex_v1::nIncomingParticles
size_t nIncomingParticles() const
Get the number of incoming particles.
xAOD::TruthVertex_v1::x
float x() const
Vertex x displacement.
xAOD::Vertex_v1::z
float z() const
Returns the z position.
xAOD::Vertex_v1::y
float y() const
Returns the y position.
xAOD::Vertex_v1::x
float x() const
Returns the x position.
Amg::distance
float distance(const Amg::Vector3D &p1, const Amg::Vector3D &p2)
calculates the distance between two point in 3D space
Definition GeoPrimitivesHelpers.h:54
Amg::Vector3D
Eigen::Matrix< double, 3, 1 > Vector3D
Definition GeoPrimitives.h:47
HepMC::is_simulation_vertex
bool is_simulation_vertex(const T &v)
Method to establish if the vertex was created during simulation (TODO migrate to be based on status).
Definition MagicNumbers.h:361
HepMC::uniqueID
int uniqueID(const T &p)
Definition MagicNumbers.h:117
HepMC::UNDEFINED_ID
constexpr int UNDEFINED_ID
Definition MagicNumbers.h:57
HepMC::is_simulation_particle
bool is_simulation_particle(const T &p)
Method to establish if a particle (or barcode) was created during the simulation (TODO update to be s...
Definition MagicNumbers.h:355
Rec
Gaudi Tools.
Definition FakeTrackBuilder.h:10
Trk::anyDirection
@ anyDirection
Definition PropDirection.h:22
Trk::Perigee
ParametersT< TrackParametersDim, Charged, PerigeeSurface > Perigee
Definition Tracking/TrkEvent/TrkParameters/TrkParameters/TrackParameters.h:33
Trk::nonInteractingMuon
@ nonInteractingMuon
Definition ParticleHypothesis.h:39
Trk::TrackParameters
ParametersBase< TrackParametersDim, Charged > TrackParameters
Definition Tracking/TrkEvent/TrkParameters/TrkParameters/TrackParameters.h:27
Trk::pixelEndCap0
@ pixelEndCap0
three pixel discs (on each side)
Definition Tracking/TrkEvent/TrkTrackSummary/TrkTrackSummary/TrackSummary.h:237
xAOD::TrackParticle
TrackParticle_v1 TrackParticle
Reference the current persistent version:
Definition Event/xAOD/xAODTracking/xAODTracking/TrackParticle.h:13
xAOD::Vertex
Vertex_v1 Vertex
Define the latest version of the vertex class.
Definition Event/xAOD/xAODTracking/xAODTracking/Vertex.h:16
xAOD::TruthParticle
TruthParticle_v1 TruthParticle
Typedef to implementation.
Definition Event/xAOD/xAODTruth/xAODTruth/TruthParticle.h:15
xAOD::expectInnermostPixelLayerHit
@ expectInnermostPixelLayerHit
Do we expect a 0th-layer barrel hit for this track?
Definition TrackingPrimitives.h:237
xAOD::numberOfNextToInnermostPixelLayerHits
@ numberOfNextToInnermostPixelLayerHits
these are the hits in the 1st pixel barrel layer
Definition TrackingPrimitives.h:249
xAOD::expectNextToInnermostPixelLayerHit
@ expectNextToInnermostPixelLayerHit
Do we expect a 1st-layer barrel hit for this track?
Definition TrackingPrimitives.h:248
xAOD::numberOfInnermostPixelLayerHits
@ numberOfInnermostPixelLayerHits
these are the hits in the 0th pixel barrel layer
Definition TrackingPrimitives.h:238
sincos.h
Helper to simultaneously calculate sin and cos of the same angle.
CxxUtils::sincos
Helper to simultaneously calculate sin and cos of the same angle.
Definition sincos.h:39
Generated on for ATLAS Offline Software by doxygen 1.14.0