GsfTools/src/MultiGaussianStateTransform.cc

0001 #include "TrackingTools/GsfTools/interface/MultiGaussianStateTransform.h"
0002 #include "TrackingTools/GsfTools/interface/GetComponents.h"
0003 #include "TrackingTools/GsfTools/interface/MultiGaussianState1D.h"
0004 #include "DataFormats/TrackReco/interface/TrackExtraFwd.h"
0005 #include "DataFormats/GsfTrackReco/interface/GsfTrackExtraFwd.h"
0006 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0007 #include "FWCore/Utilities/interface/Exception.h"
0008
0009 MultiGaussianState<MultiGaussianStateTransform::N> MultiGaussianStateTransform::outerMultiState(
0010     const reco::GsfTrack& tk) {
0011   const reco::GsfTrackExtraRef& extra(tk.gsfExtra());
0012   return multiState(extra->outerStateLocalParameters(), extra->outerStateCovariances(), extra->outerStateWeights());
0013 }
0014
0015 MultiGaussianState<MultiGaussianStateTransform::N> MultiGaussianStateTransform::innerMultiState(
0016     const reco::GsfTrack& tk) {
0017   const reco::GsfTrackExtraRef& extra(tk.gsfExtra());
0018   return multiState(extra->innerStateLocalParameters(), extra->innerStateCovariances(), extra->innerStateWeights());
0019 }
0020
0021 MultiGaussianState1D MultiGaussianStateTransform::outerMultiState1D(const reco::GsfTrack& tk, unsigned int index) {
0022   if (index >= N)
0023     throw cms::Exception("LogicError") << "MultiGaussianStateTransform: index out of range";
0024
0025   const reco::GsfTrackExtraRef& extra(tk.gsfExtra());
0026   return multiState1D(
0027       extra->outerStateLocalParameters(), extra->outerStateCovariances(), extra->outerStateWeights(), index);
0028 }
0029
0030 MultiGaussianState1D MultiGaussianStateTransform::innerMultiState1D(const reco::GsfTrack& tk, unsigned int index) {
0031   if (index >= N)
0032     throw cms::Exception("LogicError") << "MultiGaussianStateTransform: index out of range";
0033
0034   const reco::GsfTrackExtraRef& extra(tk.gsfExtra());
0035   return multiState1D(
0036       extra->innerStateLocalParameters(), extra->innerStateCovariances(), extra->innerStateWeights(), index);
0037 }
0038
0039 MultiGaussianState<MultiGaussianStateTransform::N> MultiGaussianStateTransform::multiState(
0040     const std::vector<MultiGaussianState<N>::Vector>& parameters,
0041     const std::vector<MultiGaussianState<N>::Matrix>& covariances,
0042     const std::vector<double>& weights) {
0043   unsigned int nc = parameters.size();
0044   MultiGaussianState<N>::SingleStateContainer components;
0045   components.reserve(nc);
0046   for (unsigned int i = 0; i < nc; ++i) {
0047     MultiGaussianState<N>::SingleStatePtr sgs(
0048         new MultiGaussianState<N>::SingleState(parameters[i], covariances[i], weights[i]));
0049     components.push_back(sgs);
0050   }
0051   return MultiGaussianState<N>(components);
0052 }
0053
0054 MultiGaussianState1D MultiGaussianStateTransform::multiState1D(
0055     const std::vector<MultiGaussianState<N>::Vector>& parameters,
0056     const std::vector<MultiGaussianState<N>::Matrix>& covariances,
0057     const std::vector<double>& weights,
0058     unsigned int index) {
0059   unsigned int nc = parameters.size();
0060   MultiGaussianState1D::SingleState1dContainer components;
0061   components.reserve(nc);
0062   for (unsigned int i = 0; i < nc; ++i) {
0063     components.push_back(SingleGaussianState1D(parameters[i](index), covariances[i](index, index), weights[i]));
0064   }
0065   return MultiGaussianState1D(components);
0066 }
0067
0068 MultiGaussianState<5> MultiGaussianStateTransform::multiState(const TrajectoryStateOnSurface tsos) {
0069   GetComponents comps(tsos);
0070   auto const& tsosComponents = comps();
0071   MultiGaussianState<5>::SingleStateContainer components;
0072   components.reserve(tsosComponents.size());
0073   for (auto i = tsosComponents.begin(); i != tsosComponents.end(); ++i) {
0074     MultiGaussianState<5>::SingleStatePtr sgs(
0075         new MultiGaussianState<5>::SingleState(i->localParameters().vector(), i->localError().matrix(), i->weight()));
0076     components.push_back(sgs);
0077   }
0078   return MultiGaussianState<5>(components);
0079 }
0080
0081 MultiGaussianState1D MultiGaussianStateTransform::multiState1D(const TrajectoryStateOnSurface tsos,
0082                                                                unsigned int index) {
0083   if (index >= N)
0084     throw cms::Exception("LogicError") << "MultiGaussianStateTransform: index out of range";
0085   GetComponents comps(tsos);
0086   auto const& tsosComponents = comps();
0087   MultiGaussianState1D::SingleState1dContainer components;
0088   components.reserve(tsosComponents.size());
0089   for (auto i = tsosComponents.begin(); i != tsosComponents.end(); ++i) {
0090     components.push_back(SingleGaussianState1D(
0091         i->localParameters().vector()(index), i->localError().matrix()(index, index), i->weight()));
0092   }
0093   return MultiGaussianState1D(components);
0094 }
0095
0096 TrajectoryStateOnSurface MultiGaussianStateTransform::tsosFromSingleState(const SingleGaussianState<5>& singleState,
0097                                                                           const TrajectoryStateOnSurface refTsos) {
0098   const LocalTrajectoryParameters& refPars(refTsos.localParameters());
0099   double pzSign = refPars.pzSign();
0100   bool charged = refPars.charge() != 0;
0101   LocalTrajectoryParameters pars(singleState.mean(), pzSign, charged);
0102   LocalTrajectoryError errs(singleState.covariance());
0103   // return state (doesn't use weight of the single state)
0104   return TrajectoryStateOnSurface(pars, errs, refTsos.surface(), refTsos.magneticField());
0105 }