CMS 3D CMS Logo

List of all members | Public Member Functions | Public Attributes
HGCalEEAlgo Struct Reference

Public Member Functions

void ConstructAlgo (cms::DDParsingContext &ctxt, xml_h e)
 
void ConstructLayers (const dd4hep::Volume module, cms::DDParsingContext &ctxt, xml_h e)
 
 HGCalEEAlgo ()=delete
 
 HGCalEEAlgo (cms::DDParsingContext &ctxt, xml_h e)
 
void PositionSensitive (cms::DDParsingContext &ctxt, xml_h e, const dd4hep::Volume &glog, double rin, double rout, double zpos, int layertype, int layercenter)
 

Public Attributes

int absorbMode_
 
double alpha_
 
int choiceType_
 
std::unordered_set< int > copies_
 
std::vector< int > copyNumber_
 
double cosAlpha_
 
int firstLayer_
 
double fracAreaMin_
 
HGCalGeomTools geomTools_
 
std::vector< int > layerCenter_
 
std::vector< int > layers_
 
std::vector< int > layerSense_
 
std::vector< double > layerThick_
 
std::vector< int > layerType_
 
std::vector< std::string > materials_
 
dd4hep::Volume mother_
 
std::vector< std::string > names_
 
int nCutRadPar_
 
std::vector< double > rad100to200_
 
std::vector< double > rad200to300_
 
std::vector< double > rMaxFront_
 
std::vector< double > rMinFront_
 
int sectors_
 
int sensitiveMode_
 
std::vector< double > slopeB_
 
std::vector< double > slopeT_
 
std::vector< double > thick_
 
std::vector< std::string > wafers_
 
double waferSepar_
 
double waferSize_
 
std::unique_ptr< HGCalWaferTypewaferType_
 
std::vector< double > zFrontB_
 
std::vector< double > zFrontT_
 
double zMinBlock_
 
double zMinRadPar_
 

Detailed Description

Definition at line 29 of file DDHGCalEEAlgo.cc.

Constructor & Destructor Documentation

◆ HGCalEEAlgo() [1/2]

HGCalEEAlgo::HGCalEEAlgo ( )
delete

◆ HGCalEEAlgo() [2/2]

HGCalEEAlgo::HGCalEEAlgo ( cms::DDParsingContext ctxt,
xml_h  e 
)
inline

Definition at line 68 of file DDHGCalEEAlgo.cc.

References writedatasetfile::args, cms::convert2mm(), angle_units::operators::convertRadToDeg(), funct::cos(), MillePedeFileConverter_cfg::e, mps_fire::i, cuy::ii, dqmdumpme::k, cms::DDNamespace::name(), and cms::DDNamespace::volume().

