Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-01-17 23:31:51

0001 #ifndef RECOLOCALTRACKER_SISTRIPZEROSUPPRESSION_SISTRIPAPVRESTORER_H
0002 #define RECOLOCALTRACKER_SISTRIPZEROSUPPRESSION_SISTRIPAPVRESTORER_H
0003 
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/Framework/interface/EventSetup.h"
0006 #include "DataFormats/Common/interface/Handle.h"
0007 #include "FWCore/Framework/interface/EventSetup.h"
0008 #include "FWCore/Framework/interface/ESWatcher.h"
0009 #include "CondFormats/SiStripObjects/interface/SiStripNoises.h"
0010 #include "CondFormats/SiStripObjects/interface/SiStripPedestals.h"
0011 #include "CalibFormats/SiStripObjects/interface/SiStripQuality.h"
0012 #include "CondFormats/DataRecord/interface/SiStripNoisesRcd.h"
0013 #include "CondFormats/DataRecord/interface/SiStripPedestalsRcd.h"
0014 #include "CalibTracker/Records/interface/SiStripQualityRcd.h"
0015 #include "DataFormats/Common/interface/DetSet.h"
0016 #include "DataFormats/Common/interface/DetSetVector.h"
0017 #include "DataFormats/DetId/interface/DetId.h"
0018 #include "DataFormats/SiStripDigi/interface/SiStripRawDigi.h"
0019 #include "DataFormats/SiStripDigi/interface/SiStripProcessedRawDigi.h"
0020 #include "RecoLocalTracker/SiStripZeroSuppression/interface/SiStripCommonModeNoiseSubtractor.h"
0021 
0022 #include <vector>
0023 #include <cstdint>
0024 
0025 class SiStripAPVRestorer {
0026   friend class SiStripRawProcessingFactory;
0027 
0028 protected:
0029   SiStripAPVRestorer(const edm::ParameterSet& conf, edm::ConsumesCollector);
0030 
0031 public:
0032   virtual ~SiStripAPVRestorer(){};
0033 
0034   using digi_t = int16_t;
0035   using digivector_t = std::vector<digi_t>;
0036   using digimap_t = std::map<uint16_t, digi_t>;
0037   using medians_t = std::vector<std::pair<short, float>>;
0038   using baselinemap_t = std::map<uint16_t, digivector_t>;
0039 
0040   void init(const edm::EventSetup& es);
0041 
0042   uint16_t inspect(uint32_t detId, uint16_t firstAPV, const digivector_t& digis, const medians_t& vmedians);
0043   void restore(uint16_t firstAPV, digivector_t& digis);
0044 
0045   uint16_t inspectAndRestore(uint32_t detId,
0046                              uint16_t firstAPV,
0047                              const digivector_t& rawDigisPedSubtracted,
0048                              digivector_t& processedRawDigi,
0049                              const medians_t& vmedians);
0050 
0051   void loadMeanCMMap(const edm::Event&);
0052 
0053   const baselinemap_t& getBaselineMap() const { return baselineMap_; }
0054   const std::map<uint16_t, digimap_t>& getSmoothedPoints() const { return smoothedMaps_; }
0055   const std::vector<bool>& getAPVFlags() const { return apvFlagsBool_; }
0056 
0057 private:
0058   using CMMap = std::map<uint32_t, std::vector<float>>;  //detId, Vector of MeanCM per detId
0059   constexpr static uint16_t nTotStripsPerAPV = 128;
0060 
0061   uint16_t nullInspect(uint16_t firstAPV, const digivector_t& digis);
0062   uint16_t abnormalBaselineInspect(uint16_t firstAPV, const digivector_t& digis);
0063   uint16_t baselineFollowerInspect(uint16_t firstAPV, const digivector_t& digis);
0064   uint16_t baselineAndSaturationInspect(uint16_t firstAPV, const digivector_t& digis);
0065   uint16_t forceRestoreInspect(uint16_t firstAPV, const digivector_t& digis);
0066   uint16_t hybridFormatInspect(uint16_t firstAPV, const digivector_t& digis);
0067   uint16_t hybridEmulationInspect(uint16_t firstAPV, const digivector_t& digis);
0068 
0069   void flatRestore(uint16_t apvN, uint16_t firstAPV, digivector_t& digis);
0070   bool checkBaseline(const std::vector<int16_t>& baseline) const;
0071   void baselineFollowerRestore(uint16_t apvN, uint16_t firstAPV, float median, digivector_t& digis);
0072   void derivativeFollowerRestore(uint16_t apvN, uint16_t firstAPV, digivector_t& digis);
0073 
0074   void baselineFollower(const digimap_t&, digivector_t& baseline, float median);
0075   bool flatRegionsFinder(const digivector_t& adcs, digimap_t& smoothedpoints, uint16_t apvN);
0076 
0077   void baselineCleaner(const digivector_t& adcs, digimap_t& smoothedpoints, uint16_t apvN);
0078   void cleaner_MonotonyChecker(digimap_t& smoothedpoints);
0079   void cleaner_HighSlopeChecker(digimap_t& smoothedpoints);
0080   void cleaner_LocalMinimumAdder(const digivector_t& adcs, digimap_t& smoothedpoints, uint16_t apvN);
0081 
0082   void createCMMapRealPed(const edm::DetSetVector<SiStripRawDigi>& input);
0083   void createCMMapCMstored(const edm::DetSetVector<SiStripProcessedRawDigi>& input);
0084 
0085 private:  // members
0086   edm::EDGetTokenT<edm::DetSetVector<SiStripRawDigi>> siStripRawDigiToken_;
0087   edm::EDGetTokenT<edm::DetSetVector<SiStripProcessedRawDigi>> siStripProcessedRawDigiToken_;
0088 
0089   edm::ESGetToken<SiStripQuality, SiStripQualityRcd> qualityToken_;
0090   edm::ESGetToken<SiStripNoises, SiStripNoisesRcd> noiseToken_;
0091   edm::ESGetToken<SiStripPedestals, SiStripPedestalsRcd> pedestalToken_;
0092   const SiStripQuality* qualityHandle;
0093   const SiStripNoises* noiseHandle;
0094   const SiStripPedestals* pedestalHandle;
0095   edm::ESWatcher<SiStripQualityRcd> qualityWatcher_;
0096   edm::ESWatcher<SiStripNoisesRcd> noiseWatcher_;
0097   edm::ESWatcher<SiStripPedestalsRcd> pedestalWatcher_;
0098 
0099   // event state
0100   CMMap meanCMmap_;
0101   // state
0102   uint32_t detId_;
0103   std::vector<std::string> apvFlags_;
0104   std::vector<bool> apvFlagsBool_;
0105   std::vector<bool> apvFlagsBoolOverride_;
0106   std::vector<float> median_;
0107   std::vector<bool> badAPVs_;
0108   std::map<uint16_t, digimap_t> smoothedMaps_;
0109   baselinemap_t baselineMap_;
0110 
0111   //--------------------------------------------------
0112   // Configurable Parameters of Algorithm
0113   //--------------------------------------------------
0114   bool forceNoRestore_;
0115   std::string inspectAlgo_;
0116   std::string restoreAlgo_;
0117   bool useRealMeanCM_;
0118   int32_t meanCM_;
0119   uint32_t deltaCMThreshold_;     // for BaselineFollower inspect
0120   double fraction_;               // fraction of strips deviating from nominal baseline
0121   uint32_t deviation_;            // ADC value of deviation from nominal baseline
0122   double restoreThreshold_;       // for Null inspect (fraction of adc=0 channels)
0123   uint32_t nSaturatedStrip_;      // for BaselineAndSaturation inspect
0124   uint32_t nSigmaNoiseDerTh_;     // threshold for rejecting hits strips
0125   uint32_t consecThreshold_;      // minimum length of flat region
0126   uint32_t nSmooth_;              // for smoothing and local minimum determination (odd number)
0127   uint32_t distortionThreshold_;  // (max-min) of flat regions to trigger baseline follower
0128   bool applyBaselineCleaner_;
0129   uint32_t cleaningSequence_;
0130   int32_t slopeX_;
0131   int32_t slopeY_;
0132   uint32_t hitStripThreshold_;  // height above median when strip is definitely a hit
0133   uint32_t minStripsToFit_;     // minimum strips to try spline algo (otherwise default to median)
0134   bool applyBaselineRejection_;
0135   double filteredBaselineMax_;
0136   double filteredBaselineDerivativeSumSquare_;
0137   int gradient_threshold_;
0138   int last_gradient_;
0139   int size_window_;
0140   int width_cluster_;
0141 };
0142 #endif