File indexing completed on 2024-04-06 11:59:46
0001 #include "ShallowDigisProducer.h"
0002
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/Framework/interface/EventSetup.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "DataFormats/Common/interface/DetSetVector.h"
0008 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0009 #include "DataFormats/SiStripDigi/interface/SiStripDigi.h"
0010
0011 ShallowDigisProducer::ShallowDigisProducer(const edm::ParameterSet& conf)
0012 : inputTags(conf.getParameter<std::vector<edm::InputTag> >("DigiProducersList")), noisesToken_(esConsumes()) {
0013 produces<std::vector<unsigned> >("id");
0014 produces<std::vector<unsigned> >("subdet");
0015 produces<std::vector<unsigned> >("strip");
0016 produces<std::vector<unsigned> >("adc");
0017 produces<std::vector<float> >("noise");
0018 }
0019
0020 void ShallowDigisProducer::insert(products& p, edm::Event& e) {
0021 e.put(std::move(p.id), "id");
0022 e.put(std::move(p.subdet), "subdet");
0023 e.put(std::move(p.strip), "strip");
0024 e.put(std::move(p.adc), "adc");
0025 e.put(std::move(p.noise), "noise");
0026 }
0027
0028 template <class T>
0029 inline void ShallowDigisProducer::recordDigis(const T& digiCollection, products& p, const SiStripNoises& noises) {
0030 for (auto const& set : digiCollection) {
0031 SiStripNoises::Range detNoiseRange = noises.getRange(set.detId());
0032 for (auto const& digi : set) {
0033 p.id->push_back(set.detId());
0034 p.subdet->push_back((set.detId() >> 25) & 0x7);
0035 p.strip->push_back(digi.strip());
0036 p.adc->push_back(digi.adc());
0037 p.noise->push_back(noises.getNoise(digi.strip(), detNoiseRange));
0038 }
0039 }
0040 }
0041
0042 void ShallowDigisProducer::produce(edm::Event& e, const edm::EventSetup& es) {
0043 products p;
0044 edm::Handle<edm::DetSetVector<SiStripDigi> > inputOld;
0045 edm::Handle<edmNew::DetSetVector<SiStripDigi> > inputNew;
0046 const auto& noises = es.getData(noisesToken_);
0047 if (findInput(inputOld, e))
0048 recordDigis(*inputOld, p, noises);
0049 else if (findInput(inputNew, e))
0050 recordDigis(*inputNew, p, noises);
0051 else
0052 edm::LogWarning("Input Not Found");
0053 insert(p, e);
0054 }
0055
0056 template <class T>
0057 inline bool ShallowDigisProducer::findInput(edm::Handle<T>& handle, const edm::Event& e) {
0058 for (auto const& inputTag : inputTags) {
0059 e.getByLabel(inputTag, handle);
0060 if (handle.isValid() && !handle->empty()) {
0061 LogDebug("Input") << inputTag;
0062 return true;
0063 }
0064 }
0065 return false;
0066 }