SiStripClusterizer/plugins/ClustersFromRawProducer.cc

0001 /*
0002  */
0003 #include "RecoLocalTracker/SiStripClusterizer/interface/StripClusterizerAlgorithmFactory.h"
0004 #include "RecoLocalTracker/SiStripZeroSuppression/interface/SiStripRawProcessingFactory.h"
0005
0006 #include "RecoLocalTracker/SiStripClusterizer/interface/StripClusterizerAlgorithm.h"
0007 #include "RecoLocalTracker/SiStripZeroSuppression/interface/SiStripRawProcessingAlgorithms.h"
0008
0009 #include "DataFormats/SiStripCluster/interface/SiStripCluster.h"
0010 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0011
0012 #include "DataFormats/FEDRawData/interface/FEDRawDataCollection.h"
0013 #include "EventFilter/SiStripRawToDigi/interface/SiStripFEDBuffer.h"
0014 #include "DataFormats/SiStripCommon/interface/SiStripConstants.h"
0015
0016 #include "CalibFormats/SiStripObjects/interface/SiStripDetCabling.h"
0017
0018 #include "FWCore/Framework/interface/stream/EDProducer.h"
0019 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0020 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0021 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0022 #include "FWCore/Utilities/interface/InputTag.h"
0023 #include "FWCore/Framework/interface/Event.h"
0024 #include "FWCore/Framework/interface/EventSetup.h"
0025 #include "FWCore/Framework/interface/ESHandle.h"
0026 #include "FWCore/Utilities/interface/Likely.h"
0027
0028 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0029 #include <sstream>
0030 #include <memory>
0031 #include <atomic>
0032 #include <mutex>
0033
0034 // #define VIDEBUG
0035 #ifdef VIDEBUG
0036 #include <iostream>
0037 #define COUT std::cout << "VI "
0038 #else
0039 #define COUT LogDebug("")
0040 #endif
0041
0042 namespace {
0043   std::unique_ptr<sistrip::FEDBuffer> fillBuffer(int fedId, const FEDRawDataCollection& rawColl) {
0044     std::unique_ptr<sistrip::FEDBuffer> buffer;
0045
0046     // Retrieve FED raw data for given FED
0047     const FEDRawData& rawData = rawColl.FEDData(fedId);
0048
0049     // Check on FEDRawData pointer
0050     const auto st_buffer = sistrip::preconstructCheckFEDBuffer(rawData);
0051     if UNLIKELY (sistrip::FEDBufferStatusCode::SUCCESS != st_buffer) {
0052       if (edm::isDebugEnabled()) {
0053         edm::LogWarning(sistrip::mlRawToCluster_)
0054             << "[ClustersFromRawProducer::" << __func__ << "]" << st_buffer << " for FED ID " << fedId;
0055       }
0056       return buffer;
0057     }
0058     buffer = std::make_unique<sistrip::FEDBuffer>(rawData);
0059     const auto st_chan = buffer->findChannels();
0060     if UNLIKELY (sistrip::FEDBufferStatusCode::SUCCESS != st_chan) {
0061       if (edm::isDebugEnabled()) {
0062         edm::LogWarning(sistrip::mlRawToCluster_)
0063             << "Exception caught when creating FEDBuffer object for FED " << fedId << ": " << st_chan;
0064       }
0065       buffer.reset();
0066       return buffer;
0067     }
0068     if UNLIKELY (!buffer->doChecks(false)) {
0069       if (edm::isDebugEnabled()) {
0070         edm::LogWarning(sistrip::mlRawToCluster_)
0071             << "Exception caught when creating FEDBuffer object for FED " << fedId << ": FED Buffer check fails";
0072       }
0073       buffer.reset();
0074       return buffer;
0075     }
0076
0077     /*
0078     // dump of FEDRawData to stdout
0079     if ( dump_ ) {
0080     std::stringstream ss;
0081     RawToDigiUnpacker::dumpRawData( fedId, rawData, ss );
0082     LogTrace(mlRawToDigi_)
0083     << ss.str();
0084     }
0085     */
0086
0087     return buffer;
0088   }
0089
0090   class ClusterFiller final : public StripClusterizerAlgorithm::output_t::Getter {
0091   public:
0092     ClusterFiller(const FEDRawDataCollection& irawColl,
0093                   StripClusterizerAlgorithm& iclusterizer,
0094                   SiStripRawProcessingAlgorithms& irawAlgos,
0095                   bool idoAPVEmulatorCheck,
0096                   bool legacy,
0097                   bool hybridZeroSuppressed)
0098         : rawColl(irawColl),
0099           clusterizer(iclusterizer),
0100           conditions(iclusterizer.conditions()),
0101           rawAlgos(irawAlgos),
0102           doAPVEmulatorCheck(idoAPVEmulatorCheck),
0103           legacy_(legacy),
0104           hybridZeroSuppressed_(hybridZeroSuppressed) {
0105       incTot(clusterizer.conditions().allDetIds().size());
0106       for (auto& d : done)
0107         d = nullptr;
0108     }
0109
0110     ~ClusterFiller() override { printStat(); }
0111
0112     void fill(StripClusterizerAlgorithm::output_t::TSFastFiller& record) override;
0113
0114   private:
0115     std::unique_ptr<sistrip::FEDBuffer> buffers[1024];
0116     std::atomic<sistrip::FEDBuffer*> done[1024];
0117
0118     const FEDRawDataCollection& rawColl;
0119
0120     StripClusterizerAlgorithm& clusterizer;
0121     const SiStripClusterizerConditions& conditions;
0122     SiStripRawProcessingAlgorithms& rawAlgos;
0123
0124     // March 2012: add flag for disabling APVe check in configuration
0125     bool doAPVEmulatorCheck;
0126
0127     bool legacy_;
0128     bool hybridZeroSuppressed_;
0129
0130 #ifdef VIDEBUG
0131     struct Stat {
0132       Stat() : totDet(0), detReady(0), detSet(0), detAct(0), detNoZ(0), detAbrt(0), totClus(0) {}
0133       std::atomic<int> totDet;    // all dets
0134       std::atomic<int> detReady;  // dets "updated"
0135       std::atomic<int> detSet;    // det actually set not empty
0136       std::atomic<int> detAct;    // det actually set with content
0137       std::atomic<int> detNoZ;    // det actually set with content
0138       std::atomic<int> detAbrt;   // det aborted
0139       std::atomic<int> totClus;   // total number of clusters
0140     };
0141
0142     mutable Stat stat;
0143     // void zeroStat() const { stat = std::move(Stat()); }
0144     void incTot(int n) const { stat.totDet = n; }
0145     void incReady() const { stat.detReady++; }
0146     void incSet() const { stat.detSet++; }
0147     void incAct() const { stat.detAct++; }
0148     void incNoZ() const { stat.detNoZ++; }
0149     void incAbrt() const { stat.detAbrt++; }
0150     void incClus(int n) const { stat.totClus += n; }
0151     void printStat() const {
0152       COUT << "VI clusters " << stat.totDet << ',' << stat.detReady << ',' << stat.detSet << ',' << stat.detAct << ','
0153            << stat.detNoZ << ',' << stat.detAbrt << ',' << stat.totClus << std::endl;
0154     }
0155
0156 #else
0157     static void zeroStat() {}
0158     static void incTot(int) {}
0159     static void incReady() {}
0160     static void incSet() {}
0161     static void incAct() {}
0162     static void incNoZ() {}
0163     static void incAbrt() {}
0164     static void incClus(int) {}
0165     static void printStat() {}
0166 #endif
0167   };
0168 }  // namespace
0169
0170 class SiStripClusterizerFromRaw final : public edm::stream::EDProducer<> {
0171 public:
0172   explicit SiStripClusterizerFromRaw(const edm::ParameterSet& conf)
0173       : onDemand(conf.getParameter<bool>("onDemand")),
0174         clusterizer_(StripClusterizerAlgorithmFactory::create(consumesCollector(),
0175                                                               conf.getParameter<edm::ParameterSet>("Clusterizer"))),
0176         rawAlgos_(SiStripRawProcessingFactory::create(conf.getParameter<edm::ParameterSet>("Algorithms"),
0177                                                       consumesCollector())),
0178         doAPVEmulatorCheck_(conf.getParameter<bool>("DoAPVEmulatorCheck")),
0179         legacy_(conf.getParameter<bool>("LegacyUnpacker")),
0180         hybridZeroSuppressed_(conf.getParameter<bool>("HybridZeroSuppressed")) {
0181     productToken_ = consumes<FEDRawDataCollection>(conf.getParameter<edm::InputTag>("ProductLabel"));
0182     produces<edmNew::DetSetVector<SiStripCluster> >();
0183     assert(clusterizer_.get());
0184     assert(rawAlgos_.get());
0185   }
0186
0187   void produce(edm::Event& ev, const edm::EventSetup& es) override {
0188     initialize(es);
0189
0190     // get raw data
0191     edm::Handle<FEDRawDataCollection> rawData;
0192     ev.getByToken(productToken_, rawData);
0193
0194     std::unique_ptr<edmNew::DetSetVector<SiStripCluster> > output(
0195         onDemand ? new edmNew::DetSetVector<SiStripCluster>(
0196                        std::shared_ptr<edmNew::DetSetVector<SiStripCluster>::Getter>(std::make_shared<ClusterFiller>(
0197                            *rawData, *clusterizer_, *rawAlgos_, doAPVEmulatorCheck_, legacy_, hybridZeroSuppressed_)),
0198                        clusterizer_->conditions().allDetIds())
0199                  : new edmNew::DetSetVector<SiStripCluster>());
0200
0201     if (onDemand)
0202       assert(output->onDemand());
0203
0204     output->reserve(15000, 24 * 10000);
0205
0206     if (!onDemand) {
0207       run(*rawData, *output);
0208       output->shrink_to_fit();
0209       COUT << output->dataSize() << " clusters from " << output->size() << " modules" << std::endl;
0210     }
0211
0212     ev.put(std::move(output));
0213   }
0214
0215   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0216
0217 private:
0218   void initialize(const edm::EventSetup& es);
0219
0220   void run(const FEDRawDataCollection& rawColl, edmNew::DetSetVector<SiStripCluster>& output);
0221
0222 private:
0223   bool onDemand;
0224
0225   edm::EDGetTokenT<FEDRawDataCollection> productToken_;
0226
0227   std::unique_ptr<StripClusterizerAlgorithm> clusterizer_;
0228   std::unique_ptr<SiStripRawProcessingAlgorithms> rawAlgos_;
0229
0230   // March 2012: add flag for disabling APVe check in configuration
0231   bool doAPVEmulatorCheck_;
0232
0233   bool legacy_;
0234   bool hybridZeroSuppressed_;
0235 };
0236
0237 void SiStripClusterizerFromRaw::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0238   edm::ParameterSetDescription desc;
0239
0240   desc.add("ProductLabel", edm::InputTag("rawDataCollector"));
0241   desc.add<std::string>("ConditionsLabel", "");
0242   desc.add("onDemand", true);
0243   desc.add("DoAPVEmulatorCheck", true);
0244   desc.add("LegacyUnpacker", false);
0245   desc.add("HybridZeroSuppressed", false);
0246
0247   edm::ParameterSetDescription clusterizer;
0248   StripClusterizerAlgorithmFactory::fillDescriptions(clusterizer);
0249   desc.add("Clusterizer", clusterizer);
0250
0251   edm::ParameterSetDescription algorithms;
0252   SiStripRawProcessingFactory::fillDescriptions(algorithms);
0253   desc.add("Algorithms", algorithms);
0254
0255   descriptions.addWithDefaultLabel(desc);
0256 }
0257
0258 #include "FWCore/Framework/interface/MakerMacros.h"
0259 DEFINE_FWK_MODULE(SiStripClusterizerFromRaw);
0260
0261 void SiStripClusterizerFromRaw::initialize(const edm::EventSetup& es) {
0262   (*clusterizer_).initialize(es);
0263   (*rawAlgos_).initialize(es);
0264 }
0265
0266 void SiStripClusterizerFromRaw::run(const FEDRawDataCollection& rawColl, edmNew::DetSetVector<SiStripCluster>& output) {
0267   ClusterFiller filler(rawColl, *clusterizer_, *rawAlgos_, doAPVEmulatorCheck_, legacy_, hybridZeroSuppressed_);
0268
0269   // loop over good det in cabling
0270   for (auto idet : clusterizer_->conditions().allDetIds()) {
0271     StripClusterizerAlgorithm::output_t::TSFastFiller record(output, idet);
0272
0273     filler.fill(record);
0274
0275     if (record.empty())
0276       record.abort();
0277   }  // end loop over dets
0278 }
0279
0280 namespace {
0281   class StripByStripAdder {
0282   public:
0283     typedef std::output_iterator_tag iterator_category;
0284     typedef void value_type;
0285     typedef void difference_type;
0286     typedef void pointer;
0287     typedef void reference;
0288
0289     StripByStripAdder(StripClusterizerAlgorithm& clusterizer,
0290                       StripClusterizerAlgorithm::State& state,
0291                       StripClusterizerAlgorithm::output_t::TSFastFiller& record)
0292         : clusterizer_(clusterizer), state_(state), record_(record) {}
0293
0294     StripByStripAdder& operator=(SiStripDigi digi) {
0295       clusterizer_.stripByStripAdd(state_, digi.strip(), digi.adc(), record_);
0296       return *this;
0297     }
0298
0299     StripByStripAdder& operator*() { return *this; }
0300     StripByStripAdder& operator++() { return *this; }
0301     StripByStripAdder& operator++(int) { return *this; }
0302
0303   private:
0304     StripClusterizerAlgorithm& clusterizer_;
0305     StripClusterizerAlgorithm::State& state_;
0306     StripClusterizerAlgorithm::output_t::TSFastFiller& record_;
0307   };
0308
0309   template <typename Container>
0310   class ADC_back_inserter {
0311   public:
0312     ADC_back_inserter(Container& c) : c_(c) {}
0313
0314     ADC_back_inserter& operator=(SiStripRawDigi digi) {
0315       c_.push_back(digi.adc());
0316       return *this;
0317     }
0318     ADC_back_inserter& operator*() { return *this; }
0319     ADC_back_inserter& operator++() { return *this; }
0320     ADC_back_inserter& operator++(int) { return *this; }
0321
0322   private:
0323     Container& c_;
0324   };
0325 }  // namespace
0326
0327 void ClusterFiller::fill(StripClusterizerAlgorithm::output_t::TSFastFiller& record) {
0328   try {  // edmNew::CapacityExaustedException
0329     incReady();
0330
0331     auto idet = record.id();
0332
0333     COUT << "filling " << idet << std::endl;
0334
0335     auto const& det = clusterizer.stripByStripBegin(idet);
0336     if (!det.valid())
0337       return;
0338     StripClusterizerAlgorithm::State state(det);
0339
0340     incSet();
0341     record.reserve(16);
0342     // Loop over apv-pairs of det
0343     for (auto const conn : conditions.currentConnection(det)) {
0344       if UNLIKELY (!conn)
0345         continue;
0346
0347       const uint16_t fedId = conn->fedId();
0348
0349       // If fed id is null or connection is invalid continue
0350       if UNLIKELY (!fedId || !conn->isConnected()) {
0351         continue;
0352       }
0353
0354       // If Fed hasnt already been initialised, extract data and initialise
0355       sistrip::FEDBuffer* buffer = done[fedId];
0356       if (!buffer) {
0357         buffer = fillBuffer(fedId, rawColl).release();
0358         if (!buffer) {
0359           continue;
0360         }
0361         sistrip::FEDBuffer* exp = nullptr;
0362         if (done[fedId].compare_exchange_strong(exp, buffer))
0363           buffers[fedId].reset(buffer);
0364         else {
0365           delete buffer;
0366           buffer = done[fedId];
0367         }
0368       }
0369       assert(buffer);
0370
0371       buffer->setLegacyMode(legacy_);
0372
0373       // check channel
0374       const uint8_t fedCh = conn->fedCh();
0375
0376       if UNLIKELY (!buffer->channelGood(fedCh, doAPVEmulatorCheck)) {
0377         if (edm::isDebugEnabled()) {
0378           std::ostringstream ss;
0379           ss << "Problem unpacking channel " << fedCh << " on FED " << fedId;
0380           edm::LogWarning(sistrip::mlRawToCluster_) << ss.str();
0381         }
0382         continue;
0383       }
0384
0385       // Determine APV std::pair number
0386       uint16_t ipair = conn->apvPairNumber();
0387
0388       const sistrip::FEDReadoutMode mode = buffer->readoutMode();
0389       const sistrip::FEDLegacyReadoutMode lmode =
0390           legacy_ ? buffer->legacyReadoutMode() : sistrip::READOUT_MODE_LEGACY_INVALID;
0391
0392       using namespace sistrip;
0393       if LIKELY (fedchannelunpacker::isZeroSuppressed(mode, legacy_, lmode)) {
0394         auto perStripAdder = StripByStripAdder(clusterizer, state, record);
0395         const auto isNonLite = fedchannelunpacker::isNonLiteZS(mode, legacy_, lmode);
0396         const uint8_t pCode = (isNonLite ? buffer->packetCode(legacy_, fedCh) : 0);
0397         auto st_ch = fedchannelunpacker::StatusCode::SUCCESS;
0398         if LIKELY (!hybridZeroSuppressed_) {
0399           st_ch = fedchannelunpacker::unpackZeroSuppressed(
0400               buffer->channel(fedCh), perStripAdder, ipair * 256, isNonLite, mode, legacy_, lmode, pCode);
0401         } else {
0402           const uint32_t id = conn->detId();
0403           edm::DetSet<SiStripDigi> unpDigis{id};
0404           unpDigis.reserve(256);
0405           st_ch = fedchannelunpacker::unpackZeroSuppressed(
0406               buffer->channel(fedCh), std::back_inserter(unpDigis), ipair * 256, isNonLite, mode, legacy_, lmode, pCode);
0407           if (fedchannelunpacker::StatusCode::SUCCESS == st_ch) {
0408             edm::DetSet<SiStripDigi> suppDigis{id};
0409             rawAlgos.suppressHybridData(unpDigis, suppDigis, ipair * 2);
0410             std::copy(std::begin(suppDigis), std::end(suppDigis), perStripAdder);
0411           }
0412         }
0413         if (fedchannelunpacker::StatusCode::SUCCESS != st_ch && edm::isDebugEnabled()) {
0414           edm::LogWarning(sistrip::mlRawToCluster_)
0415               << "Unordered clusters for channel " << fedCh << " on FED " << fedId << ": " << toString(st_ch);
0416           continue;
0417         }
0418       } else {
0419         auto st_ch = fedchannelunpacker::StatusCode::SUCCESS;
0420         if (fedchannelunpacker::isVirginRaw(mode, legacy_, lmode)) {
0421           std::vector<int16_t> digis;
0422           st_ch = fedchannelunpacker::unpackVirginRaw(
0423               buffer->channel(fedCh), ADC_back_inserter(digis), buffer->channel(fedCh).packetCode());
0424           if (fedchannelunpacker::StatusCode::SUCCESS == st_ch) {
0425             //process raw
0426             uint32_t id = conn->detId();
0427             edm::DetSet<SiStripDigi> zsdigis(id);
0428             //rawAlgos_->subtractorPed->subtract( id, ipair*256, digis);
0429             //rawAlgos_->subtractorCMN->subtract( id, digis);
0430             //rawAlgos_->suppressor->suppress( digis, zsdigis);
0431             uint16_t firstAPV = ipair * 2;
0432             rawAlgos.suppressVirginRawData(id, firstAPV, digis, zsdigis);
0433             for (const auto digi : zsdigis) {
0434               clusterizer.stripByStripAdd(state, digi.strip(), digi.adc(), record);
0435             }
0436           }
0437         } else if (fedchannelunpacker::isProcessedRaw(mode, legacy_, lmode)) {
0438           std::vector<int16_t> digis;
0439           st_ch = fedchannelunpacker::unpackProcessedRaw(buffer->channel(fedCh), ADC_back_inserter(digis));
0440           if (fedchannelunpacker::StatusCode::SUCCESS == st_ch) {
0441             //process raw
0442             uint32_t id = conn->detId();
0443             edm::DetSet<SiStripDigi> zsdigis(id);
0444             //rawAlgos_->subtractorCMN->subtract( id, digis);
0445             //rawAlgos_->suppressor->suppress( digis, zsdigis);
0446             uint16_t firstAPV = ipair * 2;
0447             rawAlgos.suppressProcessedRawData(id, firstAPV, digis, zsdigis);
0448             for (edm::DetSet<SiStripDigi>::const_iterator it = zsdigis.begin(); it != zsdigis.end(); it++) {
0449               clusterizer.stripByStripAdd(state, it->strip(), it->adc(), record);
0450             }
0451           }
0452         } else {
0453           edm::LogWarning(sistrip::mlRawToCluster_)
0454               << "[ClustersFromRawProducer::" << __func__ << "]"
0455               << " FEDRawData readout mode " << mode << " from FED id " << fedId << " not supported.";
0456         }
0457         if (fedchannelunpacker::StatusCode::SUCCESS != st_ch && edm::isDebugEnabled()) {
0458           edm::LogWarning(sistrip::mlRawToCluster_)
0459               << "[ClustersFromRawProducer::" << __func__ << "]" << toString(st_ch) << " from FED id " << fedId
0460               << " channel " << fedCh;
0461         }
0462       }
0463     }  // end loop over conn
0464
0465     clusterizer.stripByStripEnd(state, record);
0466
0467     incAct();
0468
0469     if (record.full()) {
0470       edm::LogError(sistrip::mlRawToCluster_) << "too many Sistrip Clusters to fit space allocated for OnDemand for "
0471                                               << record.id() << ' ' << record.size();
0472       record.abort();
0473       incAbrt();
0474     }
0475
0476     if (!record.empty())
0477       incNoZ();
0478
0479     COUT << "filled " << record.size() << std::endl;
0480     for (auto const& cl : record)
0481       COUT << cl.firstStrip() << ',' << cl.amplitudes().size() << std::endl;
0482     incClus(record.size());
0483   } catch (edmNew::CapacityExaustedException const&) {
0484     edm::LogError(sistrip::mlRawToCluster_) << "too many Sistrip Clusters to fit space allocated for OnDemand";
0485   }
0486 }