CMS 3D CMS Logo

/data/doxygen/doxygen-1.7.3/gen/CMSSW_4_2_8/src/RecoPixelVertexing/PixelVertexFinding/src/PixelVertexProducer.cc

Go to the documentation of this file.
00001 #include "RecoPixelVertexing/PixelVertexFinding/interface/PixelVertexProducer.h"
00002 #include "FWCore/Framework/interface/MakerMacros.h"
00003 #include "DataFormats/VertexReco/interface/Vertex.h"
00004 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
00005 #include "DataFormats/VertexReco/interface/VertexFwd.h"
00006 #include "RecoPixelVertexing/PixelVertexFinding/interface/DivisiveVertexFinder.h"
00007 #include "FWCore/Utilities/interface/InputTag.h"
00008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00009 #include <memory>
00010 #include <string>
00011 #include <cmath>
00012 
00013 PixelVertexProducer::PixelVertexProducer(const edm::ParameterSet& conf) 
00014   : conf_(conf), verbose_(0), dvf_(0), ptMin_(1.0)
00015 {
00016   // Register my product
00017   produces<reco::VertexCollection>();
00018 
00019   // Setup shop
00020   verbose_           = conf.getParameter<int>("Verbosity"); // 0 silent, 1 chatty, 2 loud
00021   std::string finder = conf.getParameter<std::string>("Finder"); // DivisiveVertexFinder
00022   bool useError      = conf.getParameter<bool>("UseError"); // true
00023   bool wtAverage     = conf.getParameter<bool>("WtAverage"); // true
00024   double zOffset     = conf.getParameter<double>("ZOffset"); // 5.0 sigma
00025   double zSeparation = conf.getParameter<double>("ZSeparation"); // 0.05 cm
00026   int ntrkMin        = conf.getParameter<int>("NTrkMin"); // 3
00027   // Tracking requirements before sending a track to be considered for vtx
00028   ptMin_ = conf_.getParameter<double>("PtMin"); // 1.0 GeV
00029 
00030   if (finder == "DivisiveVertexFinder") {
00031     if (verbose_ > 0) edm::LogInfo("PixelVertexProducer") << ": Using the DivisiveVertexFinder\n";
00032     dvf_ = new DivisiveVertexFinder(zOffset, ntrkMin, useError, zSeparation, wtAverage, verbose_);
00033   }
00034   else { // Finder not supported, or you made a mistake in your request
00035     // throw an exception once I figure out how CMSSW does this
00036   }
00037 }
00038 
00039 
00040 PixelVertexProducer::~PixelVertexProducer() {
00041   delete dvf_;
00042 }
00043 
00044 void PixelVertexProducer::produce(edm::Event& e, const edm::EventSetup& es) {
00045 
00046   // First fish the pixel tracks out of the event
00047   edm::Handle<reco::TrackCollection> trackCollection;
00048   edm::InputTag trackCollName = conf_.getParameter<edm::InputTag>("TrackCollection");
00049   e.getByLabel(trackCollName,trackCollection);
00050   const reco::TrackCollection tracks = *(trackCollection.product());
00051   if (verbose_ > 0) edm::LogInfo("PixelVertexProducer") << ": Found " << tracks.size() << " tracks in TrackCollection called " << trackCollName << "\n";
00052   
00053 
00054   // Second, make a collection of pointers to the tracks we want for the vertex finder
00055   reco::TrackRefVector trks;
00056   for (unsigned int i=0; i<tracks.size(); i++) {
00057     if (tracks[i].pt() > ptMin_)     
00058       trks.push_back( reco::TrackRef(trackCollection, i) );
00059   }
00060   if (verbose_ > 0) edm::LogInfo("PixelVertexProducer") << ": Selected " << trks.size() << " of these tracks for vertexing\n";
00061 
00062   edm::Handle<reco::BeamSpot> bsHandle;
00063   edm::InputTag bsName = conf_.getParameter<edm::InputTag>("beamSpot");
00064   e.getByLabel(bsName,bsHandle);
00065   math::XYZPoint myPoint(0.,0.,0.);
00066   if (bsHandle.isValid()) myPoint = math::XYZPoint(bsHandle->x0(),bsHandle->y0(), 0. ); //FIXME: fix last coordinate with vertex.z() at same time
00067 
00068   // Third, ship these tracks off to be vertexed
00069   std::auto_ptr<reco::VertexCollection> vertexes(new reco::VertexCollection);
00070   bool ok;
00071   if (conf_.getParameter<bool>("Method2")) {
00072     ok = dvf_->findVertexesAlt(trks,       // input
00073                              *vertexes,myPoint); // output
00074     if (verbose_ > 0) edm::LogInfo("PixelVertexProducer") << "Method2 returned status of " << ok;
00075   }
00076   else {
00077     ok = dvf_->findVertexes(trks,       // input
00078                             *vertexes); // output
00079     if (verbose_ > 0) edm::LogInfo("PixelVertexProducer") << "Method1 returned status of " << ok;
00080   }
00081 
00082   if (verbose_ > 0) {
00083     edm::LogInfo("PixelVertexProducer") << ": Found " << vertexes->size() << " vertexes\n";
00084     for (unsigned int i=0; i<vertexes->size(); ++i) {
00085       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) );
00086     }
00087   }
00088 
00089 
00090   if(bsHandle.isValid())
00091    {
00092     const reco::BeamSpot & bs = *bsHandle;
00093 
00094     for (unsigned int i=0; i<vertexes->size(); ++i) {
00095       double z=(*vertexes)[i].z();
00096       double x=bs.x0()+bs.dxdz()*(z-bs.z0());
00097       double y=bs.y0()+bs.dydz()*(z-bs.z0()); 
00098       reco::Vertex v( reco::Vertex::Point(x,y,z), (*vertexes)[i].error(),(*vertexes)[i].chi2() , (*vertexes)[i].ndof() , (*vertexes)[i].tracksSize());
00099        //Copy also the tracks 
00100        for (std::vector<reco::TrackBaseRef >::const_iterator it = (*vertexes)[i].tracks_begin();
00101             it !=(*vertexes)[i].tracks_end(); it++) {
00102             v.add( *it );
00103        }
00104       (*vertexes)[i]=v;
00105 
00106     }
00107   }
00108    else
00109   {
00110       edm::LogWarning("PixelVertexProducer") << "No beamspot found. Using returning vertexes with (0,0,Z) ";
00111   } 
00112 
00113   e.put(vertexes);
00114 }