#include <IPProducer.h>

Inheritance diagram for IPProducer< Container, Base, Helper >:

Public Types
typedef std::vector < reco::IPTagInfo< Container, Base > >	Product

Public Types inherited from edm::stream::EDProducer<>
using	CacheTypes = CacheContexts< T...>

using	GlobalCache = typename CacheTypes::GlobalCache

using	HasAbility = AbilityChecker< T...>

using	InputProcessBlockCache = typename CacheTypes::InputProcessBlockCache

using	LuminosityBlockCache = typename CacheTypes::LuminosityBlockCache

using	LuminosityBlockContext = LuminosityBlockContextT< LuminosityBlockCache, RunCache, GlobalCache >

using	LuminosityBlockSummaryCache = typename CacheTypes::LuminosityBlockSummaryCache

using	RunCache = typename CacheTypes::RunCache

using	RunContext = RunContextT< RunCache, GlobalCache >

using	RunSummaryCache = typename CacheTypes::RunSummaryCache

Public Member Functions
template<>
void	fillDescriptions (edm::ConfigurationDescriptions &descriptions)

template<>
void	fillDescriptions (edm::ConfigurationDescriptions &descriptions)

	IPProducer (const edm::ParameterSet &)

void	produce (edm::Event &, const edm::EventSetup &) override

	~IPProducer () override

Public Member Functions inherited from edm::stream::EDProducer<>
	EDProducer ()=default

	EDProducer (const EDProducer &)=delete

bool	hasAbilityToProduceInBeginLumis () const final

bool	hasAbilityToProduceInBeginProcessBlocks () const final

bool	hasAbilityToProduceInBeginRuns () const final

bool	hasAbilityToProduceInEndLumis () const final

bool	hasAbilityToProduceInEndProcessBlocks () const final

bool	hasAbilityToProduceInEndRuns () const final

const EDProducer &	operator= (const EDProducer &)=delete

Static Public Member Functions
static void	fillDescriptions (edm::ConfigurationDescriptions &descriptions)

Private Member Functions
void	checkEventSetup (const edm::EventSetup &iSetup)

Private Attributes
unsigned long long	m_calibrationCacheId2D

unsigned long long	m_calibrationCacheId3D

bool	m_computeGhostTrack

bool	m_computeProbabilities

const edm::ParameterSet &	m_config

double	m_cutMaxChiSquared

double	m_cutMaxLIP

double	m_cutMaxTIP

double	m_cutMinPt

int	m_cutPixelHits

int	m_cutTotalHits

bool	m_directionWithGhostTrack

bool	m_directionWithTracks

double	m_ghostTrackPriorDeltaR

Helper	m_helper

std::unique_ptr < HistogramProbabilityEstimator >	m_probabilityEstimator

bool	m_useDB

bool	m_useTrackQuality

edm::ESGetToken < TrackProbabilityCalibration, BTagTrackProbability2DRcd >	token_calib2D

edm::ESGetToken < TrackProbabilityCalibration, BTagTrackProbability3DRcd >	token_calib3D

edm::EDGetTokenT < reco::VertexCollection >	token_primaryVertex

edm::ESGetToken < TransientTrackBuilder, TransientTrackRecord >	token_trackBuilder

Detailed Description

template<class Container, class Base, class Helper>
class IPProducer< Container, Base, Helper >

Description: <one line="" class="" summary>="">

Implementation: <Notes on="" implementation>="">

Definition at line 130 of file IPProducer.h.

Member Typedef Documentation

template<class Container, class Base, class Helper>

typedef std::vector<reco::IPTagInfo<Container, Base> > IPProducer< Container, Base, Helper >::Product

Definition at line 132 of file IPProducer.h.

Constructor & Destructor Documentation

template<class Container , class Base , class Helper >

IPProducer< Container, Base, Helper >::IPProducer ( const edm::ParameterSet & iConfig )

explicit

