LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
OpDetReadoutGeometry.cxx
Go to the documentation of this file.
1 //
7 
8 #include "CLHEP/Geometry/Transform3D.h"
9 #include "CLHEP/Vector/Rotation.h"
10 #include "CLHEP/Vector/ThreeVector.h"
11 #include "Geant4/G4LogicalVolume.hh"
12 #include "Geant4/G4PVPlacement.hh"
13 #include "Geant4/G4RotationMatrix.hh"
14 #include "Geant4/G4ThreeVector.hh"
15 #include "Geant4/G4Types.hh"
16 #include "Geant4/G4VPhysicalVolume.hh"
17 
23 
24 #include <sstream>
25 
27 #include "cetlib_except/exception.h"
29 
30 namespace larg4 {
31 
33  const G4String name,
34  bool useLitePhotons /* = false */)
35  : G4VUserParallelWorld(name), fUseLitePhotons(useLitePhotons)
36  {
37  fOpDetSensitiveName = OpDetSensitiveName;
38  }
39 
41 
43  {
44  mf::LogInfo("OpDetReadoutGeometry")
45  << "constructing parallel world, looking for " << fOpDetSensitiveName;
46 
47  // Get an empty parallel world volume
48  G4VPhysicalVolume* ParallelWorld = GetWorld();
49 
50  // Start with empty vectors
51  std::vector<G4LogicalVolume*> OpDetVolumes;
52  std::vector<G4Transform3D> OpDetTransformations;
53 
54  // Get the primary world volume
55  G4VPhysicalVolume* WorldPhysical = g4b::DetectorConstruction::GetWorld();
56 
57  // Find the OpDet volumes
58  std::vector<G4Transform3D> EmptyVector;
59  EmptyVector.clear();
61  WorldPhysical, fOpDetSensitiveName, EmptyVector, OpDetVolumes, OpDetTransformations);
62 
63  fOpDetVolumes = OpDetVolumes;
64  fOpDetTransformations = OpDetTransformations;
65 
66  // Get the OpDet Lookup Table
68 
69  // Create sensitive detector
70  OpDetSensitiveDetector* TheSD =
71  new OpDetSensitiveDetector("OpDetSensitiveDetector", fUseLitePhotons);
72 
73  if (OpDetVolumes.size() > 0) {
74 
75  // Make placements
76  for (unsigned int i = 0; i != OpDetVolumes.size(); i++) {
77  std::stringstream VolumeName;
78  VolumeName.flush();
79  VolumeName.str("OpDetVolume_");
80  VolumeName << i;
81 
82  G4Transform3D TheTransform = OpDetTransformations.at(i);
83 
84  G4VSolid* TheSolid = OpDetVolumes.at(i)->GetSolid();
85  G4Material* TheMaterial = OpDetVolumes.at(i)->GetMaterial();
86  G4LogicalVolume* TheLogVolume =
87  new G4LogicalVolume(TheSolid, TheMaterial, VolumeName.str().c_str());
88 
89  TheLogVolume->SetSensitiveDetector(TheSD);
90 
91  G4PVPlacement* ThePlacement = new G4PVPlacement(
92  TheTransform, VolumeName.str().c_str(), TheLogVolume, ParallelWorld, false, 0);
93 
94  // CLHEP::Hep3Vector trans = ThePlacement->GetTranslation();
95 
96  TheOpDetLookup->AddPhysicalVolume(ThePlacement);
97  }
98  }
99 
100  // Now add any optically parameterized volumes
101 
102  std::vector<G4LogicalVolume*> OpParamVolumesFound;
103  std::vector<G4Transform3D> OpParamTransformationsFound;
104 
106 
107  std::vector<std::string> OpParamModels = lgp->OpticalParamModels();
108  std::vector<std::string> OpParamVolumes = lgp->OpticalParamVolumes();
109  std::vector<int> OpParamOrientations = lgp->OpticalParamOrientations();
110  ;
111  std::vector<std::vector<std::vector<double>>> OpParamParameters = lgp->OpticalParamParameters();
112 
113  if ((OpParamModels.size() != OpParamVolumes.size()) ||
114  (OpParamModels.size() != OpParamOrientations.size()) ||
115  (OpParamModels.size() != OpParamParameters.size())) {
116  throw cet::exception("OpDetReadoutGeometry")
117  << "sizes of OpParam specification vectors do not match\n";
118  }
119 
120  for (size_t imodel = 0; imodel != OpParamVolumes.size(); ++imodel) {
121  EmptyVector.clear();
122  FindVolumes(WorldPhysical,
123  OpParamVolumes.at(imodel),
124  EmptyVector,
125  OpParamVolumesFound,
126  OpParamTransformationsFound);
127  mf::LogInfo("OpDetReadoutGeometry")
128  << "Found " << OpParamVolumesFound.size() << " volumes of name "
129  << OpParamVolumes.at(imodel) << " to attach optical parameterization "
130  << OpParamModels.at(imodel) << std::endl;
131 
132  // Since the same named model may be instantiated more than once with
133  // different parameters, create a unique sensitive detector name for this
134  // instance
135 
136  std::stringstream SDName("");
137  SDName << OpParamModels.at(imodel) << "_" << imodel;
138 
139  OpParamSD* ParamSD = new OpParamSD(SDName.str().c_str(),
140  OpParamModels.at(imodel),
141  OpParamOrientations.at(imodel),
142  OpParamParameters.at(imodel));
143 
144  if (OpParamVolumesFound.size() > 0) {
145  for (unsigned int ivol = 0; ivol != OpParamVolumes.size(); ivol++) {
146  std::stringstream VolumeName;
147  VolumeName.flush();
148  VolumeName.str("OpParamVolume_");
149  VolumeName << ivol;
150 
151  G4Transform3D TheTransform = OpParamTransformationsFound.at(ivol);
152 
153  G4VSolid* TheSolid = OpParamVolumesFound.at(ivol)->GetSolid();
154  G4Material* TheMaterial = OpParamVolumesFound.at(ivol)->GetMaterial();
155  G4LogicalVolume* TheLogVolume =
156  new G4LogicalVolume(TheSolid, TheMaterial, VolumeName.str().c_str());
157 
158  G4PVPlacement* ThePlacement = new G4PVPlacement(
159  TheTransform, VolumeName.str().c_str(), TheLogVolume, ParallelWorld, false, 0);
160 
161  TheLogVolume->SetSensitiveDetector(ParamSD);
162 
163  // This line just suppressed a compiler warning
164  ThePlacement->GetTranslation();
165  }
166  }
167  }
168  }
169 
170  void OpDetReadoutGeometry::FindVolumes(G4VPhysicalVolume* PhysicalVolume,
171  G4String OpDetName,
172  std::vector<G4Transform3D> TransformSoFar,
173  std::vector<G4LogicalVolume*>& OpDetVolumes,
174  std::vector<G4Transform3D>& OpDetTransformations)
175  {
176 
177  // Add the next layer of transformation to the vector
178  G4ThreeVector Translation = PhysicalVolume->GetObjectTranslation();
179  G4RotationMatrix Rotation = PhysicalVolume->GetObjectRotationValue();
180  G4Transform3D NextTransform(Rotation, Translation);
181 
182  TransformSoFar.push_back(NextTransform);
183 
184  // Check if this volume is a OpDet
185  G4String OpDetNameUnderscore = OpDetName + "_";
186  G4String VolumeName = PhysicalVolume->GetName();
187  if ((VolumeName == OpDetName) ||
188  (VolumeName.find(OpDetNameUnderscore, 0, OpDetNameUnderscore.length()) == 0)) {
189 
190  // We found a OpDet! Store its volume and global transformation
191  G4ThreeVector Trans(0, 0, 0);
192  G4RotationMatrix Rot(0, 0, 0);
193  G4Transform3D TotalTransform(Rot, Trans);
194  //for ( std::vector<G4Transform3D>::reverse_iterator it = TransformSoFar.rbegin();
195  // it != TransformSoFar.rend(); ++it )
196 
197  for (std::vector<G4Transform3D>::iterator it = TransformSoFar.begin();
198  it != TransformSoFar.end();
199  ++it) {
200  CLHEP::Hep3Vector trans = (*it).getTranslation();
201  CLHEP::HepRotation rot = (*it).getRotation();
202  TotalTransform = TotalTransform * (*it);
203  }
204 
205  OpDetVolumes.push_back(PhysicalVolume->GetLogicalVolume());
206  OpDetTransformations.push_back(TotalTransform);
207  }
208  else {
209  // We did not find a OpDet. Keep looking through daughters.
210  G4LogicalVolume* LogicalVolume = PhysicalVolume->GetLogicalVolume();
211 
212  // Loop through the daughters of the volume
213  G4int NumberDaughters = LogicalVolume->GetNoDaughters();
214  for (G4int i = 0; i != NumberDaughters; ++i) {
215  // Get the ith daughter volume
216  G4VPhysicalVolume* Daughter = LogicalVolume->GetDaughter(i);
217 
218  // Recursively step into this volume
219  FindVolumes(Daughter, OpDetName, TransformSoFar, OpDetVolumes, OpDetTransformations);
220  }
221  }
222  }
223 
224 }
void FindVolumes(G4VPhysicalVolume *, G4String, std::vector< G4Transform3D >, std::vector< G4LogicalVolume * > &, std::vector< G4Transform3D > &)
Build Geant4 geometry from GDML.
intermediate_table::iterator iterator
Store parameters for running LArG4.
const std::vector< std::vector< std::vector< double > > > & OpticalParamParameters() const
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
static G4VPhysicalVolume * GetWorld()
Geant4 interface.
std::vector< G4Transform3D > fOpDetTransformations
const std::vector< std::string > & OpticalParamVolumes() const
OpDetReadoutGeometry(G4String OpDetSensitiveName, const G4String name="OpDetReadoutGeometry", bool useLitePhotons=false)
bool const fUseLitePhotons
Pass-through option for sensitive detector.
std::vector< G4LogicalVolume * > fOpDetVolumes
const std::vector< std::string > & OpticalParamModels() const
const std::vector< int > & OpticalParamOrientations() const
void AddPhysicalVolume(G4VPhysicalVolume *)
Definition: OpDetLookup.cxx:83
static OpDetLookup * Instance()
Definition: OpDetLookup.cxx:31
OpDetLookup * TheOpDetLookup
Definition: OpDetLookup.cxx:22
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33