File indexing completed on 2024-04-06 12:18:32
0001
0002
0003
0004
0005
0006 #include "FWCore/Framework/interface/Frameworkfwd.h"
0007 #include "FWCore/Framework/interface/global/EDProducer.h"
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/Framework/interface/MakerMacros.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0011 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0012 #include "DataFormats/VertexReco/interface/Vertex.h"
0013 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0014 #include "DataFormats/Scouting/interface/Run3ScoutingVertex.h"
0015 #include "DataFormats/Math/interface/libminifloat.h"
0016
0017 #include <memory>
0018 #include <utility>
0019
0020 class HLTScoutingPrimaryVertexProducer : public edm::global::EDProducer<> {
0021 public:
0022 explicit HLTScoutingPrimaryVertexProducer(const edm::ParameterSet&);
0023
0024 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0025
0026 private:
0027 void produce(edm::StreamID sid, edm::Event& iEvent, edm::EventSetup const& setup) const final;
0028
0029 edm::EDGetTokenT<reco::VertexCollection> const vertexCollToken_;
0030 int const mantissaPrecision_;
0031 };
0032
0033 HLTScoutingPrimaryVertexProducer::HLTScoutingPrimaryVertexProducer(const edm::ParameterSet& iConfig)
0034 : vertexCollToken_(consumes(iConfig.getParameter<edm::InputTag>("vertexCollection"))),
0035 mantissaPrecision_(iConfig.getParameter<int>("mantissaPrecision")) {
0036 produces<Run3ScoutingVertexCollection>("primaryVtx");
0037 }
0038
0039 void HLTScoutingPrimaryVertexProducer::produce(edm::StreamID sid,
0040 edm::Event& iEvent,
0041 edm::EventSetup const& setup) const {
0042 auto outVertices = std::make_unique<Run3ScoutingVertexCollection>();
0043
0044 auto const vertexCollHandle = iEvent.getHandle(vertexCollToken_);
0045
0046 if (vertexCollHandle.isValid()) {
0047 outVertices->reserve(vertexCollHandle->size());
0048 for (auto const& vtx : *vertexCollHandle) {
0049 outVertices->emplace_back(
0050 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.x(), mantissaPrecision_),
0051 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.y(), mantissaPrecision_),
0052 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.z(), mantissaPrecision_),
0053 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.zError(), mantissaPrecision_),
0054 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.xError(), mantissaPrecision_),
0055 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.yError(), mantissaPrecision_),
0056 vtx.tracksSize(),
0057 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.chi2(), mantissaPrecision_),
0058 vtx.ndof(),
0059 vtx.isValid(),
0060 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.covariance(0, 1), mantissaPrecision_),
0061 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.covariance(0, 2), mantissaPrecision_),
0062 MiniFloatConverter::reduceMantissaToNbitsRounding(vtx.covariance(1, 2), mantissaPrecision_));
0063 }
0064 }
0065
0066 iEvent.put(std::move(outVertices), "primaryVtx");
0067 }
0068
0069 void HLTScoutingPrimaryVertexProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0070 edm::ParameterSetDescription desc;
0071 desc.add<edm::InputTag>("vertexCollection", edm::InputTag("hltPixelVertices"))
0072 ->setComment("InputTag of input collection of primary vertices");
0073 desc.add<int>("mantissaPrecision", 10)->setComment("default float16, change to 23 for float32");
0074 descriptions.addWithDefaultLabel(desc);
0075 }
0076
0077 DEFINE_FWK_MODULE(HLTScoutingPrimaryVertexProducer);