Belle II Software  release-05-01-25
DosiCreator.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Martin Ritter, Igal Jaegle *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <beast/dosi/geometry/DosiCreator.h>
12 #include <beast/dosi/simulation/SensitiveDetector.h>
13 
14 #include <geometry/Materials.h>
15 #include <geometry/CreatorFactory.h>
16 #include <framework/gearbox/GearDir.h>
17 
18 #include <cmath>
19 #include <boost/format.hpp>
20 #include <boost/foreach.hpp>
21 #include <boost/algorithm/string.hpp>
22 
23 #include <G4LogicalVolume.hh>
24 #include <G4PVPlacement.hh>
25 
26 //Shapes
27 #include <G4Box.hh>
28 #include <G4UserLimits.hh>
29 
30 //Visualization Attributes
31 #include <G4VisAttributes.hh>
32 
33 using namespace std;
34 using namespace boost;
35 
36 namespace Belle2 {
43  namespace dosi {
44 
45  // Register the creator
47  geometry::CreatorFactory<DosiCreator> DosiFactory("DOSICreator");
48 
49  DosiCreator::DosiCreator(): m_sensitive(0)
50  {
52  }
53 
54  DosiCreator::~DosiCreator()
55  {
56  if (m_sensitive) delete m_sensitive;
57  }
58 
59  void DosiCreator::create(const GearDir& content, G4LogicalVolume& topVolume, geometry::GeometryTypes /* type */)
60  {
61  //lets get the stepsize parameter with a default value of 5 µm
62  double stepSize = content.getLength("stepSize", 5 * CLHEP::um);
63 
64  G4VisAttributes* red = new G4VisAttributes(G4Colour(1, 0, 0));
65  red->SetForceAuxEdgeVisible(true);
66  G4VisAttributes* green = new G4VisAttributes(G4Colour(0, 1, 0));
67  green->SetForceAuxEdgeVisible(true);
68  G4VisAttributes* gray = new G4VisAttributes(G4Colour(.5, .5, .5));
69  gray->SetForceAuxEdgeVisible(true);
70  G4VisAttributes* coppercolor = new G4VisAttributes(G4Colour(218. / 255., 138. / 255., 103. / 255.));
71  coppercolor->SetForceAuxEdgeVisible(true);
72 
73  //Lets loop over all the Active nodes
74  BOOST_FOREACH(const GearDir & activeParams, content.getNodes("Active")) {
75  G4double dx_dosi = activeParams.getLength("dx_dosi") / 2.*CLHEP::cm;
76  G4double dy_dosi = activeParams.getLength("dy_dosi") / 2.*CLHEP::cm;
77  G4double dz_dosi = activeParams.getLength("dz_dosi") / 2.*CLHEP::cm;
78  double thetaZ = activeParams.getAngle("ThetaZ");
79  G4VSolid* s_dosi = new G4Box("s_dosi", dx_dosi, dy_dosi, dz_dosi);
80  //G4LogicalVolume* l_dosi = new G4LogicalVolume(s_dosi, geometry::Materials::get("BGO"), "l_dosi", 0, m_sensitive);
81  G4LogicalVolume* l_dosi = new G4LogicalVolume(s_dosi, geometry::Materials::get("G4_SILICON_DIOXIDE"), "l_dosi", 0, m_sensitive);
82 
83  l_dosi->SetVisAttributes(green);
84  //Lets limit the Geant4 stepsize inside the volume
85  l_dosi->SetUserLimits(new G4UserLimits(stepSize));
86  double z_pos[100];
87  double r_pos[100];
88  int dim = 0;
89  for (double z : activeParams.getArray("z", {0})) {
90  z *= CLHEP::cm;
91  z_pos[dim] = z;
92  dim++;
93  }
94  dim = 0;
95  for (double r : activeParams.getArray("r", {0})) {
96  r *= CLHEP::cm;
97  r_pos[dim] = r + dz_dosi;
98  dim++;
99  }
100 
101  int detID = 0;
102  G4Transform3D transform;
103  for (double phi : activeParams.getArray("Phi", {M_PI / 2})) {
104  //phi *= CLHEP::deg;
105  for (int i = 0; i < dim; i++) {
106  transform = G4RotateZ3D(phi - M_PI / 2) * G4Translate3D(0, r_pos[i], z_pos[i]) * G4RotateX3D(-M_PI / 2 - thetaZ);
107  new G4PVPlacement(transform, l_dosi, TString::Format("p_dosi_%d", detID).Data() , &topVolume, false, detID);
108  detID++;
109  }
110  }
111  }
112  }
113  } // dosi namespace
115 } // Belle2 namespace
Belle2::dosi::DosiCreator::m_sensitive
SensitiveDetector * m_sensitive
SensitiveDetector DOSI.
Definition: DosiCreator.h:36
Belle2::gearbox::Interface::getAngle
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:301
Belle2::geometry::Materials::get
static G4Material * get(const std::string &name)
Find given material.
Definition: Materials.h:65
Belle2::dosi::DosiCreator::create
virtual void create(const GearDir &content, G4LogicalVolume &topVolume, geometry::GeometryTypes type)
Function to actually create the geometry, has to be overridden by derived classes.
Definition: DosiCreator.cc:59
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::GearDir
GearDir is the basic class used for accessing the parameter store.
Definition: GearDir.h:41
Belle2::gearbox::Interface::getLength
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:261
Belle2::gearbox::Interface::getArray
std::vector< double > getArray(const std::string &path) const noexcept(false)
Get the parameter path as a list of double values converted to the standard unit.
Definition: Interface.cc:133
Belle2::PXD::SensitiveDetector
VXD::SensitiveDetector< PXDSimHit, PXDTrueHit > SensitiveDetector
The PXD Sensitive Detector class.
Definition: SensitiveDetector.h:36
Belle2::dosi::DosiFactory
geometry::CreatorFactory< DosiCreator > DosiFactory("DOSICreator")
Creator creates the DOSI geometry.
Belle2::geometry::GeometryTypes
GeometryTypes
Flag indiciating the type of geometry to be used.
Definition: GeometryManager.h:39