318 typedef std::map<const Track *, TransientTrack> TransientTrackMap;
337 if (groomedFatJetsHandle->size() > fatJetsHandle->size())
339 <<
"There are more groomed (" << groomedFatJetsHandle->size() <<
") than original fat jets (" 340 << fatJetsHandle->size() <<
"). Please check that the two jet collections belong to each other.";
348 unsigned int bsCovSrc[7] = {
351 double sigmaZ = 0.0, beamWidth = 0.0;
355 bsCovSrc[3] = bsCovSrc[4] = bsCovSrc[5] = bsCovSrc[6] = 1;
362 bsCovSrc[0] = bsCovSrc[1] = 2;
363 bsCovSrc[3] = bsCovSrc[4] = bsCovSrc[5] = 1;
368 bsCovSrc[0] = bsCovSrc[1] = bsCovSrc[2] = 2;
381 std::vector<std::vector<int> > clusteredSVs(
trackIPTagInfos->size(), std::vector<int>());
384 std::vector<fastjet::PseudoJet> fjInputs;
388 std::vector<edm::Ptr<reco::Candidate> > constituents = it->getJetConstituents();
389 std::vector<edm::Ptr<reco::Candidate> >::const_iterator
m;
390 for (
m = constituents.begin();
m != constituents.end(); ++
m) {
393 edm::LogWarning(
"NullTransverseMomentum") <<
"dropping input candidate with pt=0";
396 if (it->isWeighted()) {
399 <<
"TemplatedSecondaryVertexProducer: No weights (e.g. PUPPI) given for weighted jet collection" 401 float w = (*weightsHandle)[constit];
404 fastjet::PseudoJet(constit->
px() *
w, constit->
py() *
w, constit->
pz() *
w, constit->
energy() *
w));
407 fjInputs.push_back(fastjet::PseudoJet(constit->
px(), constit->
py(), constit->
pz(), constit->
energy()));
414 std::vector<edm::Ptr<reco::Candidate> > constituents = it->jet()->getJetConstituents();
415 std::vector<edm::Ptr<reco::Candidate> >::const_iterator
m;
416 for (
m = constituents.begin();
m != constituents.end(); ++
m) {
419 edm::LogWarning(
"NullTransverseMomentum") <<
"dropping input candidate with pt=0";
422 if (it->jet()->isWeighted()) {
425 <<
"TemplatedSecondaryVertexProducer: No weights (e.g. PUPPI) given for weighted jet collection" 427 float w = (*weightsHandle)[constit];
430 fastjet::PseudoJet(constit->
px() *
w, constit->
py() *
w, constit->
pz() *
w, constit->
energy() *
w));
433 fjInputs.push_back(fastjet::PseudoJet(constit->
px(), constit->
py(), constit->
pz(), constit->
energy()));
443 fastjet::PseudoJet
p(
446 p = fastjet::PseudoJet(it->p4().px(), it->p4().py(), it->p4().pz(), it->p4().energy());
448 p.set_user_info(
new VertexInfo(it - extSecVertex->begin()));
449 fjInputs.push_back(
p);
455 std::vector<fastjet::PseudoJet> inclusiveJets = fastjet::sorted_by_pt(
fjClusterSeq->inclusive_jets(
jetPtMin));
458 if (inclusiveJets.size() < fatJetsHandle->size())
460 <<
"There are fewer reclustered (" << inclusiveJets.size() <<
") than original fat jets (" 461 << fatJetsHandle->size()
462 <<
"). Please check that the jet algorithm and jet size match those used for the original jet collection.";
465 std::vector<int> reclusteredIndices;
466 matchReclusteredJets<edm::View<reco::Jet> >(fatJetsHandle, inclusiveJets, reclusteredIndices,
"fat");
469 std::vector<int> groomedIndices;
474 std::vector<std::vector<int> > subjetIndices;
481 for (
size_t i = 0;
i < fatJetsHandle->size(); ++
i) {
482 if (reclusteredIndices.at(
i) < 0)
485 if (fatJetsHandle->at(
i).pt() == 0)
488 <<
"The original fat jet " <<
i <<
" has Pt=0. This is not expected so the jet will be skipped.";
492 if (subjetIndices.at(
i).empty())
496 if ((
std::abs(inclusiveJets.at(reclusteredIndices.at(
i)).
pt() - fatJetsHandle->at(
i).pt()) /
498 if (fatJetsHandle->at(
i).pt() < 10.)
500 <<
"The reclustered and original fat jet " <<
i <<
" have different Pt's (" 501 << inclusiveJets.at(reclusteredIndices.at(
i)).
pt() <<
" vs " << fatJetsHandle->at(
i).pt()
502 <<
" GeV, respectively).\n" 503 <<
"Please check that the jet algorithm and jet size match those used for the original fat jet " 504 "collection and also make sure the original fat jets are uncorrected. In addition, make sure you " 505 "are not using CaloJets which are presently not supported.\n" 506 <<
"Since the mismatch is at low Pt, it is ignored and only a warning is issued.\n" 507 <<
"\nIn extremely rare instances the mismatch could be caused by a difference in the machine " 508 "precision in which case make sure the original jet collection is produced and reclustering is " 509 "performed in the same job.";
512 <<
"The reclustered and original fat jet " <<
i <<
" have different Pt's (" 513 << inclusiveJets.at(reclusteredIndices.at(
i)).
pt() <<
" vs " << fatJetsHandle->at(
i).pt()
514 <<
" GeV, respectively).\n" 515 <<
"Please check that the jet algorithm and jet size match those used for the original fat jet " 516 "collection and also make sure the original fat jets are uncorrected. In addition, make sure you " 517 "are not using CaloJets which are presently not supported.\n" 518 <<
"\nIn extremely rare instances the mismatch could be caused by a difference in the machine " 519 "precision in which case make sure the original jet collection is produced and reclustering is " 520 "performed in the same job.";
524 std::vector<fastjet::PseudoJet> constituents = inclusiveJets.at(reclusteredIndices.at(
i)).constituents();
526 std::vector<int> svIndices;
528 for (std::vector<fastjet::PseudoJet>::const_iterator it = constituents.begin(); it != constituents.end();
530 if (!it->has_user_info())
533 svIndices.push_back(it->user_info<VertexInfo>().vertexIndex());
537 for (
size_t sv = 0;
sv < svIndices.size(); ++
sv) {
539 const VTX &extSV = (*extSecVertex)[svIndices.at(
sv)];
542 fastjet::PseudoJet
p(
545 p = fastjet::PseudoJet(extSV.p4().px(), extSV.p4().py(), extSV.p4().pz(), extSV.p4().energy());
547 std::vector<double> dR2toSubjets;
549 for (
size_t sj = 0; sj < subjetIndices.at(
i).size(); ++sj)
556 int closestSubjetIdx =
557 std::distance(dR2toSubjets.begin(), std::min_element(dR2toSubjets.begin(), dR2toSubjets.end()));
559 clusteredSVs.at(subjetIndices.at(
i).at(closestSubjetIdx)).
push_back(svIndices.at(
sv));
565 <<
"There are fewer reclustered (" << inclusiveJets.size() <<
") than original jets (" 567 <<
"). Please check that the jet algorithm and jet size match those used for the original jet collection.";
570 std::vector<int> reclusteredIndices;
571 matchReclusteredJets<std::vector<IPTI> >(
trackIPTagInfos, inclusiveJets, reclusteredIndices);
575 if (reclusteredIndices.at(
i) < 0)
581 <<
"The original jet " <<
i <<
" has Pt=0. This is not expected so the jet will be skipped.";
590 <<
"The reclustered and original jet " <<
i <<
" have different Pt's (" 591 << inclusiveJets.at(reclusteredIndices.at(
i)).
pt() <<
" vs " <<
trackIPTagInfos->at(
i).jet()->pt()
592 <<
" GeV, respectively).\n" 593 <<
"Please check that the jet algorithm and jet size match those used for the original jet collection " 594 "and also make sure the original jets are uncorrected. In addition, make sure you are not using " 595 "CaloJets which are presently not supported.\n" 596 <<
"Since the mismatch is at low Pt, it is ignored and only a warning is issued.\n" 597 <<
"\nIn extremely rare instances the mismatch could be caused by a difference in the machine " 598 "precision in which case make sure the original jet collection is produced and reclustering is " 599 "performed in the same job.";
602 <<
"The reclustered and original jet " <<
i <<
" have different Pt's (" 603 << inclusiveJets.at(reclusteredIndices.at(
i)).
pt() <<
" vs " <<
trackIPTagInfos->at(
i).jet()->pt()
604 <<
" GeV, respectively).\n" 605 <<
"Please check that the jet algorithm and jet size match those used for the original jet collection " 606 "and also make sure the original jets are uncorrected. In addition, make sure you are not using " 607 "CaloJets which are presently not supported.\n" 608 <<
"\nIn extremely rare instances the mismatch could be caused by a difference in the machine " 609 "precision in which case make sure the original jet collection is produced and reclustering is " 610 "performed in the same job.";
614 std::vector<fastjet::PseudoJet> constituents = inclusiveJets.at(reclusteredIndices.at(
i)).constituents();
617 for (std::vector<fastjet::PseudoJet>::const_iterator it = constituents.begin(); it != constituents.end();
619 if (!it->has_user_info())
622 clusteredSVs.at(
i).push_back(it->user_info<VertexInfo>().vertexIndex());
630 std::vector<int> groomedIndices;
635 std::vector<std::vector<int> > subjetIndices;
642 for (
size_t i = 0;
i < fatJetsHandle->size(); ++
i) {
643 if (fatJetsHandle->at(
i).pt() == 0)
646 <<
"The original fat jet " <<
i <<
" has Pt=0. This is not expected so the jet will be skipped.";
650 if (subjetIndices.at(
i).empty())
656 size_t sv = (it - extSecVertex->begin());
659 const VTX &extSV = (*extSecVertex)[
sv];
661 GlobalVector jetDir(fatJetsHandle->at(
i).px(), fatJetsHandle->at(
i).py(), fatJetsHandle->at(
i).pz());
667 fastjet::PseudoJet
p(
670 p = fastjet::PseudoJet(extSV.p4().px(), extSV.p4().py(), extSV.p4().pz(), extSV.p4().energy());
672 std::vector<double> dR2toSubjets;
674 for (
size_t sj = 0; sj < subjetIndices.at(
i).size(); ++sj)
681 int closestSubjetIdx =
682 std::distance(dR2toSubjets.begin(), std::min_element(dR2toSubjets.begin(), dR2toSubjets.end()));
684 clusteredSVs.at(subjetIndices.at(
i).at(closestSubjetIdx)).
push_back(
sv);
690 std::unique_ptr<ConfigurableVertexReconstructor>
vertexReco;
691 std::unique_ptr<GhostTrackVertexFinder> vertexRecoGT;
693 vertexRecoGT = std::make_unique<GhostTrackVertexFinder>(
702 TransientTrackMap primariesMap;
706 auto tagInfos = std::make_unique<Product>();
708 for (
typename std::vector<IPTI>::const_iterator iterJets =
trackIPTagInfos->begin();
714 const Vertex &
pv = *iterJets->primaryVertex();
716 std::set<TransientTrack> primaries;
719 TransientTrackMap::iterator
pos = primariesMap.lower_bound(iter->get());
721 if (
pos != primariesMap.end() &&
pos->first == iter->get())
722 primaries.insert(
pos->second);
725 primariesMap.insert(
pos, std::make_pair(iter->get(),
track));
726 primaries.insert(
track);
733 GlobalVector jetDir(jetRef->momentum().x(), jetRef->momentum().y(), jetRef->momentum().z());
739 const std::vector<reco::btag::TrackIPData> &ipData = iterJets->impactParameterData();
743 std::vector<TransientTrack> fitTracks;
744 std::vector<GhostTrackState> gtStates;
745 std::unique_ptr<GhostTrackPrediction> gtPred;
747 gtPred = std::make_unique<GhostTrackPrediction>(*iterJets->ghostTrack());
749 for (
unsigned int i = 0;
i <
indices.size();
i++) {
755 trackData.back().first =
indices[
i];
767 if (
pos != primariesMap.end()) {
768 primaries.erase(
pos->second);
769 fitTrack =
pos->second;
771 fitTrack = trackBuilder->
build(trackRef);
772 fitTracks.push_back(fitTrack);
779 gtState.linearize(*gtPred,
true, gtPred->lambda(
pos));
780 gtState.setWeight(ipData[
indices[
i]].ghostTrackWeight);
781 gtStates.push_back(gtState);
785 std::unique_ptr<GhostTrack> ghostTrack;
787 ghostTrack = std::make_unique<GhostTrack>(
791 GlobalVector(iterJets->ghostTrack()->px(), iterJets->ghostTrack()->py(), iterJets->ghostTrack()->pz()),
795 iterJets->ghostTrack()->chi2(),
796 iterJets->ghostTrack()->ndof());
800 std::vector<VTX> extAssoCollection;
801 std::vector<TransientVertex> fittedSVs;
802 std::vector<SecondaryVertex> SVs;
807 fittedSVs = vertexRecoGT->vertices(
pv, *ghostTrack);
814 fittedSVs = vertexRecoGT->vertices(
pv, *
beamSpot, *ghostTrack);
823 for (
unsigned int i = 0;
i < 7;
i++) {
824 unsigned int covSrc = bsCovSrc[
i];
825 for (
unsigned int j = 0;
j < 7;
j++) {
827 if (!covSrc || bsCovSrc[
j] != covSrc)
829 else if (covSrc == 1)
831 else if (
j < 3 &&
i < 3)
846 fittedSVs = vertexRecoGT->vertices(
pv,
bs, *ghostTrack);
852 std::vector<TransientTrack> primaries_(primaries.begin(), primaries.end());
854 fittedSVs = vertexRecoGT->vertices(
pv, *
beamSpot, primaries_, *ghostTrack);
861 std::remove_copy_if(boost::make_transform_iterator(fittedSVs.begin(), svBuilder),
862 boost::make_transform_iterator(fittedSVs.end(), svBuilder),
863 std::back_inserter(SVs),
870 for (
size_t iExtSv = 0; iExtSv < clusteredSVs.at(
jetIdx).size(); iExtSv++) {
871 const VTX &extVertex = (*extSecVertex)[clusteredSVs.at(
jetIdx).at(iExtSv)];
872 if (extVertex.p4().M() < 0.3)
874 extAssoCollection.push_back(extVertex);
877 for (
size_t iExtSv = 0; iExtSv < extSecVertex->size(); iExtSv++) {
878 const VTX &extVertex = (*extSecVertex)[iExtSv];
881 extVertex.p4().M() < 0.3)
883 extAssoCollection.push_back(extVertex);
888 std::remove_copy_if(boost::make_transform_iterator(extAssoCollection.begin(), svBuilder),
889 boost::make_transform_iterator(extAssoCollection.end(), svBuilder),
890 std::back_inserter(SVs),
899 extAssoCollection.clear();
905 std::vector<typename TemplatedSecondaryVertexTagInfo<IPTI, VTX>::VertexData> svData;
907 svData.resize(vtxIndices.size());
908 for (
unsigned int idx = 0;
idx < vtxIndices.size();
idx++) {
911 svData[
idx].vertex =
sv;
912 svData[
idx].dist1d =
sv.dist1d();
913 svData[
idx].dist2d =
sv.dist2d();
914 svData[
idx].dist3d =
sv.dist3d();
math::Error< dimension >::type CovarianceMatrix
reco::Vertex::Point convertPos(const GlobalPoint &p)
edm::EDGetTokenT< edm::View< reco::Jet > > token_fatJets
edm::ESGetToken< TransientTrackBuilder, TransientTrackRecord > token_trackBuilder
const math::XYZPoint & position(const reco::Vertex &sv)
T getParameter(std::string const &) const
virtual double energy() const =0
energy
TrackSelector trackSelector
virtual double pt() const =0
transverse momentum
ClusterSequencePtr fjClusterSeq
edm::EDGetTokenT< reco::BeamSpot > token_BeamSpot
virtual double pz() const =0
z coordinate of momentum vector
reco::Vertex::Error convertError(const GlobalError &ge)
constexpr bool isUninitialized() const noexcept
VertexFilter vertexFilter
reco::btag::SortCriteria sortCriterium
Log< level::Error, false > LogError
const reco::Track * toTrack(const reco::TrackBaseRef &t)
static GhostTrackVertexFinder::FitType getGhostTrackFitType(const std::string &name)
ConstraintType constraint
reco::TransientTrack build(const reco::Track *p) const
edm::EDGetTokenT< std::vector< IPTI > > token_trackIPTagInfo
IPTI::input_container::value_type input_item
bool isNull() const
Checks for null.
reco::btag::IndexedTrackData IndexedTrackData
GlobalVector flightDirection(const reco::Vertex &pv, const reco::Vertex &sv)
edm::EDGetTokenT< edm::ValueMap< float > > token_weights
Abs< T >::type abs(const T &t)
VertexSorting< SecondaryVertex > vertexSorting
virtual double py() const =0
y coordinate of momentum vector
ESHandle< T > getHandle(const ESGetToken< T, R > &iToken) const
void matchGroomedJets(const edm::Handle< edm::View< reco::Jet > > &jets, const edm::Handle< edm::View< reco::Jet > > &matchedJets, std::vector< int > &matchedIndices)
void matchSubjets(const std::vector< int > &groomedIndices, const edm::Handle< edm::View< reco::Jet > > &groomedJets, const edm::Handle< std::vector< IPTI > > &subjets, std::vector< std::vector< int > > &matchedIndices)
void markUsedTracks(TrackDataVector &trackData, const input_container &trackRefs, const SecondaryVertex &sv, size_t idx)
edm::EDGetTokenT< edm::View< VTX > > token_extSVCollection
virtual double px() const =0
x coordinate of momentum vector
edm::ParameterSet vtxRecoPSet
edm::EDGetTokenT< edm::View< reco::Jet > > token_groomedFatJets
deadvectors [0] push_back({0.0175431, 0.538005, 6.80997, 13.29})
std::vector< reco::btag::IndexedTrackData > TrackDataVector
boost::indirect_iterator< typename seq_t::const_iterator > const_iterator
Log< level::Warning, false > LogWarning
JetDefPtr fjJetDefinition
std::vector< TrackBaseRef >::const_iterator trackRef_iterator
The iteratator for the vector<TrackRef>
std::pair< unsigned int, TrackData > IndexedTrackData
Global3DVector GlobalVector
IPTI::input_container input_container
Public Types inherited from edm::stream::EDProducer<>
1.8.14