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.
10 #include <memory>
11 #include <string>
12 #include <cmath>
13 
15  : conf_(conf), verbose_(0), dvf_(0), ptMin_(1.0)
16 {
17  // Register my product
18  produces<reco::VertexCollection>();
19 
20  // Setup shop
21  verbose_ = conf.getParameter<int>("Verbosity"); // 0 silent, 1 chatty, 2 loud
22  std::string finder = conf.getParameter<std::string>("Finder"); // DivisiveVertexFinder
23  bool useError = conf.getParameter<bool>("UseError"); // true
24  bool wtAverage = conf.getParameter<bool>("WtAverage"); // true
25  double zOffset = conf.getParameter<double>("ZOffset"); // 5.0 sigma
26  double zSeparation = conf.getParameter<double>("ZSeparation"); // 0.05 cm
27  int ntrkMin = conf.getParameter<int>("NTrkMin"); // 3
28  // Tracking requirements before sending a track to be considered for vtx
29  ptMin_ = conf_.getParameter<double>("PtMin"); // 1.0 GeV
30 
31  if (finder == "DivisiveVertexFinder") {
32  if (verbose_ > 0) edm::LogInfo("PixelVertexProducer") << ": Using the DivisiveVertexFinder\n";
33  dvf_ = new DivisiveVertexFinder(zOffset, ntrkMin, useError, zSeparation, wtAverage, verbose_);
34  }
35  else { // Finder not supported, or you made a mistake in your request
36  // throw an exception once I figure out how CMSSW does this
37  }
38 }
39 
40 
42  delete dvf_;
43 }
44 
46 
47  // First fish the pixel tracks out of the event
48  edm::Handle<reco::TrackCollection> trackCollection;
49  edm::InputTag trackCollName = conf_.getParameter<edm::InputTag>("TrackCollection");
50  e.getByLabel(trackCollName,trackCollection);
51  const reco::TrackCollection tracks = *(trackCollection.product());
52  if (verbose_ > 0) edm::LogInfo("PixelVertexProducer") << ": Found " << tracks.size() << " tracks in TrackCollection called " << trackCollName << "\n";
53 
54 
55  // Second, make a collection of pointers to the tracks we want for the vertex finder
57  for (unsigned int i=0; i<tracks.size(); i++) {
58  if (tracks[i].pt() > ptMin_)
59  trks.push_back( reco::TrackRef(trackCollection, i) );
60  }
61  if (verbose_ > 0) edm::LogInfo("PixelVertexProducer") << ": Selected " << trks.size() << " of these tracks for vertexing\n";
62 
64  edm::InputTag bsName = conf_.getParameter<edm::InputTag>("beamSpot");
65  e.getByLabel(bsName,bsHandle);
66  math::XYZPoint myPoint(0.,0.,0.);
67  if (bsHandle.isValid()) myPoint = math::XYZPoint(bsHandle->x0(),bsHandle->y0(), 0. ); //FIXME: fix last coordinate with vertex.z() at same time
68 
69  // Third, ship these tracks off to be vertexed
70  std::auto_ptr<reco::VertexCollection> vertexes(new reco::VertexCollection);
71  bool ok;
72  if (conf_.getParameter<bool>("Method2")) {
73  ok = dvf_->findVertexesAlt(trks, // input
74  *vertexes,myPoint); // output
75  if (verbose_ > 0) edm::LogInfo("PixelVertexProducer") << "Method2 returned status of " << ok;
76  }
77  else {
78  ok = dvf_->findVertexes(trks, // input
79  *vertexes); // output
80  if (verbose_ > 0) edm::LogInfo("PixelVertexProducer") << "Method1 returned status of " << ok;
81  }
82 
83  if (verbose_ > 0) {
84  edm::LogInfo("PixelVertexProducer") << ": Found " << vertexes->size() << " vertexes\n";
85  for (unsigned int i=0; i<vertexes->size(); ++i) {
86  edm::LogInfo("PixelVertexProducer") << "Vertex number " << i << " has " << (*vertexes)[i].tracksSize() << " tracks with a position of " << (*vertexes)[i].z() << " +- " << std::sqrt( (*vertexes)[i].covariance(2,2) );
87  }
88  }
89 
90 
91  if(bsHandle.isValid())
92  {
93  const reco::BeamSpot & bs = *bsHandle;
94 
95  for (unsigned int i=0; i<vertexes->size(); ++i) {
96  double z=(*vertexes)[i].z();
97  double x=bs.x0()+bs.dxdz()*(z-bs.z0());
98  double y=bs.y0()+bs.dydz()*(z-bs.z0());
99  reco::Vertex v( reco::Vertex::Point(x,y,z), (*vertexes)[i].error(),(*vertexes)[i].chi2() , (*vertexes)[i].ndof() , (*vertexes)[i].tracksSize());
100  //Copy also the tracks
101  for (std::vector<reco::TrackBaseRef >::const_iterator it = (*vertexes)[i].tracks_begin();
102  it !=(*vertexes)[i].tracks_end(); it++) {
103  v.add( *it );
104  }
105  (*vertexes)[i]=v;
106 
107  }
108  }
109  else
110  {
111  edm::LogWarning("PixelVertexProducer") << "No beamspot found. Using returning vertexes with (0,0,Z) ";
112  }
113 
114  if(!vertexes->size() && bsHandle.isValid()){
115 
116  const reco::BeamSpot & bs = *bsHandle;
117 
119  if ( (bse.cxx() <= 0.) ||
120  (bse.cyy() <= 0.) ||
121  (bse.czz() <= 0.) ) {
123  we(0,0)=10000;
124  we(1,1)=10000;
125  we(2,2)=10000;
126  vertexes->push_back(reco::Vertex(bs.position(), we,0.,0.,0));
127 
128  edm::LogInfo("PixelVertexProducer") <<"No vertices found. Beamspot with invalid errors " << bse.matrix() << std::endl
129  << "Will put Vertex derived from dummy-fake BeamSpot into Event.\n"
130  << (*vertexes)[0].x() << "\n"
131  << (*vertexes)[0].y() << "\n"
132  << (*vertexes)[0].z() << "\n";
133  } else {
134  vertexes->push_back(reco::Vertex(bs.position(),
135  bs.rotatedCovariance3D(),0.,0.,0));
136 
137  edm::LogInfo("PixelVertexProducer") << "No vertices found. Will put Vertex derived from BeamSpot into Event:\n"
138  << (*vertexes)[0].x() << "\n"
139  << (*vertexes)[0].y() << "\n"
140  << (*vertexes)[0].z() << "\n";
141  }
142  }
143 
144  else if(!vertexes->size() && !bsHandle.isValid())
145  {
146  edm::LogWarning("PixelVertexProducer") << "No beamspot and no vertex found. No vertex returned.";
147  }
148 
149  e.put(vertexes);
150 
151 }
152 
153 
T getParameter(std::string const &) const
double z0() const
z coordinate
Definition: BeamSpot.h:69
int i
Definition: DBlmapReader.cc:9
virtual void produce(edm::Event &, const edm::EventSetup &)
PixelVertexProducer(const edm::ParameterSet &)
std::vector< Track > TrackCollection
collection of Tracks
Definition: TrackFwd.h:10
bool findVertexesAlt(const reco::TrackRefVector &trks, reco::VertexCollection &vertexes, const math::XYZPoint &bs)
std::vector< Vertex > VertexCollection
collection of Vertex objects
Definition: VertexFwd.h:9
edm::ParameterSet conf_
double double double z
double dydz() const
dydz slope
Definition: BeamSpot.h:85
ROOT::Math::SMatrix< double, 3, 3, ROOT::Math::MatRepSym< double, 3 > > AlgebraicSymMatrix33
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:85
T sqrt(T t)
Definition: SSEVec.h:46
math::XYZPoint Point
point in the space
Definition: Vertex.h:40
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:76
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:356
double dxdz() const
dxdz slope
Definition: BeamSpot.h:83
tuple conf
Definition: dbtoconf.py:185
void add(const TrackBaseRef &r, float w=1.0)
add a reference to a Track
tuple tracks
Definition: testEve_cfg.py:39
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:13
T const * product() const
Definition: Handle.h:74
DivisiveVertexFinder * dvf_
double y0() const
y coordinate
Definition: BeamSpot.h:67
void push_back(value_type const &ref)
Add a Ref&lt;C, T&gt; to the RefVector.
Definition: RefVector.h:64
size_type size() const
Size of the RefVector.
Definition: RefVector.h:89
const Point & position() const
position
Definition: BeamSpot.h:63
Definition: DDAxes.h:10
Covariance3DMatrix rotatedCovariance3D() const
Definition: BeamSpot.cc:79
mathSSE::Vec4< T > v
double x0() const
x coordinate
Definition: BeamSpot.h:65