Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 10:05:47

0001 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
0002 #include "DataFormats/Common/interface/Handle.h"
0003 #include "TrackingTools/TransientTrack/interface/GsfTransientTrack.h"
0004 #include "TrackingTools/TransientTrack/interface/TransientTrackFromFTS.h"
0005 #include "FWCore/Utilities/interface/isFinite.h"
0006 
0007 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0008 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
0009 
0010 #include "DataFormats/PatCandidates/interface/PackedCandidate.h"
0011 
0012 using namespace reco;
0013 using namespace std;
0014 using namespace edm;
0015 
0016 namespace {
0017   constexpr float defaultInvalidTrackReso = 0.350f;
0018 }
0019 
0020 TransientTrack TransientTrackBuilder::build(const Track* t) const {
0021   return TransientTrack(*t, theField, theTrackingGeometry);
0022 }
0023 
0024 TransientTrack TransientTrackBuilder::build(const Track& t) const {
0025   return TransientTrack(t, theField, theTrackingGeometry);
0026 }
0027 
0028 TransientTrack TransientTrackBuilder::build(const GsfTrack* t) const {
0029   return TransientTrack(new GsfTransientTrack(*t, theField, theTrackingGeometry));
0030 }
0031 
0032 TransientTrack TransientTrackBuilder::build(const GsfTrack& t) const {
0033   return TransientTrack(new GsfTransientTrack(t, theField, theTrackingGeometry));
0034 }
0035 
0036 TransientTrack TransientTrackBuilder::build(const CandidatePtr* t) const {
0037   reco::PFCandidatePtr tryPF(*t);
0038   edm::Ptr<pat::PackedCandidate> tryPacked(*t);
0039   if (tryPF.get() != nullptr && tryPF->isTimeValid()) {
0040     return TransientTrack(*t, tryPF->time(), tryPF->timeError(), theField, theTrackingGeometry);
0041   } else if (tryPacked.get() != nullptr && tryPacked->timeError() > 0.f) {
0042     return TransientTrack(*t, (double)tryPacked->time(), (double)tryPacked->timeError(), theField, theTrackingGeometry);
0043   }
0044   return TransientTrack(*t, theField, theTrackingGeometry);
0045 }
0046 
0047 TransientTrack TransientTrackBuilder::build(const CandidatePtr& t) const { return this->build(&t); }
0048 
0049 TransientTrack TransientTrackBuilder::build(const TrackRef* t) const {
0050   return TransientTrack(*t, theField, theTrackingGeometry);
0051 }
0052 
0053 TransientTrack TransientTrackBuilder::build(const TrackRef& t) const {
0054   return TransientTrack(t, theField, theTrackingGeometry);
0055 }
0056 
0057 TransientTrack TransientTrackBuilder::build(const GsfTrackRef* t) const {
0058   return TransientTrack(new GsfTransientTrack(*t, theField, theTrackingGeometry));
0059 }
0060 
0061 TransientTrack TransientTrackBuilder::build(const GsfTrackRef& t) const {
0062   return TransientTrack(new GsfTransientTrack(t, theField, theTrackingGeometry));
0063 }
0064 
0065 vector<TransientTrack> TransientTrackBuilder::build(const edm::Handle<reco::TrackCollection>& trkColl) const {
0066   vector<TransientTrack> ttVect;
0067   ttVect.reserve((*trkColl).size());
0068   for (unsigned int i = 0; i < (*trkColl).size(); i++) {
0069     ttVect.push_back(TransientTrack(TrackRef(trkColl, i), theField, theTrackingGeometry));
0070   }
0071   return ttVect;
0072 }
0073 
0074 vector<TransientTrack> TransientTrackBuilder::build(const edm::Handle<reco::GsfTrackCollection>& trkColl) const {
0075   vector<TransientTrack> ttVect;
0076   ttVect.reserve((*trkColl).size());
0077   for (unsigned int i = 0; i < (*trkColl).size(); i++) {
0078     ttVect.push_back(TransientTrack(new GsfTransientTrack(GsfTrackRef(trkColl, i), theField, theTrackingGeometry)));
0079   }
0080   return ttVect;
0081 }
0082 
0083 vector<TransientTrack> TransientTrackBuilder::build(const edm::Handle<edm::View<Track> >& trkColl) const {
0084   vector<TransientTrack> ttVect;
0085   ttVect.reserve((*trkColl).size());
0086   for (unsigned int i = 0; i < (*trkColl).size(); i++) {
0087     const Track* trk = &(*trkColl)[i];
0088     const GsfTrack* gsfTrack = dynamic_cast<const GsfTrack*>(trk);
0089     if (gsfTrack) {
0090       ttVect.push_back(TransientTrack(
0091           new GsfTransientTrack(RefToBase<Track>(trkColl, i).castTo<GsfTrackRef>(), theField, theTrackingGeometry)));
0092     } else {  // gsf
0093       ttVect.push_back(TransientTrack(RefToBase<Track>(trkColl, i).castTo<TrackRef>(), theField, theTrackingGeometry));
0094     }
0095   }
0096   return ttVect;
0097 }
0098 
0099 vector<TransientTrack> TransientTrackBuilder::build(const edm::Handle<reco::TrackCollection>& trkColl,
0100                                                     const edm::ValueMap<float>& trackTimes,
0101                                                     const edm::ValueMap<float>& trackTimeResos) const {
0102   vector<TransientTrack> ttVect;
0103   ttVect.reserve((*trkColl).size());
0104   for (unsigned int i = 0; i < (*trkColl).size(); i++) {
0105     TrackRef ref(trkColl, i);
0106     double time = trackTimes[ref];
0107     double timeReso = trackTimeResos[ref];
0108     timeReso =
0109         (timeReso > 1e-6 ? timeReso : defaultInvalidTrackReso);  // make the error much larger than the BS time width
0110     if (edm::isNotFinite(time)) {
0111       time = 0.0;
0112       timeReso = defaultInvalidTrackReso;
0113     }
0114     ttVect.push_back(TransientTrack(ref, time, timeReso, theField, theTrackingGeometry));
0115   }
0116   return ttVect;
0117 }
0118 
0119 vector<TransientTrack> TransientTrackBuilder::build(const edm::Handle<reco::GsfTrackCollection>& trkColl,
0120                                                     const edm::ValueMap<float>& trackTimes,
0121                                                     const edm::ValueMap<float>& trackTimeResos) const {
0122   vector<TransientTrack> ttVect;
0123   ttVect.reserve((*trkColl).size());
0124   for (unsigned int i = 0; i < (*trkColl).size(); i++) {
0125     GsfTrackRef ref(trkColl, i);
0126     double time = trackTimes[ref];
0127     double timeReso = trackTimeResos[ref];
0128     timeReso =
0129         (timeReso > 1e-6 ? timeReso : defaultInvalidTrackReso);  // make the error much larger than the BS time width
0130     if (edm::isNotFinite(time)) {
0131       time = 0.0;
0132       timeReso = defaultInvalidTrackReso;
0133     }
0134     ttVect.push_back(TransientTrack(new GsfTransientTrack(ref, time, timeReso, theField, theTrackingGeometry)));
0135   }
0136   return ttVect;
0137 }
0138 
0139 vector<TransientTrack> TransientTrackBuilder::build(const edm::Handle<edm::View<Track> >& trkColl,
0140                                                     const edm::ValueMap<float>& trackTimes,
0141                                                     const edm::ValueMap<float>& trackTimeResos) const {
0142   vector<TransientTrack> ttVect;
0143   ttVect.reserve((*trkColl).size());
0144   for (unsigned int i = 0; i < (*trkColl).size(); i++) {
0145     const Track* trk = &(*trkColl)[i];
0146     const GsfTrack* gsfTrack = dynamic_cast<const GsfTrack*>(trk);
0147     if (gsfTrack) {
0148       GsfTrackRef ref = RefToBase<Track>(trkColl, i).castTo<GsfTrackRef>();
0149       double time = trackTimes[ref];
0150       double timeReso = trackTimeResos[ref];
0151       timeReso =
0152           (timeReso > 1e-6 ? timeReso : defaultInvalidTrackReso);  // make the error much larger than the BS time width
0153       if (edm::isNotFinite(time)) {
0154         time = 0.0;
0155         timeReso = defaultInvalidTrackReso;
0156       }
0157       ttVect.push_back(TransientTrack(new GsfTransientTrack(
0158           RefToBase<Track>(trkColl, i).castTo<GsfTrackRef>(), time, timeReso, theField, theTrackingGeometry)));
0159     } else {  // gsf
0160       TrackRef ref = RefToBase<Track>(trkColl, i).castTo<TrackRef>();
0161       double time = trackTimes[ref];
0162       double timeReso = trackTimeResos[ref];
0163       timeReso =
0164           (timeReso > 1e-6 ? timeReso : defaultInvalidTrackReso);  // make the error much larger than the BS time width
0165       if (edm::isNotFinite(time)) {
0166         time = 0.0;
0167         timeReso = defaultInvalidTrackReso;
0168       }
0169       ttVect.push_back(TransientTrack(
0170           RefToBase<Track>(trkColl, i).castTo<TrackRef>(), time, timeReso, theField, theTrackingGeometry));
0171     }
0172   }
0173   return ttVect;
0174 }
0175 
0176 vector<TransientTrack> TransientTrackBuilder::build(const edm::Handle<reco::TrackCollection>& trkColl,
0177                                                     const reco::BeamSpot& beamSpot) const {
0178   vector<TransientTrack> ttVect = build(trkColl);
0179   for (unsigned int i = 0; i < ttVect.size(); i++) {
0180     ttVect[i].setBeamSpot(beamSpot);
0181   }
0182   return ttVect;
0183 }
0184 
0185 vector<TransientTrack> TransientTrackBuilder::build(const edm::Handle<reco::GsfTrackCollection>& trkColl,
0186                                                     const reco::BeamSpot& beamSpot) const {
0187   vector<TransientTrack> ttVect = build(trkColl);
0188   for (unsigned int i = 0; i < ttVect.size(); i++) {
0189     ttVect[i].setBeamSpot(beamSpot);
0190   }
0191   return ttVect;
0192 }
0193 
0194 vector<TransientTrack> TransientTrackBuilder::build(const edm::Handle<edm::View<Track> >& trkColl,
0195                                                     const reco::BeamSpot& beamSpot) const {
0196   vector<TransientTrack> ttVect = build(trkColl);
0197   for (unsigned int i = 0; i < ttVect.size(); i++) {
0198     ttVect[i].setBeamSpot(beamSpot);
0199   }
0200   return ttVect;
0201 }
0202 
0203 vector<TransientTrack> TransientTrackBuilder::build(const edm::Handle<reco::TrackCollection>& trkColl,
0204                                                     const reco::BeamSpot& beamSpot,
0205                                                     const edm::ValueMap<float>& trackTimes,
0206                                                     const edm::ValueMap<float>& trackTimeResos) const {
0207   vector<TransientTrack> ttVect = build(trkColl, trackTimes, trackTimeResos);
0208   for (unsigned int i = 0; i < ttVect.size(); i++) {
0209     ttVect[i].setBeamSpot(beamSpot);
0210   }
0211   return ttVect;
0212 }
0213 
0214 vector<TransientTrack> TransientTrackBuilder::build(const edm::Handle<reco::GsfTrackCollection>& trkColl,
0215                                                     const reco::BeamSpot& beamSpot,
0216                                                     const edm::ValueMap<float>& trackTimes,
0217                                                     const edm::ValueMap<float>& trackTimeResos) const {
0218   vector<TransientTrack> ttVect = build(trkColl, trackTimes, trackTimeResos);
0219   for (unsigned int i = 0; i < ttVect.size(); i++) {
0220     ttVect[i].setBeamSpot(beamSpot);
0221   }
0222   return ttVect;
0223 }
0224 
0225 vector<TransientTrack> TransientTrackBuilder::build(const edm::Handle<edm::View<Track> >& trkColl,
0226                                                     const reco::BeamSpot& beamSpot,
0227                                                     const edm::ValueMap<float>& trackTimes,
0228                                                     const edm::ValueMap<float>& trackTimeResos) const {
0229   vector<TransientTrack> ttVect = build(trkColl, trackTimes, trackTimeResos);
0230   for (unsigned int i = 0; i < ttVect.size(); i++) {
0231     ttVect[i].setBeamSpot(beamSpot);
0232   }
0233   return ttVect;
0234 }
0235 
0236 TransientTrack TransientTrackBuilder::build(const FreeTrajectoryState& fts) const {
0237   return TransientTrack(new TransientTrackFromFTS(fts));
0238 }