CMS 3D CMS Logo

HGCalCell.cc
Go to the documentation of this file.
3 #include <vector>
4 
5 HGCalCell::HGCalCell(double waferSize, int32_t nFine, int32_t nCoarse) {
6  ncell_[0] = nFine;
7  ncell_[1] = nCoarse;
8  for (int k = 0; k < 2; ++k) {
9  cellX_[k] = waferSize / (3 * ncell_[k]);
10  cellY_[k] = sqrt3By2_ * cellX_[k];
11  }
12  edm::LogVerbatim("HGCalGeom") << "HGCalCell initialized with waferSize " << waferSize << " number of cells " << nFine
13  << ":" << nCoarse;
14 }
15 
16 std::pair<double, double> HGCalCell::cellUV2XY1(int32_t u, int32_t v, int32_t placementIndex, int32_t type) {
17  if (type != 0)
18  type = 1;
19  double x(0), y(0);
20  switch (placementIndex) {
22  x = (1.5 * (v - u) + 0.5) * cellX_[type];
23  y = (v + u - 2 * ncell_[type] + 1) * cellY_[type];
24  break;
26  x = (1.5 * (v - ncell_[type]) + 1.0) * cellX_[type];
27  y = (2 * u - v - ncell_[type]) * cellY_[type];
28  break;
30  x = (1.5 * (u - ncell_[type]) + 0.5) * cellX_[type];
31  y = -(2 * v - u - ncell_[type] + 1) * cellY_[type];
32  break;
34  x = -(1.5 * (v - u) + 0.5) * cellX_[type];
35  y = -(v + u - 2 * ncell_[type] + 1) * cellY_[type];
36  break;
38  x = -(1.5 * (v - ncell_[type]) + 1) * cellX_[type];
39  y = -(2 * u - v - ncell_[type]) * cellY_[type];
40  break;
42  x = -(1.5 * (u - ncell_[type]) + 0.5) * cellX_[type];
43  y = (2 * v - u - ncell_[type] + 1) * cellY_[type];
44  break;
46  x = (1.5 * (u - v) - 0.5) * cellX_[type];
47  y = (v + u - 2 * ncell_[type] + 1) * cellY_[type];
48  break;
50  x = -(1.5 * (v - ncell_[type]) + 1.0) * cellX_[type];
51  y = (2 * u - v - ncell_[type]) * cellY_[type];
52  break;
54  x = -(1.5 * (u - ncell_[type]) + 0.5) * cellX_[type];
55  y = -(2 * v - u - ncell_[type] + 1) * cellY_[type];
56  break;
58  x = -(1.5 * (u - v) - 0.5) * cellX_[type];
59  y = -(v + u - 2 * ncell_[type] + 1) * cellY_[type];
60  break;
62  x = (1.5 * (v - ncell_[type]) + 1) * cellX_[type];
63  y = -(2 * u - v - ncell_[type]) * cellY_[type];
64  break;
65  default:
66  x = (1.5 * (u - ncell_[type]) + 0.5) * cellX_[type];
67  y = (2 * v - u - ncell_[type] + 1) * cellY_[type];
68  break;
69  }
70  return std::make_pair(x, y);
71 }
72 
73 std::pair<double, double> HGCalCell::cellUV2XY2(int32_t u, int32_t v, int32_t placementIndex, int32_t type) {
74  if (type != 0)
75  type = 1;
76  double x(0), y(0);
77  if (placementIndex < HGCalCell::cellPlacementExtra) {
78  double x0 = (1.5 * (u - v) - 0.5) * cellX_[type];
79  double y0 = (u + v - 2 * ncell_[type] + 1) * cellY_[type];
80  const std::vector<double> fcos = {1.0, 0.5, -0.5, -1.0, -0.5, 0.5};
81  const std::vector<double> fsin = {0.0, sqrt3By2_, sqrt3By2_, 0.0, -sqrt3By2_, -sqrt3By2_};
82  x = x0 * fcos[placementIndex] - y0 * fsin[placementIndex];
83  y = x0 * fsin[placementIndex] + y0 * fcos[placementIndex];
84  } else {
85  double x0 = (1.5 * (v - ncell_[type]) + 1.0) * cellX_[type];
86  double y0 = (2 * u - v - ncell_[type]) * cellY_[type];
87  const std::vector<double> fcos = {0.5, 1.0, 0.5, -0.5, -1.0, -0.5};
88  const std::vector<double> fsin = {sqrt3By2_, 0.0, -sqrt3By2_, -sqrt3By2_, 0.0, sqrt3By2_};
89  x = x0 * fcos[placementIndex - HGCalCell::cellPlacementExtra] -
90  y0 * fsin[placementIndex - HGCalCell::cellPlacementExtra];
91  y = x0 * fsin[placementIndex - HGCalCell::cellPlacementExtra] +
92  y0 * fcos[placementIndex - HGCalCell::cellPlacementExtra];
93  }
94  return std::make_pair(x, y);
95 }
96 
97 std::pair<int, int> HGCalCell::cellUV2Cell(int32_t u, int32_t v, int32_t placementIndex, int32_t type) {
98  if (type != 0)
99  type = 1;
100  int cell(0), cellx(0), cellt(HGCalCell::fullCell);
101  if (placementIndex >= HGCalCell::cellPlacementExtra) {
102  const std::vector<int> itype0 = {0, 7, 8, 9, 10, 11, 6, 3, 4, 5, 4, 5, 3};
103  const std::vector<int> itype1 = {0, 0, 1, 2, 3, 4, 5, 0, 1, 2, 0, 1, 2};
104  const std::vector<int> itype2 = {0, 11, 6, 7, 8, 9, 10, 5, 3, 4, 3, 4, 5};
105  const std::vector<int> itype3 = {0, 4, 5, 0, 1, 2, 3, 2, 0, 1, 2, 0, 1};
106  const std::vector<int> itype4 = {0, 9, 10, 11, 6, 7, 8, 4, 5, 3, 5, 3, 4};
107  const std::vector<int> itype5 = {0, 2, 3, 4, 5, 0, 1, 1, 2, 0, 1, 2, 0};
108  if (u == 0 && v == 0) {
109  cellx = 1;
110  cellt = HGCalCell::cornerCell;
111  } else if (u == 0 && (v - u) == (ncell_[type] - 1)) {
112  cellx = 2;
113  cellt = HGCalCell::cornerCell;
114  } else if ((v - u) == (ncell_[type] - 1) && v == (2 * ncell_[type] - 1)) {
115  cellx = 3;
116  cellt = HGCalCell::cornerCell;
117  } else if (u == (2 * ncell_[type] - 1) && v == (2 * ncell_[type] - 1)) {
118  cellx = 4;
119  cellt = HGCalCell::cornerCell;
120  } else if (u == (2 * ncell_[type] - 1) && (u - v) == ncell_[type]) {
121  cellx = 5;
122  cellt = HGCalCell::cornerCell;
123  } else if ((u - v) == ncell_[type] && v == 0) {
124  cellx = 6;
125  cellt = HGCalCell::cornerCell;
126  } else if (u == 0) {
127  cellx = 7;
128  cellt = HGCalCell::truncatedCell;
129  } else if ((v - u) == (ncell_[type] - 1)) {
130  cellx = 10;
131  cellt = HGCalCell::extendedCell;
132  } else if (v == (2 * ncell_[type] - 1)) {
133  cellx = 8;
134  cellt = HGCalCell::truncatedCell;
135  } else if (u == (2 * ncell_[type] - 1)) {
136  cellx = 11;
137  cellt = HGCalCell::extendedCell;
138  } else if ((u - v) == ncell_[type]) {
139  cellx = 9;
140  cellt = HGCalCell::truncatedCell;
141  } else if (v == 0) {
142  cellx = 12;
143  cellt = HGCalCell::extendedCell;
144  }
145  switch (placementIndex) {
147  cell = itype0[cellx];
148  break;
150  cell = itype1[cellx];
151  break;
153  cell = itype2[cellx];
154  break;
156  cell = itype3[cellx];
157  break;
159  cell = itype4[cellx];
160  break;
161  default:
162  cell = itype5[cellx];
163  break;
164  }
165  } else {
166  const std::vector<int> itype0 = {0, 1, 2, 3, 4, 5, 0, 1, 2, 0, 0, 1, 2};
167  const std::vector<int> itype1 = {0, 8, 9, 10, 11, 6, 7, 4, 5, 3, 4, 5, 3};
168  const std::vector<int> itype2 = {0, 3, 4, 5, 0, 1, 2, 2, 0, 1, 1, 2, 0};
169  const std::vector<int> itype3 = {0, 10, 11, 6, 7, 8, 9, 5, 3, 4, 5, 3, 4};
170  const std::vector<int> itype4 = {0, 5, 0, 1, 2, 3, 4, 0, 1, 2, 2, 0, 1};
171  const std::vector<int> itype5 = {0, 6, 7, 8, 9, 10, 11, 3, 4, 5, 3, 4, 5};
172  if (u == 0 && v == 0) {
173  cellx = 1;
174  cellt = HGCalCell::cornerCell;
175  } else if (v == 0 && (u - v) == (ncell_[type])) {
176  cellx = 2;
177  cellt = HGCalCell::cornerCell;
178  } else if ((u - v) == (ncell_[type]) && u == (2 * ncell_[type] - 1)) {
179  cellx = 3;
180  cellt = HGCalCell::cornerCell;
181  } else if (u == (2 * ncell_[type] - 1) && v == (2 * ncell_[type] - 1)) {
182  cellx = 4;
183  cellt = HGCalCell::cornerCell;
184  } else if (v == (2 * ncell_[type] - 1) && (v - u) == (ncell_[type] - 1)) {
185  cellx = 5;
186  cellt = HGCalCell::cornerCell;
187  } else if ((v - u) == (ncell_[type] - 1) && u == 0) {
188  cellx = 6;
189  cellt = HGCalCell::cornerCell;
190  } else if (v == 0) {
191  cellx = 10;
192  cellt = HGCalCell::extendedCell;
193  } else if ((u - v) == ncell_[type]) {
194  cellx = 7;
195  cellt = HGCalCell::truncatedCell;
196  } else if (u == (2 * ncell_[type] - 1)) {
197  cellx = 11;
198  cellt = HGCalCell::extendedCell;
199  } else if (v == (2 * ncell_[type] - 1)) {
200  cellx = 8;
201  cellt = HGCalCell::truncatedCell;
202  } else if ((v - u) == (ncell_[type] - 1)) {
203  cellx = 12;
204  cellt = HGCalCell::extendedCell;
205  } else if (u == 0) {
206  cellx = 9;
207  cellt = HGCalCell::truncatedCell;
208  }
209  switch (placementIndex) {
211  cell = itype0[cellx];
212  break;
214  cell = itype1[cellx];
215  break;
217  cell = itype2[cellx];
218  break;
220  cell = itype3[cellx];
221  break;
223  cell = itype4[cellx];
224  break;
225  default:
226  cell = itype5[cellx];
227  break;
228  }
229  }
230  return std::make_pair(cell, cellt);
231 }
232 
233 int HGCalCell::cellPlacementIndex(int32_t iz, int32_t frontBack, int32_t orient) {
234  int32_t indx = ((iz * frontBack) > 0) ? orient : (orient + HGCalCell::cellPlacementExtra);
235  return indx;
236 }
237 
238 std::pair<int32_t, int32_t> HGCalCell::cellOrient(int32_t placementIndex) {
239  int32_t orient = (placementIndex >= HGCalCell::cellPlacementExtra) ? (placementIndex - HGCalCell::cellPlacementExtra)
240  : placementIndex;
241  int32_t frontBackZside = (placementIndex >= HGCalCell::cellPlacementExtra) ? 1 : -1;
242  return std::make_pair(orient, frontBackZside);
243 }
244 
245 std::pair<int32_t, int32_t> HGCalCell::cellType(int32_t u, int32_t v, int32_t ncell, int32_t placementIndex) {
246  int cell(0), cellx(0), cellt(HGCalCell::fullCell);
247  if (placementIndex >= HGCalCell::cellPlacementExtra) {
248  const std::vector<int> itype0 = {HGCalCell::centralCell,
261  const std::vector<int> itype1 = {HGCalCell::centralCell,
274  const std::vector<int> itype2 = {HGCalCell::centralCell,
287  const std::vector<int> itype3 = {HGCalCell::centralCell,
300  const std::vector<int> itype4 = {HGCalCell::centralCell,
313  const std::vector<int> itype5 = {HGCalCell::centralCell,
326  if (u == 0 && v == 0) {
327  cellx = 1;
328  cellt = HGCalCell::cornerCell;
329  } else if (u == 0 && (v - u) == (ncell - 1)) {
330  cellx = 2;
331  cellt = HGCalCell::cornerCell;
332  } else if ((v - u) == (ncell - 1) && v == (2 * ncell - 1)) {
333  cellx = 3;
334  cellt = HGCalCell::cornerCell;
335  } else if (u == (2 * ncell - 1) && v == (2 * ncell - 1)) {
336  cellx = 4;
337  cellt = HGCalCell::cornerCell;
338  } else if (u == (2 * ncell - 1) && (u - v) == ncell) {
339  cellx = 5;
340  cellt = HGCalCell::cornerCell;
341  } else if ((u - v) == ncell && v == 0) {
342  cellx = 6;
343  cellt = HGCalCell::cornerCell;
344  } else if (u == 0) {
345  cellx = 7;
346  cellt = HGCalCell::truncatedCell;
347  } else if ((v - u) == (ncell - 1)) {
348  cellx = 10;
349  cellt = HGCalCell::extendedCell;
350  } else if (v == (2 * ncell - 1)) {
351  cellx = 8;
352  cellt = HGCalCell::truncatedCell;
353  } else if (u == (2 * ncell - 1)) {
354  cellx = 11;
355  cellt = HGCalCell::extendedCell;
356  } else if ((u - v) == ncell) {
357  cellx = 9;
358  cellt = HGCalCell::truncatedCell;
359  } else if (v == 0) {
360  cellx = 12;
361  cellt = HGCalCell::extendedCell;
362  }
363  switch (placementIndex) {
365  cell = itype0[cellx];
366  break;
368  cell = itype1[cellx];
369  break;
371  cell = itype2[cellx];
372  break;
374  cell = itype3[cellx];
375  break;
377  cell = itype4[cellx];
378  break;
379  default:
380  cell = itype5[cellx];
381  break;
382  }
383  } else {
384  const std::vector<int> itype0 = {HGCalCell::centralCell,
397  const std::vector<int> itype1 = {HGCalCell::centralCell,
410  const std::vector<int> itype2 = {HGCalCell::centralCell,
423  const std::vector<int> itype3 = {HGCalCell::centralCell,
436  const std::vector<int> itype4 = {HGCalCell::centralCell,
449  const std::vector<int> itype5 = {HGCalCell::centralCell,
462  if (u == 0 && v == 0) {
463  cellx = 1;
464  cellt = HGCalCell::cornerCell;
465  } else if (v == 0 && (u - v) == (ncell)) {
466  cellx = 2;
467  cellt = HGCalCell::cornerCell;
468  } else if ((u - v) == (ncell) && u == (2 * ncell - 1)) {
469  cellx = 3;
470  cellt = HGCalCell::cornerCell;
471  } else if (u == (2 * ncell - 1) && v == (2 * ncell - 1)) {
472  cellx = 4;
473  cellt = HGCalCell::cornerCell;
474  } else if (v == (2 * ncell - 1) && (v - u) == (ncell - 1)) {
475  cellx = 5;
476  cellt = HGCalCell::cornerCell;
477  } else if ((v - u) == (ncell - 1) && u == 0) {
478  cellx = 6;
479  cellt = HGCalCell::cornerCell;
480  } else if (v == 0) {
481  cellx = 10;
482  cellt = HGCalCell::extendedCell;
483  } else if ((u - v) == ncell) {
484  cellx = 7;
485  cellt = HGCalCell::truncatedCell;
486  } else if (u == (2 * ncell - 1)) {
487  cellx = 11;
488  cellt = HGCalCell::extendedCell;
489  } else if (v == (2 * ncell - 1)) {
490  cellx = 8;
491  cellt = HGCalCell::truncatedCell;
492  } else if ((v - u) == (ncell - 1)) {
493  cellx = 12;
494  cellt = HGCalCell::extendedCell;
495  } else if (u == 0) {
496  cellx = 9;
497  cellt = HGCalCell::truncatedCell;
498  }
499  switch (placementIndex) {
501  cell = itype0[cellx];
502  break;
504  cell = itype1[cellx];
505  break;
507  cell = itype2[cellx];
508  break;
510  cell = itype3[cellx];
511  break;
513  cell = itype4[cellx];
514  break;
515  default:
516  cell = itype5[cellx];
517  break;
518  }
519  }
520  return std::make_pair(cell, cellt);
521 }
Log< level::Info, true > LogVerbatim
static int32_t cellPlacementIndex(int32_t iz, int32_t frontBack, int32_t orient)
Definition: HGCalCell.cc:233
static constexpr int32_t fullCell
Definition: HGCalCell.h:28
const double sqrt3By2_
Definition: HGCalCell.h:60
static constexpr int32_t cellPlacementIndex8
Definition: HGCalCell.h:19
static std::pair< int32_t, int32_t > cellType(int32_t u, int32_t v, int32_t ncell, int32_t placementIndex)
Definition: HGCalCell.cc:245
static constexpr int32_t topCorner
Definition: HGCalCell.h:44
static constexpr int32_t cellPlacementIndex3
Definition: HGCalCell.h:14
static constexpr int32_t leftEdge
Definition: HGCalCell.h:36
int32_t ncell_[2]
Definition: HGCalCell.h:61
static constexpr int32_t cellPlacementIndex10
Definition: HGCalCell.h:21
std::pair< int32_t, int32_t > cellUV2Cell(int32_t u, int32_t v, int32_t placementIndex, int32_t type)
Definition: HGCalCell.cc:97
static constexpr int32_t cellPlacementIndex0
Definition: HGCalCell.h:11
static constexpr int32_t topLeftCorner
Definition: HGCalCell.h:43
static std::pair< int32_t, int32_t > cellOrient(int32_t placementIndex)
Definition: HGCalCell.cc:238
static constexpr int32_t bottomLeftCorner
Definition: HGCalCell.h:42
static constexpr int32_t topRightCorner
Definition: HGCalCell.h:45
static constexpr int32_t bottomRightEdge
Definition: HGCalCell.h:40
static constexpr int32_t cellPlacementIndex9
Definition: HGCalCell.h:20
HGCalCell(double waferSize, int32_t nFine, int32_t nCoarse)
Definition: HGCalCell.cc:5
static constexpr int32_t bottomLeftEdge
Definition: HGCalCell.h:35
static constexpr int32_t cellPlacementIndex7
Definition: HGCalCell.h:18
static constexpr int32_t bottomRightCorner
Definition: HGCalCell.h:46
static constexpr int32_t truncatedCell
Definition: HGCalCell.h:30
static constexpr int32_t centralCell
Definition: HGCalCell.h:34
std::pair< double, double > cellUV2XY1(int32_t u, int32_t v, int32_t placementIndex, int32_t type)
Definition: HGCalCell.cc:16
static constexpr int32_t extendedCell
Definition: HGCalCell.h:31
static constexpr int32_t topRightEdge
Definition: HGCalCell.h:38
static constexpr int32_t cellPlacementIndex11
Definition: HGCalCell.h:22
std::pair< double, double > cellUV2XY2(int32_t u, int32_t v, int32_t placementIndex, int32_t type)
Definition: HGCalCell.cc:73
static constexpr int32_t rightEdge
Definition: HGCalCell.h:39
static constexpr int32_t bottomCorner
Definition: HGCalCell.h:41
static constexpr int32_t cornerCell
Definition: HGCalCell.h:29
static constexpr int32_t cellPlacementIndex4
Definition: HGCalCell.h:15
static constexpr int32_t cellPlacementIndex1
Definition: HGCalCell.h:12
double cellX_[2]
Definition: HGCalCell.h:62
static constexpr int32_t cellPlacementIndex6
Definition: HGCalCell.h:17
double cellY_[2]
Definition: HGCalCell.h:62
static constexpr int32_t cellPlacementIndex2
Definition: HGCalCell.h:13
static constexpr int32_t topLeftEdge
Definition: HGCalCell.h:37
static constexpr int32_t cellPlacementExtra
Definition: HGCalCell.h:24