File indexing completed on 2023-03-17 11:09:42
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011 #include <memory>
0012 #include "FWCore/Framework/interface/Frameworkfwd.h"
0013 #include "FWCore/Framework/interface/global/EDProducer.h"
0014 #include "FWCore/Framework/interface/Event.h"
0015 #include "FWCore/Framework/interface/MakerMacros.h"
0016 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0017 #include "DataFormats/Common/interface/AssociationMap.h"
0018 #include "DataFormats/Common/interface/getRef.h"
0019 #include "DataFormats/Common/interface/ValueMap.h"
0020 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0021 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
0022 #include "DataFormats/TrackReco/interface/HitPattern.h"
0023 #include "DataFormats/TrackReco/interface/Track.h"
0024 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0025 #include "DataFormats/TrackReco/interface/TrackBase.h"
0026 #include "DataFormats/VertexReco/interface/Vertex.h"
0027 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0028 #include "DataFormats/Scouting/interface/Run3ScoutingMuon.h"
0029 #include "DataFormats/Scouting/interface/Run3ScoutingTrack.h"
0030 #include "DataFormats/Scouting/interface/Run3ScoutingVertex.h"
0031 #include "DataFormats/MuonReco/interface/MuonTrackLinks.h"
0032 #include "DataFormats/MuonReco/interface/MuonFwd.h"
0033 #include "DataFormats/Math/interface/libminifloat.h"
0034
0035 class HLTScoutingTrackProducer : public edm::global::EDProducer<> {
0036 public:
0037 explicit HLTScoutingTrackProducer(const edm::ParameterSet&);
0038 ~HLTScoutingTrackProducer() override;
0039
0040 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0041
0042 private:
0043 void produce(edm::StreamID sid, edm::Event& iEvent, edm::EventSetup const& setup) const final;
0044
0045 const edm::EDGetTokenT<reco::TrackCollection> otherTrackCollection_;
0046 const edm::EDGetTokenT<reco::VertexCollection> vertexCollection_;
0047
0048 const int mantissaPrecision_;
0049 const double vtxMinDist_;
0050 const double ptMin_;
0051 };
0052
0053
0054
0055
0056 HLTScoutingTrackProducer::HLTScoutingTrackProducer(const edm::ParameterSet& iConfig)
0057 : otherTrackCollection_(consumes<reco::TrackCollection>(iConfig.getParameter<edm::InputTag>("OtherTracks"))),
0058 vertexCollection_(consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("vertexCollection"))),
0059 mantissaPrecision_(iConfig.getParameter<int>("mantissaPrecision")),
0060 vtxMinDist_(iConfig.getParameter<double>("vtxMinDist")),
0061 ptMin_(iConfig.getParameter<double>("ptMin")) {
0062
0063 produces<Run3ScoutingTrackCollection>();
0064 }
0065
0066 HLTScoutingTrackProducer::~HLTScoutingTrackProducer() = default;
0067
0068
0069 void HLTScoutingTrackProducer::produce(edm::StreamID sid, edm::Event& iEvent, edm::EventSetup const& setup) const {
0070 using namespace edm;
0071
0072 std::unique_ptr<Run3ScoutingTrackCollection> outTrack(new Run3ScoutingTrackCollection());
0073
0074 Handle<reco::TrackCollection> otherTrackCollection;
0075 Handle<reco::VertexCollection> vertexCollection;
0076
0077 if (iEvent.getByToken(otherTrackCollection_, otherTrackCollection)) {
0078
0079 for (auto& trk : *otherTrackCollection) {
0080 int vtxInd = -1;
0081 double min_dist = vtxMinDist_;
0082 int vtxIt = 0;
0083
0084 if (trk.pt() < ptMin_)
0085 continue;
0086
0087 if (iEvent.getByToken(vertexCollection_, vertexCollection)) {
0088 for (auto& vrt : *vertexCollection) {
0089 double min_dist_tmp = pow(trk.dz(vrt.position()), 2);
0090
0091 if (min_dist_tmp < min_dist) {
0092 min_dist = min_dist_tmp;
0093 vtxInd = vtxIt;
0094 }
0095
0096 vtxIt++;
0097 }
0098 }
0099
0100
0101 outTrack->emplace_back(
0102 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.pt(), mantissaPrecision_),
0103 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.eta(), mantissaPrecision_),
0104 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.phi(), mantissaPrecision_),
0105 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.chi2(), mantissaPrecision_),
0106 trk.ndof(),
0107 trk.charge(),
0108 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.dxy(), mantissaPrecision_),
0109 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.dz(), mantissaPrecision_),
0110 trk.hitPattern().numberOfValidPixelHits(),
0111 trk.hitPattern().trackerLayersWithMeasurement(),
0112 trk.hitPattern().numberOfValidStripHits(),
0113 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.qoverp(), mantissaPrecision_),
0114 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.lambda(), mantissaPrecision_),
0115 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.dxyError(), mantissaPrecision_),
0116 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.dzError(), mantissaPrecision_),
0117 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.qoverpError(), mantissaPrecision_),
0118 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.lambdaError(), mantissaPrecision_),
0119 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.phiError(), mantissaPrecision_),
0120 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.dsz(), mantissaPrecision_),
0121 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.dszError(), mantissaPrecision_),
0122 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.covariance(0, 1), mantissaPrecision_),
0123 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.covariance(0, 2), mantissaPrecision_),
0124 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.covariance(0, 3), mantissaPrecision_),
0125 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.covariance(0, 4), mantissaPrecision_),
0126 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.covariance(1, 2), mantissaPrecision_),
0127 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.covariance(1, 3), mantissaPrecision_),
0128 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.covariance(1, 4), mantissaPrecision_),
0129 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.covariance(2, 3), mantissaPrecision_),
0130 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.covariance(2, 4), mantissaPrecision_),
0131 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.covariance(3, 4), mantissaPrecision_),
0132 vtxInd,
0133 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.vx(), mantissaPrecision_),
0134 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.vy(), mantissaPrecision_),
0135 MiniFloatConverter::reduceMantissaToNbitsRounding(trk.vz(), mantissaPrecision_));
0136 }
0137 }
0138
0139 iEvent.put(std::move(outTrack));
0140 }
0141
0142
0143 void HLTScoutingTrackProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0144 edm::ParameterSetDescription desc;
0145 desc.add<edm::InputTag>("OtherTracks", edm::InputTag("hltPixelTracksL3MuonNoVtx"));
0146 desc.add<edm::InputTag>("vertexCollection", edm::InputTag("hltPixelVertices"));
0147
0148 desc.add<int>("mantissaPrecision", 10)->setComment("default float16, change to 23 for float32");
0149 desc.add<double>("vtxMinDist", 0.01);
0150 desc.add<double>("ptMin", 0.3);
0151 descriptions.add("hltScoutingTrackProducer", desc);
0152 }
0153
0154
0155 #include "FWCore/Framework/interface/MakerMacros.h"
0156 DEFINE_FWK_MODULE(HLTScoutingTrackProducer);