Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-07-01 02:23:05

0001 #ifndef DataFormats_SiPixel_Cluster_SiPixelCluster_h
0002 #define DataFormats_SiPixel_Cluster_SiPixelCluster_h
0003 
0004 //---------------------------------------------------------------------------
0005 //!  \class SiPixelCluster
0006 //!  \brief Pixel cluster -- collection of neighboring pixels above threshold
0007 //!
0008 //!  Class to contain and store all the topological information of pixel clusters:
0009 //!  charge, global size, size and the barycenter in x and y
0010 //!  local directions. It builds a vector of SiPixel (which is
0011 //!  an inner class) and a container of channels.
0012 //!
0013 //!  March 2007: Edge methods moved to RectangularPixelTopology class (V.Chiochia)
0014 //!  Feb 2008: Modify the Pixel class from float to shorts
0015 //!  May   2008: Offset based packing (D.Fehling / A. Rizzi)
0016 //!  Sep 2012: added Max back, removed detId (V.I.)
0017 //!  sizeX and sizeY now clipped at 127
0018 //!  July 2017 make it compatible with PhaseII, remove errs....
0019 //---------------------------------------------------------------------------
0020 
0021 #include <vector>
0022 #include <cstdint>
0023 #include <cassert>
0024 #include <limits>
0025 
0026 class PixelDigi;
0027 
0028 class SiPixelCluster {
0029 public:
0030   // FIXME make it POD
0031   class Pixel {
0032   public:
0033     constexpr Pixel() : x(0), y(0), adc(0) {}  // for root
0034     constexpr Pixel(int pix_x, int pix_y, int pix_adc) : x(pix_x), y(pix_y), adc(pix_adc) {}
0035     uint16_t x;
0036     uint16_t y;
0037     uint16_t adc;
0038   };
0039 
0040   //--- Integer shift in x and y directions.
0041   // FIXME make    it POD
0042   class Shift {
0043   public:
0044     constexpr Shift(int dx, int dy) : dx_(dx), dy_(dy) {}
0045     constexpr Shift() : dx_(0), dy_(0) {}
0046     constexpr int dx() const { return dx_; }
0047     constexpr int dy() const { return dy_; }
0048 
0049   private:
0050     int dx_;
0051     int dy_;
0052   };
0053 
0054   //--- Position of a SiPixel
0055   // FIXME make    it POD
0056   class PixelPos {
0057   public:
0058     constexpr PixelPos() : row_(0), col_(0) {}
0059     constexpr PixelPos(int row, int col) : row_(row), col_(col) {}
0060     constexpr int row() const { return row_; }
0061     constexpr int col() const { return col_; }
0062     constexpr PixelPos operator+(const Shift& shift) const { return PixelPos(row() + shift.dx(), col() + shift.dy()); }
0063 
0064   private:
0065     int row_;
0066     int col_;
0067   };
0068 
0069   typedef std::vector<PixelDigi>::const_iterator PixelDigiIter;
0070   typedef std::pair<PixelDigiIter, PixelDigiIter> PixelDigiRange;
0071 
0072   static constexpr unsigned int MAXSPAN = 255;
0073   static constexpr unsigned int MAXPOS = 2047;
0074 
0075   static constexpr uint16_t invalidClusterId = std::numeric_limits<uint16_t>::max();
0076 
0077   /** Construct from a range of digis that form a cluster and from 
0078    *  a DetID. The range is assumed to be non-empty.
0079    */
0080 
0081   SiPixelCluster() = default;
0082   ~SiPixelCluster() = default;
0083   SiPixelCluster(SiPixelCluster const&) = default;
0084   SiPixelCluster(SiPixelCluster&&) = default;
0085   SiPixelCluster& operator=(SiPixelCluster const&) = default;
0086   SiPixelCluster& operator=(SiPixelCluster&&) = default;
0087 
0088   SiPixelCluster(unsigned int isize,
0089                  uint16_t const* adcs,
0090                  uint16_t const* xpos,
0091                  uint16_t const* ypos,
0092                  uint16_t xmin,
0093                  uint16_t ymin,
0094                  uint16_t id = invalidClusterId)
0095       : thePixelOffset(2 * isize), thePixelADC(adcs, adcs + isize), theOriginalClusterId(id) {
0096     uint16_t maxCol = 0;
0097     uint16_t maxRow = 0;
0098     for (unsigned int i = 0; i < isize; ++i) {
0099       uint16_t xoffset = xpos[i] - xmin;
0100       uint16_t yoffset = ypos[i] - ymin;
0101       thePixelOffset[i * 2] = std::min(uint16_t(MAXSPAN), xoffset);
0102       thePixelOffset[i * 2 + 1] = std::min(uint16_t(MAXSPAN), yoffset);
0103       if (xoffset > maxRow)
0104         maxRow = xoffset;
0105       if (yoffset > maxCol)
0106         maxCol = yoffset;
0107     }
0108     packRow(xmin, maxRow);
0109     packCol(ymin, maxCol);
0110   }
0111 
0112   // obsolete (only for regression tests)
0113   SiPixelCluster(const PixelPos& pix, int adc);
0114   void add(const PixelPos& pix, int adc);
0115 
0116   // Analog linear average position (barycenter)
0117   float x() const {
0118     float qm = 0.0;
0119     int isize = thePixelADC.size();
0120     for (int i = 0; i < isize; ++i)
0121       qm += float(thePixelADC[i]) * (thePixelOffset[i * 2] + minPixelRow() + 0.5f);
0122     return qm / charge();
0123   }
0124 
0125   float y() const {
0126     float qm = 0.0;
0127     int isize = thePixelADC.size();
0128     for (int i = 0; i < isize; ++i)
0129       qm += float(thePixelADC[i]) * (thePixelOffset[i * 2 + 1] + minPixelCol() + 0.5f);
0130     return qm / charge();
0131   }
0132 
0133   // Return number of pixels.
0134   int size() const { return thePixelADC.size(); }
0135 
0136   // Return cluster dimension in the x direction.
0137   int sizeX() const { return rowSpan() + 1; }
0138 
0139   // Return cluster dimension in the y direction.
0140   int sizeY() const { return colSpan() + 1; }
0141 
0142   inline int charge() const {
0143     int qm = 0;
0144     int isize = thePixelADC.size();
0145     for (int i = 0; i < isize; ++i)
0146       qm += thePixelADC[i];
0147     return qm;
0148   }  // Return total cluster charge.
0149 
0150   inline int minPixelRow() const { return theMinPixelRow; }             // The min x index.
0151   inline int maxPixelRow() const { return minPixelRow() + rowSpan(); }  // The max x index.
0152   inline int minPixelCol() const { return theMinPixelCol; }             // The min y index.
0153   inline int maxPixelCol() const { return minPixelCol() + colSpan(); }  // The max y index.
0154 
0155   const std::vector<uint8_t>& pixelOffset() const { return thePixelOffset; }
0156   const std::vector<uint16_t>& pixelADC() const { return thePixelADC; }
0157 
0158   // obsolete, use single pixel access below
0159   const std::vector<Pixel> pixels() const {
0160     std::vector<Pixel> oldPixVector;
0161     int isize = thePixelADC.size();
0162     oldPixVector.reserve(isize);
0163     for (int i = 0; i < isize; ++i) {
0164       oldPixVector.push_back(pixel(i));
0165     }
0166     return oldPixVector;
0167   }
0168 
0169   // infinite faster than above...
0170   Pixel pixel(int i) const {
0171     return Pixel(minPixelRow() + thePixelOffset[i * 2], minPixelCol() + thePixelOffset[i * 2 + 1], thePixelADC[i]);
0172   }
0173 
0174 private:
0175   static int overflow_(uint16_t span) { return span == uint16_t(MAXSPAN); }
0176 
0177 public:
0178   int colSpan() const { return thePixelColSpan; }
0179 
0180   int rowSpan() const { return thePixelRowSpan; }
0181 
0182   bool overflowCol() const { return overflow_(thePixelColSpan); }
0183 
0184   bool overflowRow() const { return overflow_(thePixelRowSpan); }
0185 
0186   bool overflow() const { return overflowCol() || overflowRow(); }
0187 
0188   void packCol(uint16_t ymin, uint16_t yspan) {
0189     theMinPixelCol = ymin;
0190     thePixelColSpan = std::min(yspan, uint16_t(MAXSPAN));
0191   }
0192   void packRow(uint16_t xmin, uint16_t xspan) {
0193     theMinPixelRow = xmin;
0194     thePixelRowSpan = std::min(xspan, uint16_t(MAXSPAN));
0195   }
0196 
0197   // ggiurgiu@fnal.gov, 01/05/12
0198   // Getters and setters for the newly added data members (err_x and err_y). See below.
0199   void setSplitClusterErrorX(float errx) { err_x = errx; }
0200   void setSplitClusterErrorY(float erry) { err_y = erry; }
0201   float getSplitClusterErrorX() const { return err_x; }
0202   float getSplitClusterErrorY() const { return err_y; }
0203 
0204   // the original id (they get sorted)
0205   auto originalId() const { return theOriginalClusterId; }
0206   void setOriginalId(uint16_t id) { theOriginalClusterId = id; }
0207 
0208 private:
0209   std::vector<uint8_t> thePixelOffset;
0210   std::vector<uint16_t> thePixelADC;
0211 
0212   uint16_t theMinPixelRow = MAXPOS;  // Minimum pixel index in the x direction (low edge).
0213   uint16_t theMinPixelCol = MAXPOS;  // Minimum pixel index in the y direction (left edge).
0214   uint8_t thePixelRowSpan = 0;       // Span pixel index in the x direction (low edge).
0215   uint8_t thePixelColSpan = 0;       // Span pixel index in the y direction (left edge).
0216 
0217   uint16_t theOriginalClusterId = invalidClusterId;
0218 
0219   float err_x = -99999.9f;
0220   float err_y = -99999.9f;
0221 };
0222 
0223 // Comparison operators  (needed by DetSetVector)
0224 inline bool operator<(const SiPixelCluster& one, const SiPixelCluster& other) {
0225   if (one.minPixelRow() < other.minPixelRow()) {
0226     return true;
0227   } else if (one.minPixelRow() > other.minPixelRow()) {
0228     return false;
0229   } else if (one.minPixelCol() < other.minPixelCol()) {
0230     return true;
0231   } else {
0232     return false;
0233   }
0234 }
0235 
0236 #include "DataFormats/Common/interface/DetSetVector.h"
0237 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0238 #include "DataFormats/Common/interface/Ref.h"
0239 #include "DataFormats/Common/interface/DetSetRefVector.h"
0240 
0241 typedef edm::DetSetVector<SiPixelCluster> SiPixelClusterCollection;
0242 typedef edm::Ref<SiPixelClusterCollection, SiPixelCluster> SiPixelClusterRef;
0243 typedef edm::DetSetRefVector<SiPixelCluster> SiPixelClusterRefVector;
0244 typedef edm::RefProd<SiPixelClusterCollection> SiPixelClusterRefProd;
0245 
0246 typedef edmNew::DetSetVector<SiPixelCluster> SiPixelClusterCollectionNew;
0247 typedef edm::Ref<SiPixelClusterCollectionNew, SiPixelCluster> SiPixelClusterRefNew;
0248 #endif