CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
PixelVertexProducer.cc
Go to the documentation of this file.
9 #include <memory>
10 #include <string>
11 #include <cmath>
12 
14  // 0 silent, 1 chatty, 2 loud
15  : verbose_(conf.getParameter<int>("Verbosity")),
16  // 1.0 GeV
17  ptMin_(conf.getParameter<double>("PtMin")),
18  method2(conf.getParameter<bool>("Method2")),
19  trackCollName(conf.getParameter<edm::InputTag>("TrackCollection")),
20  token_Tracks(consumes<reco::TrackCollection>(trackCollName)),
21  token_BeamSpot(consumes<reco::BeamSpot>(conf.getParameter<edm::InputTag>("beamSpot"))) {
22  // Register my product
23  produces<reco::VertexCollection>();
24 
25  // Setup shop
26  std::string finder = conf.getParameter<std::string>("Finder"); // DivisiveVertexFinder
27  bool useError = conf.getParameter<bool>("UseError"); // true
28  bool wtAverage = conf.getParameter<bool>("WtAverage"); // true
29  double zOffset = conf.getParameter<double>("ZOffset"); // 5.0 sigma
30  double zSeparation = conf.getParameter<double>("ZSeparation"); // 0.05 cm
31  int ntrkMin = conf.getParameter<int>("NTrkMin"); // 3
32  // Tracking requirements before sending a track to be considered for vtx
33 
34  double track_pt_min = ptMin_;
35  double track_pt_max = 10.;
36  double track_chi2_max = 9999999.;
37  double track_prob_min = -1.;
38 
39  if (conf.exists("PVcomparer")) {
40  edm::ParameterSet PVcomparerPSet = conf.getParameter<edm::ParameterSet>("PVcomparer");
41  track_pt_min = PVcomparerPSet.getParameter<double>("track_pt_min");
42  if (track_pt_min != ptMin_) {
43  if (track_pt_min < ptMin_)
44  edm::LogInfo("PixelVertexProducer")
45  << "minimum track pT setting differs between PixelVertexProducer (" << ptMin_ << ") and PVcomparer ("
46  << track_pt_min << ") [PVcomparer considers tracks w/ lower threshold than PixelVertexProducer does] !!!";
47  else
48  edm::LogInfo("PixelVertexProducer") << "minimum track pT setting differs between PixelVertexProducer ("
49  << ptMin_ << ") and PVcomparer (" << track_pt_min << ") !!!";
50  }
51  track_pt_max = PVcomparerPSet.getParameter<double>("track_pt_max");
52  track_chi2_max = PVcomparerPSet.getParameter<double>("track_chi2_max");
53  track_prob_min = PVcomparerPSet.getParameter<double>("track_prob_min");
54  }
55 
56  if (finder == "DivisiveVertexFinder") {
57  if (verbose_ > 0)
58  edm::LogInfo("PixelVertexProducer") << ": Using the DivisiveVertexFinder\n";
59  dvf_ = new DivisiveVertexFinder(track_pt_min,
60  track_pt_max,
61  track_chi2_max,
62  track_prob_min,
63  zOffset,
64  ntrkMin,
65  useError,
66  zSeparation,
67  wtAverage,
68  verbose_);
69  } else { // Finder not supported, or you made a mistake in your request
70  // throw an exception once I figure out how CMSSW does this
71  }
72 }
73 
75 
77  // First fish the pixel tracks out of the event
79  e.getByToken(token_Tracks, trackCollection);
80  const reco::TrackCollection tracks = *(trackCollection.product());
81  if (verbose_ > 0)
82  edm::LogInfo("PixelVertexProducer") << ": Found " << tracks.size() << " tracks in TrackCollection called "
83  << trackCollName << "\n";
84 
85  // Second, make a collection of pointers to the tracks we want for the vertex finder
87  for (unsigned int i = 0; i < tracks.size(); i++) {
88  if (tracks[i].pt() > ptMin_)
89  trks.push_back(reco::TrackRef(trackCollection, i));
90  }
91  if (verbose_ > 0)
92  edm::LogInfo("PixelVertexProducer") << ": Selected " << trks.size() << " of these tracks for vertexing\n";
93 
95  e.getByToken(token_BeamSpot, bsHandle);
96  math::XYZPoint myPoint(0., 0., 0.);
97  if (bsHandle.isValid())
98  //FIXME: fix last coordinate with vertex.z() at same time
99  myPoint = math::XYZPoint(bsHandle->x0(), bsHandle->y0(), 0.);
100 
101  // Third, ship these tracks off to be vertexed
102  auto vertexes = std::make_unique<reco::VertexCollection>();
103  bool ok;
104  if (method2) {
105  ok = dvf_->findVertexesAlt(trks, // input
106  *vertexes,
107  myPoint); // output
108  if (verbose_ > 0)
109  edm::LogInfo("PixelVertexProducer") << "Method2 returned status of " << ok;
110  } else {
111  ok = dvf_->findVertexes(trks, // input
112  *vertexes); // output
113  if (verbose_ > 0)
114  edm::LogInfo("PixelVertexProducer") << "Method1 returned status of " << ok;
115  }
116 
117  if (verbose_ > 0) {
118  edm::LogInfo("PixelVertexProducer") << ": Found " << vertexes->size() << " vertexes\n";
119  for (unsigned int i = 0; i < vertexes->size(); ++i) {
120  edm::LogInfo("PixelVertexProducer")
121  << "Vertex number " << i << " has " << (*vertexes)[i].tracksSize() << " tracks with a position of "
122  << (*vertexes)[i].z() << " +- " << std::sqrt((*vertexes)[i].covariance(2, 2));
123  }
124  }
125 
126  if (bsHandle.isValid()) {
127  const reco::BeamSpot& bs = *bsHandle;
128 
129  for (unsigned int i = 0; i < vertexes->size(); ++i) {
130  double z = (*vertexes)[i].z();
131  double x = bs.x0() + bs.dxdz() * (z - bs.z0());
132  double y = bs.y0() + bs.dydz() * (z - bs.z0());
134  (*vertexes)[i].error(),
135  (*vertexes)[i].chi2(),
136  (*vertexes)[i].ndof(),
137  (*vertexes)[i].tracksSize());
138  //Copy also the tracks
139  for (std::vector<reco::TrackBaseRef>::const_iterator it = (*vertexes)[i].tracks_begin();
140  it != (*vertexes)[i].tracks_end();
141  it++) {
142  v.add(*it);
143  }
144  (*vertexes)[i] = v;
145  }
146  } else {
147  edm::LogWarning("PixelVertexProducer") << "No beamspot found. Using returning vertexes with (0,0,Z) ";
148  }
149 
150  if (vertexes->empty() && bsHandle.isValid()) {
151  const reco::BeamSpot& bs = *bsHandle;
152 
154  if ((bse.cxx() <= 0.) || (bse.cyy() <= 0.) || (bse.czz() <= 0.)) {
156  we(0, 0) = 10000;
157  we(1, 1) = 10000;
158  we(2, 2) = 10000;
159  vertexes->push_back(reco::Vertex(bs.position(), we, 0., 0., 0));
160 
161  edm::LogInfo("PixelVertexProducer")
162  << "No vertices found. Beamspot with invalid errors " << bse.matrix() << std::endl
163  << "Will put Vertex derived from dummy-fake BeamSpot into Event.\n"
164  << (*vertexes)[0].x() << "\n"
165  << (*vertexes)[0].y() << "\n"
166  << (*vertexes)[0].z() << "\n";
167  } else {
168  vertexes->push_back(reco::Vertex(bs.position(), bs.rotatedCovariance3D(), 0., 0., 0));
169 
170  edm::LogInfo("PixelVertexProducer") << "No vertices found. Will put Vertex derived from BeamSpot into Event:\n"
171  << (*vertexes)[0].x() << "\n"
172  << (*vertexes)[0].y() << "\n"
173  << (*vertexes)[0].z() << "\n";
174  }
175  }
176 
177  else if (vertexes->empty() && !bsHandle.isValid()) {
178  edm::LogWarning("PixelVertexProducer") << "No beamspot and no vertex found. No vertex returned.";
179  }
180 
181  e.put(std::move(vertexes));
182 }
T getParameter(std::string const &) const
double z0() const
z coordinate
Definition: BeamSpot.h:65
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:131
const edm::EDGetTokenT< reco::TrackCollection > token_Tracks
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:525
PixelVertexProducer(const edm::ParameterSet &)
std::vector< Track > TrackCollection
collection of Tracks
Definition: TrackFwd.h:14
bool findVertexesAlt(const reco::TrackRefVector &trks, reco::VertexCollection &vertexes, const math::XYZPoint &bs)
bool exists(std::string const &parameterName) const
checks if a parameter exists
const edm::InputTag trackCollName
const edm::EDGetTokenT< reco::BeamSpot > token_BeamSpot
double dydz() const
dydz slope
Definition: BeamSpot.h:80
T sqrt(T t)
Definition: SSEVec.h:19
math::XYZPoint Point
point in the space
Definition: Vertex.h:39
bool findVertexes(const reco::TrackRefVector &trks, reco::VertexCollection &vertexes)
Run the divisive algorithm and return a vector of vertexes for the input track collection.
bool isValid() const
Definition: HandleBase.h:70
void add(const TrackBaseRef &r, float w=1.0)
add a reference to a Track
Definition: Vertex.cc:79
double dxdz() const
dxdz slope
Definition: BeamSpot.h:78
T const * product() const
Definition: Handle.h:69
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
DivisiveVertexFinder * dvf_
fixed size matrix
HLT enums.
double y0() const
y coordinate
Definition: BeamSpot.h:63
ROOT::Math::SMatrix< double, 3, 3, ROOT::Math::MatRepSym< double, 3 > > AlgebraicSymMatrix33
void push_back(value_type const &ref)
Add a Ref<C, T> to the RefVector.
Definition: RefVector.h:67
size_type size() const
Size of the RefVector.
Definition: RefVector.h:102
const Point & position() const
position
Definition: BeamSpot.h:59
void produce(edm::Event &, const edm::EventSetup &) override
Covariance3DMatrix rotatedCovariance3D() const
Definition: BeamSpot.cc:73
def move(src, dest)
Definition: eostools.py:511
double x0() const
x coordinate
Definition: BeamSpot.h:61