Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 11:23:33

0001 #include "RecoVertex/VertexPrimitives/interface/TransientVertex.h"
0002 #include "RecoVertex/VertexPrimitives/interface/VertexException.h"
0003 #include "TrackingTools/TransientTrack/interface/TrackTransientTrack.h"
0004 #include "TrackingTools/TransientTrack/interface/CandidatePtrTransientTrack.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 #include <algorithm>
0007 
0008 using namespace std;
0009 using namespace reco;
0010 
0011 TransientVertex::TransientVertex()
0012     : theVertexState(),
0013       theOriginalTracks(),
0014       theChi2(0),
0015       theNDF(0),
0016       vertexValid(false),
0017       withPrior(false),
0018       theWeightMapIsAvailable(false),
0019       theCovMapAvailable(false),
0020       withRefittedTracks(false) {}
0021 
0022 TransientVertex::TransientVertex(const GlobalPoint& pos,
0023                                  const GlobalError& posError,
0024                                  const std::vector<TransientTrack>& tracks,
0025                                  float chi2)
0026     : theVertexState(pos, posError),
0027       theOriginalTracks(tracks),
0028       theChi2(chi2),
0029       theNDF(0),
0030       vertexValid(true),
0031       withPrior(false),
0032       theWeightMapIsAvailable(false),
0033       theCovMapAvailable(false),
0034       withRefittedTracks(false) {
0035   theNDF = 2. * theOriginalTracks.size() - 3.;
0036 }
0037 
0038 TransientVertex::TransientVertex(const GlobalPoint& pos,
0039                                  const double time,
0040                                  const GlobalError& posTimeError,
0041                                  const std::vector<TransientTrack>& tracks,
0042                                  float chi2)
0043     : theVertexState(pos, time, posTimeError),
0044       theOriginalTracks(tracks),
0045       theChi2(chi2),
0046       theNDF(0),
0047       vertexValid(true),
0048       withPrior(false),
0049       theWeightMapIsAvailable(false),
0050       theCovMapAvailable(false),
0051       withRefittedTracks(false) {
0052   theNDF = 2. * theOriginalTracks.size() - 3.;
0053 }
0054 
0055 TransientVertex::TransientVertex(const GlobalPoint& pos,
0056                                  const GlobalError& posError,
0057                                  const std::vector<TransientTrack>& tracks,
0058                                  float chi2,
0059                                  float ndf)
0060     : theVertexState(pos, posError),
0061       theOriginalTracks(tracks),
0062       theChi2(chi2),
0063       theNDF(ndf),
0064       vertexValid(true),
0065       withPrior(false),
0066       theWeightMapIsAvailable(false),
0067       theCovMapAvailable(false),
0068       withRefittedTracks(false) {}
0069 
0070 TransientVertex::TransientVertex(const GlobalPoint& pos,
0071                                  const double time,
0072                                  const GlobalError& posTimeError,
0073                                  const std::vector<TransientTrack>& tracks,
0074                                  float chi2,
0075                                  float ndf)
0076     : theVertexState(pos, time, posTimeError),
0077       theOriginalTracks(tracks),
0078       theChi2(chi2),
0079       theNDF(ndf),
0080       vertexValid(true),
0081       withPrior(false),
0082       theWeightMapIsAvailable(false),
0083       theCovMapAvailable(false),
0084       withRefittedTracks(false) {}
0085 
0086 TransientVertex::TransientVertex(const GlobalPoint& priorPos,
0087                                  const GlobalError& priorErr,
0088                                  const GlobalPoint& pos,
0089                                  const GlobalError& posError,
0090                                  const std::vector<TransientTrack>& tracks,
0091                                  float chi2)
0092     : thePriorVertexState(priorPos, priorErr),
0093       theVertexState(pos, posError),
0094       theOriginalTracks(tracks),
0095       theChi2(chi2),
0096       theNDF(0),
0097       vertexValid(true),
0098       withPrior(true),
0099       theWeightMapIsAvailable(false),
0100       theCovMapAvailable(false),
0101       withRefittedTracks(false) {
0102   theNDF = 2. * theOriginalTracks.size();
0103 }
0104 
0105 TransientVertex::TransientVertex(const GlobalPoint& priorPos,
0106                                  const double priorTime,
0107                                  const GlobalError& priorErr,
0108                                  const GlobalPoint& pos,
0109                                  const double time,
0110                                  const GlobalError& posError,
0111                                  const std::vector<TransientTrack>& tracks,
0112                                  float chi2)
0113     : thePriorVertexState(priorPos, priorTime, priorErr),
0114       theVertexState(pos, time, posError),
0115       theOriginalTracks(tracks),
0116       theChi2(chi2),
0117       theNDF(0),
0118       vertexValid(true),
0119       withPrior(true),
0120       theWeightMapIsAvailable(false),
0121       theCovMapAvailable(false),
0122       withRefittedTracks(false) {
0123   theNDF = 2. * theOriginalTracks.size();
0124 }
0125 
0126 TransientVertex::TransientVertex(const GlobalPoint& priorPos,
0127                                  const GlobalError& priorErr,
0128                                  const GlobalPoint& pos,
0129                                  const GlobalError& posError,
0130                                  const std::vector<TransientTrack>& tracks,
0131                                  float chi2,
0132                                  float ndf)
0133     : thePriorVertexState(priorPos, priorErr),
0134       theVertexState(pos, posError),
0135       theOriginalTracks(tracks),
0136       theChi2(chi2),
0137       theNDF(ndf),
0138       vertexValid(true),
0139       withPrior(true),
0140       theWeightMapIsAvailable(false),
0141       theCovMapAvailable(false),
0142       withRefittedTracks(false) {}
0143 
0144 TransientVertex::TransientVertex(const GlobalPoint& priorPos,
0145                                  const double priorTime,
0146                                  const GlobalError& priorErr,
0147                                  const GlobalPoint& pos,
0148                                  const double time,
0149                                  const GlobalError& posError,
0150                                  const std::vector<TransientTrack>& tracks,
0151                                  float chi2,
0152                                  float ndf)
0153     : thePriorVertexState(priorPos, priorTime, priorErr),
0154       theVertexState(pos, time, posError),
0155       theOriginalTracks(tracks),
0156       theChi2(chi2),
0157       theNDF(ndf),
0158       vertexValid(true),
0159       withPrior(true),
0160       theWeightMapIsAvailable(false),
0161       theCovMapAvailable(false),
0162       withRefittedTracks(false) {}
0163 
0164 TransientVertex::TransientVertex(const VertexState& state, const std::vector<TransientTrack>& tracks, float chi2)
0165     : theVertexState(state),
0166       theOriginalTracks(tracks),
0167       theChi2(chi2),
0168       theNDF(0),
0169       vertexValid(true),
0170       withPrior(false),
0171       theWeightMapIsAvailable(false),
0172       theCovMapAvailable(false),
0173       withRefittedTracks(false) {
0174   theNDF = 2. * theOriginalTracks.size() - 3.;
0175 }
0176 
0177 TransientVertex::TransientVertex(const VertexState& state,
0178                                  const std::vector<TransientTrack>& tracks,
0179                                  float chi2,
0180                                  float ndf)
0181     : theVertexState(state),
0182       theOriginalTracks(tracks),
0183       theChi2(chi2),
0184       theNDF(ndf),
0185       vertexValid(true),
0186       withPrior(false),
0187       theWeightMapIsAvailable(false),
0188       theCovMapAvailable(false),
0189       withRefittedTracks(false) {}
0190 
0191 TransientVertex::TransientVertex(const VertexState& prior,
0192                                  const VertexState& state,
0193                                  const std::vector<TransientTrack>& tracks,
0194                                  float chi2)
0195     : thePriorVertexState(prior),
0196       theVertexState(state),
0197       theOriginalTracks(tracks),
0198       theChi2(chi2),
0199       theNDF(0),
0200       vertexValid(true),
0201       withPrior(true),
0202       theWeightMapIsAvailable(false),
0203       theCovMapAvailable(false),
0204       withRefittedTracks(false) {
0205   theNDF = 2. * theOriginalTracks.size();
0206 }
0207 
0208 TransientVertex::TransientVertex(const VertexState& prior,
0209                                  const VertexState& state,
0210                                  const std::vector<TransientTrack>& tracks,
0211                                  float chi2,
0212                                  float ndf)
0213     : thePriorVertexState(prior),
0214       theVertexState(state),
0215       theOriginalTracks(tracks),
0216       theChi2(chi2),
0217       theNDF(ndf),
0218       vertexValid(true),
0219       withPrior(true),
0220       theWeightMapIsAvailable(false),
0221       theCovMapAvailable(false),
0222       withRefittedTracks(false) {}
0223 
0224 void TransientVertex::weightMap(const TransientTrackToFloatMap& theMap) {
0225   theWeightMap = theMap;
0226   theWeightMapIsAvailable = true;
0227 }
0228 
0229 void TransientVertex::refittedTracks(const std::vector<reco::TransientTrack>& refittedTracks) {
0230   if (refittedTracks.empty())
0231     throw VertexException("TransientVertex::refittedTracks: No refitted tracks stored in input container");
0232   theRefittedTracks = refittedTracks;
0233   withRefittedTracks = true;
0234 }
0235 
0236 void TransientVertex::tkToTkCovariance(const TTtoTTmap& covMap) {
0237   theCovMap = covMap;
0238   theCovMapAvailable = true;
0239 }
0240 
0241 float TransientVertex::trackWeight(const TransientTrack& track) const {
0242   if (!theWeightMapIsAvailable) {
0243     std::vector<TransientTrack>::const_iterator foundTrack =
0244         find(theOriginalTracks.begin(), theOriginalTracks.end(), track);
0245     return ((foundTrack != theOriginalTracks.end()) ? 1. : 0.);
0246   }
0247   TransientTrackToFloatMap::const_iterator it = theWeightMap.find(track);
0248   if (it != theWeightMap.end()) {
0249     return (it->second);
0250   }
0251   return 0.;
0252 }
0253 
0254 AlgebraicMatrix33 TransientVertex::tkToTkCovariance(const TransientTrack& t1, const TransientTrack& t2) const {
0255   if (!theCovMapAvailable) {
0256     throw VertexException("TransientVertex::Track-to-track covariance matrices not available");
0257   }
0258   const TransientTrack* tr1;
0259   const TransientTrack* tr2;
0260   if (t1 < t2) {
0261     tr1 = &t1;
0262     tr2 = &t2;
0263   } else {
0264     tr1 = &t2;
0265     tr2 = &t1;
0266   }
0267   TTtoTTmap::const_iterator it = theCovMap.find(*tr1);
0268   if (it != theCovMap.end()) {
0269     const TTmap& tm = it->second;
0270     TTmap::const_iterator nit = tm.find(*tr2);
0271     if (nit != tm.end()) {
0272       return (nit->second);
0273     } else {
0274       throw VertexException("TransientVertex::requested Track-to-track covariance matrix does not exist");
0275     }
0276   } else {
0277     throw VertexException("TransientVertex::requested Track-to-track covariance matrix does not exist");
0278   }
0279 }
0280 
0281 TransientTrack TransientVertex::originalTrack(const TransientTrack& refTrack) const {
0282   if (theRefittedTracks.empty())
0283     throw VertexException("TransientVertex::requested No refitted tracks stored in vertex");
0284   std::vector<TransientTrack>::const_iterator it = find(theRefittedTracks.begin(), theRefittedTracks.end(), refTrack);
0285   if (it == theRefittedTracks.end())
0286     throw VertexException(
0287         "TransientVertex::requested Refitted track not found in list.\n address used for comparison: ")
0288         << refTrack.basicTransientTrack();
0289   size_t pos = it - theRefittedTracks.begin();
0290   return theOriginalTracks[pos];
0291 }
0292 
0293 TransientTrack TransientVertex::refittedTrack(const TransientTrack& track) const {
0294   if (theRefittedTracks.empty())
0295     throw VertexException("TransientVertex::requested No refitted tracks stored in vertex");
0296   std::vector<TransientTrack>::const_iterator it = find(theOriginalTracks.begin(), theOriginalTracks.end(), track);
0297   if (it == theOriginalTracks.end())
0298     throw VertexException("transientVertex::requested Track not found in list.\n address used for comparison: ")
0299         << track.basicTransientTrack();
0300   size_t pos = it - theOriginalTracks.begin();
0301   return theRefittedTracks[pos];
0302 }
0303 
0304 TransientVertex::operator reco::Vertex() const {
0305   //If the vertex is invalid, return an invalid TV !
0306   if (!isValid())
0307     return Vertex();
0308 
0309   Vertex vertex(Vertex::Point(theVertexState.position()),
0310                 //  RecoVertex::convertError(theVertexState.error()),
0311                 theVertexState.error4D().matrix4D(),
0312                 theVertexState.time(),
0313                 totalChiSquared(),
0314                 degreesOfFreedom(),
0315                 theOriginalTracks.size());
0316   for (std::vector<TransientTrack>::const_iterator i = theOriginalTracks.begin(); i != theOriginalTracks.end(); ++i) {
0317     //     const TrackTransientTrack* ttt = dynamic_cast<const TrackTransientTrack*>((*i).basicTransientTrack());
0318     //     if ((ttt!=0) && (ttt->persistentTrackRef().isNonnull()))
0319     //     {
0320     //       TrackRef tr = ttt->persistentTrackRef();
0321     //       TrackBaseRef tbr(tr);
0322     if (withRefittedTracks) {
0323       vertex.add((*i).trackBaseRef(), refittedTrack(*i).track(), trackWeight(*i));
0324     } else {
0325       vertex.add((*i).trackBaseRef(), trackWeight(*i));
0326     }
0327     //}
0328   }
0329   return vertex;
0330 }
0331 
0332 TransientVertex::operator reco::VertexCompositePtrCandidate() const {
0333   using namespace reco;
0334   if (!isValid())
0335     return VertexCompositePtrCandidate();
0336 
0337   VertexCompositePtrCandidate vtxCompPtrCand;
0338 
0339   vtxCompPtrCand.setTime(vertexState().time());
0340   vtxCompPtrCand.setCovariance(vertexState().error4D().matrix4D());
0341   vtxCompPtrCand.setChi2AndNdof(totalChiSquared(), degreesOfFreedom());
0342   vtxCompPtrCand.setVertex(Candidate::Point(position().x(), position().y(), position().z()));
0343 
0344   Candidate::LorentzVector p4;
0345   for (std::vector<reco::TransientTrack>::const_iterator tt = theOriginalTracks.begin(); tt != theOriginalTracks.end();
0346        ++tt) {
0347     if (trackWeight(*tt) < 0.5)
0348       continue;
0349 
0350     const CandidatePtrTransientTrack* cptt = dynamic_cast<const CandidatePtrTransientTrack*>(tt->basicTransientTrack());
0351     if (cptt == nullptr)
0352       edm::LogError("DynamicCastingFailed") << "Casting of TransientTrack to CandidatePtrTransientTrack failed!";
0353     else {
0354       p4 += cptt->candidate()->p4();
0355       vtxCompPtrCand.addDaughter(cptt->candidate());
0356     }
0357   }
0358 
0359   //TODO: if has refitted tracks we should scale the candidate p4 to the refitted one
0360   vtxCompPtrCand.setP4(p4);
0361   return vtxCompPtrCand;
0362 }