316 typedef std::map<const Track *, TransientTrack> TransientTrackMap;
336 if (groomedFatJetsHandle->size() > fatJetsHandle->size())
338 <<
"There are more groomed (" << groomedFatJetsHandle->size() <<
") than original fat jets ("
339 << fatJetsHandle->size() <<
"). Please check that the two jet collections belong to each other.";
347 unsigned int bsCovSrc[7] = {
350 double sigmaZ = 0.0, beamWidth = 0.0;
354 bsCovSrc[3] = bsCovSrc[4] = bsCovSrc[5] = bsCovSrc[6] = 1;
361 bsCovSrc[0] = bsCovSrc[1] = 2;
362 bsCovSrc[3] = bsCovSrc[4] = bsCovSrc[5] = 1;
367 bsCovSrc[0] = bsCovSrc[1] = bsCovSrc[2] = 2;
380 std::vector<std::vector<int> > clusteredSVs(
trackIPTagInfos->size(), std::vector<int>());
383 std::vector<fastjet::PseudoJet> fjInputs;
387 std::vector<edm::Ptr<reco::Candidate> > constituents = it->getJetConstituents();
388 std::vector<edm::Ptr<reco::Candidate> >::const_iterator
m;
389 for (
m = constituents.begin();
m != constituents.end(); ++
m) {
391 if (constit->
pt() == 0) {
392 edm::LogWarning(
"NullTransverseMomentum") <<
"dropping input candidate with pt=0";
395 if (it->isWeighted()) {
398 <<
"TemplatedSecondaryVertexProducer: No weights (e.g. PUPPI) given for weighted jet collection"
400 float w = (*weightsHandle)[constit];
402 fastjet::PseudoJet(constit->
px() *
w, constit->
py() *
w, constit->
pz() *
w, constit->
energy() *
w));
404 fjInputs.push_back(fastjet::PseudoJet(constit->
px(), constit->
py(), constit->
pz(), constit->
energy()));
411 std::vector<edm::Ptr<reco::Candidate> > constituents = it->jet()->getJetConstituents();
412 std::vector<edm::Ptr<reco::Candidate> >::const_iterator
m;
413 for (
m = constituents.begin();
m != constituents.end(); ++
m) {
415 if (constit->
pt() == 0) {
416 edm::LogWarning(
"NullTransverseMomentum") <<
"dropping input candidate with pt=0";
419 if (it->jet()->isWeighted()) {
422 <<
"TemplatedSecondaryVertexProducer: No weights (e.g. PUPPI) given for weighted jet collection"
424 float w = (*weightsHandle)[constit];
426 fastjet::PseudoJet(constit->
px() *
w, constit->
py() *
w, constit->
pz() *
w, constit->
energy() *
w));
428 fjInputs.push_back(fastjet::PseudoJet(constit->
px(), constit->
py(), constit->
pz(), constit->
energy()));
438 fastjet::PseudoJet
p(
441 p = fastjet::PseudoJet(it->p4().px(), it->p4().py(), it->p4().pz(), it->p4().energy());
443 p.set_user_info(
new VertexInfo(it - extSecVertex->begin()));
444 fjInputs.push_back(
p);
450 std::vector<fastjet::PseudoJet> inclusiveJets = fastjet::sorted_by_pt(
fjClusterSeq->inclusive_jets(
jetPtMin));
453 if (inclusiveJets.size() < fatJetsHandle->size())
455 <<
"There are fewer reclustered (" << inclusiveJets.size() <<
") than original fat jets ("
456 << fatJetsHandle->size()
457 <<
"). Please check that the jet algorithm and jet size match those used for the original jet collection.";
460 std::vector<int> reclusteredIndices;
461 matchReclusteredJets<edm::View<reco::Jet> >(fatJetsHandle, inclusiveJets, reclusteredIndices,
"fat");
464 std::vector<int> groomedIndices;
469 std::vector<std::vector<int> > subjetIndices;
476 for (
size_t i = 0;
i < fatJetsHandle->size(); ++
i) {
477 if (reclusteredIndices.at(
i) < 0)
480 if (fatJetsHandle->at(
i).pt() == 0)
483 <<
"The original fat jet " <<
i <<
" has Pt=0. This is not expected so the jet will be skipped.";
487 if (subjetIndices.at(
i).empty())
491 if ((
std::abs(inclusiveJets.at(reclusteredIndices.at(
i)).
pt() - fatJetsHandle->at(
i).pt()) /
493 if (fatJetsHandle->at(
i).pt() < 10.)
495 <<
"The reclustered and original fat jet " <<
i <<
" have different Pt's ("
496 << inclusiveJets.at(reclusteredIndices.at(
i)).
pt() <<
" vs " << fatJetsHandle->at(
i).pt()
497 <<
" GeV, respectively).\n"
498 <<
"Please check that the jet algorithm and jet size match those used for the original fat jet "
499 "collection and also make sure the original fat jets are uncorrected. In addition, make sure you "
500 "are not using CaloJets which are presently not supported.\n"
501 <<
"Since the mismatch is at low Pt, it is ignored and only a warning is issued.\n"
502 <<
"\nIn extremely rare instances the mismatch could be caused by a difference in the machine "
503 "precision in which case make sure the original jet collection is produced and reclustering is "
504 "performed in the same job.";
507 <<
"The reclustered and original fat jet " <<
i <<
" have different Pt's ("
508 << inclusiveJets.at(reclusteredIndices.at(
i)).
pt() <<
" vs " << fatJetsHandle->at(
i).pt()
509 <<
" GeV, respectively).\n"
510 <<
"Please check that the jet algorithm and jet size match those used for the original fat jet "
511 "collection and also make sure the original fat jets are uncorrected. In addition, make sure you "
512 "are not using CaloJets which are presently not supported.\n"
513 <<
"\nIn extremely rare instances the mismatch could be caused by a difference in the machine "
514 "precision in which case make sure the original jet collection is produced and reclustering is "
515 "performed in the same job.";
519 std::vector<fastjet::PseudoJet> constituents = inclusiveJets.at(reclusteredIndices.at(
i)).constituents();
521 std::vector<int> svIndices;
523 for (std::vector<fastjet::PseudoJet>::const_iterator it = constituents.begin(); it != constituents.end();
525 if (!it->has_user_info())
528 svIndices.push_back(it->user_info<VertexInfo>().vertexIndex());
532 for (
size_t sv = 0;
sv < svIndices.size(); ++
sv) {
534 const VTX &extSV = (*extSecVertex)[svIndices.at(
sv)];
537 fastjet::PseudoJet
p(
540 p = fastjet::PseudoJet(extSV.p4().px(), extSV.p4().py(), extSV.p4().pz(), extSV.p4().energy());
542 std::vector<double> dR2toSubjets;
544 for (
size_t sj = 0; sj < subjetIndices.at(
i).size(); ++sj)
551 int closestSubjetIdx =
552 std::distance(dR2toSubjets.begin(), std::min_element(dR2toSubjets.begin(), dR2toSubjets.end()));
554 clusteredSVs.at(subjetIndices.at(
i).at(closestSubjetIdx)).push_back(svIndices.at(
sv));
560 <<
"There are fewer reclustered (" << inclusiveJets.size() <<
") than original jets ("
562 <<
"). Please check that the jet algorithm and jet size match those used for the original jet collection.";
565 std::vector<int> reclusteredIndices;
566 matchReclusteredJets<std::vector<IPTI> >(
trackIPTagInfos, inclusiveJets, reclusteredIndices);
570 if (reclusteredIndices.at(
i) < 0)
576 <<
"The original jet " <<
i <<
" has Pt=0. This is not expected so the jet will be skipped.";
585 <<
"The reclustered and original jet " <<
i <<
" have different Pt's ("
586 << inclusiveJets.at(reclusteredIndices.at(
i)).
pt() <<
" vs " <<
trackIPTagInfos->at(
i).jet()->pt()
587 <<
" GeV, respectively).\n"
588 <<
"Please check that the jet algorithm and jet size match those used for the original jet collection "
589 "and also make sure the original jets are uncorrected. In addition, make sure you are not using "
590 "CaloJets which are presently not supported.\n"
591 <<
"Since the mismatch is at low Pt, it is ignored and only a warning is issued.\n"
592 <<
"\nIn extremely rare instances the mismatch could be caused by a difference in the machine "
593 "precision in which case make sure the original jet collection is produced and reclustering is "
594 "performed in the same job.";
597 <<
"The reclustered and original jet " <<
i <<
" have different Pt's ("
598 << inclusiveJets.at(reclusteredIndices.at(
i)).
pt() <<
" vs " <<
trackIPTagInfos->at(
i).jet()->pt()
599 <<
" GeV, respectively).\n"
600 <<
"Please check that the jet algorithm and jet size match those used for the original jet collection "
601 "and also make sure the original jets are uncorrected. In addition, make sure you are not using "
602 "CaloJets which are presently not supported.\n"
603 <<
"\nIn extremely rare instances the mismatch could be caused by a difference in the machine "
604 "precision in which case make sure the original jet collection is produced and reclustering is "
605 "performed in the same job.";
609 std::vector<fastjet::PseudoJet> constituents = inclusiveJets.at(reclusteredIndices.at(
i)).constituents();
612 for (std::vector<fastjet::PseudoJet>::const_iterator it = constituents.begin(); it != constituents.end();
614 if (!it->has_user_info())
617 clusteredSVs.at(
i).push_back(it->user_info<VertexInfo>().vertexIndex());
625 std::vector<int> groomedIndices;
630 std::vector<std::vector<int> > subjetIndices;
637 for (
size_t i = 0;
i < fatJetsHandle->size(); ++
i) {
638 if (fatJetsHandle->at(
i).pt() == 0)
641 <<
"The original fat jet " <<
i <<
" has Pt=0. This is not expected so the jet will be skipped.";
645 if (subjetIndices.at(
i).empty())
651 size_t sv = (it - extSecVertex->begin());
654 const VTX &extSV = (*extSecVertex)[
sv];
656 GlobalVector jetDir(fatJetsHandle->at(
i).px(), fatJetsHandle->at(
i).py(), fatJetsHandle->at(
i).pz());
662 fastjet::PseudoJet
p(
665 p = fastjet::PseudoJet(extSV.p4().px(), extSV.p4().py(), extSV.p4().pz(), extSV.p4().energy());
667 std::vector<double> dR2toSubjets;
669 for (
size_t sj = 0; sj < subjetIndices.at(
i).size(); ++sj)
676 int closestSubjetIdx =
677 std::distance(dR2toSubjets.begin(), std::min_element(dR2toSubjets.begin(), dR2toSubjets.end()));
679 clusteredSVs.at(subjetIndices.at(
i).at(closestSubjetIdx)).push_back(
sv);
685 std::unique_ptr<ConfigurableVertexReconstructor>
vertexReco;
686 std::unique_ptr<GhostTrackVertexFinder> vertexRecoGT;
697 TransientTrackMap primariesMap;
701 auto tagInfos = std::make_unique<Product>();
703 for (
typename std::vector<IPTI>::const_iterator iterJets =
trackIPTagInfos->begin();
709 const Vertex &
pv = *iterJets->primaryVertex();
711 std::set<TransientTrack> primaries;
714 TransientTrackMap::iterator
pos = primariesMap.lower_bound(iter->get());
716 if (
pos != primariesMap.end() &&
pos->first == iter->get())
717 primaries.insert(
pos->second);
720 primariesMap.insert(
pos, std::make_pair(iter->get(),
track));
721 primaries.insert(
track);
728 GlobalVector jetDir(jetRef->momentum().x(), jetRef->momentum().y(), jetRef->momentum().z());
734 const std::vector<reco::btag::TrackIPData> &ipData = iterJets->impactParameterData();
738 std::vector<TransientTrack> fitTracks;
739 std::vector<GhostTrackState> gtStates;
740 std::unique_ptr<GhostTrackPrediction> gtPred;
744 for (
unsigned int i = 0;
i <
indices.size();
i++) {
750 trackData.back().first =
indices[
i];
762 if (
pos != primariesMap.end()) {
763 primaries.erase(
pos->second);
764 fitTrack =
pos->second;
766 fitTrack = trackBuilder->
build(trackRef);
767 fitTracks.push_back(fitTrack);
774 gtState.linearize(*gtPred,
true, gtPred->lambda(
pos));
775 gtState.setWeight(ipData[
indices[
i]].ghostTrackWeight);
776 gtStates.push_back(gtState);
780 std::unique_ptr<GhostTrack> ghostTrack;
786 GlobalVector(iterJets->ghostTrack()->px(), iterJets->ghostTrack()->py(), iterJets->ghostTrack()->pz()),
790 iterJets->ghostTrack()->chi2(),
791 iterJets->ghostTrack()->ndof()));
795 std::vector<VTX> extAssoCollection;
796 std::vector<TransientVertex> fittedSVs;
797 std::vector<SecondaryVertex> SVs;
802 fittedSVs = vertexRecoGT->vertices(
pv, *ghostTrack);
809 fittedSVs = vertexRecoGT->vertices(
pv, *
beamSpot, *ghostTrack);
818 for (
unsigned int i = 0;
i < 7;
i++) {
819 unsigned int covSrc = bsCovSrc[
i];
820 for (
unsigned int j = 0;
j < 7;
j++) {
822 if (!covSrc || bsCovSrc[
j] != covSrc)
824 else if (covSrc == 1)
826 else if (
j < 3 &&
i < 3)
841 fittedSVs = vertexRecoGT->vertices(
pv,
bs, *ghostTrack);
847 std::vector<TransientTrack> primaries_(primaries.begin(), primaries.end());
849 fittedSVs = vertexRecoGT->vertices(
pv, *
beamSpot, primaries_, *ghostTrack);
856 std::remove_copy_if(boost::make_transform_iterator(fittedSVs.begin(), svBuilder),
857 boost::make_transform_iterator(fittedSVs.end(), svBuilder),
858 std::back_inserter(SVs),
865 for (
size_t iExtSv = 0; iExtSv < clusteredSVs.at(
jetIdx).size(); iExtSv++) {
866 const VTX &extVertex = (*extSecVertex)[clusteredSVs.at(
jetIdx).at(iExtSv)];
867 if (extVertex.p4().M() < 0.3)
869 extAssoCollection.push_back(extVertex);
872 for (
size_t iExtSv = 0; iExtSv < extSecVertex->size(); iExtSv++) {
873 const VTX &extVertex = (*extSecVertex)[iExtSv];
876 extVertex.p4().M() < 0.3)
878 extAssoCollection.push_back(extVertex);
883 std::remove_copy_if(boost::make_transform_iterator(extAssoCollection.begin(), svBuilder),
884 boost::make_transform_iterator(extAssoCollection.end(), svBuilder),
885 std::back_inserter(SVs),
894 extAssoCollection.clear();
900 std::vector<typename TemplatedSecondaryVertexTagInfo<IPTI, VTX>::VertexData> svData;
902 svData.resize(vtxIndices.size());
903 for (
unsigned int idx = 0;
idx < vtxIndices.size();
idx++) {
906 svData[
idx].vertex =
sv;
907 svData[
idx].dist1d =
sv.dist1d();
908 svData[
idx].dist2d =
sv.dist2d();
909 svData[
idx].dist3d =
sv.dist3d();
Public Types inherited from edm::stream::EDProducer<>
1.8.16