68  {
69  cms::DDNamespace ns(ctxt, e, true);
71 
72  mother_ = ns.volume(args.parentName());
73  wafers_ = args.value<std::vector<std::string>>("WaferNames");
74 #ifdef EDM_ML_DEBUG
75  edm::LogVerbatim("HGCalGeom") << "DDHGCalEEAlgo: " << wafers_.size() << " wafers";
76  for (unsigned int i = 0; i < wafers_.size(); ++i)
77  edm::LogVerbatim("HGCalGeom") << "Wafer[" << i << "] " << wafers_[i];
78 #endif
79 
80  materials_ = args.value<std::vector<std::string>>("MaterialNames");
81  names_ = args.value<std::vector<std::string>>("VolumeNames");
82  thick_ = args.value<std::vector<double>>("Thickness");
83  copyNumber_.resize(materials_.size(), 1);
84 #ifdef EDM_ML_DEBUG
85  edm::LogVerbatim("HGCalGeom") << "DDHGCalEEAlgo: " << materials_.size() << " types of volumes";
86  for (unsigned int i = 0; i < names_.size(); ++i)
87  edm::LogVerbatim("HGCalGeom") << "Volume [" << i << "] " << names_[i] << " of thickness "
88  << cms::convert2mm(thick_[i]) << " filled with " << materials_[i]
89  << " first copy number " << copyNumber_[i];
90 #endif
91 
92  layers_ = args.value<std::vector<int>>("Layers");
93  layerThick_ = args.value<std::vector<double>>("LayerThick");
94 #ifdef EDM_ML_DEBUG
95  edm::LogVerbatim("HGCalGeom") << "There are " << layers_.size() << " blocks";
96  for (unsigned int i = 0; i < layers_.size(); ++i)
97  edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] of thickness " << cms::convert2mm(layerThick_[i])
98  << " with " << layers_[i] << " layers";
99 #endif
100 
101  layerType_ = args.value<std::vector<int>>("LayerType");
102  layerSense_ = args.value<std::vector<int>>("LayerSense");
103  firstLayer_ = args.value<int>("FirstLayer");
104  absorbMode_ = args.value<int>("AbsorberMode");
105  sensitiveMode_ = args.value<int>("SensitiveMode");
106 #ifdef EDM_ML_DEBUG
107  edm::LogVerbatim("HGCalGeom") << "First Layer " << firstLayer_ << " and "
108  << "Absober:Sensitive mode " << absorbMode_ << ":" << sensitiveMode_;
109 #endif
110  layerCenter_ = args.value<std::vector<int>>("LayerCenter");
111 #ifdef EDM_ML_DEBUG
112  for (unsigned int i = 0; i < layerCenter_.size(); ++i)
113  edm::LogVerbatim("HGCalGeom") << "LayerCenter [" << i << "] " << layerCenter_[i];
114 #endif
115  if (firstLayer_ > 0) {
116  for (unsigned int i = 0; i < layerType_.size(); ++i) {
117  if (layerSense_[i] > 0) {
118  int ii = layerType_[i];
120 #ifdef EDM_ML_DEBUG
121  edm::LogVerbatim("HGCalGeom") << "First copy number for layer type " << i << ":" << ii << " with "
122  << materials_[ii] << " changed to " << copyNumber_[ii];
123 #endif
124  break;
125  }
126  }
127  } else {
128  firstLayer_ = 1;
129  }
130 #ifdef EDM_ML_DEBUG
131  edm::LogVerbatim("HGCalGeom") << "There are " << layerType_.size() << " layers";
132  for (unsigned int i = 0; i < layerType_.size(); ++i)
133  edm::LogVerbatim("HGCalGeom") << "Layer [" << i << "] with material type " << layerType_[i] << " sensitive class "
134  << layerSense_[i];
135 #endif
136  zMinBlock_ = args.value<double>("zMinBlock");
137 
138  rad100to200_ = args.value<std::vector<double>>("rad100to200");
139  rad200to300_ = args.value<std::vector<double>>("rad200to300");
140  zMinRadPar_ = args.value<double>("zMinForRadPar");
141  choiceType_ = args.value<int>("choiceType");
142  nCutRadPar_ = args.value<int>("nCornerCut");
143  fracAreaMin_ = args.value<double>("fracAreaMin");
144  waferSize_ = args.value<double>("waferSize");
145  waferSepar_ = args.value<double>("SensorSeparation");
146  sectors_ = args.value<int>("Sectors");
147  alpha_ = (1._pi) / sectors_;
148  cosAlpha_ = cos(alpha_);
149 #ifdef EDM_ML_DEBUG
150  edm::LogVerbatim("HGCalGeom") << "zStart " << cms::convert2mm(zMinBlock_)
151  << " radius for wafer type separation uses " << rad100to200_.size()
152  << " parameters; zmin " << cms::convert2mm(zMinRadPar_) << " cutoff " << choiceType_
153  << ":" << nCutRadPar_ << ":" << fracAreaMin_ << " wafer width "
154  << cms::convert2mm(waferSize_) << " separations " << cms::convert2mm(waferSepar_)
155  << " sectors " << sectors_ << ":" << convertRadToDeg(alpha_) << ":" << cosAlpha_;
156  for (unsigned int k = 0; k < rad100to200_.size(); ++k)
157  edm::LogVerbatim("HGCalGeom") << "[" << k << "] 100-200 " << rad100to200_[k] << " 200-300 " << rad200to300_[k];
158 #endif
159 
160  slopeB_ = args.value<std::vector<double>>("SlopeBottom");
161  zFrontB_ = args.value<std::vector<double>>("ZFrontBottom");
162  rMinFront_ = args.value<std::vector<double>>("RMinFront");
163  slopeT_ = args.value<std::vector<double>>("SlopeTop");
164  zFrontT_ = args.value<std::vector<double>>("ZFrontTop");
165  rMaxFront_ = args.value<std::vector<double>>("RMaxFront");
166 #ifdef EDM_ML_DEBUG
167  for (unsigned int i = 0; i < slopeB_.size(); ++i)
168  edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] Zmin " << cms::convert2mm(zFrontB_[i]) << " Rmin "
169  << cms::convert2mm(rMinFront_[i]) << " Slope " << slopeB_[i];
170  for (unsigned int i = 0; i < slopeT_.size(); ++i)
171  edm::LogVerbatim("HGCalGeom") << "Block [" << i << "] Zmin " << cms::convert2mm(zFrontT_[i]) << " Rmax "
172  << cms::convert2mm(rMaxFront_[i]) << " Slope " << slopeT_[i];
173 #endif
174 
175 #ifdef EDM_ML_DEBUG
176  edm::LogVerbatim("HGCalGeom") << "DDHGCalEEAlgo: NameSpace " << ns.name();
177 #endif
178 
179  waferType_ = std::make_unique<HGCalWaferType>(rad100to200_,
180  rad200to300_,
183  choiceType_,
184  nCutRadPar_,
185  fracAreaMin_);
186 
187  ConstructAlgo(ctxt, e);
188  }
std::vector< int > layerSense_
Log< level::Info, true > LogVerbatim
std::vector< double > rMaxFront_
double zMinBlock_
constexpr NumType convertRadToDeg(NumType radians)
Definition: angle_units.h:21
constexpr NumType convert2mm(NumType length)
Definition: DDutils.h:7
std::vector< double > slopeB_
std::vector< int > layerCenter_
dd4hep::Volume mother_
std::vector< double > rMinFront_
std::unique_ptr< HGCalWaferType > waferType_
double waferSepar_
double cosAlpha_
std::vector< int > layerType_
double fracAreaMin_
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
std::vector< std::string > materials_
double waferSize_
std::vector< std::string > names_
std::vector< double > thick_
double zMinRadPar_
ii
Definition: cuy.py:589
std::vector< double > rad200to300_
std::vector< double > zFrontT_
std::vector< int > layers_
void ConstructAlgo(cms::DDParsingContext &ctxt, xml_h e)
std::vector< double > rad100to200_
std::vector< double > layerThick_
std::vector< double > slopeT_
std::vector< double > zFrontB_
std::vector< std::string > wafers_
std::vector< int > copyNumber_