Definition at line 191 of file IPProducer.h.

     : m_config(iConfig), m_helper(iConfig, consumesCollector()) {
   m_calibrationCacheId3D = 0;
   m_calibrationCacheId2D = 0;
 
   token_primaryVertex = consumes<reco::VertexCollection>(m_config.getParameter<edm::InputTag>("primaryVertex"));
   token_trackBuilder =
       esConsumes<TransientTrackBuilder, TransientTrackRecord>(edm::ESInputTag("", "TransientTrackBuilder"));
   token_calib2D = esConsumes<TrackProbabilityCalibration, BTagTrackProbability2DRcd>();
   token_calib3D = esConsumes<TrackProbabilityCalibration, BTagTrackProbability3DRcd>();
 
   m_computeProbabilities = m_config.getParameter<bool>("computeProbabilities");
   m_computeGhostTrack = m_config.getParameter<bool>("computeGhostTrack");
   m_ghostTrackPriorDeltaR = m_config.getParameter<double>("ghostTrackPriorDeltaR");
   m_cutPixelHits = m_config.getParameter<int>("minimumNumberOfPixelHits");
   m_cutTotalHits = m_config.getParameter<int>("minimumNumberOfHits");
   m_cutMaxTIP = m_config.getParameter<double>("maximumTransverseImpactParameter");
   m_cutMinPt = m_config.getParameter<double>("minimumTransverseMomentum");
   m_cutMaxChiSquared = m_config.getParameter<double>("maximumChiSquared");
   m_cutMaxLIP = m_config.getParameter<double>("maximumLongitudinalImpactParameter");
   m_directionWithTracks = m_config.getParameter<bool>("jetDirectionUsingTracks");
   m_directionWithGhostTrack = m_config.getParameter<bool>("jetDirectionUsingGhostTrack");
   m_useTrackQuality = m_config.getParameter<bool>("useTrackQuality");
 
   if (m_computeGhostTrack)
     produces<reco::TrackCollection>("ghostTracks");
   produces<Product>();
 }

template<class Container , class Base , class Helper >

IPProducer< Container, Base, Helper >::~IPProducer ( )

override

Definition at line 221 of file IPProducer.h.

221 {}

Member Function Documentation

template<class Container , class Base , class Helper >

void IPProducer< Container, Base, Helper >::checkEventSetup ( const edm::EventSetup & iSetup )

private

Definition at line 420 of file IPProducer.h.

References edm::eventsetup::EventSetupRecord::cacheIdentifier(), edm::EventSetup::get(), edm::EventSetup::getHandle(), and edm::ESHandle< class >::product().

                                                                                      {
   const edm::eventsetup::EventSetupRecord& re2D = iSetup.get<BTagTrackProbability2DRcd>();
   const edm::eventsetup::EventSetupRecord& re3D = iSetup.get<BTagTrackProbability3DRcd>();
   unsigned long long cacheId2D = re2D.cacheIdentifier();
   unsigned long long cacheId3D = re3D.cacheIdentifier();
 
   if (cacheId2D != m_calibrationCacheId2D || cacheId3D != m_calibrationCacheId3D)  //Calibration changed
   {
     //iSetup.get<BTagTrackProbabilityRcd>().get(calib);
     edm::ESHandle<TrackProbabilityCalibration> calib2DHandle = iSetup.getHandle(token_calib2D);
     edm::ESHandle<TrackProbabilityCalibration> calib3DHandle = iSetup.getHandle(token_calib3D);
 
     const TrackProbabilityCalibration* ca2D = calib2DHandle.product();
     const TrackProbabilityCalibration* ca3D = calib3DHandle.product();
 
     m_probabilityEstimator = std::make_unique<HistogramProbabilityEstimator>(ca3D, ca2D);
   }
   m_calibrationCacheId3D = cacheId3D;
   m_calibrationCacheId2D = cacheId2D;
 }

template<class Container, class Base, class Helper>

static void IPProducer< Container, Base, Helper >::fillDescriptions ( edm::ConfigurationDescriptions & descriptions )

static

template<>

void IPProducer< reco::TrackRefVector, reco::JTATagInfo, IPProducerHelpers::FromJTA >::fillDescriptions ( edm::ConfigurationDescriptions & descriptions )

inline

Definition at line 444 of file IPProducer.h.

References edm::ParameterSetDescription::add(), edm::ConfigurationDescriptions::addDefault(), submitPVResolutionJobs::desc, and HLT_FULL_cff::InputTag.

                                                 {
   edm::ParameterSetDescription desc;
   desc.add<double>("maximumTransverseImpactParameter", 0.2);
   desc.add<int>("minimumNumberOfHits", 8);
   desc.add<double>("minimumTransverseMomentum", 1.0);
   desc.add<edm::InputTag>("primaryVertex", edm::InputTag("offlinePrimaryVertices"));
   desc.add<double>("maximumLongitudinalImpactParameter", 17.0);
   desc.add<bool>("computeGhostTrack", true);
   desc.add<double>("ghostTrackPriorDeltaR", 0.03);
   desc.add<edm::InputTag>("jetTracks", edm::InputTag("ak4JetTracksAssociatorAtVertexPF"));
   desc.add<bool>("jetDirectionUsingGhostTrack", false);
   desc.add<int>("minimumNumberOfPixelHits", 2);
   desc.add<bool>("jetDirectionUsingTracks", false);
   desc.add<bool>("computeProbabilities", true);
   desc.add<bool>("useTrackQuality", false);
   desc.add<double>("maximumChiSquared", 5.0);
   descriptions.addDefault(desc);
 }

