Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-08-05 23:13:47

0001 
0002 
0003 #include "FWCore/Framework/interface/MakerMacros.h"
0004 #include "FWCore/Framework/interface/Frameworkfwd.h"
0005 #include "FWCore/Framework/interface/stream/EDProducer.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "FWCore/Utilities/interface/InputTag.h"
0008 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0009 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0010 #include "DataFormats/SiStripCluster/interface/SiStripApproximateCluster.h"
0011 #include "DataFormats/SiStripCluster/interface/SiStripCluster.h"
0012 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0013 #include "DataFormats/Common/interface/DetSetVector.h"
0014 
0015 #include <vector>
0016 #include <memory>
0017 
0018 class SiStripClusters2ApproxClusters : public edm::stream::EDProducer<> {
0019 public:
0020   explicit SiStripClusters2ApproxClusters(const edm::ParameterSet& conf);
0021   void produce(edm::Event&, const edm::EventSetup&) override;
0022 
0023   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0024 
0025 private:
0026   edm::InputTag inputClusters;
0027   edm::EDGetTokenT<edmNew::DetSetVector<SiStripCluster> > clusterToken;
0028 
0029   unsigned int maxNSat;
0030 };
0031 
0032 SiStripClusters2ApproxClusters::SiStripClusters2ApproxClusters(const edm::ParameterSet& conf) {
0033   inputClusters = conf.getParameter<edm::InputTag>("inputClusters");
0034   maxNSat = conf.getParameter<unsigned int>("maxSaturatedStrips");
0035 
0036   clusterToken = consumes<edmNew::DetSetVector<SiStripCluster> >(inputClusters);
0037   produces<edmNew::DetSetVector<SiStripApproximateCluster> >();
0038 }
0039 
0040 void SiStripClusters2ApproxClusters::produce(edm::Event& event, edm::EventSetup const&) {
0041   auto result = std::make_unique<edmNew::DetSetVector<SiStripApproximateCluster> >();
0042   const auto& clusterCollection = event.get(clusterToken);
0043 
0044   for (const auto& detClusters : clusterCollection) {
0045     edmNew::DetSetVector<SiStripApproximateCluster>::FastFiller ff{*result, detClusters.id()};
0046 
0047     for (const auto& cluster : detClusters)
0048       ff.push_back(SiStripApproximateCluster(cluster, maxNSat));
0049   }
0050 
0051   event.put(std::move(result));
0052 }
0053 
0054 void SiStripClusters2ApproxClusters::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0055   edm::ParameterSetDescription desc;
0056   desc.add<edm::InputTag>("inputClusters", edm::InputTag("siStripClusters"));
0057   desc.add<unsigned int>("maxSaturatedStrips", 3);
0058   descriptions.add("SiStripClusters2ApproxClusters", desc);
0059 }
0060 
0061 DEFINE_FWK_MODULE(SiStripClusters2ApproxClusters);