Member Function Documentation

◆ ConstructAlgo()

void HGCalEEAlgo::ConstructAlgo ( cms::DDParsingContext ctxt,
xml_h  e 
)
inline

Definition at line 190 of file DDHGCalEEAlgo.cc.

References MillePedeFileConverter_cfg::e, and dqmdumpme::k.

190  {
191 #ifdef EDM_ML_DEBUG
192  edm::LogVerbatim("HGCalGeom") << "==>> Constructing DDHGCalEEAlgo...";
193  copies_.clear();
194 #endif
195  dd4hep::Volume par;
196  ConstructLayers(par, ctxt, e);
197 #ifdef EDM_ML_DEBUG
198  edm::LogVerbatim("HGCalGeom") << "DDHGCalEEAlgo: " << copies_.size() << " different wafer copy numbers";
199  int k(0);
200  for (std::unordered_set<int>::const_iterator itr = copies_.begin(); itr != copies_.end(); ++itr, ++k) {
201  edm::LogVerbatim("HGCalGeom") << "Copy [" << k << "] : " << (*itr);
202  }
203  copies_.clear();
204  edm::LogVerbatim("HGCalGeom") << "<<== End of DDHGCalEEAlgo construction...";
205 #endif
206  }
Log< level::Info, true > LogVerbatim
void ConstructLayers(const dd4hep::Volume module, cms::DDParsingContext &ctxt, xml_h e)
std::unordered_set< int > copies_
dd4hep::Volume Volume

◆ ConstructLayers()

void HGCalEEAlgo::ConstructLayers ( const dd4hep::Volume  module,
cms::DDParsingContext ctxt,
xml_h  e 
)
inline

Definition at line 208 of file DDHGCalEEAlgo.cc.

