CMS 3D CMS Logo

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