11 #include <boost/container/static_vector.hpp>
13 #include "HepPDT/ParticleDataTable.hh"
19 m_cascadeTools(cascadeTools), m_eventInfo(nullptr), m_PV_minNTracks(0),
20 m_copyAllVertices(true)
27 m_cascadeTools(cascadeTools), m_eventInfo(eventinfo), m_PV_minNTracks(0),
28 m_copyAllVertices(true)
67 const float errConst = -9999999.;
77 const std::vector<const xAOD::Vertex*> &PVlist,
78 const size_t PV_minNTracks)
const {
84 size_t size = PVlist.size();
85 double lowA0zcalc = fabs(m_cascadeTools->a0z (
mom, Obj.
vtx(), PVlist[0]));
87 if ( PVlist[
i]->nTrackParticles() >= PV_minNTracks ) {
88 double a0z = fabs(m_cascadeTools->a0z(
mom, Obj.
vtx(), PVlist[
i]));
89 if(a0z < lowA0zcalc) {
99 const std::vector<const xAOD::Vertex*> &PVlist,
100 const size_t PV_minNTracks)
const {
106 size_t size = PVlist.size();
107 double lowA0calc = m_cascadeTools->a0(
mom, Obj.
vtx(), PVlist[0]);
108 for(
size_t i =1;
i<
size;
i++) {
109 if ( PVlist[
i]->nTrackParticles() >= PV_minNTracks ) {
110 double a0 = m_cascadeTools->a0(
mom, Obj.
vtx(), PVlist[
i]);
124 std::vector<const xAOD::Vertex*> goodPrimaryVertices;
125 goodPrimaryVertices.reserve(pvContainer->
size());
129 if ( thistype == Pileupvtx || thistype == Pvtx ) {
130 goodPrimaryVertices.push_back(*
ptr);
135 return goodPrimaryVertices;
142 m_PV_minNTracks = PV_minNTracks;
147 if(m_eventInfo)
return Amg::Vector3D(m_eventInfo->beamPosX(), m_eventInfo->beamPosY(), m_eventInfo->beamPosZ());
149 static const Amg::Vector3D defaultBS(-10000.,-10000.,-10000.);
156 const std::vector<const xAOD::Vertex*> &PVlist,
157 const size_t PV_minNTracks)
const {
161 for (
size_t i = 0;
i<PVlist.size(); ++
i) {
162 if ( PVlist[
i]->nTrackParticles() >= PV_minNTracks ) {
163 double z0BA = m_cascadeTools->a0(
mom,
obj.vtx(), PVlist[
i]);
164 if (z0BA < lowZ0BAcalc) {
186 TLorentzVector totalMom;
187 unsigned int NTrk =
mom.size();
188 for(
unsigned int it=0;
it<NTrk;
it++) totalMom +=
mom[
it];
189 TVector3 totP = totalMom.Vect();
192 if (
pT.mag2() > 0 ) {
196 xDOCA = xSV - pSV*
pT.dot(xT)/
pT.mag2();
198 std::cout <<
"BPhysPVCascadeTools::DocaExtrapToBeamSpot: WARNING pT == 0."
206 auto &collection =
result->vertices();
207 for(
auto v : collection)
209 std::vector<ElementLink<DataVector<xAOD::TrackParticle> > > newLinkVector;
210 for(
unsigned int i=0;
i<
v->trackParticleLinks().
size();
i++)
214 newLinkVector.push_back( mylink );
217 v->setTrackParticleLinks( newLinkVector );
223 auto &collection =
result->vertices();
224 for(
auto v : collection) {
225 std::vector<ElementLink<DataVector<xAOD::TrackParticle> > > newLinkVector;
226 for(
auto mylink :
v->trackParticleLinks()) {
229 for(
auto col : trackCols){
235 if(foundcontainer ==
nullptr){
236 throw std::runtime_error(
"Could not find track in original containers");
238 mylink.setStorableObject(*foundcontainer,
true);
239 newLinkVector.push_back( mylink );
242 v->setTrackParticleLinks( newLinkVector );
248 std::vector<const xAOD::TrackParticle*> exclTrk;
249 for(
size_t jt=0; jt<cascadeVertices.size(); jt++) {
250 for(
size_t it=0;
it<cascadeVertices[jt]->vxTrackAtVertex().
size();
it++) {
251 if(cascadeVertices[jt]->trackParticle(
it)->
charge() != 0) exclTrk.push_back(cascadeVertices[jt]->trackParticle(
it));
268 const std::vector<xAOD::Vertex*> &cascadeVertices = casc->
vertices();
269 const bool doPt = (DoVertexType & 1) != 0;
270 const bool doA0 = (DoVertexType & 2) != 0;
271 const bool doZ0 = (DoVertexType & 4) != 0;
272 const bool doZ0BA = (DoVertexType & 8) != 0;
275 std::vector<const xAOD::TrackParticle*> exclTrk = CollectAllChargedTracks(cascadeVertices);
278 const std::vector<const xAOD::Vertex*> GoodPVs = GetGoodPV(pvContainer);
280 if (GoodPVs.empty() ==
false) {
282 size_t pVmax =
std::min((
size_t)in_PV_max, GoodPVs.size());
283 std::vector<const xAOD::Vertex*> refPVvertexes;
284 std::vector<xAOD::Vertex*> refPVvertexes_toDelete;
285 std::vector<int> exitCode;
286 refPVvertexes.reserve(pVmax);
287 refPVvertexes_toDelete.reserve(pVmax);
288 exitCode.reserve(pVmax);
292 for (
size_t i =0;
i < pVmax ;
i++) {
301 if (refPV ==
nullptr) {
302 refPVvertexes.push_back(oldPV);
303 refPVvertexes_toDelete.push_back(
nullptr);
305 refPVvertexes.push_back(refPV);
306 refPVvertexes_toDelete.push_back(refPV);
309 boost::container::static_vector<size_t, 4> indexesUsed;
310 boost::container::static_vector<std::pair<size_t, xAOD::BPhysHelper::pv_type>, 4> indexestoProcess;
313 indexestoProcess.push_back(std::make_pair
317 indexestoProcess.push_back(std::make_pair( FindLowA0Index(
mom, vtx, refPVvertexes, m_PV_minNTracks),
321 indexestoProcess.push_back(std::make_pair(FindLowZIndex(
mom, vtx, refPVvertexes, m_PV_minNTracks),
325 size_t lowZBA = FindLowZ0BAIndex(
mom, vtx, refPVvertexes, m_PV_minNTracks);
326 if( lowZBA < pVmax ) {
332 for(
size_t i =0 ;
i<indexestoProcess.size();
i++){
334 auto index = indexestoProcess[
i].first;
335 auto pvtype = indexestoProcess[
i].second;
337 (refPVvertexes_toDelete.at(
index)) ? refPvContainer : pvContainer;
338 if(ParentContainer == refPvContainer &&
std::find(indexesUsed.begin(),
339 indexesUsed.end(),
index) == indexesUsed.end()) {
342 indexesUsed.push_back(
index);
344 FillBPhysHelper(
mom,
cov, vtx, refPVvertexes[
index],
345 ParentContainer, pvtype, exitCode[
index]);
350 for(
size_t x : indexesUsed) refPVvertexes_toDelete[
x] =
nullptr;
354 refPVvertexes.clear();
355 refPVvertexes_toDelete.clear();
366 size_t lowA0 = FindLowA0Index(
mom, vtx, GoodPVs, m_PV_minNTracks);
371 size_t lowZ = FindLowZIndex(
mom, vtx, GoodPVs, m_PV_minNTracks);
376 size_t lowZBA = FindLowZ0BAIndex(
mom, vtx, GoodPVs, m_PV_minNTracks);
377 if ( lowZBA < GoodPVs.size() ) {
387 if(pvContainer->
empty())
return StatusCode::FAILURE;
419 return StatusCode::SUCCESS;
426 const std::vector< std::vector<TLorentzVector> > &moms = casc->
getParticleMoms();
427 const std::vector<xAOD::Vertex*> &cascadeVertices = casc->
vertices();
429 std::vector<float>
px;
430 std::vector<float>
py;
431 std::vector<float>
pz;
432 for(
size_t jt=0; jt<moms.size(); jt++) {
433 for(
size_t it=0;
it<cascadeVertices[jt]->vxTrackAtVertex().
size();
it++) {
434 px.push_back( moms[jt][
it].Px() );
435 py.push_back( moms[jt][
it].Py() );
436 pz.push_back( moms[jt][
it].Pz() );
466 auto precedingVerticesItr = vertices.begin();
467 for(; precedingVerticesItr!=vertices.end(); ++precedingVerticesItr) {
469 if( !(*precedingVerticesItr) )
482 precedingVertexLinks.push_back( vertexLink );
487 decor(*vert) = precedingVertexLinks;
492 auto ptr = pdt->particle( pdgcode );
493 return ptr ?
ptr->mass() : 0.;