File indexing completed on 2023-03-17 11:23:31
0001 #include "RecoVertex/TrimmedKalmanVertexFinder/interface/ConfigurableTrimmedVertexFinder.h"
0002 #include "CommonTools/Statistics/interface/ChiSquaredProbability.h"
0003
0004 using namespace reco;
0005
0006 ConfigurableTrimmedVertexFinder::ConfigurableTrimmedVertexFinder(const VertexFitter<5>* vf,
0007 const VertexUpdator<5>* vu,
0008 const VertexTrackCompatibilityEstimator<5>* ve)
0009 : theClusterFinder(vf, vu, ve),
0010 theVtxFitProbCut(0.01),
0011 theTrackCompatibilityToPV(0.05),
0012 theTrackCompatibilityToSV(0.01),
0013 theMaxNbOfVertices(0) {
0014
0015 theFilter.setPtCut(1.5);
0016 }
0017
0018 void ConfigurableTrimmedVertexFinder::setParameters(const edm::ParameterSet& s) {
0019 theFilter.setPtCut(s.getParameter<double>("ptCut"));
0020 theTrackCompatibilityToPV = s.getParameter<double>("trackCompatibilityToPVcut");
0021 theTrackCompatibilityToSV = s.getParameter<double>("trackCompatibilityToSVcut");
0022 theVtxFitProbCut = s.getParameter<double>("vtxFitProbCut");
0023 theMaxNbOfVertices = s.getParameter<int>("maxNbOfVertices");
0024 }
0025
0026 std::vector<TransientVertex> ConfigurableTrimmedVertexFinder::vertices(const std::vector<TransientTrack>& tracks) const {
0027 std::vector<TransientTrack> remaining;
0028
0029 return vertices(tracks, remaining, reco::BeamSpot(), false);
0030 }
0031
0032 std::vector<TransientVertex> ConfigurableTrimmedVertexFinder::vertices(const std::vector<TransientTrack>& tracks,
0033 const reco::BeamSpot& spot) const {
0034 std::vector<TransientTrack> remaining;
0035 return vertices(tracks, remaining, spot, true);
0036 }
0037
0038 std::vector<TransientVertex> ConfigurableTrimmedVertexFinder::vertices(const std::vector<TransientTrack>& tracks,
0039 std::vector<TransientTrack>& unused,
0040 const reco::BeamSpot& spot,
0041 bool use_spot) const {
0042 resetEvent(tracks);
0043 analyseInputTracks(tracks);
0044
0045 std::vector<TransientTrack> filtered;
0046 for (std::vector<TransientTrack>::const_iterator it = tracks.begin(); it != tracks.end(); it++) {
0047 if (theFilter(*it)) {
0048 filtered.push_back(*it);
0049 } else {
0050 unused.push_back(*it);
0051 }
0052 }
0053
0054 std::vector<TransientVertex> all = vertexCandidates(filtered, unused, spot, use_spot);
0055
0056 analyseVertexCandidates(all);
0057
0058 std::vector<TransientVertex> sel = clean(all);
0059
0060 analyseFoundVertices(sel);
0061
0062 return sel;
0063 }
0064
0065 std::vector<TransientVertex> ConfigurableTrimmedVertexFinder::vertexCandidates(const std::vector<TransientTrack>& tracks,
0066 std::vector<TransientTrack>& unused,
0067 const reco::BeamSpot& spot,
0068 bool use_spot) const {
0069 std::vector<TransientVertex> cand;
0070
0071 std::vector<TransientTrack> remain = tracks;
0072
0073 while (true) {
0074 float tkCompCut = (cand.empty() ? theTrackCompatibilityToPV : theTrackCompatibilityToSV);
0075
0076
0077 theClusterFinder.setTrackCompatibilityCut(tkCompCut);
0078
0079
0080
0081 std::vector<TransientVertex> newVertices;
0082 if (cand.empty() && use_spot) {
0083 newVertices = theClusterFinder.vertices(remain, spot);
0084 } else {
0085 newVertices = theClusterFinder.vertices(remain);
0086 }
0087 if (newVertices.empty())
0088 break;
0089
0090 analyseClusterFinder(newVertices, remain);
0091
0092 for (std::vector<TransientVertex>::const_iterator iv = newVertices.begin(); iv != newVertices.end(); iv++) {
0093 if (iv->originalTracks().size() > 1) {
0094 cand.push_back(*iv);
0095 } else {
0096
0097 for (std::vector<TransientTrack>::const_iterator trk = iv->originalTracks().begin();
0098 trk != iv->originalTracks().end();
0099 ++trk) {
0100 unused.push_back(*trk);
0101 }
0102 }
0103 }
0104
0105
0106 if (theMaxNbOfVertices != 0) {
0107 if (cand.size() >= (unsigned int)theMaxNbOfVertices)
0108 break;
0109 }
0110 }
0111
0112 for (std::vector<TransientTrack>::const_iterator it = remain.begin(); it != remain.end(); it++) {
0113 unused.push_back(*it);
0114 }
0115
0116 return cand;
0117 }
0118
0119 std::vector<TransientVertex> ConfigurableTrimmedVertexFinder::clean(
0120 const std::vector<TransientVertex>& candidates) const {
0121 std::vector<TransientVertex> sel;
0122 for (std::vector<TransientVertex>::const_iterator i = candidates.begin(); i != candidates.end(); i++) {
0123 if (ChiSquaredProbability((*i).totalChiSquared(), (*i).degreesOfFreedom()) > theVtxFitProbCut) {
0124 sel.push_back(*i);
0125 }
0126 }
0127
0128 return sel;
0129 }