CMS 3D CMS Logo

PixelBrickedDigitizerAlgorithm.cc
Go to the documentation of this file.
1 #include <iostream>
2 #include <cmath>
3 
6 
12 
17 
18 // Geometry
22 
23 using namespace edm;
24 using namespace sipixelobjects;
25 
27  : PixelDigitizerAlgorithm(conf, iC)
28 /* odd_row_interchannelCoupling_next_row_(conf.getParameter<ParameterSet>("PixelBrickedDigitizerAlgorithm")
29  .getParameter<double>("Odd_row_interchannelCoupling_next_row"))
30  even_row_interchannelCoupling_next_row_(conf.getParameter<ParameterSet>("PixelBrickedDigitizerAlgorithm")
31  .getParameter<double>("Even_row_interchannelCoupling_next_row")),
32  odd_column_interchannelCoupling_next_column_(
33  conf.getParameter<ParameterSet>("PixelBrickedDigitizerAlgorithm")
34  .getParameter<double>("Odd_column_interchannelCoupling_next_column")),
35  even_column_interchannelCoupling_next_column_(
36  conf.getParameter<ParameterSet>("PixelBrickedDigitizerAlgorithm")
37  .getParameter<double>("Even_column_interchannelCoupling_next_column"))*/
38 {
39  even_row_interchannelCoupling_next_row_ = conf.getParameter<ParameterSet>("PixelBrickedDigitizerAlgorithm")
40  .getParameter<double>("Even_row_interchannelCoupling_next_row");
41  pixelFlag_ = true;
42  LogDebug("PixelBrickedDigitizerAlgorithm")
43  << "Algorithm constructed "
44  << "Configuration parameters:"
45  << "Threshold/Gain = "
46  << "threshold in electron Endcap = " << theThresholdInE_Endcap_
47  << "threshold in electron Barrel = " << theThresholdInE_Barrel_ << " " << theElectronPerADC_ << " "
48  << theAdcFullScale_ << " The delta cut-off is set to " << tMax_ << " pix-inefficiency " << addPixelInefficiency_;
49 }
51  LogDebug("PixelBrickedDigitizerAlgorithm") << "Algorithm deleted";
52 }
54  const size_t hitIndex,
55  const uint32_t tofBin,
56  const Phase2TrackerGeomDetUnit* pixdet,
57  const std::vector<DigitizerUtility::SignalPoint>& collection_points) {
58  // X - Rows, Left-Right, 160, (1.6cm) for barrel
59  // Y - Columns, Down-Up, 416, (6.4cm)
60  const Phase2TrackerTopology* topol = &pixdet->specificTopology();
61  uint32_t detID = pixdet->geographicalId().rawId();
62  signal_map_type& theSignal = _signal[detID];
63 
64  // local map to store pixels hit by 1 Hit.
65  using hit_map_type = std::map<int, float, std::less<int> >;
66  hit_map_type hit_signal;
67 
68  // Assign signals to readout channels and store sorted by channel number
69  // Iterate over collection points on the collection plane
70  for (auto const& v : collection_points) {
71  float CloudCenterX = v.position().x(); // Charge position in x
72  float CloudCenterY = v.position().y(); // in y
73  float SigmaX = v.sigma_x(); // Charge spread in x
74  float SigmaY = v.sigma_y(); // in y
75  float Charge = v.amplitude(); // Charge amplitude
76 
77  // Find the maximum cloud spread in 2D plane , assume 3*sigma
78  float CloudRight = CloudCenterX + clusterWidth_ * SigmaX;
79  float CloudLeft = CloudCenterX - clusterWidth_ * SigmaX;
80  float CloudUp = CloudCenterY + clusterWidth_ * SigmaY;
81  float CloudDown = CloudCenterY - clusterWidth_ * SigmaY;
82 
83  // Define 2D cloud limit points
84  LocalPoint PointRightUp = LocalPoint(CloudRight, CloudUp);
85  LocalPoint PointLeftDown = LocalPoint(CloudLeft, CloudDown);
86 
87  // This points can be located outside the sensor area.
88  // The conversion to measurement point does not check for that
89  // so the returned pixel index might be wrong (outside range).
90  // We rely on the limits check below to fix this.
91  // But remember whatever we do here THE CHARGE OUTSIDE THE ACTIVE
92  // PIXEL AREA IS LOST, it should not be collected.
93 
94  // Convert the 2D points to pixel indices
95  MeasurementPoint mp = topol->measurementPosition(PointRightUp);
96  //MeasurementPoint mp_bricked = topol->measurementPosition(PointRightUp);
97  int IPixRightUpX = static_cast<int>(std::floor(mp.x())); // cast reqd.
98  //int IPixRightUpY = static_cast<int>(std::floor(mp.y()));
99 
100  int numColumns = topol->ncolumns(); // det module number of cols&rows
101  int numRows = topol->nrows();
102  IPixRightUpX = numRows > IPixRightUpX ? IPixRightUpX : numRows - 1;
103 
104  //Specific to bricked geometry
105  int IPixRightUpY = static_cast<int>(mp.y() - 0.5 * (IPixRightUpX % 2));
106 
107  mp = topol->measurementPosition(PointLeftDown);
108 
109  int IPixLeftDownX = static_cast<int>(std::floor(mp.x()));
110 
111  IPixLeftDownX = 0 < IPixLeftDownX ? IPixLeftDownX : 0;
112 
113  //Specific to bricked geometry
114  int IPixLeftDownY = static_cast<int>(mp.y() - 0.5 * (IPixLeftDownX % 2)); //changed in case negative value
115 
116  IPixRightUpY = numColumns > IPixRightUpY ? IPixRightUpY : numColumns - 1;
117  IPixLeftDownY = 0 < IPixLeftDownY ? IPixLeftDownY : 0;
118 
119  // First integrate charge strips in x
120  hit_map_type x;
121  for (int ix = IPixLeftDownX; ix <= IPixRightUpX; ++ix) { // loop over x index
122  float xLB, LowerBound;
123  // Why is set to 0 if ix=0, does it meen that we accept charge
124  // outside the sensor?
125  if (ix == 0 || SigmaX == 0.) { // skip for surface segemnts
126  LowerBound = 0.;
127  } else {
128  mp = MeasurementPoint(ix, 0.0);
129  xLB = topol->localPosition(mp).x();
130  LowerBound = 1 - calcQ((xLB - CloudCenterX) / SigmaX);
131  }
132 
133  float xUB, UpperBound;
134  if (ix == numRows - 1 || SigmaX == 0.) {
135  UpperBound = 1.;
136  } else {
137  mp = MeasurementPoint(ix + 1, 0.0);
138  xUB = topol->localPosition(mp).x();
139  UpperBound = 1. - calcQ((xUB - CloudCenterX) / SigmaX);
140  }
141  float TotalIntegrationRange = UpperBound - LowerBound; // get strip
142  x.emplace(ix, TotalIntegrationRange); // save strip integral
143  }
144 
145  // Now integrate strips in y. Two maps will be filled: y and y_bricked which will both be used for the induced signal.
146  int IPixLeftDownY_bricked = IPixLeftDownY;
147  int IPixRightUpY_bricked = IPixRightUpY;
148 
149  //Specific to bricked geometry
150  IPixRightUpY = std::min(IPixRightUpY + int((IPixRightUpX % 2)), numColumns - 1);
151 
152  //This map will be twice as large as the non-bricked hit map in y to harbor both the integrated charge from the bricked and non-bricked columns.
153  hit_map_type y;
154  for (int iy = IPixLeftDownY; iy <= IPixRightUpY; ++iy) { // loop over y index
155  float yLB, LowerBound;
156  if (iy == 0 || SigmaY == 0.) {
157  LowerBound = 0.;
158  } else {
159  mp = MeasurementPoint(0.0, iy);
160  yLB = topol->localPosition(mp).y();
161  LowerBound = 1. - calcQ((yLB - CloudCenterY) / SigmaY);
162  }
163 
164  float yUB, UpperBound;
165  if (iy >= numColumns - 1 || SigmaY == 0.) {
166  UpperBound = 1.;
167  } else {
168  mp = MeasurementPoint(0.0, iy + 1);
169  yUB = topol->localPosition(mp).y();
170  UpperBound = 1. - calcQ((yUB - CloudCenterY) / SigmaY);
171  }
172 
173  float TotalIntegrationRange = UpperBound - LowerBound;
174 
175  //Even indices correspond to the non-bricked columns
176  y.emplace(2 * iy, TotalIntegrationRange); // save strip integral
177  }
178 
179  IPixLeftDownY_bricked = std::max(IPixLeftDownY_bricked - int((!(IPixLeftDownX % 2))), 0);
180  for (int iy = IPixLeftDownY_bricked; iy <= IPixRightUpY_bricked; ++iy) { // loop over y index
181  float yLB, LowerBound;
182  if (iy == 0 || SigmaY == 0.) {
183  LowerBound = 0.;
184  } else {
185  mp = MeasurementPoint(0.0, iy + 0.5);
186  yLB = topol->localPosition(mp).y();
187  LowerBound = 1. - calcQ((yLB - CloudCenterY) / SigmaY);
188  }
189 
190  float yUB, UpperBound;
191  if (iy >= numColumns || SigmaY == 0.) { // This was changed for bricked pixels
192  UpperBound = 1.;
193  } else {
194  mp = MeasurementPoint(0.0, iy + 1.5);
195  yUB = topol->localPosition(mp).y();
196  UpperBound = 1. - calcQ((yUB - CloudCenterY) / SigmaY);
197  }
198 
199  float TotalIntegrationRange = UpperBound - LowerBound;
200  //Odd indices correspond to bricked columns
201  y.emplace(2 * iy + 1, TotalIntegrationRange); // save strip integral
202  } //loop over y index
203 
204  // Get the 2D charge integrals by folding x and y strips
205  for (int ix = IPixLeftDownX; ix <= IPixRightUpX; ++ix) { // loop over x index
206  for (int iy = std::max(0, IPixLeftDownY - int((ix % 2))); iy <= IPixRightUpY; ++iy) { // loop over y index
207  int iy_considered = iy * 2 + ix % 2;
208  float ChargeFraction = Charge * x[ix] * y[iy_considered];
209 
210  int chanFired = -1;
211  if (ChargeFraction > 0.) {
212  chanFired =
214  // Load the ampl
215  hit_signal[chanFired] += ChargeFraction;
216  }
217  }
218  } //x loop
219  } //collection loop
220  // Fill the global map with all hit pixels from this event
221  float corr_time = hit.tof() - pixdet->surface().toGlobal(hit.localPosition()).mag() * c_inv;
222  for (auto const& hit_s : hit_signal) {
223  int chan = hit_s.first;
224  theSignal[chan] +=
225  (makeDigiSimLinks_ ? DigitizerUtility::Amplitude(hit_s.second, &hit, hit_s.second, corr_time, hitIndex, tofBin)
226  : DigitizerUtility::Amplitude(hit_s.second, nullptr, hit_s.second));
227  }
228 }
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
Point3DBase< Scalar, LocalTag > LocalPoint
Definition: Definitions.h:30
virtual LocalPoint localPosition(const MeasurementPoint &) const =0
virtual int ncolumns() const =0
virtual int nrows() const =0
PixelBrickedDigitizerAlgorithm(const edm::ParameterSet &conf, edm::ConsumesCollector iC)
T x() const
Definition: PV2DBase.h:43
void induce_signal(const PSimHit &hit, const size_t hitIndex, const unsigned int tofBin, const Phase2TrackerGeomDetUnit *pixdet, const std::vector< DigitizerUtility::SignalPoint > &collection_points) override
static int pixelToChannel(int row, int col)
Definition: PixelDigi.h:75
T y() const
Definition: PV2DBase.h:44
T x() const
Definition: PV3DBase.h:59
T y() const
Definition: PV3DBase.h:60
Measurement2DPoint MeasurementPoint
Measurement points are two-dimensional by default.
static PackedDigiType pixelToChannel(unsigned int row, unsigned int col)
virtual MeasurementPoint measurementPosition(const LocalPoint &) const =0
DetId geographicalId() const
The label of this GeomDet.
Definition: GeomDet.h:64
std::map< int, DigitizerUtility::Amplitude, std::less< int > > signal_map_type
constexpr double c_inv
GlobalPoint toGlobal(const Point2DBase< Scalar, LocalTag > lp) const
Definition: Surface.h:79
const Plane & surface() const
The nominal surface of the GeomDet.
Definition: GeomDet.h:37
T mag() const
The vector magnitude. Equivalent to sqrt(vec.mag2())
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
chan
lumi = TPaveText(lowX+0.38, lowY+0.061, lowX+0.45, lowY+0.161, "NDC") lumi.SetBorderSize( 0 ) lumi...
HLT enums.
virtual const PixelTopology & specificTopology() const
Returns a reference to the pixel proxy topology.
#define LogDebug(id)