File indexing completed on 2025-02-05 23:51:40
0001 #include "PhysicsTools/PatAlgos/interface/PATPrimaryVertexSelector.h"
0002 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0003 #include "DataFormats/TrackReco/interface/Track.h"
0004 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0005
0006 PATPrimaryVertexSelector::PATPrimaryVertexSelector(const edm::ParameterSet& cfg, edm::ConsumesCollector&& iC)
0007 : multiplicityCut_(cfg.getParameter<unsigned int>("minMultiplicity")),
0008 ptSumCut_(cfg.getParameter<double>("minPtSum")),
0009 trackEtaCut_(cfg.getParameter<double>("maxTrackEta")),
0010 chi2Cut_(cfg.getParameter<double>("maxNormChi2")),
0011 dr2Cut_(cfg.getParameter<double>("maxDeltaR")),
0012 dzCut_(cfg.getParameter<double>("maxDeltaZ")) {
0013
0014 dr2Cut_ *= dr2Cut_;
0015 }
0016
0017 void PATPrimaryVertexSelector::select(const edm::Handle<collection>& handle,
0018 const edm::Event& event,
0019 const edm::EventSetup& setup) {
0020 selected_.clear();
0021
0022 const math::XYZPoint beamSpot(0., 0., 0.);
0023 unsigned int multiplicity;
0024 double ptSum;
0025 for (collection::const_iterator iv = handle->begin(); iv != handle->end(); ++iv) {
0026 math::XYZVector displacement(iv->position() - beamSpot);
0027 if (iv->normalizedChi2() < chi2Cut_ && fabs(displacement.z()) < dzCut_ && displacement.perp2() < dr2Cut_) {
0028 getVertexVariables(*iv, multiplicity, ptSum);
0029 if (multiplicity >= multiplicityCut_ && ptSum > ptSumCut_)
0030 selected_.push_back(&*iv);
0031 }
0032 }
0033 sort(selected_.begin(), selected_.end(), *this);
0034 }
0035
0036 void PATPrimaryVertexSelector::fillPSetDescription(edm::ParameterSetDescription& desc) {
0037 desc.add<unsigned int>("minMultiplicity", 1);
0038 desc.add<double>("minPtSum", 0.);
0039 desc.add<double>("maxTrackEta", 9999.);
0040 desc.add<double>("maxNormChi2", 9999.);
0041 desc.add<double>("maxDeltaR", 9999.);
0042 desc.add<double>("maxDeltaZ", 9999.);
0043 }
0044
0045 bool PATPrimaryVertexSelector::operator()(const reco::Vertex* v1, const reco::Vertex* v2) const {
0046 unsigned int mult1;
0047 double ptSum1;
0048 getVertexVariables(*v1, mult1, ptSum1);
0049 unsigned int mult2;
0050 double ptSum2;
0051 getVertexVariables(*v2, mult2, ptSum2);
0052 return ptSum1 > ptSum2;
0053 }
0054
0055 void PATPrimaryVertexSelector::getVertexVariables(const reco::Vertex& vertex,
0056 unsigned int& multiplicity,
0057 double& ptSum) const {
0058 multiplicity = 0;
0059 ptSum = 0.;
0060 for (reco::Vertex::trackRef_iterator it = vertex.tracks_begin(); it != vertex.tracks_end(); ++it) {
0061 if (fabs((**it).eta()) < trackEtaCut_) {
0062 ++multiplicity;
0063 ptSum += (**it).pt();
0064 }
0065 }
0066 }