319 typedef std::map<const Track *, TransientTrack> TransientTrackMap;
340 if( groomedFatJetsHandle->size() > fatJetsHandle->size() )
341 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.";
345 unsigned int bsCovSrc[7] = { 0, };
346 double sigmaZ = 0.0, beamWidth = 0.0;
350 bsCovSrc[3] = bsCovSrc[4] = bsCovSrc[5] = bsCovSrc[6] = 1;
351 sigmaZ = beamSpot->sigmaZ();
352 beamWidth = beamSpot->BeamWidthX();
357 bsCovSrc[0] = bsCovSrc[1] = 2;
358 bsCovSrc[3] = bsCovSrc[4] = bsCovSrc[5] = 1;
359 sigmaZ = beamSpot->sigmaZ();
363 bsCovSrc[0] = bsCovSrc[1] = bsCovSrc[2] = 2;
376 std::vector<std::vector<int> > clusteredSVs(trackIPTagInfos->size(),std::vector<int>());
380 std::vector<fastjet::PseudoJet> fjInputs;
386 std::vector<edm::Ptr<reco::Candidate> > constituents = it->getJetConstituents();
387 std::vector<edm::Ptr<reco::Candidate> >::const_iterator
m;
388 for( m = constituents.begin(); m != constituents.end(); ++
m )
391 if(constit->pt() == 0)
393 edm::LogWarning(
"NullTransverseMomentum") <<
"dropping input candidate with pt=0";
396 fjInputs.push_back(fastjet::PseudoJet(constit->px(),constit->py(),constit->pz(),constit->energy()));
402 for(
typename std::vector<IPTI>::const_iterator it = trackIPTagInfos->begin(); it != trackIPTagInfos->end(); ++it)
404 std::vector<edm::Ptr<reco::Candidate> > constituents = it->jet()->getJetConstituents();
405 std::vector<edm::Ptr<reco::Candidate> >::const_iterator
m;
406 for( m = constituents.begin(); m != constituents.end(); ++
m )
409 if(constit->pt() == 0)
411 edm::LogWarning(
"NullTransverseMomentum") <<
"dropping input candidate with pt=0";
414 fjInputs.push_back(fastjet::PseudoJet(constit->px(),constit->py(),constit->pz(),constit->energy()));
421 const reco::Vertex &
pv = *(trackIPTagInfos->front().primaryVertex());
424 fastjet::PseudoJet
p(dir.
x(),dir.
y(),dir.
z(),dir.
mag());
426 p = fastjet::PseudoJet(it->p4().px(),it->p4().py(),it->p4().pz(),it->p4().energy());
428 p.set_user_info(
new VertexInfo( it - extSecVertex->begin() ));
429 fjInputs.push_back(p);
435 std::vector<fastjet::PseudoJet> inclusiveJets = fastjet::sorted_by_pt(
fjClusterSeq->inclusive_jets(
jetPtMin) );
439 if( inclusiveJets.size() < fatJetsHandle->size() )
440 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.";
443 std::vector<int> reclusteredIndices;
444 matchReclusteredJets<edm::View<reco::Jet> >(fatJetsHandle,inclusiveJets,reclusteredIndices,
"fat");
447 std::vector<int> groomedIndices;
451 std::vector<std::vector<int> > subjetIndices;
452 matchSubjets(groomedIndices,groomedFatJetsHandle,trackIPTagInfos,subjetIndices);
455 for(
size_t i=0;
i<fatJetsHandle->size(); ++
i)
457 if( reclusteredIndices.at(
i) < 0 )
continue;
459 if( fatJetsHandle->at(
i).pt() == 0 )
461 edm::LogWarning(
"NullTransverseMomentum") <<
"The original fat jet " <<
i <<
" has Pt=0. This is not expected so the jet will be skipped.";
465 if( subjetIndices.at(
i).size()==0 )
continue;
468 if( (
std::abs( inclusiveJets.at(reclusteredIndices.at(
i)).
pt() - fatJetsHandle->at(
i).pt() ) / fatJetsHandle->at(
i).pt() ) >
relPtTolerance )
470 if( fatJetsHandle->at(
i).pt() < 10. )
471 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"
472 <<
"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"
473 <<
"Since the mismatch is at low Pt, it is ignored and only a warning is issued.\n"
474 <<
"\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.";
476 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"
477 <<
"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"
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.";
482 std::vector<fastjet::PseudoJet> constituents = inclusiveJets.at(reclusteredIndices.at(
i)).constituents();
484 std::vector<int> svIndices;
486 for(std::vector<fastjet::PseudoJet>::const_iterator it = constituents.begin(); it != constituents.end(); ++it)
488 if( !it->has_user_info() )
continue;
490 svIndices.push_back( it->user_info<VertexInfo>().vertexIndex() );
494 for(
size_t sv=0; sv<svIndices.size(); ++sv)
496 const reco::Vertex &pv = *(trackIPTagInfos->front().primaryVertex());
497 const VTX &extSV = (*extSecVertex)[ svIndices.at(sv) ];
500 fastjet::PseudoJet
p(dir.
x(),dir.
y(),dir.
z(),dir.
mag());
502 p = fastjet::PseudoJet(extSV.p4().px(),extSV.p4().py(),extSV.p4().pz(),extSV.p4().energy());
504 std::vector<double> dR2toSubjets;
506 for(
size_t sj=0; sj<subjetIndices.at(
i).size(); ++sj)
507 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() ) );
510 int closestSubjetIdx = std::distance( dR2toSubjets.begin(), std::min_element(dR2toSubjets.begin(), dR2toSubjets.end()) );
512 clusteredSVs.at(subjetIndices.at(
i).at(closestSubjetIdx)).push_back( svIndices.at(sv) );
518 if( inclusiveJets.size() < trackIPTagInfos->size() )
519 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.";
522 std::vector<int> reclusteredIndices;
523 matchReclusteredJets<std::vector<IPTI> >(trackIPTagInfos,inclusiveJets,reclusteredIndices);
526 for(
size_t i=0;
i<trackIPTagInfos->size(); ++
i)
528 if( reclusteredIndices.at(
i) < 0 )
continue;
530 if( trackIPTagInfos->at(
i).jet()->pt() == 0 )
532 edm::LogWarning(
"NullTransverseMomentum") <<
"The original jet " <<
i <<
" has Pt=0. This is not expected so the jet will be skipped.";
537 if( (
std::abs( inclusiveJets.at(reclusteredIndices.at(
i)).
pt() - trackIPTagInfos->at(
i).jet()->pt() ) / trackIPTagInfos->at(
i).jet()->pt() ) >
relPtTolerance )
539 if( trackIPTagInfos->at(
i).jet()->pt() < 10. )
540 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"
541 <<
"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"
542 <<
"Since the mismatch is at low Pt, it is ignored and only a warning is issued.\n"
543 <<
"\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.";
545 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"
546 <<
"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"
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.";
551 std::vector<fastjet::PseudoJet> constituents = inclusiveJets.at(reclusteredIndices.at(
i)).constituents();
554 for(std::vector<fastjet::PseudoJet>::const_iterator it = constituents.begin(); it != constituents.end(); ++it)
556 if( !it->has_user_info() )
continue;
558 clusteredSVs.at(
i).push_back( it->user_info<VertexInfo>().vertexIndex() );
565 std::auto_ptr<ConfigurableVertexReconstructor> vertexReco;
566 std::auto_ptr<GhostTrackVertexFinder> vertexRecoGT;
578 TransientTrackMap primariesMap;
582 std::auto_ptr<Product> tagInfos(
new Product);
584 for(
typename std::vector<IPTI>::const_iterator iterJets =
585 trackIPTagInfos->begin(); iterJets != trackIPTagInfos->end();
590 const Vertex &pv = *iterJets->primaryVertex();
592 std::set<TransientTrack> primaries;
596 TransientTrackMap::iterator pos =
597 primariesMap.lower_bound(
iter->get());
599 if (pos != primariesMap.end() &&
600 pos->first ==
iter->get())
601 primaries.insert(pos->second);
606 primariesMap.insert(pos,
607 std::make_pair(
iter->get(), track));
608 primaries.insert(track);
619 std::vector<std::size_t> indices =
624 const std::vector<reco::btag::TrackIPData> &ipData =
625 iterJets->impactParameterData();
629 std::vector<TransientTrack> fitTracks;
630 std::vector<GhostTrackState> gtStates;
631 std::auto_ptr<GhostTrackPrediction> gtPred;
634 *iterJets->ghostTrack()));
636 for(
unsigned int i = 0;
i < indices.size();
i++) {
642 trackData.back().first = indices[
i];
649 trackData.back().second.svStatus =
654 TransientTrackMap::const_iterator pos =
657 if (pos != primariesMap.end()) {
658 primaries.erase(pos->second);
659 fitTrack = pos->second;
661 fitTrack = trackBuilder->build(trackRef);
662 fitTracks.push_back(fitTrack);
664 trackData.back().second.svStatus =
670 ipData[indices[
i]].closestToGhostTrack;
671 gtState.linearize(*gtPred,
true,
672 gtPred->lambda(pos));
673 gtState.setWeight(ipData[indices[i]].ghostTrackWeight);
674 gtStates.push_back(gtState);
678 std::auto_ptr<GhostTrack> ghostTrack;
685 iterJets->ghostTrack()->px(),
686 iterJets->ghostTrack()->py(),
687 iterJets->ghostTrack()->pz()),
690 iterJets->ghostTrack()->chi2(),
691 iterJets->ghostTrack()->ndof()));
696 std::vector<VTX> extAssoCollection;
697 std::vector<TransientVertex> fittedSVs;
698 std::vector<SecondaryVertex> SVs;
703 fittedSVs = vertexRecoGT->vertices(
706 fittedSVs = vertexReco->vertices(fitTracks);
711 fittedSVs = vertexRecoGT->vertices(
712 pv, *beamSpot, *ghostTrack);
714 fittedSVs = vertexReco->vertices(fitTracks,
722 for(
unsigned int i = 0; i < 7; i++) {
723 unsigned int covSrc = bsCovSrc[
i];
724 for(
unsigned int j = 0;
j < 7;
j++) {
726 if (!covSrc || bsCovSrc[
j] != covSrc)
728 else if (covSrc == 1)
729 v = beamSpot->covariance(i,
j);
738 beamSpot.
isValid() ? beamSpot->dxdz() : 0.,
739 beamSpot.
isValid() ? beamSpot->dydz() : 0.,
743 fittedSVs = vertexRecoGT->vertices(
744 pv, bs, *ghostTrack);
746 fittedSVs = vertexReco->vertices(fitTracks, bs);
750 std::vector<TransientTrack> primaries_(
751 primaries.begin(), primaries.end());
753 fittedSVs = vertexRecoGT->vertices(
754 pv, *beamSpot, primaries_,
757 fittedSVs = vertexReco->vertices(
758 primaries_, fitTracks,
764 std::remove_copy_if(boost::make_transform_iterator(
765 fittedSVs.begin(), svBuilder),
766 boost::make_transform_iterator(
767 fittedSVs.end(), svBuilder),
768 std::back_inserter(SVs),
773 for(
size_t iExtSv = 0; iExtSv < extSecVertex->size(); iExtSv++){
774 const VTX & extVertex = (*extSecVertex)[iExtSv];
777 extAssoCollection.push_back( extVertex );
782 size_t jetIdx = ( iterJets - trackIPTagInfos->begin() );
784 for(
size_t iExtSv = 0; iExtSv < clusteredSVs.at(jetIdx).size(); iExtSv++){
785 const VTX & extVertex = (*extSecVertex)[ clusteredSVs.at(jetIdx).at(iExtSv) ];
786 if( extVertex.p4().M() < 0.3 )
788 extAssoCollection.push_back( extVertex );
793 std::remove_copy_if(boost::make_transform_iterator( extAssoCollection.begin(), svBuilder),
794 boost::make_transform_iterator(extAssoCollection.end(), svBuilder),
795 std::back_inserter(SVs),
804 extAssoCollection.clear();
810 std::vector<typename TemplatedSecondaryVertexTagInfo<IPTI,VTX>::VertexData> svData;
812 svData.resize(vtxIndices.size());
813 for(
unsigned int idx = 0;
idx < vtxIndices.size();
idx++) {
816 svData[
idx].vertex = sv;
828 trackData, svData, SVs.size(),
830 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
Static Public Member Functions inherited from edm::stream::EDProducerBase
1.8.5
