LArSoft  v10_04_05
Liquid Argon Software toolkit - https://larsoft.org/
AuxDetReadoutGeometry.cxx
Go to the documentation of this file.
1 
8 
9 // Framework includes
11 
12 // LArSoft includes
16 
17 // G4 includes
18 #include "Geant4/G4LogicalVolume.hh"
19 #include "Geant4/G4Point3D.hh"
20 #include "Geant4/G4SDManager.hh"
21 
22 namespace larg4 {
23 
24  // Constructor and destructor.
26  G4String const name)
27  : G4VUserParallelWorld(name), fAuxDetGeom{auxDetGeom}, fNumSensitiveVol(0)
28  {}
29 
32  {
33  // We want to find all of the AuxDets that the Geometry service would find and make
34  // each one a G4 sensitive detector.
35 
36  // Call initial case of a function that will rucursively run through the volume tree
37  // down to max depth. Start at 0 depth with World, where the initial translation and
38  // rotation should be 0 as well
39  unsigned int MaxDepth = 8; // should be plenty
40  std::vector<const G4VPhysicalVolume*> path(MaxDepth);
42  G4Transform3D InitTransform(path[0]->GetObjectRotationValue(), path[0]->GetObjectTranslation());
43 
44  // first try to make sensitive volumes, if those are not in the gdml file (ie file was
45  // made before the introduction of sensitive volumes) fall back to looking for the aux
46  // dets
47  FindAndMakeAuxDetSensitive(path, 0, InitTransform);
48  if (fNumSensitiveVol < 1) FindAndMakeAuxDet(path, 0, InitTransform);
49  }
50 
51  //---------------------------------------------------------------
53  std::vector<const G4VPhysicalVolume*>& path,
54  unsigned int depth,
55  G4Transform3D DepthToWorld)
56  {
57  G4LogicalVolume* LogicalVolumeAtDepth = path[depth]->GetLogicalVolume();
58 
59  std::string volName(path[depth]->GetName());
60  if (volName.find("volAuxDet") != std::string::npos &&
61  volName.find("Sensitive") != std::string::npos) {
62 
63  // find world coordinate of the AuxDet origin in cm
64  G4Point3D local(0., 0., 0.);
65  G4Point3D world = DepthToWorld * local; // G4 works in mm
66  geo::Point_t const worldPos{
67  world.x() / CLHEP::cm, world.y() / CLHEP::cm, world.z() / CLHEP::cm};
68 
69  size_t adNum = 0;
70  size_t svNum = 0;
71  fAuxDetGeom->FindAuxDetSensitiveAtPosition(worldPos, adNum, svNum);
72  // N.B. This name is expected by code in LArG4:
73  std::string SDName = "AuxDetSD_AuxDet" + std::to_string(adNum) + "_" + std::to_string(svNum);
74  AuxDetReadout* adReadout = new larg4::AuxDetReadout(SDName, adNum, svNum);
75 
76  MF_LOG_DEBUG("AuxDetReadoutGeometry")
77  << "found" << path[depth]->GetName() << ", number " << adNum << ":" << svNum;
78 
79  // Tell Geant4's sensitive-detector manager about the AuxDetReadout class
80  (G4SDManager::GetSDMpointer())->AddNewDetector(adReadout);
81  LogicalVolumeAtDepth->SetSensitiveDetector(adReadout);
83  return;
84  }
85 
86  // Explore the next layer down -- unless it is a very deep geometry,
87  // recursion should end before exception is thrown.
88  unsigned int deeper = depth + 1;
89  if (deeper >= path.size()) {
90  throw cet::exception("AuxDetReadoutGeometry") << "exceeded maximum TGeoNode depth\n";
91  }
92 
93  // Note that there will be nd different branches off of path[depth]
94  G4int nd = LogicalVolumeAtDepth->GetNoDaughters();
95  for (int d = 0; d < nd; ++d) {
96 
97  // get the physvol daughter in the logicalvol
98  path[deeper] = LogicalVolumeAtDepth->GetDaughter(d);
99 
100  // keep track of the transform to world coordinates for PositionToAuxDet
101  G4Transform3D DeeperToMother(path[deeper]->GetObjectRotationValue(),
102  path[deeper]->GetObjectTranslation());
103  G4Transform3D DeeperToWorld = DepthToWorld * DeeperToMother;
104  FindAndMakeAuxDetSensitive(path, deeper, DeeperToWorld);
105  }
106  }
107 
108  //---------------------------------------------------------------
109  void AuxDetReadoutGeometry::FindAndMakeAuxDet(std::vector<const G4VPhysicalVolume*>& path,
110  unsigned int depth,
111  G4Transform3D DepthToWorld)
112  {
113  G4LogicalVolume* LogicalVolumeAtDepth = path[depth]->GetLogicalVolume();
114 
115  std::string volName(path[depth]->GetName());
116  if (volName.find("volAuxDet") != std::string::npos) {
117 
118  // find world coordinate of the AuxDet origin in cm
119  G4Point3D local(0., 0., 0.);
120  G4Point3D world = DepthToWorld * local; // G4 works in mm
121  geo::Point_t const worldPos{
122  world.x() / CLHEP::cm, world.y() / CLHEP::cm, world.z() / CLHEP::cm};
123 
124  auto const adNum = fAuxDetGeom->FindAuxDetAtPosition(worldPos);
125  // N.B. This name is expected by code in LArG4:
126  std::string SDName = "AuxDetSD_AuxDet" + std::to_string(adNum) + "_0";
127  AuxDetReadout* adReadout = new larg4::AuxDetReadout(SDName, adNum, 0);
128 
129  MF_LOG_DEBUG("AuxDetReadoutGeometry")
130  << "found" << path[depth]->GetName() << ", number " << adNum << ":0";
131 
132  // Tell Geant4's sensitive-detector manager about the AuxDetReadout class
133  (G4SDManager::GetSDMpointer())->AddNewDetector(adReadout);
134  LogicalVolumeAtDepth->SetSensitiveDetector(adReadout);
136  return;
137  }
138 
139  // Explore the next layer down -- unless it is a very deep geometry,
140  // recursion should end before exception is thrown.
141  unsigned int deeper = depth + 1;
142  if (deeper >= path.size()) {
143  throw cet::exception("AuxDetReadoutGeometry") << "exceeded maximum TGeoNode depth\n";
144  }
145 
146  // Note that there will be nd different branches off of path[depth]
147  G4int nd = LogicalVolumeAtDepth->GetNoDaughters();
148  for (int d = 0; d < nd; ++d) {
149 
150  // get the physvol daughter in the logicalvol
151  path[deeper] = LogicalVolumeAtDepth->GetDaughter(d);
152 
153  // keep track of the transform to world coordinates for PositionToAuxDet
154  G4Transform3D DeeperToMother(path[deeper]->GetObjectRotationValue(),
155  path[deeper]->GetObjectTranslation());
156  G4Transform3D DeeperToWorld = DepthToWorld * DeeperToMother;
157  FindAndMakeAuxDet(path, deeper, DeeperToWorld);
158  }
159  }
160 
161 } // namespace larg4
Build Geant4 geometry from GDML.
std::size_t FindAuxDetAtPosition(Point_t const &point, double tolerance=0) const
Returns the index of the auxiliary detector at specified location.
static G4VPhysicalVolume * GetWorld()
Geant4 interface.
geo::AuxDetGeometryCore const * fAuxDetGeom
Handle to the geometry.
void FindAndMakeAuxDet(std::vector< const G4VPhysicalVolume * > &path, unsigned int depth, G4Transform3D DepthToWorld)
AuxDetReadoutGeometry(geo::AuxDetGeometryCore const *auxDetGeom, G4String name="AuxDetReadoutGeometry")
Description of physical geometry of one set of auxiliary detectors.
uint32_t fNumSensitiveVol
number of sensitive volumes
Access the description of auxiliary detector geometry.
decltype(auto) constexpr to_string(T &&obj)
ADL-aware version of std::to_string.
Float_t d
Definition: plot.C:235
void FindAuxDetSensitiveAtPosition(Point_t const &point, std::size_t &adg, std::size_t &sv, double tolerance=0) const
Fills the indices of the sensitive auxiliary detector at location.
ROOT::Math::PositionVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Point_t
Type for representation of position in physical 3D space.
Definition: geo_vectors.h:180
Define the "parallel" geometry that&#39;s seen by the AuxDet.
#define MF_LOG_DEBUG(id)
A Geant4 sensitive detector that accumulates information.
void FindAndMakeAuxDetSensitive(std::vector< const G4VPhysicalVolume * > &path, unsigned int depth, G4Transform3D DepthToWorld)
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33