Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-03-17 11:19:34

0001 #ifndef RecoLocalTracker_SiPixelClusterizer_PixelThresholdClusterizer_H
0002 #define RecoLocalTracker_SiPixelClusterizer_PixelThresholdClusterizer_H
0003 
0004 //-----------------------------------------------------------------------
0005 //! \class PixelThresholdClusterizer
0006 //! \brief An explicit threshold-based clustering algorithm.
0007 //!
0008 //! A threshold-based clustering algorithm which clusters SiPixelDigis
0009 //! into SiPixelClusters for each DetUnit.  The algorithm is straightforward
0010 //! and purely topological: the clustering process starts with seed pixels
0011 //! and continues by adding adjacent pixels above the pixel threshold.
0012 //! Once the cluster is made, it has to be above the cluster threshold as
0013 //! well.
0014 //!
0015 //! The clusterization is performed on a matrix with size
0016 //! equal to the size of the pixel detector, each cell containing
0017 //! the ADC count of the corresponding pixel.
0018 //! The matrix is reset after each clusterization.
0019 //!
0020 //! The search starts from seed pixels, i.e. pixels with sufficiently
0021 //! large amplitudes, found at the time of filling of the matrix
0022 //! and stored in a
0023 //!
0024 //! At this point the noise and dead channels are ignored, but soon they
0025 //! won't be.
0026 //!
0027 //! SiPixelCluster contains a barrycenter, but it should be noted that that
0028 //! information is largely useless.  One must use a PositionEstimator
0029 //! class to compute the RecHit position and its error for every given
0030 //! cluster.
0031 //!
0032 //! \author Largely copied from NewPixelClusterizer in ORCA written by
0033 //!     Danek Kotlinski (PSI).   Ported to CMSSW by Petar Maksimovic (JHU).
0034 //!     DetSetVector data container implemented by V.Chiochia (Uni Zurich)
0035 //!
0036 //! Sets the PixelArrayBuffer dimensions and pixel thresholds.
0037 //! Makes clusters and stores them in theCache if the option
0038 //! useCache has been set.
0039 //-----------------------------------------------------------------------
0040 
0041 // Base class, defines SiPixelDigi and SiPixelCluster.  The latter includes
0042 // Pixel, PixelPos and Shift as inner classes.
0043 //
0044 #include "DataFormats/Common/interface/DetSetVector.h"
0045 #include "PixelClusterizerBase.h"
0046 
0047 // The private pixel buffer
0048 #include "SiPixelArrayBuffer.h"
0049 
0050 // Parameter Set:
0051 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0052 
0053 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0054 
0055 #include <vector>
0056 
0057 class dso_hidden PixelThresholdClusterizer : public PixelClusterizerBase {
0058 public:
0059   PixelThresholdClusterizer(edm::ParameterSet const& conf);
0060   ~PixelThresholdClusterizer() override;
0061 
0062   // Full I/O in DetSet
0063   void clusterizeDetUnit(const edm::DetSet<PixelDigi>& input,
0064                          const PixelGeomDetUnit* pixDet,
0065                          const TrackerTopology* tTopo,
0066                          const std::vector<short>& badChannels,
0067                          edmNew::DetSetVector<SiPixelCluster>::FastFiller& output) override {
0068     clusterizeDetUnitT(input, pixDet, tTopo, badChannels, output);
0069   }
0070   void clusterizeDetUnit(const edmNew::DetSet<SiPixelCluster>& input,
0071                          const PixelGeomDetUnit* pixDet,
0072                          const TrackerTopology* tTopo,
0073                          const std::vector<short>& badChannels,
0074                          edmNew::DetSetVector<SiPixelCluster>::FastFiller& output) override {
0075     clusterizeDetUnitT(input, pixDet, tTopo, badChannels, output);
0076   }
0077 
0078   static void fillPSetDescription(edm::ParameterSetDescription& desc);
0079 
0080 protected:
0081   template <typename T>
0082   void clusterizeDetUnitT(const T& input,
0083                           const PixelGeomDetUnit* pixDet,
0084                           const TrackerTopology* tTopo,
0085                           const std::vector<short>& badChannels,
0086                           edmNew::DetSetVector<SiPixelCluster>::FastFiller& output);
0087 
0088   //! Data storage
0089   SiPixelArrayBuffer theBuffer;                    // internal nrow * ncol matrix
0090   std::vector<SiPixelCluster::PixelPos> theSeeds;  // cached seed pixels
0091   std::vector<SiPixelCluster> theClusters;         // resulting clusters
0092 
0093   std::vector<bool> theFakePixels;  // fake pixels introduced to guide clustering
0094 
0095   std::vector<uint8_t> thePixelOccurrence;  // the number of times each pixel occurs (for tracking duplicate pixels)
0096 
0097   //! Clustering-related quantities:
0098   float thePixelThresholdInNoiseUnits;    // Pixel threshold in units of noise
0099   float theSeedThresholdInNoiseUnits;     // Pixel cluster seed in units of noise
0100   float theClusterThresholdInNoiseUnits;  // Cluster threshold in units of noise
0101 
0102   const int thePixelThreshold;       // Pixel threshold in electrons
0103   const int theSeedThreshold;        // Seed threshold in electrons
0104   const int theClusterThreshold;     // Cluster threshold in electrons
0105   const int theClusterThreshold_L1;  // Cluster threshold in electrons for Layer 1
0106   const int theConversionFactor;     // adc to electron conversion factor
0107   const int theConversionFactor_L1;  // adc to electron conversion factor for Layer 1
0108   const int theOffset;               // adc to electron conversion offset
0109   const int theOffset_L1;            // adc to electron conversion offset for Layer 1
0110 
0111   const double theElectronPerADCGain;  //  ADC to electrons conversion
0112 
0113   const bool doPhase2Calibration;  // The ADC --> electrons calibration is for phase-2 tracker
0114 
0115   const bool dropDuplicates;  // Enabling dropping duplicate pixels
0116 
0117   const int thePhase2ReadoutMode;      // Readout mode of the phase-2 IT digitizer
0118   const double thePhase2DigiBaseline;  // Threshold above which digis are measured in the phase-2 IT
0119   const int thePhase2KinkADC;          // ADC count at which the kink in the dual slop kicks in
0120 
0121   //! Geometry-related information
0122   int theNumOfRows;
0123   int theNumOfCols;
0124   uint32_t theDetid;
0125   int theLayer;
0126   const bool doMissCalibrate;  // Use calibration or not
0127   const bool doSplitClusters;
0128   //! Private helper methods:
0129   bool setup(const PixelGeomDetUnit* pixDet);
0130   void copy_to_buffer(DigiIterator begin, DigiIterator end);
0131   void copy_to_buffer(ClusterIterator begin, ClusterIterator end);
0132   void clear_buffer(DigiIterator begin, DigiIterator end);
0133   void clear_buffer(ClusterIterator begin, ClusterIterator end);
0134   SiPixelCluster make_cluster(const SiPixelCluster::PixelPos& pix,
0135                               edmNew::DetSetVector<SiPixelCluster>::FastFiller& output);
0136   // Calibrate the ADC charge to electrons
0137   int calibrate(int adc, int col, int row);
0138 };
0139 
0140 #endif