template<>

void IPProducer< std::vector< reco::CandidatePtr >, reco::JetTagInfo, IPProducerHelpers::FromJetAndCands >::fillDescriptions ( edm::ConfigurationDescriptions & descriptions )

inline

Definition at line 466 of file IPProducer.h.

References edm::ParameterSetDescription::add(), edm::ConfigurationDescriptions::addDefault(), edm::ParameterSetDescription::addOptional(), submitPVResolutionJobs::desc, and HLT_FULL_cff::InputTag.

                                                 {
   edm::ParameterSetDescription desc;
   desc.add<double>("maximumTransverseImpactParameter", 0.2);
   desc.add<int>("minimumNumberOfHits", 8);
   desc.add<double>("minimumTransverseMomentum", 1.0);
   desc.add<edm::InputTag>("primaryVertex", edm::InputTag("offlinePrimaryVertices"));
   desc.add<double>("maximumLongitudinalImpactParameter", 17.0);
   desc.add<bool>("computeGhostTrack", true);
   desc.add<double>("maxDeltaR", 0.4);
   desc.add<edm::InputTag>("candidates", edm::InputTag("particleFlow"));
   desc.add<bool>("jetDirectionUsingGhostTrack", false);
   desc.add<int>("minimumNumberOfPixelHits", 2);
   desc.add<bool>("jetDirectionUsingTracks", false);
   desc.add<bool>("computeProbabilities", true);
   desc.add<bool>("useTrackQuality", false);
   desc.add<edm::InputTag>("jets", edm::InputTag("ak4PFJetsCHS"));
   desc.add<double>("ghostTrackPriorDeltaR", 0.03);
   desc.add<double>("maximumChiSquared", 5.0);
   desc.addOptional<bool>("explicitJTA", false);
   descriptions.addDefault(desc);
 }

template<class Container , class Base , class Helper >

void IPProducer< Container, Base, Helper >::produce	(	edm::Event &	iEvent,
		const edm::EventSetup &	iSetup
	)

override

Definition at line 228 of file IPProducer.h.

