170 {
173 return StatusCode::FAILURE;
174 }
175
176 constexpr int topoN = 3;
179 return StatusCode::FAILURE;
180 }
181 std::array<SG::WriteHandle<xAOD::VertexContainer>, topoN> VtxWriteHandles; int ikey(0);
182 for(
const SG::WriteHandleKey<xAOD::VertexContainer>& key :
m_outputsKeys) {
183 VtxWriteHandles[ikey] = SG::WriteHandle<xAOD::VertexContainer>(key, ctx);
184 ATH_CHECK( VtxWriteHandles[ikey].record(std::make_unique<xAOD::VertexContainer>(), std::make_unique<xAOD::VertexAuxContainer>()) );
185 ikey++;
186 }
187
188
189
190
193 if (pvContainer.cptr()->size()==0) {
194 ATH_MSG_WARNING(
"You have no primary vertices: " << pvContainer.cptr()->size());
195 return StatusCode::RECOVERABLE;
196 }
197
198
201
202 std::vector<const xAOD::TrackParticle*> tracksJpsi1;
203 std::vector<const xAOD::TrackParticle*> tracksDiTrk1;
204 std::vector<const xAOD::TrackParticle*> tracksPsi1;
205 std::vector<const xAOD::TrackParticle*> tracksJpsi2;
206 std::vector<const xAOD::TrackParticle*> tracksDiTrk2;
207 std::vector<const xAOD::TrackParticle*> tracksPsi2;
208 std::vector<const xAOD::TrackParticle*> inputTracks;
209 std::vector<double> massesPsi1;
214 std::vector<double> massesPsi2;
219 std::vector<double> massesInputTracks;
220 massesInputTracks.reserve(massesPsi1.size() + massesPsi2.size());
221 for(auto mass : massesPsi1) massesInputTracks.push_back(mass);
222 for(auto mass : massesPsi2) massesInputTracks.push_back(mass);
223
224
229
230
231 std::vector<const xAOD::Vertex*> selectedPsi1Candidates;
232 for(auto vxcItr=psi1Container->cbegin(); vxcItr!=psi1Container->cend(); ++vxcItr) {
233
238 if(flagAcc.isAvailable(*vtx) && flagAcc(*vtx)) {
240 }
241 }
243
244
246 double mass_psi1 =
m_V0Tools->invariantMass(*vxcItr, massesPsi1);
247 if (mass_psi1 < m_psi1MassLower || mass_psi1 >
m_psi1MassUpper)
continue;
248 }
249
250 TLorentzVector p4_mu1, p4_mu2;
251 p4_mu1.SetPtEtaPhiM( (*vxcItr)->trackParticle(0)->pt(),
252 (*vxcItr)->trackParticle(0)->eta(),
254 p4_mu2.SetPtEtaPhiM( (*vxcItr)->trackParticle(1)->pt(),
255 (*vxcItr)->trackParticle(1)->eta(),
257 double mass_jpsi1 = (p4_mu1 + p4_mu2).M();
258 if (mass_jpsi1 < m_jpsi1MassLower || mass_jpsi1 >
m_jpsi1MassUpper)
continue;
259
261 TLorentzVector p4_trk1, p4_trk2;
262 p4_trk1.SetPtEtaPhiM( (*vxcItr)->trackParticle(2)->pt(),
263 (*vxcItr)->trackParticle(2)->eta(),
265 p4_trk2.SetPtEtaPhiM( (*vxcItr)->trackParticle(3)->pt(),
266 (*vxcItr)->trackParticle(3)->eta(),
268 double mass_diTrk1 = (p4_trk1 + p4_trk2).M();
270 }
271
272 selectedPsi1Candidates.push_back(*vxcItr);
273 }
274
275
276 std::vector<const xAOD::Vertex*> selectedPsi2Candidates;
277 for(auto vxcItr=psi2Container->cbegin(); vxcItr!=psi2Container->cend(); ++vxcItr) {
278
283 if(flagAcc.isAvailable(*vtx) && flagAcc(*vtx)) {
285 }
286 }
288
289
291 double mass_psi2 =
m_V0Tools->invariantMass(*vxcItr,massesPsi2);
292 if(mass_psi2 < m_psi2MassLower || mass_psi2 >
m_psi2MassUpper)
continue;
293 }
294
295 TLorentzVector p4_mu1, p4_mu2;
296 p4_mu1.SetPtEtaPhiM( (*vxcItr)->trackParticle(0)->pt(),
297 (*vxcItr)->trackParticle(0)->eta(),
299 p4_mu2.SetPtEtaPhiM( (*vxcItr)->trackParticle(1)->pt(),
300 (*vxcItr)->trackParticle(1)->eta(),
302 double mass_jpsi2 = (p4_mu1 + p4_mu2).M();
303 if (mass_jpsi2 < m_jpsi2MassLower || mass_jpsi2 >
m_jpsi2MassUpper)
continue;
304
306 TLorentzVector p4_trk1, p4_trk2;
307 p4_trk1.SetPtEtaPhiM( (*vxcItr)->trackParticle(2)->pt(),
308 (*vxcItr)->trackParticle(2)->eta(),
310 p4_trk2.SetPtEtaPhiM( (*vxcItr)->trackParticle(3)->pt(),
311 (*vxcItr)->trackParticle(3)->eta(),
313 double mass_diTrk2 = (p4_trk1 + p4_trk2).M();
315 }
316 selectedPsi2Candidates.push_back(*vxcItr);
317 }
318
319 std::vector<std::pair<const xAOD::Vertex*, const xAOD::Vertex*> > candidatePairs;
320 for(auto psi1Itr=selectedPsi1Candidates.cbegin(); psi1Itr!=selectedPsi1Candidates.cend(); ++psi1Itr) {
321 tracksPsi1.clear();
322 tracksPsi1.reserve((*psi1Itr)->nTrackParticles());
323 for(
size_t i=0;
i<(*psi1Itr)->nTrackParticles();
i++) tracksPsi1.push_back((*psi1Itr)->trackParticle(i));
324 for(auto psi2Itr=selectedPsi2Candidates.cbegin(); psi2Itr!=selectedPsi2Candidates.cend(); ++psi2Itr) {
326 for(size_t j=0; j<(*psi2Itr)->nTrackParticles(); j++) {
327 if(std::find(tracksPsi1.cbegin(), tracksPsi1.cend(), (*psi2Itr)->trackParticle(j)) != tracksPsi1.cend()) {
skip =
true;
break; }
328 }
331 for(
size_t ic=0;
ic<candidatePairs.size();
ic++) {
334 if((psi1Vertex == *psi1Itr && psi2Vertex == *psi2Itr) || (psi1Vertex == *psi2Itr && psi2Vertex == *psi1Itr)) {
skip =
true;
break; }
335 }
336 }
338 candidatePairs.push_back(std::pair<const xAOD::Vertex*, const xAOD::Vertex*>(*psi1Itr,*psi2Itr));
339 }
340 }
341
342 std::sort( candidatePairs.begin(), candidatePairs.end(), [](std::pair<const xAOD::Vertex*, const xAOD::Vertex*>
a, std::pair<const xAOD::Vertex*, const xAOD::Vertex*> b) { return a.first->chiSquared()/a.first->numberDoF()+a.second->chiSquared()/a.second->numberDoF() < b.first->chiSquared()/b.first->numberDoF()+b.second->chiSquared()/b.second->numberDoF(); } );
344 candidatePairs.erase(candidatePairs.begin()+
m_maxCandidates, candidatePairs.end());
345 }
346
347 for(
size_t ic=0;
ic<candidatePairs.size();
ic++) {
350
351 tracksPsi1.clear();
354 if (tracksPsi1.size() != massesPsi1.size()) {
355 ATH_MSG_ERROR(
"Problems with Psi1 input: number of tracks or track mass inputs is not correct!");
356 }
357 tracksPsi2.clear();
360 if (tracksPsi2.size() != massesPsi2.size()) {
361 ATH_MSG_ERROR(
"Problems with Psi2 input: number of tracks or track mass inputs is not correct!");
362 }
363
364 tracksJpsi1.clear();
367 tracksDiTrk1.clear();
371 }
372 tracksJpsi2.clear();
375 tracksDiTrk2.clear();
379 }
380
381 TLorentzVector p4_moth;
386 }
390 }
392
393 inputTracks.clear();
402
403
408 }
411 }
414 }
417 }
420 }
423 }
424
425
426 Amg::Vector3D startingPoint((psi1Vertex->
x()+psi2Vertex->
x())/2,(psi1Vertex->
y()+psi2Vertex->
y())/2,(psi1Vertex->
z()+psi2Vertex->
z())/2);
427
428
429 std::unique_ptr<xAOD::Vertex> theResult(
m_iVertexFitter->fit(inputTracks, startingPoint, *state) );
430
431 if(theResult != nullptr){
432
433 double chi2DOF = theResult->chiSquared()/theResult->numberDoF();
435 if(chi2CutPassed) {
437 for(
size_t i=0;
i<theResult->trackParticleLinks().
size();
i++) {
440 tpLinkVector.push_back( mylink );
441 }
442 theResult->clearTracks();
443 theResult->setTrackParticleLinks( tpLinkVector );
444
446 xAOD::BPhysHelper psi1_helper(psi1Vertex_);
447 psi1_helper.setRefTrks();
451 mylink.setStorableObject(*trackContainer.ptr(), true);
452 tpLinkVector_psi1.push_back( mylink );
453 }
456
458 xAOD::BPhysHelper psi2_helper(psi2Vertex_);
459 psi2_helper.setRefTrks();
463 mylink.setStorableObject(*trackContainer.ptr(), true);
464 tpLinkVector_psi2.push_back( mylink );
465 }
468
469 VtxWriteHandles[0].ptr()->push_back(psi1Vertex_);
470 VtxWriteHandles[1].ptr()->push_back(psi2Vertex_);
471
472
477 if( vertexLink1.
isValid() ) precedingVertexLinks.push_back( vertexLink1 );
481 if( vertexLink2.
isValid() ) precedingVertexLinks.push_back( vertexLink2 );
482
484 PrecedingLinksDecor(*theResult.get()) = precedingVertexLinks;
485
486 xAOD::BPhysHypoHelper vtx(
m_hypoName, theResult.get());
487 vtx.setRefTrks();
488 vtx.setPass(true);
489
490 VtxWriteHandles[2].ptr()->push_back(theResult.release());
491 }
492 }
493 }
494
499
502 ATH_CHECK( refPvContainer.record(std::make_unique<xAOD::VertexContainer>(), std::make_unique<xAOD::VertexAuxContainer>()) );
503
504 if(VtxWriteHandles[2]->
size()>0) {
505 for(
int i=0;
i<topoN;
i++) {
508 ATH_MSG_FATAL(
"FillCandwithRefittedVertices failed - check the vertices you passed");
510 }
511 }
512 }
513 }
514 else {
515 if(VtxWriteHandles[2]->
size()>0) {
516 for(
int i=0;
i<topoN;
i++) {
519 ATH_MSG_FATAL(
"FillCandExistingVertices failed - check the vertices you passed");
521 }
522 }
523 }
524 }
525
526 return StatusCode::SUCCESS;
527 }
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_WARNING(x)
bool setElement(ElementType element)
Set to point to an element.
bool setStorableObject(BaseConstReference data, bool replace=false, IProxyDict *sg=0)
Set link to point to a new container (storable).
bool isValid() const
Test to see if the link can be dereferenced.
SG::Decorator< T, ALLOC > Decorator
SG::Accessor< T, ALLOC > Accessor
virtual double phi() const override final
The azimuthal angle ( ) of the particle (has range to .)
virtual double pt() const override final
The transverse momentum ( ) of the particle.
virtual double eta() const override final
The pseudorapidity ( ) of the particle.
float z() const
Returns the z position.
void setTrackParticleLinks(const TrackParticleLinks_t &trackParticles)
Set all track particle links at once.
size_t nTrackParticles() const
Get the number of tracks associated with this vertex.
void clearTracks()
Remove all tracks from the vertex.
const TrackParticleLinks_t & trackParticleLinks() const
Get all the particles associated with the vertex.
const TrackParticle * trackParticle(size_t i) const
Get the pointer to a given track that was used in vertex reco.
float y() const
Returns the y position.
float x() const
Returns the x position.
Eigen::Matrix< double, 3, 1 > Vector3D
ElementLink< xAOD::TrackParticleContainer > TrackParticleLink
ElementLink< xAOD::VertexContainer > VertexLink
std::vector< VertexLink > VertexLinkVector
std::vector< TrackParticleLink > TrackParticleLinkVector
::StatusCode StatusCode
StatusCode definition for legacy code.
const float SC[NC]
Cross sections for Carbon.
void sort(typename DataModel_detail::iterator< DVL > beg, typename DataModel_detail::iterator< DVL > end)
Specialization of sort for DataVector/List.
Vertex_v1 Vertex
Define the latest version of the vertex class.