References funct::abs(), cms::DDNamespace::addSolidNS(), cms::DDNamespace::addVolumeNS(), cms::convert2mm(), angle_units::operators::convertRadToDeg(), filterCSVwithJSON::copy, MillePedeFileConverter_cfg::e, mps_fire::i, cuy::ii, dqmdumpme::k, cms::DDNamespace::material(), g4SimHits_cfi::Material, callgraph::module, Skims_PA_cff::name, HGCalGeometryMode::Polyhedra, PixelTestBeamValidation_cfi::Position, cms::DDNamespace::prepend(), HGCalGeomTools::radius(), AlCaHLTBitMon_QueryRunRegistry::string, cond::impl::to_string(), and geometryCSVtoXML::zz.

208  {
209  static constexpr double tol1 = 0.01 * dd4hep::mm;
210  static constexpr double tol2 = 0.00001 * dd4hep::mm;
211  cms::DDNamespace ns(ctxt, e, true);
212 
213 #ifdef EDM_ML_DEBUG
214  edm::LogVerbatim("HGCalGeom") << "DDHGCalEEAlgo: \t\tInside Layers";
215 #endif
216 
217  double zi(zMinBlock_);
218  int laymin(0);
219  for (unsigned int i = 0; i < layers_.size(); i++) {
220  double zo = zi + layerThick_[i];
221  double routF = HGCalGeomTools::radius(zi, zFrontT_, rMaxFront_, slopeT_);
222  int laymax = laymin + layers_[i];
223  double zz = zi;
224  double thickTot(0);
225  for (int ly = laymin; ly < laymax; ++ly) {
226  int ii = layerType_[ly];
227  int copy = copyNumber_[ii];
228  double hthick = 0.5 * thick_[ii];
229  double rinB = HGCalGeomTools::radius(zo - tol1, zFrontB_, rMinFront_, slopeB_);
230  zz += hthick;
231  thickTot += thick_[ii];
232 
233  std::string name = ns.prepend(names_[ii]) + std::to_string(copy);
234 #ifdef EDM_ML_DEBUG
235  edm::LogVerbatim("HGCalGeom") << "DDHGCalEEAlgo: Layer " << ly << ":" << ii << " Front " << cms::convert2mm(zi)
236  << ", " << cms::convert2mm(routF) << " Back " << cms::convert2mm(zo) << ", "
237  << cms::convert2mm(rinB) << " superlayer thickness "
239 #endif
240 
241  std::string matName = materials_[ii];
242  dd4hep::Material matter = ns.material(matName);
243  dd4hep::Volume glog;
244  if (layerSense_[ly] < 1) {
245  std::vector<double> pgonZ, pgonRin, pgonRout;
246  if (layerSense_[ly] == 0 || absorbMode_ == 0) {
247  double rmax = routF * cosAlpha_ - tol1;
248  pgonZ.emplace_back(-hthick);
249  pgonZ.emplace_back(hthick);
250  pgonRin.emplace_back(rinB);
251  pgonRin.emplace_back(rinB);
252  pgonRout.emplace_back(rmax);
253  pgonRout.emplace_back(rmax);
254  } else {
255  HGCalGeomTools::radius(zz - hthick,
256  zz + hthick,
257  zFrontB_,
258  rMinFront_,
259  slopeB_,
260  zFrontT_,
261  rMaxFront_,
262  slopeT_,
263  -layerSense_[ly],
264  pgonZ,
265  pgonRin,
266  pgonRout);
267 #ifdef EDM_ML_DEBUG
268  edm::LogVerbatim("HGCalGeom") << "DDHGCalEEAlgo: z " << cms::convert2mm((zz - hthick)) << ":"
269  << cms::convert2mm((zz + hthick)) << " with " << pgonZ.size() << " palnes";
270  for (unsigned int isec = 0; isec < pgonZ.size(); ++isec)
271  edm::LogVerbatim("HGCalGeom") << "[" << isec << "] z " << cms::convert2mm(pgonZ[isec]) << " R "
272  << cms::convert2mm(pgonRin[isec]) << ":" << cms::convert2mm(pgonRout[isec]);
273 #endif
274  for (unsigned int isec = 0; isec < pgonZ.size(); ++isec) {
275  pgonZ[isec] -= zz;
276  pgonRout[isec] = pgonRout[isec] * cosAlpha_ - tol1;
277  }
278  }
279 
280  dd4hep::Solid solid = dd4hep::Polyhedra(sectors_, -alpha_, 2._pi, pgonZ, pgonRin, pgonRout);
281  ns.addSolidNS(ns.prepend(name), solid);
282  glog = dd4hep::Volume(solid.name(), solid, matter);
283  ns.addVolumeNS(glog);
284 
285 #ifdef EDM_ML_DEBUG
286  edm::LogVerbatim("HGCalGeom") << "DDHGCalEEAlgo: " << solid.name() << " polyhedra of " << sectors_
287  << " sectors covering " << convertRadToDeg(-alpha_) << ":"
288  << convertRadToDeg(-alpha_ + 2._pi) << " with " << pgonZ.size()
289  << " sections and filled with " << matName;
290 
291  for (unsigned int k = 0; k < pgonZ.size(); ++k)
292  edm::LogVerbatim("HGCalGeom") << "[" << k << "] z " << cms::convert2mm(pgonZ[k]) << " R "
293  << cms::convert2mm(pgonRin[k]) << ":" << cms::convert2mm(pgonRout[k]);
294 #endif
295  } else {
296  double rins =
297  (sensitiveMode_ < 1) ? rinB : HGCalGeomTools::radius(zz + hthick - tol1, zFrontB_, rMinFront_, slopeB_);
298  double routs =
299  (sensitiveMode_ < 1) ? routF : HGCalGeomTools::radius(zz - hthick, zFrontT_, rMaxFront_, slopeT_);
300  dd4hep::Solid solid = dd4hep::Tube(rins, routs, hthick, 0.0, 2._pi);
301  ns.addSolidNS(ns.prepend(name), solid);
302  glog = dd4hep::Volume(solid.name(), solid, matter);
303  ns.addVolumeNS(glog);
304 
305 #ifdef EDM_ML_DEBUG
306  edm::LogVerbatim("HGCalGeom") << "DDHGCalEEFileAlgo: " << solid.name() << " Tubs made of " << matter.name()
307  << " of dimensions " << cms::convert2mm(rinB) << ":" << cms::convert2mm(rins)
308  << ", " << cms::convert2mm(routF) << ":" << cms::convert2mm(routs) << ", "
309  << cms::convert2mm(hthick) << ", 0.0, 360.0 and position " << glog.name()
310  << " number " << copy << ":" << layerCenter_[copy - firstLayer_];
311 #endif
313  ctxt, e, glog, rins, routs, zz, layerSense_[ly], layerCenter_[copy - firstLayer_]); //, cpv);
314  }
315 
316  dd4hep::Position r1(0, 0, zz);
317  mother_.placeVolume(glog, copy, r1);
318  ++copyNumber_[ii];
319 
320 #ifdef EDM_ML_DEBUG
321  edm::LogVerbatim("HGCalGeom") << "DDHGCalEEAlgo: " << glog.name() << " number " << copy << " positioned in "
322  << module.name() << " at (0,0," << cms::convert2mm(zz) << ") with no rotation";
323 #endif
324  zz += hthick;
325  } // End of loop over layers in a block
326  zi = zo;
327  laymin = laymax;
328  if (std::abs(thickTot - layerThick_[i]) >= tol2) {
329  if (thickTot > layerThick_[i]) {
330  edm::LogError("HGCalGeom") << "Thickness of the partition " << cms::convert2mm(layerThick_[i])
331  << " is smaller than " << cms::convert2mm(thickTot)
332  << ": thickness of all its components **** ERROR ****";
333  } else {
334  edm::LogWarning("HGCalGeom") << "Thickness of the partition " << cms::convert2mm(layerThick_[i])
335  << " does not match with " << cms::convert2mm(thickTot) << " of the components";
336  }
337  }
338 
339  } // End of loop over layers in a block
340  }
std::vector< int > layerSense_
Log< level::Info, true > LogVerbatim
static void radius(double zf, double zb, std::vector< double > const &zFront1, std::vector< double > const &rFront1, std::vector< double > const &slope1, std::vector< double > const &zFront2, std::vector< double > const &rFront2, std::vector< double > const &slope2, int flag, std::vector< double > &zz, std::vector< double > &rin, std::vector< double > &rout)
std::vector< double > rMaxFront_
double zMinBlock_
constexpr NumType convertRadToDeg(NumType radians)
Definition: angle_units.h:21
constexpr NumType convert2mm(NumType length)
Definition: DDutils.h:7
std::vector< double > slopeB_
std::string to_string(const V &value)
Definition: OMSAccess.h:71
std::vector< int > layerCenter_
Log< level::Error, false > LogError
dd4hep::Volume mother_
std::vector< double > rMinFront_
double cosAlpha_
std::vector< int > layerType_
std::vector< std::string > materials_
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
void PositionSensitive(cms::DDParsingContext &ctxt, xml_h e, const dd4hep::Volume &glog, double rin, double rout, double zpos, int layertype, int layercenter)
std::vector< std::string > names_
dd4hep::Volume Volume
std::vector< double > thick_
ii
Definition: cuy.py:589
std::vector< double > zFrontT_
std::vector< int > layers_
std::vector< double > layerThick_
std::vector< double > slopeT_
std::vector< double > zFrontB_
Log< level::Warning, false > LogWarning
std::vector< int > copyNumber_