References funct::abs(), IPTools::closestApproachToJet(), reco::btag::TrackIPData::closestToGhostTrack, reco::btag::TrackIPData::closestToJetAxis, RecoVertex::convertError(), RecoVertex::convertPos(), reco::GhostTrackPrediction::direction(), VertexDistance3D::distance(), reco::btag::TrackIPData::distanceToGhostTrack, reco::btag::TrackIPData::distanceToJetAxis, reco::TrackBase::dxy(), reco::TrackBase::dz(), alignCSCRings::e, relativeConstraints::error, reco::Vertex::error(), reco::TransientTrack::field(), reco::GhostTrackFitter::fit(), edm::Event::getByToken(), edm::EventSetup::getHandle(), edm::Event::getRefBeforePut(), reco::btag::TrackIPData::ghostTrackWeight, TrajectoryStateOnSurface::globalPosition(), reco::TrackBase::hitPattern(), reco::TransientTrack::impactPointState(), reco::btag::TrackIPData::ip2d, reco::btag::TrackIPData::ip3d, TrajectoryStateOnSurface::isValid(), IPTools::jetTrackDistance(), reco::GhostTrackPrediction::lambda(), eostools::move(), reco::TrackBase::normalizedChi2(), reco::HitPattern::numberOfValidHits(), reco::TrackBase::numberOfValidHits(), reco::HitPattern::numberOfValidPixelHits(), AlCaHLTBitMon_ParallelJobs::p, fireworks::p1, fireworks::p2, reco::GhostTrackPrediction::position(), reco::Vertex::position(), reco::GhostTrackPrediction::positionError(), reco::GhostTrackPrediction::prediction(), beam_dqm_sourceclient-live_cfg::primaryVertex, reco::TrackBase::pt(), edm::Event::put(), MetAnalyzer::pv(), mps_fire::result, TrackCollections2monitor_cff::selectedTracks, IPTools::signedImpactParameter3D(), IPTools::signedTransverseImpactParameter(), reco::btag::toTrack(), reco::TransientTrack::track(), HLT_FULL_cff::track, and tracks.

                                                                                                {
   // Update probability estimator if event setup is changed
   if (m_computeProbabilities)
     checkEventSetup(iSetup);
 
   edm::Handle<reco::VertexCollection> primaryVertex;
   iEvent.getByToken(token_primaryVertex, primaryVertex);
 
   edm::ESHandle<TransientTrackBuilder> builder = iSetup.getHandle(token_trackBuilder);
   // m_algo.setTransientTrackBuilder(builder.product());
 
   // output collections
   auto result = std::make_unique<Product>();
 
   std::unique_ptr<reco::TrackCollection> ghostTracks;
   reco::TrackRefProd ghostTrackRefProd;
   if (m_computeGhostTrack) {
     ghostTracks = std::make_unique<reco::TrackCollection>();
     ghostTrackRefProd = iEvent.getRefBeforePut<reco::TrackCollection>("ghostTracks");
   }
 
   // use first pv of the collection
   reco::Vertex dummy;
   const reco::Vertex* pv = &dummy;
   edm::Ref<reco::VertexCollection> pvRef;
   if (!primaryVertex->empty()) {
     pv = &*primaryVertex->begin();
     // we always use the first vertex (at the moment)
     pvRef = edm::Ref<reco::VertexCollection>(primaryVertex, 0);
   } else {  // create a dummy PV
     reco::Vertex::Error e;
     e(0, 0) = 0.0015 * 0.0015;
     e(1, 1) = 0.0015 * 0.0015;
     e(2, 2) = 15. * 15.;
     reco::Vertex::Point p(0, 0, 0);
     dummy = reco::Vertex(p, e, 0, 0, 0);
   }
 
   std::vector<Base> baseTagInfos = m_helper.makeBaseVector(iEvent);
   for (typename std::vector<Base>::const_iterator it = baseTagInfos.begin(); it != baseTagInfos.end(); it++) {
     Container tracks = m_helper.tracks(*it);
     math::XYZVector jetMomentum = it->jet()->momentum();
 
     if (m_directionWithTracks) {
       jetMomentum *= 0.5;
       for (typename Container::const_iterator itTrack = tracks.begin(); itTrack != tracks.end(); ++itTrack)
         if (reco::btag::toTrack(*itTrack)->numberOfValidHits() >= m_cutTotalHits)  //minimal quality cuts
           jetMomentum += (*itTrack)->momentum();
     }
 
     Container selectedTracks;
     std::vector<reco::TransientTrack> transientTracks;
 
     for (typename Container::const_iterator itTrack = tracks.begin(); itTrack != tracks.end(); ++itTrack) {
       reco::TransientTrack transientTrack = builder->build(*itTrack);
       const reco::Track& track = transientTrack.track();  //**itTrack;
                                                           /*    cout << " pt " <<  track.pt() <<
                " d0 " <<  fabs(track.d0()) <<
                " #hit " <<    track.hitPattern().numberOfValidHits()<<
                " ipZ " <<   fabs(track.dz()-pv->z())<<
                " chi2 " <<  track.normalizedChi2()<<
                " #pixel " <<    track.hitPattern().numberOfValidPixelHits()<< endl;
 */
       if (track.pt() > m_cutMinPt &&
           track.hitPattern().numberOfValidHits() >= m_cutTotalHits &&  // min num tracker hits
           track.hitPattern().numberOfValidPixelHits() >= m_cutPixelHits &&
           track.normalizedChi2() < m_cutMaxChiSquared && std::abs(track.dxy(pv->position())) < m_cutMaxTIP &&
           std::abs(track.dz(pv->position())) < m_cutMaxLIP) {
         //   std::cout << "selected" << std::endl;
         selectedTracks.push_back(*itTrack);
         transientTracks.push_back(transientTrack);
       }
     }
     //  std::cout <<"SIZE: " << transientTracks.size() << std::endl;
     GlobalVector direction(jetMomentum.x(), jetMomentum.y(), jetMomentum.z());
 
     std::unique_ptr<reco::GhostTrack> ghostTrack;
     reco::TrackRef ghostTrackRef;
     if (m_computeGhostTrack) {
       reco::GhostTrackFitter fitter;
       GlobalPoint origin = RecoVertex::convertPos(pv->position());
       GlobalError error = RecoVertex::convertError(pv->error());
       ghostTrack = std::make_unique<reco::GhostTrack>(
           fitter.fit(origin, error, direction, m_ghostTrackPriorDeltaR, transientTracks));
 
       /*
     if (std::sqrt(jetMomentum.Perp2()) > 30) {
         double offset = ghostTrack->prediction().lambda(origin);
         std::cout << "------------------ jet pt " << std::sqrt(jetMomentum.Perp2()) << std::endl;
         const std::vector<GhostTrackState> *states = &ghostTrack->states();
         for(std::vector<GhostTrackState>::const_iterator state = states->begin();
             state != states->end(); ++state) {
             double dist = state->lambda() - offset;
             double err = state->lambdaError(ghostTrack->prediction(), error);
             double ipSig = IPTools::signedImpactParameter3D(state->track(), direction, *pv).second.significance();
             double axisDist = state->axisDistance(ghostTrack->prediction());
             std::cout << state->track().impactPointState().freeState()->momentum().perp()
                       << ": " << dist << "/" << err << " [" << (dist / err) << "], ipsig = " << ipSig << ", dist = " << axisDist << ", w = " << state->weight() << std::endl;
         }
     }
 */
       ghostTrackRef = reco::TrackRef(ghostTrackRefProd, ghostTracks->size());
       ghostTracks->push_back(*ghostTrack);
 
       if (m_directionWithGhostTrack) {
         const reco::GhostTrackPrediction& pred = ghostTrack->prediction();
         double lambda = pred.lambda(origin);
         dummy = reco::Vertex(RecoVertex::convertPos(pred.position(lambda)),
                              RecoVertex::convertError(pred.positionError(lambda)),
                              0,
                              0,
                              0);
         pv = &dummy;
         direction = pred.direction();
       }
     }
 
     std::vector<float> prob2D, prob3D;
     std::vector<reco::btag::TrackIPData> ipData;
 
     for (unsigned int ind = 0; ind < transientTracks.size(); ind++) {
       const reco::TransientTrack& transientTrack = transientTracks[ind];
       const reco::Track& track = transientTrack.track();
 
       reco::btag::TrackIPData trackIP;
       trackIP.ip3d = IPTools::signedImpactParameter3D(transientTrack, direction, *pv).second;
       trackIP.ip2d = IPTools::signedTransverseImpactParameter(transientTrack, direction, *pv).second;
 
       TrajectoryStateOnSurface closest =
           IPTools::closestApproachToJet(transientTrack.impactPointState(), *pv, direction, transientTrack.field());
       if (closest.isValid())
         trackIP.closestToJetAxis = closest.globalPosition();
 
       // TODO: cross check if it is the same using other methods
       trackIP.distanceToJetAxis = IPTools::jetTrackDistance(transientTrack, direction, *pv).second;
 
       if (ghostTrack.get()) {
         const std::vector<reco::GhostTrackState>& states = ghostTrack->states();
         std::vector<reco::GhostTrackState>::const_iterator pos =
             std::find_if(states.begin(), states.end(), [&](auto& arg) { return arg.track() == transientTrack; });
 
         if (pos != states.end() && pos->isValid()) {
           VertexDistance3D dist;
           const reco::GhostTrackPrediction& pred = ghostTrack->prediction();
           GlobalPoint p1 = pos->tsos().globalPosition();
           GlobalError e1 = pos->tsos().cartesianError().position();
           GlobalPoint p2 = pred.position(pos->lambda());
           GlobalError e2 = pred.positionError(pos->lambda());
           trackIP.closestToGhostTrack = p1;
           trackIP.distanceToGhostTrack = dist.distance(VertexState(p1, e1), VertexState(p2, e2));
           trackIP.ghostTrackWeight = pos->weight();
         } else {
           trackIP.distanceToGhostTrack = Measurement1D(0., -1.);
           trackIP.ghostTrackWeight = 0.;
         }
       } else {
         trackIP.distanceToGhostTrack = Measurement1D(0., -1.);
         trackIP.ghostTrackWeight = 1.;
       }
 
       ipData.push_back(trackIP);
 
       if (m_computeProbabilities) {
         //probability with 3D ip
         std::pair<bool, double> probability = m_probabilityEstimator->probability(
             m_useTrackQuality, 0, ipData.back().ip3d.significance(), track, *(it->jet()), *pv);
         prob3D.push_back(probability.first ? probability.second : -1.);
 
         //probability with 2D ip
         probability = m_probabilityEstimator->probability(
             m_useTrackQuality, 1, ipData.back().ip2d.significance(), track, *(it->jet()), *pv);
         prob2D.push_back(probability.first ? probability.second : -1.);
       }
     }
 
     result->push_back(
         typename Product::value_type(ipData, prob2D, prob3D, selectedTracks, *it, pvRef, direction, ghostTrackRef));
   }
 
   if (m_computeGhostTrack)
     iEvent.put(std::move(ghostTracks), "ghostTracks");
   iEvent.put(std::move(result));
 }