ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4STRead.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4STRead.cc
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // class G4STRead Implementation
28 //
29 // History:
30 // - Created. Zoltan Torzsok, November 2007
31 // -------------------------------------------------------------------------
32 
33 #include <fstream>
34 
35 #include "G4STRead.hh"
36 
37 #include "G4Material.hh"
38 #include "G4Box.hh"
39 #include "G4QuadrangularFacet.hh"
40 #include "G4TriangularFacet.hh"
41 #include "G4TessellatedSolid.hh"
42 #include "G4LogicalVolume.hh"
43 #include "G4PVPlacement.hh"
44 #include "G4AffineTransform.hh"
45 #include "G4VoxelLimits.hh"
46 
47 void G4STRead::TessellatedRead(const std::string& line)
48 {
49  if (tessellatedList.size()>0)
50  {
51  tessellatedList.back()->SetSolidClosed(true);
52  // Finish the previous solid at first!
53  }
54 
55  std::istringstream stream(line.substr(2));
56 
57  G4String name;
58  stream >> name;
59 
60  G4TessellatedSolid* tessellated = new G4TessellatedSolid(name);
61  volumeMap[tessellated] =
62  new G4LogicalVolume(tessellated, solid_material, name+"_LV" , 0, 0, 0);
63  tessellatedList.push_back(tessellated);
64 
65 #ifdef G4VERBOSE
66  G4cout << "G4STRead: Reading solid: " << name << G4endl;
67 #endif
68 }
69 
70 void G4STRead::FacetRead(const std::string& line)
71 {
72  if (tessellatedList.size()==0)
73  {
74  G4Exception("G4STRead::FacetRead()", "ReadError", FatalException,
75  "A solid must be defined before defining a facet!");
76  }
77 
78  if (line[2]=='3') // Triangular facet
79  {
80  G4double x1,y1,z1;
81  G4double x2,y2,z2;
82  G4double x3,y3,z3;
83 
84  std::istringstream stream(line.substr(4));
85  stream >> x1 >> y1 >> z1 >> x2 >> y2 >> z2 >> x3 >> y3 >> z3;
86  tessellatedList.back()->
87  AddFacet(new G4TriangularFacet(G4ThreeVector(x1,y1,z1),
88  G4ThreeVector(x2,y2,z2),
89  G4ThreeVector(x3,y3,z3), ABSOLUTE));
90  }
91  else if (line[2]=='4') // Quadrangular facet
92  {
93  G4double x1,y1,z1;
94  G4double x2,y2,z2;
95  G4double x3,y3,z3;
96  G4double x4,y4,z4;
97 
98  std::istringstream stream(line.substr(4));
99  stream >> x1 >> y1 >> z1 >> x2 >> y2 >> z2
100  >> x3 >> y3 >> z3 >> x4 >> y4 >> z4;
101  tessellatedList.back()->
102  AddFacet(new G4QuadrangularFacet(G4ThreeVector(x1,y1,z1),
103  G4ThreeVector(x2,y2,z2),
104  G4ThreeVector(x3,y3,z3),
105  G4ThreeVector(x4,y4,z4), ABSOLUTE));
106  }
107  else
108  {
109  G4Exception("G4STRead::FacetRead()", "ReadError", FatalException,
110  "Number of vertices per facet should be either 3 or 4!");
111  }
112 }
113 
114 void G4STRead::PhysvolRead(const std::string& line)
115 {
116  G4int level;
117  G4String name;
118  G4double r1,r2,r3;
119  G4double r4,r5,r6;
120  G4double r7,r8,r9;
121  G4double pX,pY,pZ;
122  G4double n1,n2,n3,n4,n5;
123 
124  std::istringstream stream(line.substr(2));
125  stream >> level >> name >> r1 >> r2 >> r3 >> n1 >> r4 >> r5 >> r6
126  >> n2 >> r7 >> r8 >> r9 >> n3 >> pX >> pY >> pZ >> n4 >> n5;
127  std::string::size_type idx = name.rfind("_");
128  if (idx!=std::string::npos)
129  {
130  name.resize(idx);
131  }
132  else
133  {
134  G4Exception("G4STRead::PhysvolRead()", "ReadError",
135  FatalException, "Invalid input stream!");
136  return;
137  }
138 
139  G4cout << "G4STRead: Placing tessellated solid: " << name << G4endl;
140 
141  G4TessellatedSolid* tessellated = 0;
142 
143  for (size_t i=0; i<tessellatedList.size(); i++)
144  { // Find the volume for this physvol!
145  if (tessellatedList[i]->GetName() == G4String(name))
146  {
147  tessellated = tessellatedList[i];
148  break;
149  }
150  }
151 
152  if (tessellated == 0)
153  {
154  G4String error_msg = "Referenced solid '" + name + "' not found!";
155  G4Exception("G4STRead::PhysvolRead()", "ReadError",
156  FatalException, error_msg);
157  }
158  if (volumeMap.find(tessellated) == volumeMap.end())
159  {
160  G4String error_msg = "Referenced solid '" + name
161  + "' is not associated with a logical volume!";
162  G4Exception("G4STRead::PhysvolRead()", "InvalidSetup",
163  FatalException, error_msg);
164  }
165  const G4RotationMatrix rot(G4ThreeVector(r1,r2,r3),
166  G4ThreeVector(r4,r5,r6),
167  G4ThreeVector(r7,r8,r9));
168  const G4ThreeVector pos(pX,pY,pZ);
169 
171  volumeMap[tessellated], name+"_PV", world_volume, 0, 0);
172  // Note: INVERSE of rotation is needed!!!
173 
174  G4double minx,miny,minz;
175  G4double maxx,maxy,maxz;
176  const G4VoxelLimits limits;
177 
178  tessellated->CalculateExtent(kXAxis,limits,
179  G4AffineTransform(rot,pos),minx,maxx);
180  tessellated->CalculateExtent(kYAxis,limits,
181  G4AffineTransform(rot,pos),miny,maxy);
182  tessellated->CalculateExtent(kZAxis,limits,
183  G4AffineTransform(rot,pos),minz,maxz);
184 
185  if (world_extent.x() < std::fabs(minx))
186  { world_extent.setX(std::fabs(minx)); }
187  if (world_extent.y() < std::fabs(miny))
188  { world_extent.setY(std::fabs(miny)); }
189  if (world_extent.z() < std::fabs(minz))
190  { world_extent.setZ(std::fabs(minz)); }
191  if (world_extent.x() < std::fabs(maxx))
192  { world_extent.setX(std::fabs(maxx)); }
193  if (world_extent.y() < std::fabs(maxy))
194  { world_extent.setY(std::fabs(maxy)); }
195  if (world_extent.z() < std::fabs(maxz))
196  { world_extent.setZ(std::fabs(maxz)); }
197 }
198 
200 {
201 #ifdef G4VERBOSE
202  G4cout << "G4STRead: Reading '" << name << "'..." << G4endl;
203 #endif
204  std::ifstream GeomFile(name);
205 
206  if (!GeomFile)
207  {
208  G4String error_msg = "Cannot open file: " + name;
209  G4Exception("G4STRead::ReadGeom()", "ReadError",
210  FatalException, error_msg);
211  }
212 
213  tessellatedList.clear();
214  volumeMap.clear();
215  std::string line;
216 
217  while (getline(GeomFile,line))
218  {
219  if (line[0] == 'f') { TessellatedRead(line); } else
220  if (line[0] == 'p') { FacetRead(line); }
221  }
222 
223  if (tessellatedList.size()>0) // Finish the last solid!
224  {
225  tessellatedList.back()->SetSolidClosed(true);
226  }
227 
228  G4cout << "G4STRead: Reading '" << name << "' done." << G4endl;
229 }
230 
232 {
233 #ifdef G4VERBOSE
234  G4cout << "G4STRead: Reading '" << name << "'..." << G4endl;
235 #endif
236  std::ifstream TreeFile(name);
237 
238  if (!TreeFile)
239  {
240  G4String error_msg = "Cannot open file: " + name;
241  G4Exception("G4STRead::ReadTree()", "ReadError",
242  FatalException, error_msg);
243  }
244 
245  std::string line;
246 
247  while (getline(TreeFile,line))
248  {
249  if (line[0] == 'g') { PhysvolRead(line); }
250  }
251 
252  G4cout << "G4STRead: Reading '" << name << "' done." << G4endl;
253 }
254 
256 G4STRead::Read(const G4String& name, G4Material* mediumMaterial,
257  G4Material* solidMaterial)
258 {
259  if (mediumMaterial == 0)
260  {
261  G4Exception("G4STRead::Read()", "InvalidSetup", FatalException,
262  "Pointer to medium material is not valid!");
263  }
264  if (solidMaterial == 0)
265  {
266  G4Exception("G4STRead::Read()", "InvalidSetup", FatalException,
267  "Pointer to solid material is not valid!");
268  }
269 
270  solid_material = solidMaterial;
271 
272  world_box = new G4Box("TessellatedWorldBox",kInfinity,kInfinity,kInfinity);
273  // We don't know the extent of the world yet!
274  world_volume = new G4LogicalVolume(world_box, mediumMaterial,
275  "TessellatedWorldLV", 0, 0, 0);
276  world_extent = G4ThreeVector(0,0,0);
277 
278  ReadGeom(name+".geom");
279  ReadTree(name+".tree");
280 
281  // Now setting the world extent ...
282  //
289 
290  return world_volume;
291 }