◆ PositionSensitive()

void HGCalEEAlgo::PositionSensitive ( cms::DDParsingContext ctxt,
xml_h  e,
const dd4hep::Volume &  glog,
double  rin,
double  rout,
double  zpos,
int  layertype,
int  layercenter 
)
inline

Definition at line 342 of file DDHGCalEEAlgo.cc.

References funct::abs(), cms::convert2mm(), filterCSVwithJSON::copy, TCMET_cfi::corner, PVValHelper::dy, MillePedeFileConverter_cfg::e, createfilelist::int, gpuVertexFinder::iv, HGCalParameters::k_CornerSize, N, EgHLTOffHistBins_cfi::nr, gpuPixelDoublets::ntot, HGCalTypes::packTypeUV(), PixelTestBeamValidation_cfi::Position, dttmaxenums::R, HGCalGeomTools::shiftXY(), mathSSE::sqrt(), findQualityFiles::v, cms::DDNamespace::volume(), and HGCalGeomTools::waferCorner().

349  {
350  cms::DDNamespace ns(ctxt, e, true);
351  static const double sqrt3 = std::sqrt(3.0);
352  double r = 0.5 * (waferSize_ + waferSepar_);
353  double R = 2.0 * r / sqrt3;
354  double dy = 0.75 * R;
355  int N = (int)(0.5 * rout / r) + 2;
356  const auto& xyoff = geomTools_.shiftXY(layercenter, (waferSize_ + waferSepar_));
357 #ifdef EDM_ML_DEBUG
358  int ium(0), ivm(0), iumAll(0), ivmAll(0), kount(0), ntot(0), nin(0);
359  std::vector<int> ntype(6, 0);
360  edm::LogVerbatim("HGCalGeom") << "DDHGCalEEAlgo: " << glog.name() << " rout " << cms::convert2mm(rout) << " N " << N
361  << " for maximum u, v; r " << cms::convert2mm(r) << " R " << cms::convert2mm(R)
362  << " dy " << cms::convert2mm(dy) << " Shift " << cms::convert2mm(xyoff.first) << ":"
363  << cms::convert2mm(xyoff.second) << " WaferSize "
365 #endif
366 
367  for (int u = -N; u <= N; ++u) {
368  for (int v = -N; v <= N; ++v) {
369  int nr = 2 * v;
370  int nc = -2 * u + v;
371  double xpos = xyoff.first + nc * r;
372  double ypos = xyoff.second + nr * dy;
373  const auto& corner = HGCalGeomTools::waferCorner(xpos, ypos, r, R, rin, rout, false);
374 #ifdef EDM_ML_DEBUG
375  int iu = std::abs(u);
376  int iv = std::abs(v);
377  ++ntot;
378  if (((corner.first <= 0) && std::abs(u) < 5 && std::abs(v) < 5) || (std::abs(u) < 2 && std::abs(v) < 2)) {
379  edm::LogVerbatim("HGCalGeom") << "DDHGCalEEAlgo: " << glog.name() << " R " << cms::convert2mm(rin) << ":"
380  << cms::convert2mm(rout) << "\n Z " << cms::convert2mm(zpos) << " LayerType "
381  << layertype << " u " << u << " v " << v << " with " << corner.first
382  << " corners";
383  }
384 #endif
385  if (corner.first > 0) {
386  int type = waferType_->getType(cms::convert2mm(xpos), cms::convert2mm(ypos), cms::convert2mm(zpos));
387  int copy = HGCalTypes::packTypeUV(type, u, v);
388 #ifdef EDM_ML_DEBUG
389  if (iu > ium)
390  ium = iu;
391  if (iv > ivm)
392  ivm = iv;
393  kount++;
394  if (copies_.count(copy) == 0)
395  copies_.insert(copy);
396 #endif
397  if (corner.first == (int)(HGCalParameters::k_CornerSize)) {
398 #ifdef EDM_ML_DEBUG
399  if (iu > iumAll)
400  iumAll = iu;
401  if (iv > ivmAll)
402  ivmAll = iv;
403  ++nin;
404 #endif
405 
406  dd4hep::Position tran(xpos, ypos, 0.0);
407  if (layertype > 1)
408  type += 3;
409  glog.placeVolume(ns.volume(wafers_[type]), copy, tran);
410 #ifdef EDM_ML_DEBUG
411  ++ntype[type];
412  edm::LogVerbatim("HGCalGeom")
413  << " DDHGCalEEAlgo: " << wafers_[type] << " number " << copy << " positioned in " << glog.name()
414  << " at (" << cms::convert2mm(xpos) << ", " << cms::convert2mm(ypos) << ",0) with no rotation";
415 #endif
416  }
417  }
418  }
419  }
420 
421 #ifdef EDM_ML_DEBUG
422  edm::LogVerbatim("HGCalGeom") << " DDHGCalEEAlgo: Maximum # of u " << ium << ":" << iumAll << " # of v " << ivm
423  << ":" << ivmAll << " and " << nin << ":" << kount << ":" << ntot << " wafers ("
424  << ntype[0] << ":" << ntype[1] << ":" << ntype[2] << ":" << ntype[3] << ":"
425  << ntype[4] << ":" << ntype[5] << ") for " << glog.name() << " R "
426  << cms::convert2mm(rin) << ":" << cms::convert2mm(rout);
427 #endif
428  }
Log< level::Info, true > LogVerbatim
int32_t *__restrict__ iv
HGCalGeomTools geomTools_
constexpr NumType convert2mm(NumType length)
Definition: DDutils.h:7
std::unique_ptr< HGCalWaferType > waferType_
std::pair< double, double > shiftXY(int waferPosition, double waferSize) const
static constexpr uint32_t k_CornerSize
double waferSepar_
static std::pair< int32_t, int32_t > waferCorner(double xpos, double ypos, double r, double R, double rMin, double rMax, bool oldBug=false)
std::unordered_set< int > copies_
T sqrt(T t)
Definition: SSEVec.h:19
double waferSize_
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
#define N
Definition: blowfish.cc:9
static int32_t packTypeUV(int type, int u, int v)
Definition: HGCalTypes.cc:3
std::vector< std::string > wafers_
__shared__ uint32_t ntot

