Belle II Software development
GeoPXDCreator.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9#include <pxd/geometry/GeoPXDCreator.h>
10#include <vxd/geometry/GeoCache.h>
11#include <pxd/geometry/SensorInfo.h>
12#include <pxd/simulation/SensitiveDetector.h>
13
14#include <simulation/background/BkgSensitiveDetector.h>
15
16#include <geometry/Materials.h>
17#include <geometry/CreatorFactory.h>
18#include <geometry/utilities.h>
19#include <framework/gearbox/GearDir.h>
20#include <framework/gearbox/Unit.h>
21#include <framework/logging/Logger.h>
22
23#include <cmath>
24#include <boost/format.hpp>
25
26#include <G4LogicalVolume.hh>
27#include <G4PVPlacement.hh>
28
29//Shapes
30#include <G4Box.hh>
31#include <G4Tubs.hh>
32#include <G4Polycone.hh>
33#include <G4SubtractionSolid.hh>
34#include <G4Region.hh>
35
36using namespace std;
37
38namespace Belle2 {
43 using namespace geometry;
45 namespace PXD {
46
49
51 {
52 for (SensorInfo* sensorInfo : m_SensorInfo) delete sensorInfo;
53 m_SensorInfo.clear();
54 }
55
57 {
58
59 const PXDSensorInfoPar& infoPar = dynamic_cast<const PXDSensorInfoPar&>(*sensor.getSensorInfo());
60
61 SensorInfo* info = new SensorInfo(
62 VxdID(0, 0, 0),
63 infoPar.getWidth(),
64 infoPar.getLength(),
65 infoPar.getThickness(),
66 infoPar.getUCells(),
67 infoPar.getVCells1(),
68 infoPar.getSplitLength(),
69 infoPar.getVCells2()
70 );
71 info->setDEPFETParams(
72 infoPar.getBulkDoping() / (Unit::um * Unit::um * Unit::um),
73 infoPar.getBackVoltage(),
74 infoPar.getTopVoltage(),
81 infoPar.getGateDepth(),
82 infoPar.getDoublePixel(),
83 infoPar.getChargeThreshold(),
84 infoPar.getNoiseFraction()
85 );
86 info->setIntegrationWindow(
87 infoPar.getIntegrationStart(),
88 infoPar.getIntegrationEnd()
89 );
90
91 m_SensorInfo.push_back(info);
92 return info;
93 }
94
96 {
98 VxdID(0, 0, 0),
99 sensor.getLength("width"),
100 sensor.getLength("length"),
101 sensor.getLength("height"),
102 sensor.getInt("pixelsU"),
103 sensor.getInt("pixelsV[1]"),
104 sensor.getLength("splitLength", 0),
105 sensor.getInt("pixelsV[2]", 0)
106 );
107
108 info->setDEPFETParams(
109 sensor.getDouble("BulkDoping"),
110 sensor.getWithUnit("BackVoltage"),
111 sensor.getWithUnit("TopVoltage"),
112 sensor.getLength("SourceBorderSmallPixel"),
113 sensor.getLength("ClearBorderSmallPixel"),
114 sensor.getLength("DrainBorderSmallPixel"),
115 sensor.getLength("SourceBorderLargePixel"),
116 sensor.getLength("ClearBorderLargePixel"),
117 sensor.getLength("DrainBorderLargePixel"),
118 sensor.getLength("GateDepth"),
119 sensor.getBool("DoublePixel"),
120 sensor.getDouble("ChargeThreshold"),
121 sensor.getDouble("NoiseFraction")
122 );
123 info->setIntegrationWindow(
124 sensor.getTime("IntegrationStart"),
125 sensor.getTime("IntegrationEnd")
126 );
127 return info;
128 }
129
130
132 {
133 // Create an empty payload
134 PXDGeometryPar pxdGeometryPar;
135
136 //Read prefix ('SVD' or 'PXD')
137 pxdGeometryPar.setPrefix(m_prefix);
138
139 //Read some global parameters
140 VXDGlobalPar globals((float)content.getDouble("ElectronTolerance", 100),
141 (float)content.getDouble("MinimumElectrons", 10),
142 content.getLength("ActiveStepSize", 0.005),
143 content.getBool("ActiveChips", false),
144 content.getBool("SeeNeutrons", false),
145 content.getBool("OnlyPrimaryTrueHits", false),
146 content.getBool("OnlyActiveMaterial", false),
147 (float)content.getLength("DistanceTolerance", 0.005),
148 content.getString("DefaultMaterial", "Air")
149 );
150 pxdGeometryPar.setGlobalParams(globals);
151
152 //Read envelope parameters
153 GearDir envelopeParams(content, "Envelope/");
154 VXDEnvelopePar envelope(envelopeParams.getString("Name", ""),
155 envelopeParams.getString("Material", "Air"),
156 envelopeParams.getString("Color", ""),
157 envelopeParams.getAngle("minPhi", 0),
158 envelopeParams.getAngle("maxPhi", 2 * M_PI),
159 (envelopeParams.getNodes("InnerPoints/point").size() > 0)
160 );
161
162 for (const GearDir& point : envelopeParams.getNodes("InnerPoints/point")) {
163 pair<double, double> ZXPoint(point.getLength("z"), point.getLength("x"));
164 envelope.getInnerPoints().push_back(ZXPoint);
165 }
166 for (const GearDir& point : envelopeParams.getNodes("OuterPoints/point")) {
167 pair<double, double> ZXPoint(point.getLength("z"), point.getLength("x"));
168 envelope.getOuterPoints().push_back(ZXPoint);
169 }
170 pxdGeometryPar.setEnvelope(envelope);
171
172 // Read alignment for detector m_prefix ('PXD' or 'SVD')
173 string pathAlign = (boost::format("Align[@component='%1%']/") % m_prefix).str();
174 GearDir paramsAlign(GearDir(content, "Alignment/"), pathAlign);
175 if (!paramsAlign) {
176 B2WARNING("Could not find alignment parameters for component " << m_prefix);
177 return pxdGeometryPar;
178 }
179 pxdGeometryPar.getAlignmentMap()[m_prefix] = VXDAlignmentPar(paramsAlign.getLength("du"),
180 paramsAlign.getLength("dv"),
181 paramsAlign.getLength("dw"),
182 paramsAlign.getAngle("alpha"),
183 paramsAlign.getAngle("beta"),
184 paramsAlign.getAngle("gamma")
185 );
186
187 //Read the definition of all sensor types
188 GearDir components(content, "Components/");
189 for (const GearDir& paramsSensor : components.getNodes("Sensor")) {
190 string sensorTypeID = paramsSensor.getString("@type");
191
192 VXDGeoSensorPar sensor(paramsSensor.getString("Material"),
193 paramsSensor.getString("Color", ""),
194 paramsSensor.getLength("width"),
195 paramsSensor.getLength("width2", 0),
196 paramsSensor.getLength("length"),
197 paramsSensor.getLength("height"),
198 paramsSensor.getAngle("angle", 0),
199 paramsSensor.getBool("@slanted", false)
200 );
201 sensor.setActive(VXDGeoComponentPar(
202 paramsSensor.getString("Material"),
203 paramsSensor.getString("Active/Color", "#f00"),
204 paramsSensor.getLength("Active/width"),
205 paramsSensor.getLength("Active/width2", 0),
206 paramsSensor.getLength("Active/length"),
207 paramsSensor.getLength("Active/height")
209 "Active",
210 paramsSensor.getLength("Active/u"),
211 paramsSensor.getLength("Active/v"),
212 paramsSensor.getString("Active/w", "center"),
213 paramsSensor.getLength("Active/woffset", 0)
214 ));
215
216 PXDSensorInfoPar* pxdInfo = readSensorInfo(GearDir(paramsSensor, "Active"));
217 sensor.setSensorInfo(pxdInfo);
218 sensor.setComponents(getSubComponents(paramsSensor));
219 pxdGeometryPar.getSensorMap()[sensorTypeID] = sensor;
220 pxdGeometryPar.getSensorInfos().push_back(pxdInfo);
221 }
222
223 //Build all ladders including Sensors
224 GearDir support(content, "Support/");
225 readHalfShellSupport(support, pxdGeometryPar);
226
227 for (const GearDir& shell : content.getNodes("HalfShell")) {
228
229 string shellName = m_prefix + "." + shell.getString("@name");
230 string pathShell = (boost::format("Align[@component='%1%']/") % shellName).str();
231 GearDir paramsShell(GearDir(content, "Alignment/"), pathShell);
232 if (!paramsShell) {
233 B2WARNING("Could not find alignment parameters for component " << shellName);
234 return pxdGeometryPar;
235 }
236 pxdGeometryPar.getAlignmentMap()[shellName] = VXDAlignmentPar(paramsShell.getLength("du"),
237 paramsShell.getLength("dv"),
238 paramsShell.getLength("dw"),
239 paramsShell.getAngle("alpha"),
240 paramsShell.getAngle("beta"),
241 paramsShell.getAngle("gamma")
242 );
243
244 VXDHalfShellPar halfShell(shell.getString("@name"), shell.getAngle("shellAngle", 0));
245
246 for (const GearDir& layer : shell.getNodes("Layer")) {
247 int layerID = layer.getInt("@id");
248
249 readLadder(layerID, components, pxdGeometryPar);
250
251 //Loop over defined ladders
252 for (const GearDir& ladder : layer.getNodes("Ladder")) {
253 int ladderID = ladder.getInt("@id");
254 double phi = ladder.getAngle("phi", 0);
255 readLadderComponents(layerID, ladderID, content, pxdGeometryPar);
256 halfShell.addLadder(layerID, ladderID, phi);
257 }
258 }
259 pxdGeometryPar.getHalfShells().push_back(halfShell);
260 }
261
262 //Create diamond radiation sensors if defined and in background mode
263 GearDir radiationDir(content, "RadiationSensors");
264 if (pxdGeometryPar.getGlobalParams().getActiveChips() && radiationDir) {
265 VXDGeoRadiationSensorsPar radiationSensors(
266 m_prefix,
267 radiationDir.getBool("insideEnvelope"),
268 radiationDir.getLength("width"),
269 radiationDir.getLength("length"),
270 radiationDir.getLength("height"),
271 radiationDir.getString("material")
272 );
273
274 //Add radiation sensor positions
275 for (GearDir& position : radiationDir.getNodes("position")) {
276 VXDGeoRadiationSensorsPositionPar diamonds(position.getLength("z"),
277 position.getLength("radius"),
278 position.getAngle("theta")
279 );
280
281 //Loop over all phi positions
282 for (GearDir& sensor : position.getNodes("phi")) {
283 //Add sensor with angle and id
284 diamonds.addSensor(sensor.getInt("@id"), sensor.getAngle());
285 }
286 radiationSensors.addPosition(diamonds);
287 }
288 pxdGeometryPar.setRadiationSensors(radiationSensors);
289 }
290
291 return pxdGeometryPar;
292 }
293
294 void GeoPXDCreator::readHalfShellSupport(const GearDir& support, PXDGeometryPar& pxdGeometryPar)
295 {
296 for (const GearDir& endflange : support.getNodes("Endflange")) {
297 VXDPolyConePar endflangePar(
298 endflange.getString("@name"),
299 endflange.getString("Material", "Air"),
300 endflange.getAngle("minPhi", 0),
301 endflange.getAngle("maxPhi", 2 * M_PI),
302 (endflange.getNodes("Cutout").size() > 0),
303 endflange.getLength("Cutout/width1", 0.),
304 endflange.getLength("Cutout/width2", 0.),
305 endflange.getLength("Cutout/height", 0.),
306 endflange.getLength("Cutout/depth", 0.)
307 );
308
309 for (const GearDir& plane : endflange.getNodes("Plane")) {
310 VXDPolyConePlanePar planePar(
311 plane.getLength("posZ"),
312 plane.getLength("innerRadius"),
313 plane.getLength("outerRadius")
314 );
315 endflangePar.getPlanes().push_back(planePar);
316 }
317 pxdGeometryPar.getEndflanges().push_back(endflangePar);
318 }
319
320 // Cout outs for endflanges
321 pxdGeometryPar.setNCutOuts(support.getInt("Cutout/count"));
322 pxdGeometryPar.setCutOutWidth(support.getLength("Cutout/width"));
323 pxdGeometryPar.setCutOutHeight(support.getLength("Cutout/height"));
324 pxdGeometryPar.setCutOutShift(support.getLength("Cutout/shift"));
325 pxdGeometryPar.setCutOutRPhi(support.getLength("Cutout/rphi"));
326 pxdGeometryPar.setCutOutStartPhi(support.getAngle("Cutout/startPhi"));
327 pxdGeometryPar.setCutOutDeltaPhi(support.getAngle("Cutout/deltaPhi"));
328
329 //Create Carbon cooling tubes
330 pxdGeometryPar.setNTubes(support.getInt("CarbonTubes/count"));
331 pxdGeometryPar.setTubesMinZ(support.getLength("CarbonTubes/minZ"));
332 pxdGeometryPar.setTubesMaxZ(support.getLength("CarbonTubes/maxZ"));
333 pxdGeometryPar.setTubesMinR(support.getLength("CarbonTubes/innerRadius"));
334 pxdGeometryPar.setTubesMaxR(support.getLength("CarbonTubes/outerRadius"));
335 pxdGeometryPar.setTubesRPhi(support.getLength("CarbonTubes/rphi"));
336 pxdGeometryPar.setTubesStartPhi(support.getAngle("CarbonTubes/startPhi"));
337 pxdGeometryPar.setTubesDeltaPhi(support.getAngle("CarbonTubes/deltaPhi"));
338 pxdGeometryPar.setTubesMaterial(support.getString("CarbonTubes/Material", "Carbon"));
339
340 return;
341 }
342
343 void GeoPXDCreator::createGeometry(const PXDGeometryPar& parameters, G4LogicalVolume& topVolume, geometry::GeometryTypes)
344 {
345
346 m_activeStepSize = parameters.getGlobalParams().getActiveStepSize() / Unit::mm;
347 m_activeChips = parameters.getGlobalParams().getActiveChips();
348 m_seeNeutrons = parameters.getGlobalParams().getSeeNeutrons();
349 m_onlyPrimaryTrueHits = parameters.getGlobalParams().getOnlyPrimaryTrueHits();
350 m_distanceTolerance = parameters.getGlobalParams().getDistanceTolerance();
351 m_electronTolerance = parameters.getGlobalParams().getElectronTolerance();
352 m_minimumElectrons = parameters.getGlobalParams().getMinimumElectrons();
353 m_onlyActiveMaterial = parameters.getGlobalParams().getOnlyActiveMaterial();
354 m_defaultMaterial = parameters.getGlobalParams().getDefaultMaterial();
355
356 G4Material* material = Materials::get(m_defaultMaterial);
357 if (!material) B2FATAL("Default Material of VXD, '" << m_defaultMaterial << "', could not be found");
358
359
360 //Build envelope
361 G4LogicalVolume* envelope(0);
362 G4VPhysicalVolume* physEnvelope{nullptr};
363 if (!parameters.getEnvelope().getExists()) {
364 B2INFO("Could not find definition for " + m_prefix + " Envelope, placing directly in top volume");
365 envelope = &topVolume;
366 } else {
367 double minZ(0), maxZ(0);
368 G4Polycone* envelopeCone = geometry::createRotationSolid("Envelope",
369 parameters.getEnvelope().getInnerPoints(),
370 parameters.getEnvelope().getOuterPoints(),
371 parameters.getEnvelope().getMinPhi(),
372 parameters.getEnvelope().getMaxPhi(),
373 minZ, maxZ
374 );
375 envelope = new G4LogicalVolume(envelopeCone, material, m_prefix + ".Envelope");
376 setVisibility(*envelope, false);
377 physEnvelope = new G4PVPlacement(getAlignment(parameters.getAlignment(m_prefix)), envelope, m_prefix + ".Envelope",
378 &topVolume, false, 1);
379
380 // Set up region for production cuts
381 G4Region* aRegion = new G4Region("PXDEnvelope");
382 envelope->SetRegion(aRegion);
383 aRegion->AddRootLogicalVolume(envelope);
384 }
385
386 //Read the definition of all sensor types
387 for (const pair<const string, VXDGeoSensorPar>& typeAndSensor : parameters.getSensorMap()) {
388 const string& sensorTypeID = typeAndSensor.first;
389 const VXDGeoSensorPar& paramsSensor = typeAndSensor.second;
390 VXDGeoSensor sensor(
391 paramsSensor.getMaterial(),
392 paramsSensor.getColor(),
393 paramsSensor.getWidth() / Unit::mm,
394 paramsSensor.getWidth2() / Unit::mm,
395 paramsSensor.getLength() / Unit::mm,
396 paramsSensor.getHeight() / Unit::mm,
397 paramsSensor.getSlanted()
398 );
399 sensor.setActive(VXDGeoComponent(
400 paramsSensor.getMaterial(),
401 paramsSensor.getActiveArea().getColor(),
402 paramsSensor.getActiveArea().getWidth() / Unit::mm,
403 paramsSensor.getActiveArea().getWidth2() / Unit::mm,
404 paramsSensor.getActiveArea().getLength() / Unit::mm,
405 paramsSensor.getActiveArea().getHeight() / Unit::mm
407 "Active",
408 paramsSensor.getActivePlacement().getU() / Unit::mm,
409 paramsSensor.getActivePlacement().getV() / Unit::mm,
410 paramsSensor.getActivePlacement().getW(),
411 paramsSensor.getActivePlacement().getWOffset() / Unit::mm
412 ));
413 sensor.setSensorInfo(createSensorInfo(paramsSensor));
414
415 vector<VXDGeoPlacement> subcomponents;
416 const auto& components = paramsSensor.getComponents();
417 subcomponents.reserve(components.size());
418 std::transform(components.begin(), components.end(), std::back_inserter(subcomponents),
419 [](auto const & component) {
420 return VXDGeoPlacement(component.getName(),
421 component.getU() / Unit::mm,
422 component.getV() / Unit::mm,
423 component.getW(),
424 component.getWOffset() / Unit::mm
425 );
426 });
427 sensor.setComponents(subcomponents);
428 m_sensorMap[sensorTypeID] = sensor;
429 }
430
431 //Read the component cache from DB
432 for (const string& name : parameters.getComponentInsertOder()) {
433 if (m_componentCache.find(name) != m_componentCache.end()) {
434 // already created due to being a sub component of a previous
435 // component. Seems fishy since the information of this component
436 // is in the db at least twice so we could run into
437 // inconsistencies.
438 B2WARNING("Component " << name << " already created from previous subcomponents, should not be here");
439 continue;
440 }
441 const VXDGeoComponentPar& paramsComponent = parameters.getComponent(name);
443 paramsComponent.getMaterial(),
444 paramsComponent.getColor(),
445 paramsComponent.getWidth() / Unit::mm,
446 paramsComponent.getWidth2() / Unit::mm,
447 paramsComponent.getLength() / Unit::mm,
448 paramsComponent.getHeight() / Unit::mm
449 );
450 double angle = paramsComponent.getAngle();
451
452
453 if (c.getWidth() <= 0 || c.getLength() <= 0 || c.getHeight() <= 0) {
454 B2DEBUG(100, "One dimension empty, using auto resize for component");
455 } else {
456 G4VSolid* solid = createTrapezoidal(m_prefix + "." + name, c.getWidth(), c.getWidth2(), c.getLength(), c.getHeight(), angle);
457 c.setVolume(new G4LogicalVolume(solid, Materials::get(c.getMaterial()), m_prefix + "." + name));
458 }
459
460 vector<VXDGeoPlacement> subComponents;
461 const auto& paramsSubComponents = paramsComponent.getSubComponents();
462 subComponents.reserve(paramsSubComponents.size());
463 std::transform(paramsSubComponents.begin(), paramsSubComponents.end(), std::back_inserter(subComponents),
464 [](auto const & paramsSubComponent) {
465 return VXDGeoPlacement(paramsSubComponent.getName(),
466 paramsSubComponent.getU() / Unit::mm,
467 paramsSubComponent.getV() / Unit::mm,
468 paramsSubComponent.getW(),
469 paramsSubComponent.getWOffset() / Unit::mm
470 );
471 });
472 createSubComponents(m_prefix + "." + name, c, subComponents);
473 if (m_activeChips && parameters.getSensitiveChipID(name) >= 0) {
474 int chipID = parameters.getSensitiveChipID(name);
475 B2DEBUG(50, "Creating BkgSensitiveDetector for component " << name << " with chipID " << chipID);
476 BkgSensitiveDetector* sensitive = new BkgSensitiveDetector(m_prefix.c_str(), chipID);
477 c.getVolume()->SetSensitiveDetector(sensitive);
478 m_sensitive.push_back(sensitive);
479 }
480
481 m_componentCache[name] = c;
482 }
483
484 //Build all ladders including Sensors
485 VXD::GeoVXDAssembly shellSupport = createHalfShellSupport(parameters);
486
487 //const std::vector<VXDHalfShellPar>& HalfShells = parameters.getHalfShells();
488 for (const VXDHalfShellPar& shell : parameters.getHalfShells()) {
489 string shellName = shell.getName();
490 m_currentHalfShell = m_prefix + "." + shellName;
491 G4Transform3D shellAlignment = getAlignment(parameters.getAlignment(m_currentHalfShell));
492
493 // Remember shell coordinate system (into which ladders are inserted)
495
496 //Place shell support
497 double shellAngle = shell.getShellAngle(); // Only used to move support, not active volumes!
498 if (!m_onlyActiveMaterial) shellSupport.place(envelope, shellAlignment * G4RotateZ3D(shellAngle));
499
500 //const std::map< int, std::vector<std::pair<int, double>> >& Layers = shell.getLayers();
501 for (const std::pair<const int, std::vector<std::pair<int, double>> >& layer : shell.getLayers()) {
502 int layerID = layer.first;
503 const std::vector<std::pair<int, double>>& Ladders = layer.second;
504
505
506 setCurrentLayer(layerID, parameters);
507
508 //Place Layer support
510 if (!m_onlyActiveMaterial) layerSupport.place(envelope, shellAlignment * G4RotateZ3D(shellAngle));
512
513 //Loop over defined ladders
514 for (const std::pair<int, double>& ladder : Ladders) {
515 int ladderID = ladder.first;
516 double phi = ladder.second;
517
518 G4Transform3D ladderPlacement = placeLadder(ladderID, phi, envelope, shellAlignment, parameters);
519 if (!m_onlyActiveMaterial) ladderSupport.place(envelope, ladderPlacement);
520 }
521
522 }
523 }
524
525 //Now build cache with all transformations
526 if (physEnvelope) {
528 } else {
529 //create a temporary placement of the top volume.
530 G4PVPlacement topPlacement(nullptr, G4ThreeVector(0, 0, 0), &topVolume,
531 "temp_Top", nullptr, false, 1, false);
532 //and search for all VXD sensitive sensors within
534 }
535
536 //Create diamond radiation sensors if defined and in background mode
537 if (m_activeChips) {
538 if (parameters.getRadiationSensors().getSubDetector() == "") {
539 B2DEBUG(10, "Apparently no radiation sensors defined, skipping");
540 } else {
541 createDiamonds(parameters.getRadiationSensors(), topVolume, *envelope);
542 }
543 }
544 }
545
547 const VXDGeoSensorPlacement& placement)
548 {
549 SensorInfo* sensorInfo = new SensorInfo(dynamic_cast<const SensorInfo&>(*sensor.getSensorInfo()));
550 sensorInfo->setID(sensorID);
551 if (placement.getFlipV()) sensorInfo->flipVSegmentation();
552 SensitiveDetector* sensitive = new SensitiveDetector(sensorInfo);
553 return sensitive;
554 }
555
557
559
561 {
562 VXD::GeoVXDAssembly supportAssembly;
563
564 if (!parameters.getBuildSupport()) return supportAssembly;
565
566
567 // Create the Endlanges
568 const std::vector<VXDPolyConePar> Endflanges = parameters.getEndflanges();
569 for (const VXDPolyConePar& endflange : Endflanges) {
570
571 double minZ(0), maxZ(0);
572 string name = endflange.getName();
573
574 // Create a polycone
575 double minPhi = endflange.getMinPhi();
576 double dPhi = endflange.getMaxPhi() - minPhi;
577 int nPlanes = endflange.getPlanes().size();
578 if (nPlanes < 2) {
579 B2ERROR("Polycone needs at least two planes");
580 return supportAssembly;
581 }
582 std::vector<double> z(nPlanes, 0);
583 std::vector<double> rMin(nPlanes, 0);
584 std::vector<double> rMax(nPlanes, 0);
585 int index(0);
586 minZ = numeric_limits<double>::infinity();
587 maxZ = -numeric_limits<double>::infinity();
588
589 const std::vector<VXDPolyConePlanePar> Planes = endflange.getPlanes();
590 for (const VXDPolyConePlanePar& plane : Planes) {
591 z[index] = plane.getPosZ() / Unit::mm;
592 minZ = min(minZ, z[index]);
593 maxZ = max(maxZ, z[index]);
594 rMin[index] = plane.getInnerRadius() / Unit::mm;
595 rMax[index] = plane.getOuterRadius() / Unit::mm;
596 ++index;
597 }
598
599 G4VSolid* supportCone = new G4Polycone(name, minPhi, dPhi, nPlanes, z.data(), rMin.data(), rMax.data());
600
601
602 //Cutout boxes to make place for modules
603
604 //We have the z dimensions of the polycon. Let's
605 //add 1mm on each side to make sure we don't run into problems when the
606 //surfaces match
607 minZ -= 1. / Unit::mm;
608 maxZ += 1. / Unit::mm;
609
610
611 //Now get the number of cutouts and their size/position/angle
612 int nCutouts = parameters.getNCutOuts();
613 double sizeX = parameters.getCutOutWidth() / Unit::mm / 2.;
614 double sizeY = parameters.getCutOutHeight() / Unit::mm / 2.;
615 double sizeZ = (maxZ - minZ) / 2.;
616 G4ThreeVector origin(
617 parameters.getCutOutShift() / Unit::mm,
618 parameters.getCutOutRPhi() / Unit::mm,
619 minZ + sizeZ
620 );
621
622 double phi0 = parameters.getCutOutStartPhi();
623 double dphi = parameters.getCutOutDeltaPhi();
624 for (int i = 0; i < nCutouts; ++i) {
625 G4Box* box = new G4Box("Cutout", sizeX, sizeY, sizeZ);
626 G4Transform3D placement = G4RotateZ3D(phi0 + i * dphi) * G4Translate3D(origin);
627 G4VSolid* supportConeOld = supportCone;
628 supportCone = new G4SubtractionSolid("PXD Support endflange", supportConeOld, box, placement);
629 }
630
631
632 string materialName = endflange.getMaterial();
633 G4Material* material = geometry::Materials::get(materialName);
634 if (!material) B2FATAL("Material '" << materialName << "', required by PXD component " << name << ", could not be found");
635
636 G4LogicalVolume* volume = new G4LogicalVolume(supportCone, material, name);
637 geometry::setColor(*volume, "#ccc4");
638 supportAssembly.add(volume);
639
640 }
641
642
643 //Create Carbon cooling tubes
644 {
645 int nTubes = parameters.getNTubes();
646 double minZ = parameters.getTubesMinZ() / Unit::mm;
647 double maxZ = parameters.getTubesMaxZ() / Unit::mm;
648 double minR = parameters.getTubesMinR() / Unit::mm;
649 double maxR = parameters.getTubesMaxR() / Unit::mm;
650 double sizeZ = (maxZ - minZ) / 2.;
651 double shiftX = parameters.getTubesRPhi() / Unit::mm;
652 double shiftY = 0;
653 double shiftZ = minZ + sizeZ;
654 double phi0 = parameters.getTubesStartPhi();
655 double dphi = parameters.getTubesDeltaPhi();
656 string material = parameters.getTubesMaterial();
657
658 G4Tubs* tube = new G4Tubs("CarbonTube", minR, maxR, sizeZ, 0, 2 * M_PI);
659 G4LogicalVolume* tubeVol = new G4LogicalVolume(tube, geometry::Materials::get(material), "CarbonTube");
660 geometry::setColor(*tubeVol, "#000");
661 for (int i = 0; i < nTubes; ++i) {
662 G4Transform3D placement = G4RotateZ3D(phi0 + i * dphi) * G4Translate3D(shiftX, shiftY, shiftZ);
663 supportAssembly.add(tubeVol, placement);
664 }
665 }
666
667 return supportAssembly;
668 }
669
670
671 }
673}
The Class for BeamBackground Sensitive Detector.
GearDir is the basic class used for accessing the parameter store.
Definition: GearDir.h:31
virtual std::string getString(const std::string &path="") const noexcept(false) override
Get the parameter path as a string.
Definition: GearDir.h:69
The Class for VXD geometry.
void setNTubes(int nTubes)
set number of tubes
void setTubesRPhi(double rphi)
set tubes tubes RPhi
void setTubesMaxR(double maxR)
set tubes maxR
const std::vector< VXDPolyConePar > & getEndflanges() const
get Endflanges
void setNCutOuts(int nCutouts)
set number of cutouts
void setCutOutShift(double shift)
set shift of cutouts
std::vector< PXDSensorInfoPar * > & getSensorInfos()
get sensorInfos
void setCutOutWidth(double width)
set width of cutouts
void setTubesMinZ(double minZ)
set tubes minZ
void setCutOutHeight(double height)
set height of cutouts
void setCutOutDeltaPhi(double delta)
set deltaphi of cutouts
void setTubesDeltaPhi(double delta)
set tubes tubes DeltaPhi
void setTubesMaxZ(double maxZ)
set tubes maxZ
void setTubesMinR(double minR)
set tubes minR
void setCutOutStartPhi(double start)
set start phi of cutouts
void setCutOutRPhi(double rphi)
set rphi of cutouts
void setTubesStartPhi(double start)
set tubes tubes StartPhi
void setTubesMaterial(const std::string &material)
set tubes tubes material
The Class for VXD geometry.
double getTopVoltage() const
Return the voltage at the top of the sensor.
double getClearBorderSmallPitch() const
Return the distance along u between the clear side of a small pixel and the start of the internal gat...
double getChargeThreshold() const
Get the charge threshold in ADU for the sensor.
double getGateDepth() const
Return the gate depth for the sensor.
double getDrainBorderLargePitch() const
Return the distance along v between the drain side of a large pixel and the start of the internal gat...
double getClearBorderLargePitch() const
Return the distance along u between the clear side of a large pixel and the start of the internal gat...
double getNoiseFraction() const
Get the noise fraction for the sensor.
double getBulkDoping() const
Return the bulk doping of the Silicon sensor.
bool getDoublePixel() const
Return true if the Sensor is a double pixel structure: every other pixel is mirrored along v.
double getSourceBorderSmallPitch() const
Return the distance along v between the source side of a small pixel and the start of the internal ga...
double getSourceBorderLargePitch() const
Return the distance along v between the source side of a large pixel and the start of the internal ga...
double getIntegrationEnd() const
Return the end of the integration window, the timeframe the PXD is sensitive.
double getBackVoltage() const
Return the voltage at the backside of the sensor.
double getIntegrationStart() const
Return the start of the integration window, the timeframe the PXD is sensitive.
double getDrainBorderSmallPitch() const
Return the distance along v between the drain side of a small pixel and the start of the internal gat...
virtual VXD::GeoVXDAssembly createLayerSupport()
Create support structure for a PXD Layer.
void createGeometry(const PXDGeometryPar &parameters, G4LogicalVolume &topVolume, geometry::GeometryTypes type)
Create the geometry from a parameter object.
virtual VXD::SensitiveDetectorBase * createSensitiveDetector(VxdID sensorID, const VXDGeoSensor &sensor, const VXDGeoSensorPlacement &placement) override
Return a SensitiveDetector implementation for a given sensor.
std::vector< SensorInfo * > m_SensorInfo
Vector of points to SensorInfo objects.
void readHalfShellSupport(const GearDir &support, PXDGeometryPar &pxdGeometryPar)
Create support structure for VXD Half Shell, that means everything that does not depend on layer or s...
PXDGeometryPar createConfiguration(const GearDir &param)
Create a parameter object from the Gearbox XML parameters.
PXDSensorInfoPar * readSensorInfo(const GearDir &sensor)
Read the sensor definitions from the gearbox.
virtual ~GeoPXDCreator()
The destructor of the GeoPXDCreator class.
virtual VXD::GeoVXDAssembly createLadderSupport()
Create support structure for a PXD Ladder.
virtual VXD::GeoVXDAssembly createHalfShellSupport(const PXDGeometryPar &parameters)
Create support structure for PXD Half Shell, that means everything that does not depend on layer or s...
virtual VXD::SensorInfoBase * createSensorInfo(const VXDGeoSensorPar &sensor) override
Read the sensor definitions from the database.
Specific implementation of SensorInfo for PXD Sensors which provides additional pixel specific inform...
Definition: SensorInfo.h:23
void flipVSegmentation()
Flip the Pitch segmentation along v.
Definition: SensorInfo.h:43
void setID(VxdID id)
Change the SensorID, useful to copy the SensorInfo from one sensor and use it for another.
Definition: SensorInfo.h:37
static const double mm
[millimeters]
Definition: Unit.h:70
static const double um
[micrometers]
Definition: Unit.h:71
The Class for VXD Alignment payload.
The Class for VXD Envelope parameters.
The Class for VXD geometry component.
double getWidth() const
get the width of the component
double getHeight() const
get the height of the component
const std::vector< VXDGeoPlacementPar > & getSubComponents() const
get sub components
const std::string & getColor() const
get the name of the color for the component
double getAngle() const
get the angle of the component
double getWidth2() const
get the forward width of the component, 0 for rectangular
const std::string & getMaterial() const
get the name of the Material for the component
double getLength() const
get the length of the component
Class holding all parameters for an VXD geometry component.
The Class for VXD placement payload.
double getU() const
get local u coordinate where to place the component
double getV() const
get local v coordinate where to place the component
const std::string & getW() const
get local w position where to place the component
double getWOffset() const
get offset to local w position where to place the component
Class holding all parameters to place a VXD geometry subcomponent.
The Class for VXD Radiation Sensor parameters.
void addPosition(const VXDGeoRadiationSensorsPositionPar &position)
add radiation sensor position
The Class for VXD Radiation Sensor Position parameters.
void addSensor(int id, double phi)
add sensor with individual id
The Class for VXD Sensor payload.
const VXDGeoPlacementPar & getActivePlacement() const
get the placement description for the active area
const std::vector< VXDGeoPlacementPar > & getComponents() const
get the list of sub components
const VXDGeoComponentPar & getActiveArea() const
get the component description for the active area
bool getSlanted() const
return wether or not the sensor is slanted (usually only the first sensor in layers 4-6)
Struct holding the information where a sensor should be placed inside the ladder.
bool getFlipV() const
check whether or not the sensor should be flipped around the V coordinate
Struct holding all parameters for a completeVXD Sensor.
const std::vector< VXDHalfShellPar > & getHalfShells() const
get half-shell
void setGlobalParams(const VXDGlobalPar &globals)
set global parameters
void setRadiationSensors(const VXDGeoRadiationSensorsPar &diamonds)
set radiation sensor parameters
const std::map< std::string, VXDGeoSensorPar > & getSensorMap() const
get sensor map
void setEnvelope(const VXDEnvelopePar &envelope)
set envelope parameters
const VXDGlobalPar & getGlobalParams() const
get global parameters
void setPrefix(const std::string &prefix)
set prefix
std::map< std::string, VXDAlignmentPar > & getAlignmentMap()
get alignmant map
The Class for VXD global paramter payload.
Definition: VXDGlobalPar.h:24
bool getActiveChips() const
Get whether chips are sensitive
Definition: VXDGlobalPar.h:48
The Class for VXD half shell payload.
void addLadder(int layerID, int ladderID, double phi)
add ladder
The Class for VXD PolyCone, possibly with coutouts.
std::vector< VXDPolyConePlanePar > & getPlanes(void)
Get planes.
The Class for VXD Polycone Plane.
const std::list< std::pair< double, double > > & getInnerPoints() const
get inner XZ points
const std::list< std::pair< double, double > > & getOuterPoints() const
get outer XZ points
double getWidth() const
Return the (backward) width of the sensor.
int getUCells() const
Return number of pixel/strips in u direction.
double getThickness() const
Return the thickness of the sensor.
int getVCells2() const
Return number of pixel/strips in v direction for second segment.
double getSplitLength() const
Return the splitLength of the sensor.
int getVCells1() const
Return number of pixel/strips in v direction for first segment.
double getLength() const
Return the length of the sensor.
void findVolumes(G4VPhysicalVolume *envelope)
Search a given Geometry for Sensors.
Definition: GeoCache.cc:78
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition: GeoCache.cc:214
void addHalfShellPlacement(VxdID halfShell, const G4Transform3D &placement)
Remember how half-shell is placed into world volume.
Definition: GeoCache.cc:232
Class to group some Geant4 volumes and place them all at once with a given transformation matrix.
void place(G4LogicalVolume *mother, const G4Transform3D &transform)
Place all the volumes already added to the assembly in the given mother.
void add(G4LogicalVolume *volume, const G4Transform3D &transform=G4Transform3D())
Add a volume to the assembly.
std::map< std::string, VXDGeoSensor > m_sensorMap
Map containing Information about all defined sensor types.
float m_minimumElectrons
minimum number of electrons to be deposited by a particle to be saved
G4VSolid * createTrapezoidal(const std::string &name, double width, double width2, double length, double &height, double angle=0)
Create a trapezoidal solid.
virtual void readLadderComponents(int layerID, int ladderID, GearDir content, VXDGeometryPar &vxdGeometryPar)
Read parameters for ladder components and their alignment corresponding to the given ladder id.
void createDiamonds(const VXDGeoRadiationSensorsPar &params, G4LogicalVolume &topVolume, G4LogicalVolume &envelopeVolume)
Create diamond radiation sensors.
std::vector< Simulation::SensitiveDetectorBase * > m_sensitive
List to all created sensitive detector instances.
std::string m_prefix
Prefix to prepend to all volume names.
bool m_onlyActiveMaterial
If this is true, only active Materials will be placed for tracking studies.
double m_activeStepSize
Stepsize to be used inside active volumes.
float m_distanceTolerance
tolerance for Geant4 steps to be merged to a single step
virtual void readLadder(int layer, GearDir components, VXDGeometryPar &geoparameters)
Read parameters for a ladder in layer with given ID from gearbox and layer store them in payload.
std::map< std::string, Belle2::VxdID > m_halfShellVxdIDs
Used for translation of half-shell name into a VxdID to consitently handle it in hierarchy.
bool m_onlyPrimaryTrueHits
If true only create TrueHits from primary particles and ignore secondaries.
float m_electronTolerance
tolerance for the energy deposition in electrons to be merged in a single step
G4Transform3D getAlignment(const VXDAlignmentPar &params)
Get Alignment from paylead object.
GeoVXDAssembly createSubComponents(const std::string &name, VXDGeoComponent &component, std::vector< VXDGeoPlacement > placements, bool originCenter=true, bool allowOutside=false)
Place a list of subcomponents into an component.
std::map< std::string, VXDGeoComponent > m_componentCache
Cache of all previously created components.
G4Transform3D placeLadder(int ladderID, double phi, G4LogicalVolume *volume, const G4Transform3D &placement, const VXDGeometryPar &parameters)
Place ladder corresponding to the given ladder id into volume setLayer has to be called first to set ...
virtual void setCurrentLayer(int layer, const VXDGeometryPar &parameters)
Read parameters for given layer and store in m_ladder.
std::vector< VXDGeoPlacementPar > getSubComponents(const GearDir &path)
Return vector of VXDGeoPlacements with all the components defined inside a given path.
bool m_seeNeutrons
Make sensitive detectors also see neutrons.
std::string m_defaultMaterial
Name of the Material to be used for Air.
bool m_activeChips
Make also chips sensitive.
std::string m_currentHalfShell
Current half-shell being processed (need to know ladder parent for hierarchy)
Base class for Sensitive Detector implementation of PXD and SVD.
Sensitive Detector implementation of PXD and SVD.
Base class to provide Sensor Information for PXD and SVD.
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:33
double getAngle(const std::string &path="") const noexcept(false)
Get the parameter path as a double converted to the standard angle unit.
Definition: Interface.h:299
double getLength(const std::string &path="") const noexcept(false)
Get the parameter path as a double converted to the standard length unit.
Definition: Interface.h:259
std::vector< GearDir > getNodes(const std::string &path="") const
Get vector of GearDirs which point to all the nodes the given path evaluates to.
Definition: Interface.cc:21
bool getBool(const std::string &path="") const noexcept(false)
Get the parameter path as a bool.
Definition: Interface.cc:80
int getInt(const std::string &path="") const noexcept(false)
Get the parameter path as a int.
Definition: Interface.cc:60
static G4Material * get(const std::string &name)
Find given material.
Definition: Materials.h:63
VXD::SensitiveDetector< PXDSimHit, PXDTrueHit > SensitiveDetector
The PXD Sensitive Detector class.
geometry::CreatorFactory< GeoPXDCreator > GeoPXDFactory("PXDCreator")
Register the creator.
G4Polycone * createRotationSolid(const std::string &name, const GearDir &params, double &minZ, double &maxZ)
Create a solid by roating two polylines around the Z-Axis.
Definition: utilities.cc:203
void setVisibility(G4LogicalVolume &volume, bool visible)
Helper function to quickly set the visibility of a given volume.
Definition: utilities.cc:108
void setColor(G4LogicalVolume &volume, const std::string &color)
Set the color of a logical volume.
Definition: utilities.cc:100
GeometryTypes
Flag indiciating the type of geometry to be used.
Abstract base class for different kinds of events.
STL namespace.
Very simple class to provide an easy way to register creators with the CreatorManager.