CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_1_8_patch9/src/PhysicsTools/RecoUtils/src/CandCommonVertexFitter.cc

Go to the documentation of this file.
00001 #include "PhysicsTools/RecoUtils/interface/CandCommonVertexFitter.h"
00002 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
00003 #include "DataFormats/Candidate/interface/VertexCompositeCandidate.h"
00004 #include "FWCore/Utilities/interface/EDMException.h"
00005 #include <sstream>
00006 using namespace reco;
00007 using namespace std;
00008 
00009 void CandCommonVertexFitterBase::set(VertexCompositeCandidate & c) const {
00010   if(bField_ == 0)
00011     throw edm::Exception(edm::errors::InvalidReference)
00012       << "B-Field was not set up CandCommonVertexFitter.\n"
00013       << "the following method must be called before fitting a candidate:\n"
00014       << " CandCommonVertexFitter:.set( const MagneticField * )" << endl;
00015   vector<TransientTrack> tracks;
00016   vector<Candidate *> daughters;
00017   vector<RecoCandidate::TrackType> trackTypes;
00018   fill(tracks, daughters, trackTypes, c);
00019   assert(tracks.size() == daughters.size());
00020   TransientVertex vertex;
00021   if(fit(vertex, tracks)) {
00022     tracks = vertex.refittedTracks();    
00023     Candidate::Point vtx(vertex.position());
00024     c.setVertex(vtx);
00025     vector<TransientTrack>::const_iterator trackIt = tracks.begin(), tracksEnd = tracks.end();
00026     vector<Candidate *>::const_iterator daughterIt = daughters.begin();
00027     vector<RecoCandidate::TrackType>::const_iterator trackTypeIt = trackTypes.begin();
00028     Candidate::LorentzVector mp4(0, 0, 0, 0);
00029     for(; trackIt != tracksEnd; ++ trackIt, ++ daughterIt, ++trackTypeIt) {
00030       const Track & track = trackIt->track();
00031       Candidate & daughter = * * daughterIt;
00032       double px = track.px(), py = track.py(), pz = track.pz(), p = track.p();
00033       double energy;
00034       daughter.setVertex( vtx );
00035       if(*trackTypeIt == RecoCandidate::recoTrackType) {
00036         double mass = daughter.mass();
00037         energy = sqrt(p*p + mass*mass);
00038       } else {
00039         energy = daughter.energy();
00040         double scale = energy / p;
00041         px *= scale; py *= scale; pz *= scale; 
00042       }
00043       Candidate::LorentzVector dp4(px, py, pz, energy);
00044       daughter.setP4(dp4);
00045       mp4 += dp4;
00046     }
00047     c.setP4(mp4);
00048     Vertex v = vertex;
00049     c.setChi2AndNdof(chi2_ = v.chi2(), ndof_ = v.ndof());
00050     v.fill(cov_);
00051     c.setCovariance(cov_);
00052   } else {
00053     c.setChi2AndNdof(chi2_ = -1, ndof_ = 0);
00054     c.setCovariance(cov_ = CovarianceMatrix(ROOT::Math::SMatrixIdentity())); 
00055   }
00056 }
00057 
00058 void CandCommonVertexFitterBase::fill(vector<TransientTrack> & tracks, 
00059                                       vector<Candidate *> & daughters,
00060                                       vector<RecoCandidate::TrackType> & trackTypes,
00061                                       Candidate & c) const {
00062   size_t nDau = c.numberOfDaughters();
00063   for(unsigned int j = 0; j < nDau ; ++j) {
00064     Candidate * d = c.daughter(j);
00065     if(d == 0) {
00066       ostringstream message;
00067       message << "Can't access in write mode candidate daughters. "
00068               << "pdgId = " << c.pdgId() << ".\n";
00069       const Candidate * d1 = c.daughter(j);
00070       if(d1 == 0)
00071         message << "Null daughter also found in read-only mode\n";
00072       else
00073         message << "Daughter found in read-only mode with id: " << d1->pdgId() << "\n";
00074       throw edm::Exception(edm::errors::InvalidReference) << message.str();
00075     }
00076     if(d->numberOfDaughters() > 0)
00077       fill(tracks, daughters, trackTypes, * d);
00078     else {
00079       const Track * trk = d->get<const Track *>();
00080       RecoCandidate::TrackType type = d->get<RecoCandidate::TrackType>();
00081       if(trk != 0) {
00082         tracks.push_back(TransientTrack(* trk, bField_));
00083         daughters.push_back(d);
00084         trackTypes.push_back(type);
00085       } else {
00086         cerr << ">>> warning: candidate of type " << d->pdgId() 
00087              << " has no track reference." << endl;
00088       }
00089     }
00090   }
00091 }