Member Data Documentation

◆ absorbMode_

int HGCalEEAlgo::absorbMode_

Definition at line 45 of file DDHGCalEEAlgo.cc.

◆ alpha_

double HGCalEEAlgo::alpha_

Definition at line 64 of file DDHGCalEEAlgo.cc.

◆ choiceType_

int HGCalEEAlgo::choiceType_

Definition at line 51 of file DDHGCalEEAlgo.cc.

◆ copies_

std::unordered_set<int> HGCalEEAlgo::copies_

Definition at line 63 of file DDHGCalEEAlgo.cc.

◆ copyNumber_

std::vector<int> HGCalEEAlgo::copyNumber_

Definition at line 38 of file DDHGCalEEAlgo.cc.

◆ cosAlpha_

double HGCalEEAlgo::cosAlpha_

Definition at line 64 of file DDHGCalEEAlgo.cc.

◆ firstLayer_

int HGCalEEAlgo::firstLayer_

Definition at line 44 of file DDHGCalEEAlgo.cc.

◆ fracAreaMin_

double HGCalEEAlgo::fracAreaMin_

Definition at line 53 of file DDHGCalEEAlgo.cc.

◆ geomTools_

HGCalGeomTools HGCalEEAlgo::geomTools_

Definition at line 30 of file DDHGCalEEAlgo.cc.

