323 typedef std::map<const Track *, TransientTrack> TransientTrackMap;
344 if( groomedFatJetsHandle->size() > fatJetsHandle->size() )
345 edm::LogError(
"TooManyGroomedJets") <<
"There are more groomed (" << groomedFatJetsHandle->size() <<
") than original fat jets (" << fatJetsHandle->size() <<
"). Please check that the two jet collections belong to each other.";
349 unsigned int bsCovSrc[7] = { 0, };
350 double sigmaZ = 0.0, beamWidth = 0.0;
354 bsCovSrc[3] = bsCovSrc[4] = bsCovSrc[5] = bsCovSrc[6] = 1;
355 sigmaZ = beamSpot->sigmaZ();
356 beamWidth = beamSpot->BeamWidthX();
361 bsCovSrc[0] = bsCovSrc[1] = 2;
362 bsCovSrc[3] = bsCovSrc[4] = bsCovSrc[5] = 1;
363 sigmaZ = beamSpot->sigmaZ();
367 bsCovSrc[0] = bsCovSrc[1] = bsCovSrc[2] = 2;
380 std::vector<std::vector<int> > clusteredSVs(trackIPTagInfos->size(),std::vector<int>());
384 std::vector<fastjet::PseudoJet> fjInputs;
390 std::vector<edm::Ptr<reco::Candidate> > constituents = it->getJetConstituents();
391 std::vector<edm::Ptr<reco::Candidate> >::const_iterator
m;
392 for( m = constituents.begin(); m != constituents.end(); ++
m )
395 if(constit->pt() == 0)
397 edm::LogWarning(
"NullTransverseMomentum") <<
"dropping input candidate with pt=0";
400 fjInputs.push_back(fastjet::PseudoJet(constit->px(),constit->py(),constit->pz(),constit->energy()));
406 for(
typename std::vector<IPTI>::const_iterator it = trackIPTagInfos->begin(); it != trackIPTagInfos->end(); ++it)
408 std::vector<edm::Ptr<reco::Candidate> > constituents = it->jet()->getJetConstituents();
409 std::vector<edm::Ptr<reco::Candidate> >::const_iterator
m;
410 for( m = constituents.begin(); m != constituents.end(); ++
m )
413 if(constit->pt() == 0)
415 edm::LogWarning(
"NullTransverseMomentum") <<
"dropping input candidate with pt=0";
418 fjInputs.push_back(fastjet::PseudoJet(constit->px(),constit->py(),constit->pz(),constit->energy()));
425 const reco::Vertex &
pv = *(trackIPTagInfos->front().primaryVertex());
428 fastjet::PseudoJet
p(dir.
x(),dir.
y(),dir.
z(),dir.
mag());
430 p = fastjet::PseudoJet(it->p4().px(),it->p4().py(),it->p4().pz(),it->p4().energy());
432 p.set_user_info(
new VertexInfo( it - extSecVertex->begin() ));
433 fjInputs.push_back(p);
439 std::vector<fastjet::PseudoJet> inclusiveJets = fastjet::sorted_by_pt(
fjClusterSeq->inclusive_jets(
jetPtMin) );
443 if( inclusiveJets.size() < fatJetsHandle->size() )
444 edm::LogError(
"TooFewReclusteredJets") <<
"There are fewer reclustered (" << inclusiveJets.size() <<
") than original fat jets (" << fatJetsHandle->size() <<
"). Please check that the jet algorithm and jet size match those used for the original jet collection.";
447 std::vector<int> reclusteredIndices;
448 matchReclusteredJets<edm::View<reco::Jet> >(fatJetsHandle,inclusiveJets,reclusteredIndices,
"fat");
451 std::vector<int> groomedIndices;
455 std::vector<std::vector<int> > subjetIndices;
456 matchSubjets(groomedIndices,groomedFatJetsHandle,trackIPTagInfos,subjetIndices);
459 for(
size_t i=0;
i<fatJetsHandle->size(); ++
i)
461 if( reclusteredIndices.at(
i) < 0 )
continue;
463 if( fatJetsHandle->at(
i).pt() == 0 )
465 edm::LogWarning(
"NullTransverseMomentum") <<
"The original fat jet " <<
i <<
" has Pt=0. This is not expected so the jet will be skipped.";
469 if( subjetIndices.at(
i).size()==0 )
continue;
472 if( (
std::abs( inclusiveJets.at(reclusteredIndices.at(
i)).
pt() - fatJetsHandle->at(
i).pt() ) / fatJetsHandle->at(
i).pt() ) >
relPtTolerance )
474 if( fatJetsHandle->at(
i).pt() < 10. )
475 edm::LogWarning(
"JetPtMismatchAtLowPt") <<
"The reclustered and original fat jet " <<
i <<
" have different Pt's (" << inclusiveJets.at(reclusteredIndices.at(
i)).
pt() <<
" vs " << fatJetsHandle->at(
i).pt() <<
" GeV, respectively).\n"
476 <<
"Please check that the jet algorithm and jet size match those used for the original fat jet collection and also make sure the original fat jets are uncorrected. In addition, make sure you are not using CaloJets which are presently not supported.\n"
477 <<
"Since the mismatch is at low Pt, it is ignored and only a warning is issued.\n"
478 <<
"\nIn extremely rare instances the mismatch could be caused by a difference in the machine precision in which case make sure the original jet collection is produced and reclustering is performed in the same job.";
480 edm::LogError(
"JetPtMismatch") <<
"The reclustered and original fat jet " <<
i <<
" have different Pt's (" << inclusiveJets.at(reclusteredIndices.at(
i)).
pt() <<
" vs " << fatJetsHandle->at(
i).pt() <<
" GeV, respectively).\n"
481 <<
"Please check that the jet algorithm and jet size match those used for the original fat jet collection and also make sure the original fat jets are uncorrected. In addition, make sure you are not using CaloJets which are presently not supported.\n"
482 <<
"\nIn extremely rare instances the mismatch could be caused by a difference in the machine precision in which case make sure the original jet collection is produced and reclustering is performed in the same job.";
486 std::vector<fastjet::PseudoJet> constituents = inclusiveJets.at(reclusteredIndices.at(
i)).constituents();
488 std::vector<int> svIndices;
490 for(std::vector<fastjet::PseudoJet>::const_iterator it = constituents.begin(); it != constituents.end(); ++it)
492 if( !it->has_user_info() )
continue;
494 svIndices.push_back( it->user_info<VertexInfo>().vertexIndex() );
498 for(
size_t sv=0; sv<svIndices.size(); ++sv)
500 const reco::Vertex &pv = *(trackIPTagInfos->front().primaryVertex());
501 const VTX &extSV = (*extSecVertex)[ svIndices.at(sv) ];
504 fastjet::PseudoJet
p(dir.
x(),dir.
y(),dir.
z(),dir.
mag());
506 p = fastjet::PseudoJet(extSV.p4().px(),extSV.p4().py(),extSV.p4().pz(),extSV.p4().energy());
508 std::vector<double> dR2toSubjets;
510 for(
size_t sj=0; sj<subjetIndices.at(
i).size(); ++sj)
511 dR2toSubjets.push_back(
Geom::deltaR2( p.rapidity(), p.phi_std(), trackIPTagInfos->at(subjetIndices.at(
i).at(sj)).
jet()->rapidity(), trackIPTagInfos->at(subjetIndices.at(
i).at(sj)).
jet()->phi() ) );
514 int closestSubjetIdx =
std::distance( dR2toSubjets.begin(), std::min_element(dR2toSubjets.begin(), dR2toSubjets.end()) );
516 clusteredSVs.at(subjetIndices.at(
i).at(closestSubjetIdx)).push_back( svIndices.at(sv) );
522 if( inclusiveJets.size() < trackIPTagInfos->size() )
523 edm::LogError(
"TooFewReclusteredJets") <<
"There are fewer reclustered (" << inclusiveJets.size() <<
") than original jets (" << trackIPTagInfos->size() <<
"). Please check that the jet algorithm and jet size match those used for the original jet collection.";
526 std::vector<int> reclusteredIndices;
527 matchReclusteredJets<std::vector<IPTI> >(
trackIPTagInfos,inclusiveJets,reclusteredIndices);
530 for(
size_t i=0;
i<trackIPTagInfos->size(); ++
i)
532 if( reclusteredIndices.at(
i) < 0 )
continue;
534 if( trackIPTagInfos->at(
i).jet()->pt() == 0 )
536 edm::LogWarning(
"NullTransverseMomentum") <<
"The original jet " <<
i <<
" has Pt=0. This is not expected so the jet will be skipped.";
541 if( (
std::abs( inclusiveJets.at(reclusteredIndices.at(
i)).
pt() - trackIPTagInfos->at(
i).jet()->pt() ) / trackIPTagInfos->at(
i).jet()->pt() ) >
relPtTolerance )
543 if( trackIPTagInfos->at(
i).jet()->pt() < 10. )
544 edm::LogWarning(
"JetPtMismatchAtLowPt") <<
"The reclustered and original jet " <<
i <<
" have different Pt's (" << inclusiveJets.at(reclusteredIndices.at(
i)).
pt() <<
" vs " << trackIPTagInfos->at(
i).jet()->pt() <<
" GeV, respectively).\n"
545 <<
"Please check that the jet algorithm and jet size match those used for the original jet collection and also make sure the original jets are uncorrected. In addition, make sure you are not using CaloJets which are presently not supported.\n"
546 <<
"Since the mismatch is at low Pt, it is ignored and only a warning is issued.\n"
547 <<
"\nIn extremely rare instances the mismatch could be caused by a difference in the machine precision in which case make sure the original jet collection is produced and reclustering is performed in the same job.";
549 edm::LogError(
"JetPtMismatch") <<
"The reclustered and original jet " <<
i <<
" have different Pt's (" << inclusiveJets.at(reclusteredIndices.at(
i)).
pt() <<
" vs " << trackIPTagInfos->at(
i).jet()->pt() <<
" GeV, respectively).\n"
550 <<
"Please check that the jet algorithm and jet size match those used for the original jet collection and also make sure the original jets are uncorrected. In addition, make sure you are not using CaloJets which are presently not supported.\n"
551 <<
"\nIn extremely rare instances the mismatch could be caused by a difference in the machine precision in which case make sure the original jet collection is produced and reclustering is performed in the same job.";
555 std::vector<fastjet::PseudoJet> constituents = inclusiveJets.at(reclusteredIndices.at(
i)).constituents();
558 for(std::vector<fastjet::PseudoJet>::const_iterator it = constituents.begin(); it != constituents.end(); ++it)
560 if( !it->has_user_info() )
continue;
562 clusteredSVs.at(
i).push_back( it->user_info<VertexInfo>().vertexIndex() );
569 std::auto_ptr<ConfigurableVertexReconstructor>
vertexReco;
570 std::auto_ptr<GhostTrackVertexFinder> vertexRecoGT;
582 TransientTrackMap primariesMap;
588 for(
typename std::vector<IPTI>::const_iterator iterJets =
589 trackIPTagInfos->begin(); iterJets != trackIPTagInfos->end();
594 const Vertex &pv = *iterJets->primaryVertex();
596 std::set<TransientTrack> primaries;
600 TransientTrackMap::iterator pos =
601 primariesMap.lower_bound(iter->get());
603 if (pos != primariesMap.end() &&
604 pos->first == iter->get())
605 primaries.insert(pos->second);
610 primariesMap.insert(pos,
611 std::make_pair(iter->get(), track));
612 primaries.insert(track);
623 std::vector<std::size_t> indices =
628 const std::vector<reco::btag::TrackIPData> &ipData =
629 iterJets->impactParameterData();
633 std::vector<TransientTrack> fitTracks;
634 std::vector<GhostTrackState> gtStates;
635 std::auto_ptr<GhostTrackPrediction> gtPred;
638 *iterJets->ghostTrack()));
640 for(
unsigned int i = 0;
i < indices.size();
i++) {
646 trackData.back().first = indices[
i];
653 trackData.back().second.svStatus =
658 TransientTrackMap::const_iterator pos =
661 if (pos != primariesMap.end()) {
662 primaries.erase(pos->second);
663 fitTrack = pos->second;
665 fitTrack = trackBuilder->build(trackRef);
666 fitTracks.push_back(fitTrack);
668 trackData.back().second.svStatus =
674 ipData[indices[
i]].closestToGhostTrack;
675 gtState.linearize(*gtPred,
true,
676 gtPred->lambda(pos));
677 gtState.setWeight(ipData[indices[i]].ghostTrackWeight);
678 gtStates.push_back(gtState);
682 std::auto_ptr<GhostTrack> ghostTrack;
689 iterJets->ghostTrack()->px(),
690 iterJets->ghostTrack()->py(),
691 iterJets->ghostTrack()->pz()),
694 iterJets->ghostTrack()->chi2(),
695 iterJets->ghostTrack()->ndof()));
700 std::vector<VTX> extAssoCollection;
701 std::vector<TransientVertex> fittedSVs;
702 std::vector<SecondaryVertex> SVs;
707 fittedSVs = vertexRecoGT->vertices(
710 fittedSVs = vertexReco->vertices(fitTracks);
715 fittedSVs = vertexRecoGT->vertices(
716 pv, *beamSpot, *ghostTrack);
718 fittedSVs = vertexReco->vertices(fitTracks,
726 for(
unsigned int i = 0; i < 7; i++) {
727 unsigned int covSrc = bsCovSrc[
i];
728 for(
unsigned int j = 0;
j < 7;
j++) {
730 if (!covSrc || bsCovSrc[
j] != covSrc)
732 else if (covSrc == 1)
733 v = beamSpot->covariance(i,
j);
742 beamSpot.
isValid() ? beamSpot->dxdz() : 0.,
743 beamSpot.
isValid() ? beamSpot->dydz() : 0.,
747 fittedSVs = vertexRecoGT->vertices(
748 pv, bs, *ghostTrack);
750 fittedSVs = vertexReco->vertices(fitTracks, bs);
754 std::vector<TransientTrack> primaries_(
755 primaries.begin(), primaries.end());
757 fittedSVs = vertexRecoGT->vertices(
758 pv, *beamSpot, primaries_,
761 fittedSVs = vertexReco->vertices(
762 primaries_, fitTracks,
768 std::remove_copy_if(boost::make_transform_iterator(
769 fittedSVs.begin(), svBuilder),
770 boost::make_transform_iterator(
771 fittedSVs.end(), svBuilder),
772 std::back_inserter(SVs),
777 for(
size_t iExtSv = 0; iExtSv < extSecVertex->size(); iExtSv++){
778 const VTX & extVertex = (*extSecVertex)[iExtSv];
781 extAssoCollection.push_back( extVertex );
786 size_t jetIdx = ( iterJets - trackIPTagInfos->begin() );
788 for(
size_t iExtSv = 0; iExtSv < clusteredSVs.at(jetIdx).size(); iExtSv++){
789 const VTX & extVertex = (*extSecVertex)[ clusteredSVs.at(jetIdx).at(iExtSv) ];
790 if( extVertex.p4().M() < 0.3 )
792 extAssoCollection.push_back( extVertex );
797 std::remove_copy_if(boost::make_transform_iterator( extAssoCollection.begin(), svBuilder),
798 boost::make_transform_iterator(extAssoCollection.end(), svBuilder),
799 std::back_inserter(SVs),
808 extAssoCollection.clear();
814 std::vector<typename TemplatedSecondaryVertexTagInfo<IPTI,VTX>::VertexData> svData;
816 svData.resize(vtxIndices.size());
817 for(
unsigned int idx = 0;
idx < vtxIndices.size();
idx++) {
820 svData[
idx].vertex = sv;
832 trackData, svData, SVs.size(),
834 iterJets - trackIPTagInfos->begin())));
math::Error< dimension >::type CovarianceMatrix
reco::Vertex::Point convertPos(const GlobalPoint &p)
T getParameter(std::string const &) const
edm::EDGetTokenT< edm::View< reco::Jet > > token_fatJets
TrackSelector trackSelector
trackRef_iterator tracks_end() const
last iterator over tracks
ClusterSequencePtr fjClusterSeq
boost::shared_ptr< fastjet::ClusterSequence > ClusterSequencePtr
edm::EDGetTokenT< reco::BeamSpot > token_BeamSpot
reco::Vertex::Error convertError(const GlobalError &ge)
virtual Vector momentum() const
spatial momentum vector
Measurement1D dist3d() const
Measurement1D dist2d() const
double covariance(int i, int j) const
(i, j)-th element of error matrix, i, j = 0, ... 2
VertexFilter vertexFilter
reco::btag::SortCriteria sortCriterium
const reco::Track * toTrack(const reco::TrackBaseRef &t)
const Point & position() const
position
static GhostTrackVertexFinder::FitType getGhostTrackFitType(const std::string &name)
ConstraintType constraint
edm::EDGetTokenT< std::vector< IPTI > > token_trackIPTagInfo
IPTI::input_container::value_type input_item
reco::btag::IndexedTrackData IndexedTrackData
GlobalVector flightDirection(const reco::Vertex &pv, const reco::Vertex &sv)
Abs< T >::type abs(const T &t)
double deltaR2(const T1 &t1, const T2 &t2)
VertexSorting< SecondaryVertex > vertexSorting
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
Vector3DBase unit() const
edm::ParameterSet vtxRecoPSet
edm::EDGetTokenT< edm::View< reco::Jet > > token_groomedFatJets
tuple idx
DEBUGGING if hasattr(process,"trackMonIterativeTracking2012"): print "trackMonIterativeTracking2012 D...
Error error() const
return SMatrix
std::vector< reco::btag::IndexedTrackData > TrackDataVector
boost::indirect_iterator< typename seq_t::const_iterator > const_iterator
static int position[264][3]
std::vector< TemplatedSecondaryVertexTagInfo< IPTI, VTX > > Product
std::vector< TrackBaseRef >::const_iterator trackRef_iterator
The iteratator for the vector<TrackRef>
trackRef_iterator tracks_begin() const
first iterator over tracks
JetDefPtr fjJetDefinition
std::pair< unsigned int, TrackData > IndexedTrackData
Global3DVector GlobalVector
IPTI::input_container input_container
Public Member Functions inherited from edm::stream::EDProducer<>
1.8.5
