Belle II Software development
GeoSVDCreator.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 <svd/geometry/GeoSVDCreator.h>
10#include <vxd/geometry/GeoCache.h>
11#include <svd/geometry/SensorInfo.h>
12#include <svd/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 <G4Torus.hh>
33#include <G4Polycone.hh>
34#include <G4UnionSolid.hh>
35#include <G4SubtractionSolid.hh>
36#include <G4Region.hh>
37
38using namespace std;
39
40namespace Belle2 {
45 using namespace geometry;
47 namespace SVD {
48
51
53 {
54 for (SensorInfo* sensorInfo : m_SensorInfo) delete sensorInfo;
55 m_SensorInfo.clear();
56 }
57
59 {
60 const SVDSensorInfoPar& infoPar = dynamic_cast<const SVDSensorInfoPar&>(*sensor.getSensorInfo());
61
62 SensorInfo* info = new SensorInfo(
63 VxdID(0, 0, 0),
64 infoPar.getWidth(),
65 infoPar.getLength(),
66 infoPar.getThickness(),
67 infoPar.getUCells(),
68 infoPar.getVCells(),
69 infoPar.getWidth2()
70 );
71 info->setSensorParams(
72 infoPar.getStripEdgeU(),
73 infoPar.getStripEdgeV(),
74 infoPar.getDepletionVoltage(),
75 infoPar.getBiasVoltage(),
82 infoPar.getAduEquivalentU(),
83 infoPar.getAduEquivalentV(),
84 infoPar.getElectronicNoiseU(),
85 infoPar.getElectronicNoiseV(),
86 infoPar.getAduEquivalentSbwU(),
87 infoPar.getAduEquivalentSbwV(),
88 infoPar.getElectronicNoiseSbwU(),
90 );
91 m_SensorInfo.push_back(info);
92 return info;
93 }
94
96 {
97
98 const double unit_pFcm = 1;
99 // This was 1000 * Unit::fC / Unit::V / Unit::cm; // pF/cm.
100 // We only use ratios of capacities, and this gives nicer numbers.
102 VxdID(0, 0, 0),
103 sensor.getLength("width"),
104 sensor.getLength("length"),
105 sensor.getLength("height"),
106 sensor.getInt("stripsU"),
107 sensor.getInt("stripsV"),
108 sensor.getLength("width2", 0)
109 );
110
111 info->setSensorParams(
112 sensor.getWithUnit("stripEdgeU"),
113 sensor.getWithUnit("stripEdgeV"),
114 sensor.getWithUnit("DepletionVoltage"),
115 sensor.getWithUnit("BiasVoltage"),
116 sensor.getDouble("BackplaneCapacitanceU") * unit_pFcm,
117 sensor.getDouble("InterstripCapacitanceU") * unit_pFcm,
118 sensor.getDouble("CouplingCapacitanceU") * unit_pFcm,
119 sensor.getDouble("BackplaneCapacitanceV") * unit_pFcm,
120 sensor.getDouble("InterstripCapacitanceV") * unit_pFcm,
121 sensor.getDouble("CouplingCapacitanceV") * unit_pFcm,
122 sensor.getWithUnit("ADUEquivalentU"),
123 sensor.getWithUnit("ADUEquivalentV"),
124 sensor.getWithUnit("ElectronicNoiseU"),
125 sensor.getWithUnit("ElectronicNoiseV"),
126 sensor.getWithUnit("ADUEquivalentSbwU", 0),
127 sensor.getWithUnit("ADUEquivalentSbwV", 0),
128 sensor.getWithUnit("ElectronicNoiseSbwU", 0),
129 sensor.getWithUnit("ElectronicNoiseSbwV", 0)
130 );
131
132 return info;
133 }
134
136 VxdID sensorID, const VXDGeoSensor& sensor, const VXDGeoSensorPlacement&)
137 {
138 SensorInfo* sensorInfo = new SensorInfo(dynamic_cast<const SensorInfo&>(*sensor.getSensorInfo()));
139 sensorInfo->setID(sensorID);
140 SensitiveDetector* sensitive = new SensitiveDetector(sensorInfo);
141 return sensitive;
142 }
143
145 {
146 SVDGeometryPar svdGeometryPar;
147
148 //Read prefix ('SVD' or 'PXD')
149 svdGeometryPar.setPrefix(m_prefix);
150
151 //Read some global parameters
152 VXDGlobalPar globals((float)content.getDouble("ElectronTolerance", 100),
153 (float)content.getDouble("MinimumElectrons", 10),
154 content.getLength("ActiveStepSize", 0.005),
155 content.getBool("ActiveChips", false),
156 content.getBool("SeeNeutrons", false),
157 content.getBool("OnlyPrimaryTrueHits", false),
158 content.getBool("OnlyActiveMaterial", false),
159 (float)content.getLength("DistanceTolerance", 0.005),
160 content.getString("DefaultMaterial", "Air")
161 );
162 svdGeometryPar.setGlobalParams(globals);
163
164 //Read envelope parameters
165 GearDir envelopeParams(content, "Envelope/");
166 VXDEnvelopePar envelope(envelopeParams.getString("Name", ""),
167 envelopeParams.getString("Material", "Air"),
168 envelopeParams.getString("Color", ""),
169 envelopeParams.getAngle("minPhi", 0),
170 envelopeParams.getAngle("maxPhi", 2 * M_PI),
171 (envelopeParams.getNodes("InnerPoints/point").size() > 0)
172 );
173
174 for (const GearDir& point : envelopeParams.getNodes("InnerPoints/point")) {
175 pair<double, double> ZXPoint(point.getLength("z"), point.getLength("x"));
176 envelope.getInnerPoints().push_back(ZXPoint);
177 }
178 for (const GearDir& point : envelopeParams.getNodes("OuterPoints/point")) {
179 pair<double, double> ZXPoint(point.getLength("z"), point.getLength("x"));
180 envelope.getOuterPoints().push_back(ZXPoint);
181 }
182 svdGeometryPar.setEnvelope(envelope);
183
184 // Read alignment for detector m_prefix ('PXD' or 'SVD')
185 string pathAlign = (boost::format("Align[@component='%1%']/") % m_prefix).str();
186 GearDir paramsAlign(GearDir(content, "Alignment/"), pathAlign);
187 if (!paramsAlign) {
188 B2WARNING("Could not find alignment parameters for component " << m_prefix);
189 return svdGeometryPar;
190 }
191 svdGeometryPar.getAlignmentMap()[m_prefix] = VXDAlignmentPar(paramsAlign.getLength("du"),
192 paramsAlign.getLength("dv"),
193 paramsAlign.getLength("dw"),
194 paramsAlign.getAngle("alpha"),
195 paramsAlign.getAngle("beta"),
196 paramsAlign.getAngle("gamma")
197 );
198
199 //Read the definition of all sensor types
200 GearDir components(content, "Components/");
201 for (const GearDir& paramsSensor : components.getNodes("Sensor")) {
202 string sensorTypeID = paramsSensor.getString("@type");
203
204 VXDGeoSensorPar sensor(paramsSensor.getString("Material"),
205 paramsSensor.getString("Color", ""),
206 paramsSensor.getLength("width"),
207 paramsSensor.getLength("width2", 0),
208 paramsSensor.getLength("length"),
209 paramsSensor.getLength("height"),
210 paramsSensor.getAngle("angle", 0),
211 paramsSensor.getBool("@slanted", false)
212 );
213 sensor.setActive(VXDGeoComponentPar(
214 paramsSensor.getString("Material"),
215 paramsSensor.getString("Active/Color", "#f00"),
216 paramsSensor.getLength("Active/width"),
217 paramsSensor.getLength("Active/width2", 0),
218 paramsSensor.getLength("Active/length"),
219 paramsSensor.getLength("Active/height")
221 "Active",
222 paramsSensor.getLength("Active/u"),
223 paramsSensor.getLength("Active/v"),
224 paramsSensor.getString("Active/w", "center"),
225 paramsSensor.getLength("Active/woffset", 0)
226 ));
227
228 SVDSensorInfoPar* svdInfo = readSensorInfo(GearDir(paramsSensor, "Active"));
229 sensor.setSensorInfo(svdInfo);
230 sensor.setComponents(getSubComponents(paramsSensor));
231 svdGeometryPar.getSensorMap()[sensorTypeID] = sensor;
232 svdGeometryPar.getSensorInfos().push_back(svdInfo);
233 }
234
235 //Build all ladders including Sensors
236 GearDir support(content, "Support/");
237 readHalfShellSupport(support, svdGeometryPar);
238
239 for (const GearDir& shell : content.getNodes("HalfShell")) {
240
241 string shellName = m_prefix + "." + shell.getString("@name");
242 string pathShell = (boost::format("Align[@component='%1%']/") % shellName).str();
243 GearDir paramsShell(GearDir(content, "Alignment/"), pathShell);
244 if (!paramsShell) {
245 B2WARNING("Could not find alignment parameters for component " << shellName);
246 return svdGeometryPar;
247 }
248 svdGeometryPar.getAlignmentMap()[shellName] = VXDAlignmentPar(paramsShell.getLength("du"),
249 paramsShell.getLength("dv"),
250 paramsShell.getLength("dw"),
251 paramsShell.getAngle("alpha"),
252 paramsShell.getAngle("beta"),
253 paramsShell.getAngle("gamma")
254 );
255
256 VXDHalfShellPar halfShell(shell.getString("@name"), shell.getAngle("shellAngle", 0));
257
258 for (const GearDir& layer : shell.getNodes("Layer")) {
259 int layerID = layer.getInt("@id");
260
261 readLadder(layerID, components, svdGeometryPar);
262 readLayerSupport(layerID, support, svdGeometryPar);
263 readLadderSupport(layerID, support, svdGeometryPar);
264
265 //Loop over defined ladders
266 for (const GearDir& ladder : layer.getNodes("Ladder")) {
267 int ladderID = ladder.getInt("@id");
268 double phi = ladder.getAngle("phi", 0);
269 readLadderComponents(layerID, ladderID, content, svdGeometryPar);
270 halfShell.addLadder(layerID, ladderID, phi);
271 }
272 }
273 svdGeometryPar.getHalfShells().push_back(halfShell);
274 }
275
276 //Create diamond radiation sensors if defined and in background mode
277 GearDir radiationDir(content, "RadiationSensors");
278 if (svdGeometryPar.getGlobalParams().getActiveChips() && radiationDir) {
279 VXDGeoRadiationSensorsPar radiationSensors(
280 m_prefix,
281 radiationDir.getBool("insideEnvelope"),
282 radiationDir.getLength("width"),
283 radiationDir.getLength("length"),
284 radiationDir.getLength("height"),
285 radiationDir.getString("material")
286 );
287
288 //Add radiation sensor positions
289 for (GearDir& position : radiationDir.getNodes("position")) {
290 VXDGeoRadiationSensorsPositionPar diamonds(position.getLength("z"),
291 position.getLength("radius"),
292 position.getAngle("theta")
293 );
294
295 //Loop over all phi positions
296 for (GearDir& sensor : position.getNodes("phi")) {
297 //Add sensor with angle and id
298 diamonds.addSensor(sensor.getInt("@id"), sensor.getAngle());
299 }
300 radiationSensors.addPosition(diamonds);
301 }
302 svdGeometryPar.setRadiationSensors(radiationSensors);
303 }
304 return svdGeometryPar;
305 }
306
307 void GeoSVDCreator::createGeometry(const SVDGeometryPar& parameters, G4LogicalVolume& topVolume, geometry::GeometryTypes)
308 {
309
310 m_activeStepSize = parameters.getGlobalParams().getActiveStepSize() / Unit::mm;
311 m_activeChips = parameters.getGlobalParams().getActiveChips();
312 m_seeNeutrons = parameters.getGlobalParams().getSeeNeutrons();
313 m_onlyPrimaryTrueHits = parameters.getGlobalParams().getOnlyPrimaryTrueHits();
314 m_distanceTolerance = parameters.getGlobalParams().getDistanceTolerance();
315 m_electronTolerance = parameters.getGlobalParams().getElectronTolerance();
316 m_minimumElectrons = parameters.getGlobalParams().getMinimumElectrons();
317 m_onlyActiveMaterial = parameters.getGlobalParams().getOnlyActiveMaterial();
318 m_defaultMaterial = parameters.getGlobalParams().getDefaultMaterial();
319
320 G4Material* material = Materials::get(m_defaultMaterial);
321 if (!material) B2FATAL("Default Material of VXD, '" << m_defaultMaterial << "', could not be found");
322
323
324 //Build envelope
325 G4LogicalVolume* envelope(0);
326 G4VPhysicalVolume* physEnvelope{nullptr};
327 if (!parameters.getEnvelope().getExists()) {
328 B2INFO("Could not find definition for " + m_prefix + " Envelope, placing directly in top volume");
329 envelope = &topVolume;
330 } else {
331 double minZ(0), maxZ(0);
332 G4Polycone* envelopeCone = geometry::createRotationSolid("Envelope",
333 parameters.getEnvelope().getInnerPoints(),
334 parameters.getEnvelope().getOuterPoints(),
335 parameters.getEnvelope().getMinPhi(),
336 parameters.getEnvelope().getMaxPhi(),
337 minZ, maxZ
338 );
339 envelope = new G4LogicalVolume(envelopeCone, material, m_prefix + ".Envelope");
340 setVisibility(*envelope, false);
341 physEnvelope = new G4PVPlacement(getAlignment(parameters.getAlignment(m_prefix)), envelope, m_prefix + ".Envelope",
342 &topVolume, false, 1);
343
344 // Set up region for production cuts
345 G4Region* aRegion = new G4Region("SVDEnvelope");
346 envelope->SetRegion(aRegion);
347 aRegion->AddRootLogicalVolume(envelope);
348 }
349
350 //Read the definition of all sensor types
351 for (const pair<const string, VXDGeoSensorPar>& typeAndSensor : parameters.getSensorMap()) {
352 const string& sensorTypeID = typeAndSensor.first;
353 const VXDGeoSensorPar& paramsSensor = typeAndSensor.second;
354 VXDGeoSensor sensor(
355 paramsSensor.getMaterial(),
356 paramsSensor.getColor(),
357 paramsSensor.getWidth() / Unit::mm,
358 paramsSensor.getWidth2() / Unit::mm,
359 paramsSensor.getLength() / Unit::mm,
360 paramsSensor.getHeight() / Unit::mm,
361 paramsSensor.getSlanted()
362 );
363 sensor.setActive(VXDGeoComponent(
364 paramsSensor.getMaterial(),
365 paramsSensor.getActiveArea().getColor(),
366 paramsSensor.getActiveArea().getWidth() / Unit::mm,
367 paramsSensor.getActiveArea().getWidth2() / Unit::mm,
368 paramsSensor.getActiveArea().getLength() / Unit::mm,
369 paramsSensor.getActiveArea().getHeight() / Unit::mm
371 "Active",
372 paramsSensor.getActivePlacement().getU() / Unit::mm,
373 paramsSensor.getActivePlacement().getV() / Unit::mm,
374 paramsSensor.getActivePlacement().getW(),
375 paramsSensor.getActivePlacement().getWOffset() / Unit::mm
376 ));
377 sensor.setSensorInfo(createSensorInfo(paramsSensor));
378
379 vector<VXDGeoPlacement> subcomponents;
380 for (const VXDGeoPlacementPar& component : paramsSensor.getComponents()) {
381 subcomponents.push_back(VXDGeoPlacement(
382 component.getName(),
383 component.getU() / Unit::mm,
384 component.getV() / Unit::mm,
385 component.getW(),
386 component.getWOffset() / Unit::mm
387 ));
388 }
389 sensor.setComponents(subcomponents);
390 m_sensorMap[sensorTypeID] = sensor;
391 }
392
393 //Read the component cache from DB
394 for (const string& name : parameters.getComponentInsertOder()) {
395 if (m_componentCache.find(name) != m_componentCache.end()) {
396 // already created due to being a sub component of a previous
397 // component. Seems fishy since the information of this component
398 // is in the db at least twice so we could run into
399 // inconsistencies.
400 B2WARNING("Component " << name << " already created from previous subcomponents, should not be here");
401 continue;
402 }
403 const VXDGeoComponentPar& paramsComponent = parameters.getComponent(name);
405 paramsComponent.getMaterial(),
406 paramsComponent.getColor(),
407 paramsComponent.getWidth() / Unit::mm,
408 paramsComponent.getWidth2() / Unit::mm,
409 paramsComponent.getLength() / Unit::mm,
410 paramsComponent.getHeight() / Unit::mm
411 );
412 double angle = paramsComponent.getAngle();
413
414
415 if (c.getWidth() <= 0 || c.getLength() <= 0 || c.getHeight() <= 0) {
416 B2DEBUG(100, "One dimension empty, using auto resize for component");
417 } else {
418 G4VSolid* solid = createTrapezoidal(m_prefix + "." + name, c.getWidth(), c.getWidth2(), c.getLength(), c.getHeight(), angle);
419 c.setVolume(new G4LogicalVolume(solid, Materials::get(c.getMaterial()), m_prefix + "." + name));
420 }
421
422 vector<VXDGeoPlacement> subComponents;
423 for (const VXDGeoPlacementPar& paramsSubComponent : paramsComponent.getSubComponents()) {
424 subComponents.push_back(VXDGeoPlacement(
425 paramsSubComponent.getName(),
426 paramsSubComponent.getU() / Unit::mm,
427 paramsSubComponent.getV() / Unit::mm,
428 paramsSubComponent.getW(),
429 paramsSubComponent.getWOffset() / Unit::mm
430 ));
431
432 }
433
434 createSubComponents(m_prefix + "." + name, c, subComponents);
435 if (m_activeChips && parameters.getSensitiveChipID(name) >= 0) {
436 int chipID = parameters.getSensitiveChipID(name);
437 B2DEBUG(50, "Creating BkgSensitiveDetector for component " << name << " with chipID " << chipID);
438 BkgSensitiveDetector* sensitive = new BkgSensitiveDetector(m_prefix.c_str(), chipID);
439 c.getVolume()->SetSensitiveDetector(sensitive);
440 m_sensitive.push_back(sensitive);
441 }
442
443 m_componentCache[name] = c;
444 }
445
446 //Build all ladders including Sensors
447 VXD::GeoVXDAssembly shellSupport = createHalfShellSupport(parameters);
448
449 //const std::vector<VXDHalfShellPar>& HalfShells = parameters.getHalfShells();
450 for (const VXDHalfShellPar& shell : parameters.getHalfShells()) {
451 string shellName = shell.getName();
452 m_currentHalfShell = m_prefix + "." + shellName;
453 G4Transform3D shellAlignment = getAlignment(parameters.getAlignment(m_currentHalfShell));
454
455 // Remember shell coordinate system (into which ladders are inserted)
457
458 //Place shell support
459 double shellAngle = shell.getShellAngle();
460 if (!m_onlyActiveMaterial) shellSupport.place(envelope, shellAlignment * G4RotateZ3D(shellAngle));
461
462 //const std::map< int, std::vector<std::pair<int, double>> >& Layers = shell.getLayers();
463 for (const std::pair<const int, std::vector<std::pair<int, double>> >& layer : shell.getLayers()) {
464 int layerID = layer.first;
465 const std::vector<std::pair<int, double>>& Ladders = layer.second;
466
467
468 setCurrentLayer(layerID, parameters);
469
470 //Place Layer support
471 VXD::GeoVXDAssembly layerSupport = createLayerSupport(layerID, parameters);
472 if (!m_onlyActiveMaterial) layerSupport.place(envelope, shellAlignment * G4RotateZ3D(shellAngle));
473 VXD::GeoVXDAssembly ladderSupport = createLadderSupport(layerID, parameters);
474
475 //Loop over defined ladders
476 for (const std::pair<int, double>& ladder : Ladders) {
477 int ladderID = ladder.first;
478 double phi = ladder.second;
479
480 G4Transform3D ladderPlacement = placeLadder(ladderID, phi, envelope, shellAlignment, parameters);
481 if (!m_onlyActiveMaterial) ladderSupport.place(envelope, ladderPlacement);
482 }
483
484 }
485 }
486
487 //Now build cache with all transformations
488 if (physEnvelope) {
490 } else {
491 //create a temporary placement of the top volume.
492 G4PVPlacement topPlacement(nullptr, G4ThreeVector(0, 0, 0), &topVolume,
493 "temp_Top", nullptr, false, 1, false);
494 //and search for all VXD sensitive sensors within
496 }
497
498 //Create diamond radiation sensors if defined and in background mode
499 if (m_activeChips) {
500 if (parameters.getRadiationSensors().getSubDetector() == "") {
501 B2DEBUG(10, "Apparently no radiation sensors defined, skipping");
502 } else {
503 createDiamonds(parameters.getRadiationSensors(), topVolume, *envelope);
504 }
505 }
506 }
507
508 void GeoSVDCreator::readHalfShellSupport(const GearDir& support, SVDGeometryPar& svdGeometryPar)
509 {
510 if (!support) return;
511
512 for (const GearDir& params : support.getNodes("HalfShell/RotationSolid")) {
513
514 VXDRotationSolidPar rotationSolidPar(params.getString("Name", ""),
515 params.getString("Material", "Air"),
516 params.getString("Color", ""),
517 params.getAngle("minPhi", 0),
518 params.getAngle("maxPhi", 2 * M_PI),
519 (params.getNodes("InnerPoints/point").size() > 0)
520 );
521
522 for (const GearDir& point : params.getNodes("InnerPoints/point")) {
523 pair<double, double> ZXPoint(point.getLength("z"), point.getLength("x"));
524 rotationSolidPar.getInnerPoints().push_back(ZXPoint);
525 }
526 for (const GearDir& point : params.getNodes("OuterPoints/point")) {
527 pair<double, double> ZXPoint(point.getLength("z"), point.getLength("x"));
528 rotationSolidPar.getOuterPoints().push_back(ZXPoint);
529 }
530 svdGeometryPar.getRotationSolids().push_back(rotationSolidPar);
531 }
532 return;
533 }
534
535 void GeoSVDCreator::readLayerSupport(int layer, const GearDir& support, SVDGeometryPar& svdGeometryPar)
536 {
537 if (!support) return;
538
539 //Check if there are any endrings defined for this layer. If not we don't create any
540 GearDir endrings(support, (boost::format("Endrings/Layer[@id='%1%']") % layer).str());
541 if (endrings) {
542 svdGeometryPar.getEndrings()[layer] = SVDEndringsPar(support.getString("Endrings/Material"),
543 support.getLength("Endrings/length"),
544 support.getLength("Endrings/gapWidth"),
545 support.getLength("Endrings/baseThickness")
546 );
547
548 //Create the endrings
549 for (const GearDir& endring : endrings.getNodes("Endring")) {
550 SVDEndringsTypePar endringPar(endring.getString("@name"),
551 endring.getLength("z"),
552 endring.getLength("baseRadius"),
553 endring.getLength("innerRadius"),
554 endring.getLength("outerRadius"),
555 endring.getLength("horizontalBar"),
556 endring.getLength("verticalBar")
557 );
558 svdGeometryPar.getEndrings()[layer].getTypes().push_back(endringPar);
559 }
560 }
561
562 // Now let's add the cooling pipes to the Support
563 GearDir pipes(support, (boost::format("CoolingPipes/Layer[@id='%1%']") % layer).str());
564 if (pipes) {
565 svdGeometryPar.getCoolingPipes()[layer] = SVDCoolingPipesPar(support.getString("CoolingPipes/Material"),
566 support.getLength("CoolingPipes/outerDiameter"),
567 support.getLength("CoolingPipes/wallThickness"),
568 pipes.getInt("nPipes"),
569 pipes.getAngle("startPhi"),
570 pipes.getAngle("deltaPhi"),
571 pipes.getLength("radius"),
572 pipes.getLength("zstart"),
573 pipes.getLength("zend")
574 );
575 if (pipes.exists("deltaL")) svdGeometryPar.getCoolingPipes()[layer].setDeltaL(pipes.getLength("deltaL"));
576 }
577 return;
578 }
579
580 void GeoSVDCreator::readLadderSupport(int layer, const GearDir& support, SVDGeometryPar& svdGeometryPar)
581 {
582 if (!support) return;
583
584 // Check if there are any support ribs defined for this layer. If not return empty assembly
585 GearDir params(support, (boost::format("SupportRibs/Layer[@id='%1%']") % layer).str());
586 if (params) {
587 svdGeometryPar.getSupportRibs()[layer] = SVDSupportRibsPar(support.getLength("SupportRibs/spacing"),
588 support.getLength("SupportRibs/height"),
589 support.getLength("SupportRibs/inner/width"),
590 support.getLength("SupportRibs/outer/width"),
591 support.getLength("SupportRibs/inner/tabLength"),
592 support.getString("SupportRibs/outer/Material"),
593 support.getString("SupportRibs/inner/Material"),
594 support.getString("SupportRibs/outer/Color"),
595 support.getString("SupportRibs/inner/Color"),
596 support.getString("SupportRibs/endmount/Material")
597 );
598
599 // Get values for the layer if available
600 if (params.exists("spacing")) svdGeometryPar.getSupportRibs()[layer].setSpacing(params.getLength("spacing"));
601 if (params.exists("height")) svdGeometryPar.getSupportRibs()[layer].setHeight(params.getLength("height"));
602
603 for (const GearDir& box : params.getNodes("box")) {
604 SVDSupportBoxPar boxPar(box.getAngle("theta"),
605 box.getLength("z"),
606 box.getLength("r"),
607 box.getLength("length")
608 );
609 svdGeometryPar.getSupportRibs()[layer].getBoxes().push_back(boxPar);
610 }
611
612 for (const GearDir& tab : params.getNodes("tab")) {
613 SVDSupportTabPar tabPar(tab.getAngle("theta"),
614 tab.getLength("z"),
615 tab.getLength("r")
616 );
617 svdGeometryPar.getSupportRibs()[layer].getTabs().push_back(tabPar);
618 }
619
620 for (const GearDir& endmount : params.getNodes("Endmount")) {
621 SVDEndmountPar mountPar(endmount.getString("@name"),
622 endmount.getLength("height"),
623 endmount.getLength("width"),
624 endmount.getLength("length"),
625 endmount.getLength("z"),
626 endmount.getLength("r")
627 );
628 svdGeometryPar.getSupportRibs()[layer].getEndmounts().push_back(mountPar);
629 }
630 }
631 return;
632 }
633
635 {
636 VXD::GeoVXDAssembly supportAssembly;
637
638 //Half shell support is easy as we just add all the defined RotationSolids
639 double minZ(0), maxZ(0);
640
641 const std::vector<VXDRotationSolidPar>& RotationSolids = parameters.getRotationSolids();
642
643 for (const VXDRotationSolidPar& component : RotationSolids) {
644
645 string name = component.getName();
646 string material = component.getMaterial();
647
648 G4Polycone* solid = geometry::createRotationSolid(name,
649 component.getInnerPoints(),
650 component.getOuterPoints(),
651 component.getMinPhi(),
652 component.getMaxPhi(),
653 minZ, maxZ
654 );
655
656 G4LogicalVolume* volume = new G4LogicalVolume(
657 solid, geometry::Materials::get(material), m_prefix + ". " + name);
658 geometry::setColor(*volume, component.getColor());
659 supportAssembly.add(volume);
660 }
661 return supportAssembly;
662 }
663
664
665
667 {
668 VXD::GeoVXDAssembly supportAssembly;
669
670 //Check if there are any endrings defined for this layer. If not we don't create any
671 if (parameters.getEndringsExist(layer)) {
672 const SVDEndringsPar& support = parameters.getEndrings(layer);
673
674 string material = support.getMaterial();
675 double length = support.getLength() / Unit::mm / 2.0;
676 double gapWidth = support.getGapWidth() / Unit::mm;
677 double baseThickness = support.getBaseThickness() / Unit::mm / 2.0;
678
679 //Create the endrings
680 const std::vector<SVDEndringsTypePar>& Endrings = support.getTypes();
681 for (const SVDEndringsTypePar& endring : Endrings) {
682 double z = endring.getZ() / Unit::mm;
683 double baseRadius = endring.getBaseRadius() / Unit::mm;
684 double innerRadius = endring.getInnerRadius() / Unit::mm;
685 double outerRadius = endring.getOuterRadius() / Unit::mm;
686 double horiBarWidth = endring.getHorizontalBarWidth() / Unit::mm / 2.0;
687 double vertBarWidth = endring.getVerticalBarWidth() / Unit::mm / 2.0;
688
689 double angle = asin(gapWidth / innerRadius);
690 G4VSolid* endringSolid = new G4Tubs("OuterEndring", innerRadius, outerRadius, length, -M_PI / 2 + angle, M_PI - 2 * angle);
691 angle = asin(gapWidth / baseRadius);
692 G4VSolid* endringBase = new G4Tubs("InnerEndring", baseRadius, baseRadius + baseThickness, length, -M_PI / 2 + angle,
693 M_PI - 2 * angle);
694 endringSolid = new G4UnionSolid("Endring", endringSolid, endringBase);
695
696 //Now we need the bars which connect the two rings
697 double height = (innerRadius - baseRadius) / 2.0;
698 double x = vertBarWidth + gapWidth;
699 G4Box* verticalBar = new G4Box("VerticalBar", vertBarWidth, height, length);
700 G4Box* horizontalBar = new G4Box("HorizontalBar", height, horiBarWidth, length);
701 endringSolid = new G4UnionSolid("Endring", endringSolid, verticalBar, G4Translate3D(x, baseRadius + height, 0));
702 endringSolid = new G4UnionSolid("Endring", endringSolid, verticalBar, G4Translate3D(x, -(baseRadius + height), 0));
703 endringSolid = new G4UnionSolid("Endring", endringSolid, horizontalBar, G4Translate3D((baseRadius + height), 0, 0));
704
705 //Finally create the volume and add it to the assembly at the correct z position
706 G4LogicalVolume* endringVolume = new G4LogicalVolume(
707 endringSolid, geometry::Materials::get(material),
708 (boost::format("%1%.Layer%2%.%3%") % m_prefix % layer % endring.getName()).str());
709 supportAssembly.add(endringVolume, G4TranslateZ3D(z));
710 }
711 }
712
713 //Check if there are any coling pipes defined for this layer. If not we don't create any
714 if (parameters.getCoolingPipesExist(layer)) {
715 const SVDCoolingPipesPar& pipes = parameters.getCoolingPipes(layer);
716
717 string material = pipes.getMaterial();
718 double outerRadius = pipes.getOuterDiameter() / Unit::mm / 2.0;
719 double innerRadius = outerRadius - pipes.getWallThickness() / Unit::mm;
720 int nPipes = pipes.getNPipes();
721 double startPhi = pipes.getStartPhi();
722 double deltaPhi = pipes.getDeltaPhi();
723 double radius = pipes.getRadius() / Unit::mm;
724 double zstart = pipes.getZStart() / Unit::mm;
725 double zend = pipes.getZEnd() / Unit::mm;
726 double zlength = (zend - zstart) / 2.0;
727
728
729 // There are two parts: the straight pipes and the bendings. So we only need two different volumes
730 // which we place multiple times
731 G4Tubs* pipeSolid = new G4Tubs("CoolingPipe", innerRadius, outerRadius, zlength, 0, 2 * M_PI);
732 G4LogicalVolume* pipeVolume = new G4LogicalVolume(
733 pipeSolid, geometry::Materials::get(material),
734 (boost::format("%1%.Layer%2%.CoolingPipe") % m_prefix % layer).str());
735 geometry::setColor(*pipeVolume, "#ccc");
736
737#ifndef __clang_analyzer__
738 G4Torus* bendSolid = new G4Torus("CoolingBend", innerRadius, outerRadius, sin(deltaPhi / 2.0)*radius, -M_PI / 2, M_PI);
739 G4LogicalVolume* bendVolume = new G4LogicalVolume(
740 bendSolid, geometry::Materials::get(material),
741 (boost::format("%1%.Layer%2%.CoolingBend") % m_prefix % layer).str());
742#endif
743 // Last pipe may be closer, thus we need additional bending
744 if (pipes.getDeltaL() > 0) {
745 double deltaL = pipes.getDeltaL() / Unit::mm;
746 G4Torus* bendSolidLast = new G4Torus("CoolingBendLast", innerRadius, outerRadius, sin(deltaPhi / 2.0) * radius - deltaL / 2.0,
747 -M_PI / 2, M_PI);
748 G4LogicalVolume* bendVolumeLast = new G4LogicalVolume(bendSolidLast, geometry::Materials::get(material),
749 (boost::format("%1%.Layer%2%.CoolingBendLast") % m_prefix % layer).str());
750 --nPipes;
751
752 // Place the last straight pipe
753 G4Transform3D placement_pipe = G4RotateZ3D(startPhi + (nPipes - 0.5) * deltaPhi) * G4Translate3D(cos(deltaPhi / 2.0) * radius,
754 sin(deltaPhi / 2.0) * radius - deltaL, zstart + zlength);
755 supportAssembly.add(pipeVolume, placement_pipe);
756
757 // Place forward or backward bend
758 double zpos = nPipes % 2 > 0 ? zend : zstart;
759 // Calculate transformation
760 G4Transform3D placement = G4RotateZ3D(startPhi + (nPipes - 0.5) * deltaPhi) * G4Translate3D(cos(deltaPhi / 2.0) * radius,
761 -deltaL / 2.0, zpos) * G4RotateY3D(M_PI / 2);
762 // If we are at the forward side we rotate the bend by 180 degree
763 if (nPipes % 2 > 0) {
764 placement = placement * G4RotateZ3D(M_PI);
765 }
766 // And place the bend
767 supportAssembly.add(bendVolumeLast, placement);
768 }
769
770 for (int i = 0; i < nPipes; ++i) {
771 // Place the straight pipes
772 G4Transform3D placement_pipe = G4RotateZ3D(startPhi + i * deltaPhi) * G4Translate3D(radius, 0, zstart + zlength);
773 supportAssembly.add(pipeVolume, placement_pipe);
774
775 // This was the easy part, now lets add the connection between the pipes. We only need n-1 bendings
776 if (i > 0) {
777 // Place forward or backward bend
778 double zpos = i % 2 > 0 ? zend : zstart;
779 // Calculate transformation
780 G4Transform3D placement = G4RotateZ3D(startPhi + (i - 0.5) * deltaPhi) * G4Translate3D(cos(deltaPhi / 2.0) * radius, 0,
781 zpos) * G4RotateY3D(M_PI / 2);
782 // If we are at the forward side we rotate the bend by 180 degree
783 if (i % 2 > 0) {
784 placement = placement * G4RotateZ3D(M_PI);
785 }
786#ifndef __clang_analyzer__
787 // And place the bend
788 supportAssembly.add(bendVolume, placement);
789#endif
790 }
791 }
792
793 }
794
795 return supportAssembly;
796 }
797
798
799
801 {
802 VXD::GeoVXDAssembly supportAssembly;
803
804 if (!parameters.getSupportRibsExist(layer)) return supportAssembly;
805 const SVDSupportRibsPar& support = parameters.getSupportRibs(layer);
806
807 // Get the common values for all layers
808 double spacing = support.getSpacing() / Unit::mm / 2.0;
809 double height = support.getHeight() / Unit::mm / 2.0;
810 double innerWidth = support.getInnerWidth() / Unit::mm / 2.0;
811 double outerWidth = support.getOuterWidth() / Unit::mm / 2.0;
812 double tabLength = support.getTabLength() / Unit::mm / 2.0;
813 G4VSolid* inner(0);
814 G4VSolid* outer(0);
815 G4Transform3D placement;
816
817
818 // Now lets create the ribs by adding all boxes to form one union solid
819 const std::vector<SVDSupportBoxPar>& Boxes = support.getBoxes();
820 for (const SVDSupportBoxPar& box : Boxes) {
821 double theta = box.getTheta();
822 double zpos = box.getZ() / Unit::mm;
823 double rpos = box.getR() / Unit::mm;
824 double length = box.getLength() / Unit::mm / 2.0;
825 G4Box* innerBox = new G4Box("innerBox", height, innerWidth, length);
826 G4Box* outerBox = new G4Box("outerBox", height, outerWidth, length);
827 if (!inner) {
828 inner = innerBox;
829 outer = outerBox;
830 placement = G4Translate3D(rpos, 0, zpos) * G4RotateY3D(theta);
831 } else {
832 G4Transform3D relative = placement.inverse() * G4Translate3D(rpos, 0, zpos) * G4RotateY3D(theta);
833 inner = new G4UnionSolid("innerBox", inner, innerBox, relative);
834 outer = new G4UnionSolid("outerBox", outer, outerBox, relative);
835 }
836 }
837
838 // Now lets add the tabs
839 const std::vector<SVDSupportTabPar>& Tabs = support.getTabs();
840 for (const SVDSupportTabPar& tab : Tabs) {
841 double theta = tab.getTheta();
842 double zpos = tab.getZ() / Unit::mm;
843 double rpos = tab.getR() / Unit::mm;
844 G4Box* innerBox = new G4Box("innerBox", height, innerWidth, tabLength);
845 if (!inner) {
846 inner = innerBox;
847 placement = G4Translate3D(rpos, 0, zpos) * G4RotateY3D(theta);
848 } else {
849 G4Transform3D relative = placement.inverse() * G4Translate3D(rpos, 0, zpos) * G4RotateY3D(theta);
850 inner = new G4UnionSolid("innerBox", inner, innerBox, relative);
851 }
852 }
853
854 // Now lets create forward and backward endmounts for the ribs
855 const std::vector<SVDEndmountPar>& Endmounts = support.getEndmounts();
856 for (const SVDEndmountPar& endmount : Endmounts) {
857 double endMountHeight = endmount.getHeight() / Unit::mm / 2.0;
858 double endMountWidth = endmount.getWidth() / Unit::mm / 2.0;
859 double endMountLength = endmount.getLength() / Unit::mm / 2.0;
860 double zpos = endmount.getZ() / Unit::mm;
861 double rpos = endmount.getR() / Unit::mm;
862 G4VSolid* endmountBox = new G4Box("endmountBox", endMountHeight, endMountWidth, endMountLength);
863 if (outer) { // holes for the ribs
864 endmountBox = new G4SubtractionSolid("endmountBox", endmountBox, outer, G4TranslateY3D(-spacing)*placement * G4Translate3D(-rpos, 0,
865 -zpos));
866 endmountBox = new G4SubtractionSolid("endmountBox", endmountBox, outer, G4TranslateY3D(spacing)*placement * G4Translate3D(-rpos, 0,
867 -zpos));
868 }
869 G4LogicalVolume* endmountVolume = new G4LogicalVolume(
870 endmountBox, geometry::Materials::get(support.getEndmountMaterial()),
871 (boost::format("%1%.Layer%2%.%3%Endmount") % m_prefix % layer % endmount.getName()).str());
872 supportAssembly.add(endmountVolume, G4Translate3D(rpos, 0, zpos));
873 }
874
875 // If there has been at least one Box, create the volumes and add them to the assembly
876 if (inner) {
877 outer = new G4SubtractionSolid("outerBox", outer, inner);
878 G4LogicalVolume* outerVolume = new G4LogicalVolume(
880 (boost::format("%1%.Layer%2%.SupportRib") % m_prefix % layer).str());
881 G4LogicalVolume* innerVolume = new G4LogicalVolume(
883 (boost::format("%1%.Layer%2%.SupportRib.Airex") % m_prefix % layer).str());
884 geometry::setColor(*outerVolume, support.getOuterColor());
885 geometry::setColor(*innerVolume, support.getInnerColor());
886 supportAssembly.add(innerVolume, G4TranslateY3D(-spacing)*placement);
887 supportAssembly.add(innerVolume, G4TranslateY3D(spacing)*placement);
888 supportAssembly.add(outerVolume, G4TranslateY3D(-spacing)*placement);
889 supportAssembly.add(outerVolume, G4TranslateY3D(spacing)*placement);
890 }
891
892 // Done, return the finished assembly
893 return supportAssembly;
894 }
895
896
897
898 }
900}
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 SVD Cooling Pipes.
double getStartPhi() const
get start phi
double getOuterDiameter() const
get outer diameter
double getDeltaPhi() const
get delta phi
int getNPipes() const
get nPipes
double getDeltaL() const
get deltal
double getZEnd() const
get zend
double getZStart() const
get zstart
const std::string & getMaterial() const
get material
double getWallThickness() const
get wall thickness
double getRadius() const
get radius
The Class for SVD Support Rib Endmounts.
The Class for SVD Endring.
double getGapWidth() const
get gap width
double getBaseThickness() const
get base thickness
const std::string & getMaterial() const
get material
const std::vector< SVDEndringsTypePar > & getTypes() const
get types (forward/backward)
double getLength() const
get length
The Class for SVD Endring Type.
The Class for VXD geometry.
const std::map< int, SVDSupportRibsPar > & getSupportRibs() const
get support ribs
const std::map< int, SVDCoolingPipesPar > & getCoolingPipes() const
get cooling pipes
std::vector< SVDSensorInfoPar * > & getSensorInfos()
get sensorInfos
const std::vector< VXDRotationSolidPar > & getRotationSolids() const
get SVD halfshell Rotation Solids
const std::map< int, SVDEndringsPar > & getEndrings() const
get endrings
The Class for VXD geometry.
double getBiasVoltage() const
Return the bias voltage on the sensor.
double getAduEquivalentSbwU() const
Return ADU equivalent for U strips in Sbw barrel sensor.
double getElectronicNoiseSbwU() const
Return electronic noise in e- for u strips in bw barrel sensors.
double getBackplaneCapacitanceV() const
Return the backplane capacitance/cm for V-side strips.
double getAduEquivalentU() const
Return ADU equivalent for U strips.
double getElectronicNoiseV() const
Return electronic noise in e- for v strips.
double getStripEdgeV() const
Return the distance between end of strip and edge of active area.
double getAduEquivalentSbwV() const
Return ADU equivalent for V strips in Sbw barrel sensor.
double getBackplaneCapacitanceU() const
Return the backplane capacitance/cm for U-side strips.
double getInterstripCapacitanceU() const
Return the interstrip capacitance/cm for U-side strips.
double getElectronicNoiseU() const
Return electronic noise in e- for u strips.
double getStripEdgeU() const
Return the distance between end of strip and edge of active area.
double getCouplingCapacitanceU() const
Return the coupling capacitance/cm for U-side strips.
double getCouplingCapacitanceV() const
Return the coupling capacitance/cm for V-side strips.
double getAduEquivalentV() const
Return ADU equivalent for V strips.
double getDepletionVoltage() const
Return the depletion voltage of the sensor.
double getInterstripCapacitanceV() const
Return the interstrip capacitance/cm for V-side strips.
double getElectronicNoiseSbwV() const
Return electronic noise in e- for v strips in bw barrel sensors.
The Class for SVD Support Box.
The Class for SVD Support Ribs (one layer)
double getHeight() const
get height
double getInnerWidth() const
get inner width
const std::string & getOuterMaterial() const
get the name of outer material
double getTabLength() const
get tabLength
double getOuterWidth() const
get outer width
const std::string & getOuterColor() const
get the name of outer color
const std::string & getInnerMaterial() const
get the name of inner material
double getSpacing() const
get spacing
const std::string & getInnerColor() const
get the name of inner color
const std::vector< SVDSupportTabPar > & getTabs() const
get tabs
const std::vector< SVDEndmountPar > & getEndmounts() const
get endmounts
const std::string & getEndmountMaterial() const
get the name of endmount material
const std::vector< SVDSupportBoxPar > & getBoxes() const
get boxes
The Class for SVD Support Rib Tab.
virtual VXD::SensitiveDetectorBase * createSensitiveDetector(VxdID sensorID, const VXDGeoSensor &sensor, const VXDGeoSensorPlacement &placement) override
Return a SensitiveDetector implementation for a given sensor.
void readLayerSupport(int layer, const GearDir &support, SVDGeometryPar &svdGeometryPar)
Create support structure for a SVD Layer.
void readLadderSupport(int layer, const GearDir &support, SVDGeometryPar &svdGeometryPar)
Create support structure for a SVD Ladder.
SVDGeometryPar createConfiguration(const GearDir &param)
Create a parameter object from the Gearbox XML parameters.
void readHalfShellSupport(const GearDir &support, SVDGeometryPar &svdGeometryPar)
Create support structure for SVD Half Shell, that means everything that does not depend on layer or s...
SVDSensorInfoPar * readSensorInfo(const GearDir &sensor)
Read the sensor definitions from gearbox.
std::vector< SensorInfo * > m_SensorInfo
Vector of pointers to SensorInfo objects.
virtual VXD::GeoVXDAssembly createLayerSupport(int, const SVDGeometryPar &parameters)
Create support structure for a SVD Layer.
virtual ~GeoSVDCreator()
The destructor of the GeoSVDCreator class.
virtual VXD::GeoVXDAssembly createLadderSupport(int, const SVDGeometryPar &parameters)
Create support structure for a SVD Ladder.
virtual VXD::GeoVXDAssembly createHalfShellSupport(const SVDGeometryPar &parameters)
Create support structure for SVD Half Shell, that means everything that does not depend on layer or s...
void createGeometry(const SVDGeometryPar &parameters, G4LogicalVolume &topVolume, geometry::GeometryTypes type)
Create the geometry from a parameter object.
virtual VXD::SensorInfoBase * createSensorInfo(const VXDGeoSensorPar &sensor) override
Read the sensor definitions from the database.
Specific implementation of SensorInfo for SVD Sensors which provides additional sensor specific infor...
Definition: SensorInfo.h:25
void setID(VxdID id)
Change the SensorID.
Definition: SensorInfo.h:85
static const double mm
[millimeters]
Definition: Unit.h:70
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.
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 Envelope parameters.
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 getVCells() const
Return number of pixel/strips in v direction.
int getUCells() const
Return number of pixel/strips in u direction.
double getThickness() const
Return the thickness of the sensor.
double getWidth2() const
Return forward width for a slanted sensor.
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
bool exists(const std::string &path="") const
Check if a given parameter path exists.
Definition: Interface.h:54
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< SVDSimHit, SVDTrueHit > SensitiveDetector
The SVD Sensitive Detector class.
geometry::CreatorFactory< GeoSVDCreator > GeoSVDFactory("SVDCreator")
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.