◆ layerCenter_

std::vector<int> HGCalEEAlgo::layerCenter_

Definition at line 43 of file DDHGCalEEAlgo.cc.

◆ layers_

std::vector<int> HGCalEEAlgo::layers_

Definition at line 39 of file DDHGCalEEAlgo.cc.

◆ layerSense_

std::vector<int> HGCalEEAlgo::layerSense_

Definition at line 42 of file DDHGCalEEAlgo.cc.

◆ layerThick_

std::vector<double> HGCalEEAlgo::layerThick_

Definition at line 40 of file DDHGCalEEAlgo.cc.

◆ layerType_

std::vector<int> HGCalEEAlgo::layerType_

Definition at line 41 of file DDHGCalEEAlgo.cc.

◆ materials_

std::vector<std::string> HGCalEEAlgo::materials_

Definition at line 35 of file DDHGCalEEAlgo.cc.

◆ mother_

dd4hep::Volume HGCalEEAlgo::mother_

Definition at line 32 of file DDHGCalEEAlgo.cc.

◆ names_

std::vector<std::string> HGCalEEAlgo::names_

Definition at line 36 of file DDHGCalEEAlgo.cc.

◆ nCutRadPar_

int HGCalEEAlgo::nCutRadPar_

Definition at line 52 of file DDHGCalEEAlgo.cc.

