Belle II Software  release-05-01-25
SddCreator.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/srsensor/geometry/SddCreator.h>
12 #include <beast/srsensor/simulation/SensitiveDetector.h>
13 
14 #include <geometry/Materials.h>
15 #include <geometry/CreatorFactory.h>
16 #include <framework/gearbox/GearDir.h>
17 #include <framework/logging/Logger.h>
18 
19 #include <cmath>
20 #include <boost/format.hpp>
21 #include <boost/foreach.hpp>
22 #include <boost/algorithm/string.hpp>
23 
24 #include <G4LogicalVolume.hh>
25 #include <G4PVPlacement.hh>
26 
27 //Shapes
28 #include <G4UserLimits.hh>
29 #include "G4Tubs.hh"
30 
31 using namespace std;
32 using namespace boost;
33 
34 namespace Belle2 {
41  namespace srsensor {
42 
43  // Register the creator
45  geometry::CreatorFactory<SddCreator> SddFactory("SDDCreator");
46 
47  SddCreator::SddCreator(): m_sensitive(0)
48  {
50  }
51 
52  SddCreator::~SddCreator()
53  {
54  if (m_sensitive) delete m_sensitive;
55  }
56 
57  void SddCreator::create(const GearDir& content, G4LogicalVolume& topVolume, geometry::GeometryTypes /* type */)
58  {
59  //lets get the stepsize parameter with a default value of 5 µm
60  double stepSize = content.getLength("stepSize", 5 * CLHEP::um);
61 
62  //no get the array. Notice that the default framework unit is cm, so the
63  //values will be automatically converted
64  vector<double> bar_sdd = content.getArray("bar_sdd");
65  B2INFO("Contents of bar_sdd: ");
66  BOOST_FOREACH(double value, bar_sdd) {
67  B2INFO("value: " << value);
68  }
69  int sddNb = 0;
70  //Lets loop over all the Active nodes
71  BOOST_FOREACH(const GearDir & activeParams, content.getNodes("Active")) {
72 
73  //create sdd volume
74  G4double startAngle = 0.*CLHEP::deg;
75  G4double spanningAngle = 360.*CLHEP::deg;
76  G4Tubs* s_SDD = new G4Tubs("s_SDD",
77  activeParams.getLength("sdd_innerRadius")*CLHEP::cm,
78  activeParams.getLength("sdd_outerRadius")*CLHEP::cm,
79  activeParams.getLength("sdd_hz")*CLHEP::cm,
80  startAngle, spanningAngle);
81 
82  string matSDD = activeParams.getString("MaterialSDD");
83  G4LogicalVolume* l_SDD = new G4LogicalVolume(s_SDD, geometry::Materials::get(matSDD), "l_SDD", 0, m_sensitive);
84 
85  //Lets limit the Geant4 stepsize inside the volume
86  l_SDD->SetUserLimits(new G4UserLimits(stepSize));
87 
88  //position sdd volume
89  G4ThreeVector SDDpos = G4ThreeVector(
90  activeParams.getLength("x_sdd") * CLHEP::cm,
91  activeParams.getLength("y_sdd") * CLHEP::cm,
92  activeParams.getLength("z_sdd") * CLHEP::cm
93  );
94 
95  G4RotationMatrix* rot_sdd = new G4RotationMatrix();
96  rot_sdd->rotateX(activeParams.getAngle("AngleX"));
97  rot_sdd->rotateY(activeParams.getAngle("AngleY"));
98  rot_sdd->rotateZ(activeParams.getAngle("AngleZ"));
99  //geometry::setColor(*l_SDD, "#006699");
100 
101  new G4PVPlacement(rot_sdd, SDDpos, l_SDD, "p_SDD", &topVolume, false, sddNb);
102 
103  sddNb++;
104  }
105  }
106  } // srsensor namespace
108 } // Belle2 namespace
Belle2::srsensor::SddFactory
geometry::CreatorFactory< SddCreator > SddFactory("SDDCreator")
Creator creates the SDD geometry.
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::srsensor::SddCreator::m_sensitive
SensitiveDetector * m_sensitive
SensitiveDetector SDD.
Definition: SddCreator.h:36
Belle2::geometry::Materials::get
static G4Material * get(const std::string &name)
Find given material.
Definition: Materials.h:65
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::srsensor::SddCreator::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: SddCreator.cc:57
Belle2::GearDir::getString
virtual std::string getString(const std::string &path="") const noexcept(false) override
Get the parameter path as a string.
Definition: GearDir.h:79
Belle2::srsensor::Fei4Creator::m_sensitive
SensitiveDetector * m_sensitive
SensitiveDetector FEI4.
Definition: Fei4Creator.h:36
Belle2::PXD::SensitiveDetector
VXD::SensitiveDetector< PXDSimHit, PXDTrueHit > SensitiveDetector
The PXD Sensitive Detector class.
Definition: SensitiveDetector.h:36
Belle2::geometry::GeometryTypes
GeometryTypes
Flag indiciating the type of geometry to be used.
Definition: GeometryManager.h:39