◆ rad100to200_

std::vector<double> HGCalEEAlgo::rad100to200_

Definition at line 48 of file DDHGCalEEAlgo.cc.

◆ rad200to300_

std::vector<double> HGCalEEAlgo::rad200to300_

Definition at line 49 of file DDHGCalEEAlgo.cc.

◆ rMaxFront_

std::vector<double> HGCalEEAlgo::rMaxFront_

Definition at line 62 of file DDHGCalEEAlgo.cc.

◆ rMinFront_

std::vector<double> HGCalEEAlgo::rMinFront_

Definition at line 59 of file DDHGCalEEAlgo.cc.

◆ sectors_

int HGCalEEAlgo::sectors_

Definition at line 56 of file DDHGCalEEAlgo.cc.

◆ sensitiveMode_

int HGCalEEAlgo::sensitiveMode_

Definition at line 46 of file DDHGCalEEAlgo.cc.

◆ slopeB_

std::vector<double> HGCalEEAlgo::slopeB_

Definition at line 57 of file DDHGCalEEAlgo.cc.

◆ slopeT_

std::vector<double> HGCalEEAlgo::slopeT_

Definition at line 60 of file DDHGCalEEAlgo.cc.

◆ thick_

std::vector<double> HGCalEEAlgo::thick_

Definition at line 37 of file DDHGCalEEAlgo.cc.

◆ wafers_

std::vector<std::string> HGCalEEAlgo::wafers_

Definition at line 34 of file DDHGCalEEAlgo.cc.

◆ waferSepar_

double HGCalEEAlgo::waferSepar_

Definition at line 55 of file DDHGCalEEAlgo.cc.

◆ waferSize_

double HGCalEEAlgo::waferSize_

Definition at line 54 of file DDHGCalEEAlgo.cc.

◆ waferType_

std::unique_ptr<HGCalWaferType> HGCalEEAlgo::waferType_

Definition at line 31 of file DDHGCalEEAlgo.cc.

◆ zFrontB_

std::vector<double> HGCalEEAlgo::zFrontB_

Definition at line 58 of file DDHGCalEEAlgo.cc.

◆ zFrontT_

std::vector<double> HGCalEEAlgo::zFrontT_

Definition at line 61 of file DDHGCalEEAlgo.cc.

◆ zMinBlock_

double HGCalEEAlgo::zMinBlock_

Definition at line 47 of file DDHGCalEEAlgo.cc.

◆ zMinRadPar_

double HGCalEEAlgo::zMinRadPar_

Definition at line 50 of file DDHGCalEEAlgo.cc.