source: trunk/source/persistency/ascii/src/G4tgbVolume.cc @ 1228

Last change on this file since 1228 was 1228, checked in by garnier, 14 years ago

update geant4.9.3 tag
File size: 50.5 KB
Line 
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//
28// class G4tgbVolume
29
30// History:
31// - Created.                                 P.Arce, CIEMAT (November 2007)
32// -------------------------------------------------------------------------
33
34#include "G4tgbVolume.hh"
35#include "G4tgbVolumeMgr.hh"
36#include "G4tgbMaterialMgr.hh"
37#include "G4tgbRotationMatrixMgr.hh"
38#include "G4tgbPlaceParamLinear.hh"
39#include "G4tgbPlaceParamSquare.hh"
40#include "G4tgbPlaceParamCircle.hh"
41
42#include "G4tgrSolid.hh"
43#include "G4tgrSolidBoolean.hh"
44#include "G4tgrVolume.hh"
45#include "G4tgrVolumeDivision.hh"
46#include "G4tgrVolumeAssembly.hh"
47#include "G4tgrVolumeMgr.hh"
48#include "G4tgrPlace.hh"
49#include "G4tgrPlaceSimple.hh"
50#include "G4tgrPlaceDivRep.hh"
51#include "G4tgrPlaceParameterisation.hh"
52#include "G4tgrUtils.hh"
53
54#include "G4VSolid.hh"
55#include "G4UnionSolid.hh"
56#include "G4SubtractionSolid.hh"
57#include "G4IntersectionSolid.hh"
58#include "G4LogicalVolume.hh"
59#include "G4VPhysicalVolume.hh"
60#include "G4PVPlacement.hh"
61#include "G4PVDivision.hh"
62#include "G4PVReplica.hh"
63#include "G4PVParameterised.hh"
64#include "G4Box.hh"
65#include "G4Tubs.hh"
66#include "G4Cons.hh"
67#include "G4Trap.hh"
68#include "G4Sphere.hh"
69#include "G4Orb.hh"
70#include "G4Trd.hh"
71#include "G4Para.hh"
72#include "G4Torus.hh"
73#include "G4Hype.hh"
74#include "G4Polycone.hh"
75#include "G4Polyhedra.hh"
76#include "G4EllipticalTube.hh"
77#include "G4Ellipsoid.hh"
78#include "G4EllipticalCone.hh"
79#include "G4Hype.hh"
80#include "G4Tet.hh"
81#include "G4TwistedBox.hh"
82#include "G4TwistedTrap.hh"
83#include "G4TwistedTrd.hh"
84#include "G4TwistedTubs.hh"
85#include "G4AssemblyVolume.hh"
86#include "G4BREPSolidBox.hh"
87#include "G4BREPSolidCylinder.hh"
88#include "G4BREPSolidCone.hh"
89#include "G4BREPSolidSphere.hh"
90#include "G4BREPSolidTorus.hh"
91#include "G4BREPSolidPCone.hh"
92#include "G4BREPSolidPolyhedra.hh"
93#include "G4BREPSolidOpenPCone.hh"
94#include "G4TessellatedSolid.hh"
95#include "G4TriangularFacet.hh"
96#include "G4QuadrangularFacet.hh"
97#include "G4ExtrudedSolid.hh"
98
99#include "G4VisExtent.hh"
100#include "G4Material.hh"
101#include "G4RotationMatrix.hh"
102#include "G4ReflectionFactory.hh"
103
104#include "G4VisAttributes.hh"
105#include "G4RegionStore.hh"
106#include "G4tgrMessenger.hh"
107#include "G4UIcommand.hh"
108#include "G4GeometryTolerance.hh"
109
110//-------------------------------------------------------------------
111G4tgbVolume::G4tgbVolume()
112{
113}
114
115
116//-------------------------------------------------------------------
117G4tgbVolume::~G4tgbVolume()
118{
119}
120
121
122//-------------------------------------------------------------------
123G4tgbVolume::G4tgbVolume( G4tgrVolume* vol)
124{
125  theTgrVolume = vol;
126  theG4AssemblyVolume = 0;
127}
128
129
130//-------------------------------------------------------------------
131void G4tgbVolume::ConstructG4Volumes( const G4tgrPlace* place,
132                                      const G4LogicalVolume* parentLV )
133{
134#ifdef G4VERBOSE
135  if( G4tgrMessenger::GetVerboseLevel() >= 2 )
136  {
137    G4cout << G4endl <<  "@@@ G4tgbVolume::ConstructG4Volumes - " << GetName() << G4endl;
138    if( place && parentLV ) G4cout << "   place in LV " << parentLV->GetName() << G4endl;
139  }
140#endif
141  G4tgbVolumeMgr* g4vmgr = G4tgbVolumeMgr::GetInstance();
142  G4LogicalVolume* logvol = g4vmgr->FindG4LogVol( GetName() );
143  G4bool bFirstCopy = false;
144  if( (logvol == 0) ) 
145  {
146    bFirstCopy = true;
147    if( theTgrVolume->GetType() != "VOLDivision" )
148    {
149      //--- If first time build solid and LogVol
150      G4VSolid* solid = FindOrConstructG4Solid( theTgrVolume->GetSolid() ); 
151      if( solid != 0 )   // for G4AssemblyVolume it is 0
152      {
153        g4vmgr->RegisterMe( solid );
154        logvol = ConstructG4LogVol( solid );
155        g4vmgr->RegisterMe( logvol );
156        g4vmgr->RegisterChildParentLVs( logvol, parentLV ); 
157       
158      }
159    } 
160  } 
161  //--- Construct PhysVol
162  G4VPhysicalVolume* physvol = ConstructG4PhysVol( place, logvol, parentLV );
163  if( physvol != 0 )  // 0 for G4AssemblyVolumes
164  {
165    g4vmgr->RegisterMe( physvol );
166
167    if( logvol == 0 ) // case of divisions
168    {
169      logvol = physvol->GetLogicalVolume();
170    }
171  }
172  else 
173  {
174    return;
175  }
176
177  //--- If first copy build children placements in this LogVol
178  if(bFirstCopy)
179  {
180    std::pair<G4mmapspl::iterator, G4mmapspl::iterator> children
181      = G4tgrVolumeMgr::GetInstance()->GetChildren( GetName() );
182    G4mmapspl::iterator cite; 
183    for( cite = children.first; cite != children.second; cite++ )
184    {
185      //----- Call G4tgrPlace ->constructG4Volumes
186      //---- find G4tgbVolume corresponding to the G4tgrVolume
187      //     pointed by G4tgrPlace
188      G4tgrPlace* pl = const_cast<G4tgrPlace*>((*cite).second);
189      G4tgbVolume* svol = g4vmgr->FindVolume( pl->GetVolume()->GetName() );
190      //--- find copyNo
191#ifdef G4VERBOSE
192      if( G4tgrMessenger::GetVerboseLevel() >= 2 )
193      {
194        G4cout << " G4tgbVolume::ConstructG4Volumes - construct daughter " <<  pl->GetVolume()->GetName() << " # " << pl->GetCopyNo() << G4endl;
195      }
196#endif
197      svol->ConstructG4Volumes( pl, logvol );
198    }
199  }
200
201}
202
203
204
205//-------------------------------------------------------------------
206G4VSolid* G4tgbVolume::FindOrConstructG4Solid( const G4tgrSolid* sol ) 
207{
208  G4double angularTolerance = G4GeometryTolerance::GetInstance()
209                            ->GetAngularTolerance();
210
211  if( sol == 0 ) { return 0; }
212
213#ifdef G4VERBOSE
214  if( G4tgrMessenger::GetVerboseLevel() >= 2 )
215  {
216    G4cout << " G4tgbVolume::FindOrConstructG4Solid():" << G4endl
217           << "   SOLID = " << sol << G4endl
218           << "   " << sol->GetName() << " of type " << sol->GetType()
219           << G4endl;
220  }
221#endif
222
223  //----- Check if solid exists already
224  G4VSolid* solid = G4tgbVolumeMgr::GetInstance()
225                    ->FindG4Solid( sol->GetName() );
226  if( solid ) { return solid; }
227
228  // Give 'sol' as Boolean solids needs to call this method twice
229
230#ifdef G4VERBOSE
231  if( G4tgrMessenger::GetVerboseLevel() >= 2 )
232  {
233    G4cout << " G4tgbVolume::FindOrConstructG4Solid() - "
234           << sol->GetSolidParams().size() << G4endl;
235  }
236#endif
237 
238  std::vector<G4double> solParam;
239
240  // In case of BOOLEAN solids, solidParams are taken from components
241
242  if( sol->GetSolidParams().size() == 1)
243  { 
244    solParam = * sol->GetSolidParams()[ 0 ];
245  }
246
247  //----------- instantiate the appropiate G4VSolid type
248  G4String stype = sol->GetType();
249  G4String sname = sol->GetName();
250
251  if( stype == "BOX" )
252  {
253    CheckNoSolidParams( stype, 3, solParam.size() );
254    solid = new G4Box( sname, solParam[0], solParam[1], solParam[2] ); 
255
256  }
257  else if( stype == "TUBE" )
258  {
259    CheckNoSolidParams( stype, 3, solParam.size() );
260    solid = new G4Tubs( sname, solParam[0], solParam[1], solParam[2],
261                        0.*deg, 360.*deg );
262  }
263  else if( stype == "TUBS" )
264  {
265    CheckNoSolidParams( stype, 5, solParam.size() );
266    G4double phiDelta = solParam[4];
267    if( std::fabs(phiDelta - twopi) < angularTolerance ) { phiDelta = twopi; }
268    solid = new G4Tubs( sname, solParam[0], solParam[1], solParam[2],
269                        solParam[3], phiDelta );
270  }
271  else if( stype == "TRAP" )
272  {
273    if( solParam.size() == 11 )
274    {
275      solid = new G4Trap( sname, solParam[0], solParam[1], solParam[2],
276                          solParam[3], solParam[4], solParam[5], solParam[6],
277                          solParam[7], solParam[8], solParam[9], solParam[10] );
278    }
279    else if( solParam.size() == 4 )
280    {
281      solid = new G4Trap( sname, solParam[0], solParam[1]/deg,
282                          solParam[2]/deg, solParam[3]);
283    }
284    else
285    {
286      G4String ErrMessage1 = "Solid type " + stype;
287      G4String ErrMessage2 = " should have 11 or 4 parameters,\n";
288      G4String ErrMessage3 = "and it has "
289                         + G4UIcommand::ConvertToString(G4int(solParam.size()));
290      G4String ErrMessage = ErrMessage1 + ErrMessage2 + ErrMessage3 + " !";
291      G4Exception("G4tgbVolume::FindOrConstructG4Solid()",
292                  "InvalidSetup", FatalException, ErrMessage);
293    }
294   
295  }
296  else if( stype == "TRD" )
297  {
298    CheckNoSolidParams( stype, 5, solParam.size() );
299    solid = new G4Trd( sname, solParam[0], solParam[1], solParam[2],
300                       solParam[3], solParam[4] );
301  }
302  else if( stype == "PARA" )
303  {
304    CheckNoSolidParams( stype, 6, solParam.size() );
305    solid = new G4Para( sname, solParam[0], solParam[1], solParam[2],
306                        solParam[3], solParam[4], solParam[5] );
307  }
308  else if( stype == "CONE" )
309  {
310    CheckNoSolidParams( stype, 5, solParam.size() );
311    solid = new G4Cons( sname, solParam[0], solParam[1], solParam[2],
312                        solParam[3], solParam[4], 0., 360.*deg);
313  }
314  else if( stype == "CONS" )
315  {
316    CheckNoSolidParams( stype, 7, solParam.size() );
317    G4double phiDelta = solParam[6];
318    if( std::fabs(phiDelta - twopi) < angularTolerance ) { phiDelta = twopi; }
319    solid = new G4Cons( sname, solParam[0], solParam[1], solParam[2],
320                        solParam[3], solParam[4], solParam[5], phiDelta);
321  }
322  else if( stype == "SPHERE" )
323  {
324    CheckNoSolidParams( stype, 6, solParam.size() );
325    G4double phiDelta = solParam[3];
326    if( std::fabs(phiDelta - twopi) < angularTolerance ) { phiDelta = twopi; }
327    G4double thetaDelta = solParam[5];
328    if( std::fabs(thetaDelta - pi) < angularTolerance ) { thetaDelta = pi; }
329    solid = new G4Sphere( sname, solParam[0], solParam[1], solParam[2],
330                          phiDelta, solParam[4], thetaDelta);
331  }
332  else if( stype == "ORB" )
333  {
334    CheckNoSolidParams( stype, 1, solParam.size() );
335    solid = new G4Orb( sname, solParam[0] );
336  }
337  else if( stype == "TORUS" )
338  {
339    CheckNoSolidParams( stype, 5, solParam.size() );
340    G4double phiDelta = solParam[4];
341    if( std::fabs(phiDelta - twopi) < angularTolerance ) { phiDelta = twopi; }
342    solid = new G4Torus( sname, solParam[0], solParam[1], solParam[2],
343                         solParam[3], phiDelta );
344  }
345  else if( stype == "POLYCONE" )
346  {
347    size_t nplanes = size_t(solParam[2]);
348    G4bool genericPoly = false;
349    if( solParam.size() == 3+nplanes*3 )
350    { 
351      genericPoly = true;
352    }
353    else if( solParam.size() == 3+nplanes*2 )
354    { 
355      genericPoly = false;
356    }
357    else
358    {
359      G4String Err1 = "Solid type " + stype + " should have ";
360      G4String Err2 = G4UIcommand::ConvertToString(G4int(3+nplanes*3))
361                    + " (Z,Rmin,Rmax)\n";
362      G4String Err3 = "or " + G4UIcommand::ConvertToString(G4int(3+nplanes*2));
363      G4String Err4 = " (RZ corners) parameters,\n";
364      G4String Err5 = "and it has "
365                    +  G4UIcommand::ConvertToString(G4int(solParam.size()));
366      G4String ErrMessage = Err1 + Err2 + Err3 + Err4 + Err5 + " !";
367      G4Exception("G4tgbVolume::FindOrConstructG4Solid()",
368                  "InvalidSetup", FatalException, ErrMessage);
369    }
370
371    if( genericPoly )
372    {
373      std::vector<G4double>* z_p = new std::vector<G4double>;
374      std::vector<G4double>* rmin_p = new std::vector<G4double>;
375      std::vector<G4double>* rmax_p = new std::vector<G4double>;
376      for( size_t ii = 0; ii < nplanes; ii++ )
377      {
378        (*z_p).push_back( solParam[3+3*ii] );
379        (*rmin_p).push_back( solParam[3+3*ii+1] );
380        (*rmax_p).push_back(  solParam[3+3*ii+2] );
381      }
382      G4double phiTotal = solParam[1];
383      if( std::fabs(phiTotal - twopi) < angularTolerance ) { phiTotal = twopi; }
384      solid = new G4Polycone( sname, solParam[0], phiTotal, // start,delta-phi
385                              nplanes, // sections
386                              &((*z_p)[0]), &((*rmin_p)[0]), &((*rmax_p)[0]));
387    }
388    else
389    {
390      std::vector<G4double>* R_c = new std::vector<G4double>;
391      std::vector<G4double>* Z_c = new std::vector<G4double>;
392      for( size_t ii = 0; ii < nplanes; ii++ )
393      {
394        (*R_c).push_back( solParam[3+2*ii] );
395        (*Z_c).push_back( solParam[3+2*ii+1] );
396      }
397      G4double phiTotal = solParam[1];
398      if( std::fabs(phiTotal - twopi) < angularTolerance ) { phiTotal = twopi; }
399      solid = new G4Polycone( sname, solParam[0], phiTotal, // start,delta-phi
400                              nplanes, // sections
401                              &((*R_c)[0]), &((*Z_c)[0]));
402    }
403
404  }
405  else if( stype == "POLYHEDRA" )
406  {
407    size_t nplanes = size_t(solParam[3]);
408    G4bool genericPoly = false;
409    if( solParam.size() == 4+nplanes*3 )
410    { 
411      genericPoly = true;
412    }
413    else if( solParam.size() == 4+nplanes*2 )
414    { 
415      genericPoly = false;
416    }
417    else
418    {
419      G4String Err1 = "Solid type " + stype + " should have ";
420      G4String Err2 = G4UIcommand::ConvertToString(G4int(4+nplanes*3))
421                    + " (Z,Rmin,Rmax)\n";
422      G4String Err3 = "or " + G4UIcommand::ConvertToString(G4int(4+nplanes*2));
423      G4String Err4 = " (RZ corners) parameters,\n";
424      G4String Err5 = "and it has "
425                    + G4UIcommand::ConvertToString(G4int(solParam.size()));
426      G4String ErrMessage = Err1 + Err2 + Err3 + Err4 + Err5 + " !";
427      G4Exception("G4tgbVolume::FindOrConstructG4Solid()",
428                  "InvalidSetup", FatalException, ErrMessage);
429    }
430   
431    if( genericPoly )
432    {
433      std::vector<G4double>* z_p = new std::vector<G4double>;
434      std::vector<G4double>* rmin_p = new std::vector<G4double>;
435      std::vector<G4double>* rmax_p = new std::vector<G4double>;
436      for( size_t ii = 0; ii < nplanes; ii++ )
437      {
438        (*z_p).push_back( solParam[4+3*ii] );
439        (*rmin_p).push_back( solParam[4+3*ii+1] );
440        (*rmax_p).push_back(  solParam[4+3*ii+2] );
441      }
442      G4double phiTotal = solParam[1];
443      if( std::fabs(phiTotal - twopi) < angularTolerance ) { phiTotal = twopi; }
444      solid = new G4Polyhedra( sname, solParam[0], phiTotal,
445                               G4int(solParam[2]), nplanes,
446                               &((*z_p)[0]), &((*rmin_p)[0]), &((*rmax_p)[0]));
447    }
448    else
449    {
450      std::vector<G4double>* R_c = new std::vector<G4double>;
451      std::vector<G4double>* Z_c = new std::vector<G4double>;
452      for( size_t ii = 0; ii < nplanes; ii++ )
453      {
454        (*R_c).push_back( solParam[4+2*ii] );
455        (*Z_c).push_back( solParam[4+2*ii+1] );
456      }
457      G4double phiTotal = solParam[1];
458      if( std::fabs(phiTotal - twopi) < angularTolerance ) { phiTotal = twopi; }
459      solid = new G4Polyhedra( sname, solParam[0], phiTotal,
460                               G4int(solParam[2]), nplanes,
461                               &((*R_c)[0]), &((*Z_c)[0]));
462    }
463  }
464  else if( stype == "ELLIPTICALTUBE" )
465  {
466    CheckNoSolidParams( stype, 3, solParam.size() );
467    solid = new G4EllipticalTube( sname, solParam[0], solParam[1], solParam[2]);
468  }
469  else if( stype == "ELLIPSOID" )
470  {
471    CheckNoSolidParams( stype, 5, solParam.size() );
472    solid = new G4Ellipsoid( sname, solParam[0], solParam[1], solParam[2],
473                             solParam[3], solParam[4] );
474  }
475  else if( stype == "ELLIPTICALCONE" )
476  {
477    CheckNoSolidParams( stype, 4, solParam.size() );
478    solid = new G4EllipticalCone( sname, solParam[0], solParam[1],
479                                  solParam[2], solParam[3] );
480  }
481  else if( stype == "HYPE" )
482  {
483    CheckNoSolidParams( stype, 5, solParam.size() );
484    solid = new G4Hype( sname, solParam[0], solParam[1], solParam[2],
485                        solParam[3], solParam[4] );
486  }
487  else if( stype == "TET" )
488  {
489    CheckNoSolidParams( stype, 12, solParam.size() );
490    G4ThreeVector anchor(solParam[0], solParam[1], solParam[2]);
491    G4ThreeVector p2(solParam[3], solParam[4], solParam[5]);
492    G4ThreeVector p3(solParam[6], solParam[7], solParam[8]);
493    G4ThreeVector p4(solParam[9], solParam[10], solParam[11]);
494    solid = new G4Tet( sname, anchor, p2, p3, p4 );
495  }
496  else if( stype == "TWISTEDBOX" )
497  {
498    CheckNoSolidParams( stype, 4, solParam.size() );
499    solid = new G4TwistedBox( sname, solParam[0], solParam[1],
500                              solParam[2], solParam[3]);
501  }
502  else if( stype == "TWISTEDTRAP" )
503  {
504    CheckNoSolidParams( stype, 11, solParam.size() );
505    solid = new G4TwistedTrap( sname, solParam[0], solParam[1], solParam[2],
506                        solParam[3], solParam[4], solParam[5], solParam[6],
507                        solParam[7], solParam[8], solParam[9], solParam[10] );
508  }
509  else if( stype == "TWISTEDTRD" )
510  {
511    CheckNoSolidParams( stype, 6, solParam.size() );
512    solid = new G4TwistedTrd( sname, solParam[0], solParam[1], solParam[2],
513                              solParam[3], solParam[4], solParam[5]);
514  }
515  else if( stype == "TWISTEDTUBS" )
516  {
517    CheckNoSolidParams( stype, 5, solParam.size() );
518    G4double phiTotal = solParam[4];
519    if( std::fabs(phiTotal - twopi) < angularTolerance ) { phiTotal = twopi; }
520    solid = new G4TwistedTubs( sname, solParam[0], solParam[1], solParam[2],
521                               solParam[3], phiTotal);
522  }
523  else if( stype == "BREPBOX" )   // EntityType is = "Closed_Shell"
524  {
525    CheckNoSolidParams( stype, 24, solParam.size() );
526    std::vector<G4Point3D> points;
527    for( size_t ii = 0; ii < 8; ii++ )
528    {
529      points.push_back( G4Point3D(solParam[ii*3+0],
530                                  solParam[ii*3+1],
531                                  solParam[ii*3+2]) );
532    }
533    solid = new G4BREPSolidBox( sname, points[0], points[1], points[2],
534                                points[3], points[4], points[5], points[6],
535                                points[7] );
536  }
537  else if( stype == "BREPCYLINDER" )   // EntityType is = "Closed_Shell"
538  {
539    CheckNoSolidParams( stype, 11, solParam.size() );
540    solid = new G4BREPSolidCylinder( sname,
541                  G4ThreeVector( solParam[0], solParam[1], solParam[2] ),
542                  G4ThreeVector( solParam[3], solParam[4], solParam[5] ),
543                  G4ThreeVector( solParam[6], solParam[7], solParam[8] ),
544                  solParam[9], solParam[10] );
545  }
546  else if( stype == "BREPCONE" )   // EntityType is = "Closed_Shell"
547  {
548    CheckNoSolidParams( stype, 12, solParam.size() );
549    solid = new G4BREPSolidCone( sname,
550                  G4ThreeVector( solParam[0], solParam[1], solParam[2] ),
551                  G4ThreeVector( solParam[3], solParam[4], solParam[5] ),
552                  G4ThreeVector( solParam[6], solParam[7], solParam[8] ),
553                  solParam[9], solParam[10], solParam[11] );
554  }
555  else if( stype == "BREPSPHERE" )   // EntityType is = "Closed_Shell"
556  {
557    CheckNoSolidParams( stype, 10, solParam.size() );
558    solid = new G4BREPSolidSphere( sname,
559                  G4ThreeVector( solParam[0], solParam[1], solParam[2] ),
560                  G4ThreeVector( solParam[3], solParam[4], solParam[5] ),
561                  G4ThreeVector( solParam[6], solParam[7], solParam[8] ),
562                  solParam[9] );
563
564  }
565  else if( stype == "BREPTORUS" )   // EntityType is = "Closed_Shell"
566  {
567    CheckNoSolidParams( stype, 11, solParam.size() );
568    solid = new G4BREPSolidTorus( sname,
569                  G4ThreeVector( solParam[0], solParam[1], solParam[2] ),
570                  G4ThreeVector( solParam[3], solParam[4], solParam[5] ),
571                  G4ThreeVector( solParam[6], solParam[7], solParam[8] ),
572                  solParam[9], solParam[10] );
573  }
574  else if( stype == "BREPPCONE" )   // EntityType is = "Closed_Shell"
575  {
576    size_t nplanes = size_t(solParam[2]);
577    CheckNoSolidParams( stype, 4+3*nplanes, solParam.size() );
578    std::vector<G4double>* z_p = new std::vector<G4double>;
579    std::vector<G4double>* rmin_p = new std::vector<G4double>;
580    std::vector<G4double>* rmax_p = new std::vector<G4double>;
581    for( size_t ii = 0; ii < nplanes; ii++ )
582    {
583      (*z_p).push_back( solParam[4+3*ii] );
584      (*rmin_p).push_back( solParam[4+3*ii+1] );
585      (*rmax_p).push_back(  solParam[4+3*ii+2] );
586    }
587    G4double phiTotal = solParam[1];
588    if( std::fabs(phiTotal - twopi) < angularTolerance ) { phiTotal = twopi; }
589    CheckNoSolidParams( stype, 12, solParam.size() );
590    solid = new G4BREPSolidPCone( sname, solParam[0], phiTotal, // start,dph
591                                  nplanes,                      // sections
592                                  solParam[3],                  // z_start
593                                  &((*z_p)[0]), &((*rmin_p)[0]),
594                                  &((*rmax_p)[0]));
595  }
596  else if( stype == "BREPPOLYHEDRA" )   // EntityType is = "Closed_Shell"
597  {
598    size_t nplanes = size_t(solParam[3]);
599    CheckNoSolidParams( stype, 5+3*nplanes, solParam.size() );
600    std::vector<G4double>* z_p = new std::vector<G4double>;
601    std::vector<G4double>* rmin_p = new std::vector<G4double>;
602    std::vector<G4double>* rmax_p = new std::vector<G4double>;
603    for( size_t ii = 0; ii < nplanes; ii++ )
604    {
605      (*z_p).push_back( solParam[5+3*ii] );
606      (*rmin_p).push_back( solParam[5+3*ii+1] );
607      (*rmax_p).push_back(  solParam[5+3*ii+2] );
608    }
609    G4double phiTotal = solParam[1];
610    if( std::fabs(phiTotal - twopi) < angularTolerance ) { phiTotal = twopi; }
611    CheckNoSolidParams( stype, 12, solParam.size() );
612    solid = new G4BREPSolidPolyhedra( sname, solParam[0], phiTotal, // start,dph
613                                      G4int(solParam[2]), // sides
614                                      nplanes,            // sections
615                                      solParam[4],        // z_start
616                                      &((*z_p)[0]), &((*rmin_p)[0]),
617                                      &((*rmax_p)[0]));
618  }
619  else if( stype == "BREPOPENPCONE" )   // EntityType is = "Closed_Shell"
620  {
621    size_t nplanes = size_t(solParam[2]);
622    std::vector<G4double>* z_p = new std::vector<G4double>;
623    std::vector<G4double>* rmin_p = new std::vector<G4double>;
624    std::vector<G4double>* rmax_p = new std::vector<G4double>;
625    for( size_t ii = 0; ii < nplanes; ii++ )
626    {
627      (*z_p).push_back( solParam[4+3*ii] );
628      (*rmin_p).push_back( solParam[4+3*ii+1] );
629      (*rmax_p).push_back(  solParam[4+3*ii+2] );
630    }
631    G4double phiTotal = solParam[1];
632    if( std::fabs(phiTotal - twopi) < angularTolerance ) { phiTotal = twopi; }
633    CheckNoSolidParams( stype, 12, solParam.size() );
634    solid = new G4BREPSolidOpenPCone( sname, solParam[0], phiTotal, // start,dph
635                                      nplanes,                      // sections
636                                      solParam[3],                  // z_start
637                                      &((*z_p)[0]), &((*rmin_p)[0]),
638                                      &((*rmax_p)[0]));
639  }
640  else if( stype == "TESSELLATED" )
641  {
642    G4int nFacets = G4int(solParam[0]);
643    G4int jj = 0;
644    solid = new G4TessellatedSolid(sname);
645    G4TessellatedSolid* solidTS = (G4TessellatedSolid*)(solid);
646    G4VFacet* facet=0;
647   
648    for( G4int ii = 0; ii < nFacets; ii++){
649      G4int nPoints = G4int(solParam[jj+1]);
650      if( G4int(solParam.size()) < jj + nPoints*3 + 2 ) {
651        G4String Err1 = "Too small number of parameters in tesselated solid, it should be at least "  + G4UIcommand::ConvertToString(jj + nPoints*3 + 2 );
652        G4String Err2 = " facet number " + G4UIcommand::ConvertToString(ii);
653        G4String Err3 = " number of parameters is " + G4UIcommand::ConvertToString(G4int(solParam.size()));
654        G4String ErrMessage = Err1 + Err2 + Err3 + " !";
655        G4Exception("G4tgbVolume::FindOrConstructG4Solid()",
656                    "InvalidSetup", FatalException, ErrMessage);
657      }
658     
659      if( nPoints == 3 ) 
660      {
661        G4ThreeVector pt0(solParam[jj+2],solParam[jj+3],solParam[jj+4]);
662        G4ThreeVector vt1(solParam[jj+5],solParam[jj+6],solParam[jj+7]);
663        G4ThreeVector vt2(solParam[jj+8],solParam[jj+9],solParam[jj+10]);
664        G4FacetVertexType vertexType = ABSOLUTE;
665        if( solParam[jj+11] == 0 ) 
666        {
667          vertexType = ABSOLUTE;
668        } 
669        else if( solParam[jj+11] == 0 ) 
670        {
671          vertexType = RELATIVE;
672        } 
673        else 
674        {
675          G4String Err1 = "Wrong number of vertex type in tesselated solid, it should be 0 =ABSOLUTE) or 1 (=RELATIVE)";
676          G4String Err2 = " facet number " + G4UIcommand::ConvertToString(G4int(ii));
677          G4String Err3 = " vertex type is " + G4UIcommand::ConvertToString(solParam[jj+11]);
678          G4String ErrMessage = Err1 + Err2 + Err3 + " !";
679          G4Exception("G4tgbVolume::FindOrConstructG4Solid()",
680                      "InvalidSetup", FatalException, ErrMessage);
681        }
682        facet = new G4TriangularFacet( pt0, vt1, vt2, vertexType );
683      } 
684      else if( nPoints == 4 ) 
685      {
686        G4ThreeVector pt0(solParam[jj+2],solParam[jj+3],solParam[jj+4]);
687        G4ThreeVector vt1(solParam[jj+5],solParam[jj+6],solParam[jj+7]);
688        G4ThreeVector vt2(solParam[jj+8],solParam[jj+9],solParam[jj+10]);
689        G4ThreeVector vt3(solParam[jj+11],solParam[jj+12],solParam[jj+13]);
690        G4FacetVertexType vertexType = ABSOLUTE;
691        if( solParam[jj+14] == 0 ) 
692        {
693          vertexType = ABSOLUTE;
694        }
695        else if( solParam[jj+14] == 0 ) 
696        {
697          vertexType = RELATIVE;
698        } 
699        else 
700        {
701          G4String Err1 = "Wrong number of vertex type in tesselated solid, it should be 0 =ABSOLUTE) or 1 (=RELATIVE)";
702          G4String Err2 = " facet number " + G4UIcommand::ConvertToString(G4int(ii));
703          G4String Err3 = " vertex type is " + G4UIcommand::ConvertToString(solParam[jj+14]);
704          G4String ErrMessage = Err1 + Err2 + Err3 + " !";
705          G4Exception("G4tgbVolume::FindOrConstructG4Solid()",
706                      "InvalidSetup", FatalException, ErrMessage);
707        }
708        facet = new G4QuadrangularFacet( pt0, vt1, vt2, vt3, vertexType );
709      } 
710      else 
711      {
712        G4String Err1 = "Wrong number of points in tesselated solid, it should be 3 or 4";
713        G4String Err2 = " facet number " + G4UIcommand::ConvertToString(G4int(ii));
714        G4String Err3 = " number of points is " + G4UIcommand::ConvertToString(G4int(nPoints));
715        G4String ErrMessage = Err1 + Err2 + Err3 + " !";
716        G4Exception("G4tgbVolume::FindOrConstructG4Solid()",
717                    "InvalidSetup", FatalException, ErrMessage);
718      }
719     
720      solidTS->AddFacet( facet );
721      jj += nPoints*3 + 2;
722    }
723   
724  } 
725  else if( stype == "EXTRUDED" ) 
726  {
727    std::vector<G4TwoVector> polygonList;
728    std::vector<G4ExtrudedSolid::ZSection> zsectionList;
729    G4int nPolygons = G4int(solParam[0]);
730    G4int ii = 1;
731    G4int nMax = nPolygons*2+1;
732    for( ;ii < nMax; ii+=2 )
733    {
734      polygonList.push_back( G4TwoVector(solParam[ii],solParam[ii+1]) );
735    }
736    G4int nZSections = G4int(solParam[ii]);
737    nMax = nPolygons*2 + nZSections*4 + 2; 
738    ii++;
739    for( ; ii < nMax; ii+=4 )
740    {
741      G4TwoVector offset(solParam[ii+1],solParam[ii+2]);
742      zsectionList.push_back( G4ExtrudedSolid::ZSection(solParam[ii],offset,solParam[ii+3]) );
743    }
744    solid = new G4ExtrudedSolid( sname, polygonList, zsectionList );
745   
746  }
747  else if( stype.substr(0,7) == "Boolean" )
748  {
749    const G4tgrSolidBoolean* solb = dynamic_cast<const G4tgrSolidBoolean*>(sol);
750    G4VSolid* sol1 = FindOrConstructG4Solid( solb->GetSolid(0));
751    G4VSolid* sol2 = FindOrConstructG4Solid( solb->GetSolid(1));
752    G4RotationMatrix* relRotMat = G4tgbRotationMatrixMgr::GetInstance()
753      ->FindOrBuildG4RotMatrix( sol->GetRelativeRotMatName() );
754    G4ThreeVector relPlace = solb->GetRelativePlace();
755
756    if( stype == "Boolean_UNION" )
757    {
758      solid = new G4UnionSolid( sname, sol1, sol2, relRotMat, relPlace );
759    }
760    else if( stype == "Boolean_SUBTRACTION" )
761    {
762      solid = new G4SubtractionSolid( sname, sol1, sol2, relRotMat, relPlace );
763    }
764    else if( stype == "Boolean_INTERSECTION" )
765    {
766      solid = new G4IntersectionSolid( sname, sol1, sol2, relRotMat, relPlace );
767    }
768    else
769    {
770      G4String ErrMessage = "Unknown Boolean type " + stype;
771      G4Exception("G4tgbVolume::FindOrConstructG4Solid()",
772                  "InvalidSetup", FatalException, ErrMessage);
773    }
774  }
775  else
776  {
777    G4String ErrMessage = "Solids of type " + stype
778                        + " not implemented yet, sorry...";
779    G4Exception("G4tgbVolume::FindOrConstructG4Solid()", "NotImplemented",
780                FatalException, ErrMessage);
781  } 
782 
783#ifdef G4VERBOSE
784  if( G4tgrMessenger::GetVerboseLevel() >= 2 )
785  {
786    G4cout << " G4tgbVolume::FindOrConstructG4Solid()" << G4endl
787           << "   Created solid " << sname
788           << " of type " << solid->GetEntityType() << G4endl; 
789  }
790#endif
791
792#ifdef G4VERBOSE
793    if( G4tgrMessenger::GetVerboseLevel() >= 1 )
794    {
795      G4cout << " Constructing new G4Solid: " 
796             << *solid << G4endl;
797    }
798#endif
799 
800  return solid;
801}
802
803//-------------------------------------------------------------------
804void G4tgbVolume::CheckNoSolidParams( const G4String& solidType,
805                                      const unsigned int NoParamExpected,
806                                      const unsigned int NoParam )
807{
808  if( NoParamExpected != NoParam )
809  {
810    G4String Err1 = "Solid type " + solidType + " should have ";
811    G4String Err2 = G4UIcommand::ConvertToString(G4int(NoParamExpected))
812                  + " parameters,\n";
813    G4String Err3 = "and it has "
814                  + G4UIcommand::ConvertToString(G4int(NoParam));
815    G4String ErrMessage = Err1 + Err2 + Err3 + " !";
816    G4Exception("G4tgbVolume::CheckNoSolidParams()", "InvalidSetup",
817                FatalException, ErrMessage);
818  }
819}
820
821
822//-------------------------------------------------------------------
823G4LogicalVolume* G4tgbVolume::ConstructG4LogVol( const G4VSolid* solid )
824{
825  G4LogicalVolume* logvol;
826
827#ifdef G4VERBOSE
828  if( G4tgrMessenger::GetVerboseLevel() >= 2 )
829  {
830    G4cout << " G4tgbVolume::ConstructG4LogVol() - " << GetName() << G4endl;
831  }
832#endif
833
834  //----------- Get the material first
835  G4Material* mate = G4tgbMaterialMgr::GetInstance()
836              ->FindOrBuildG4Material( theTgrVolume->GetMaterialName() );
837  if( mate == 0 )
838  {
839    G4String ErrMessage = "Material not found "
840                        + theTgrVolume->GetMaterialName()
841                        + " for volume " + GetName() + ".";
842    G4Exception("G4tgbVolume::ConstructG4LogVol()", "InvalidSetup",
843                FatalException, ErrMessage);
844  }
845#ifdef G4VERBOSE
846  if( G4tgrMessenger::GetVerboseLevel() >= 2 )
847  {
848    G4cout << " G4tgbVolume::ConstructG4LogVol() -"
849           << " Material constructed: " << mate->GetName() << G4endl; 
850  }
851#endif
852 
853  //---------- Construct the LV
854  logvol = new G4LogicalVolume( const_cast<G4VSolid*>(solid),
855                                const_cast<G4Material*>(mate), GetName() );
856
857#ifdef G4VERBOSE
858    if( G4tgrMessenger::GetVerboseLevel() >= 1 )
859    {
860      G4cout << " Constructing new G4LogicalVolume: " 
861             << logvol->GetName() << " mate " << mate->GetName() << G4endl;
862    }
863#endif
864 
865  //---------- Set visibility and colour
866  if( !GetVisibility() || GetColour()[0] != -1 )
867  {
868    G4VisAttributes* visAtt = new G4VisAttributes();
869#ifdef G4VERBOSE
870    if( G4tgrMessenger::GetVerboseLevel() >= 1 )
871    {
872      G4cout << " Constructing new G4VisAttributes: " 
873             << *visAtt << G4endl;
874    }
875#endif
876 
877    if( !GetVisibility() )
878    {
879      visAtt->SetVisibility( GetVisibility() );
880    }
881    else if( GetColour()[0] != -1 )
882    {
883      // this else should not be necessary, because if the visibility
884      // is set to off, colour should have no effect. But it does not
885      // work: if you set colour and vis off, it is visualized!?!?!?
886
887      const G4double* col = GetColour();
888      if( col[3] == -1. )
889      {
890        visAtt->SetColour( G4Colour(col[0],col[1],col[2]));
891      }
892      else
893      {
894        visAtt->SetColour( G4Colour(col[0],col[1],col[2],col[3]));
895      }
896    }
897    logvol->SetVisAttributes(visAtt);
898  }
899
900#ifdef G4VERBOSE
901  if( G4tgrMessenger::GetVerboseLevel() >= 2 )
902  {
903    G4cout << " G4tgbVolume::ConstructG4LogVol() -"
904           << " Created logical volume: " << GetName() << G4endl;
905  }
906#endif
907
908  return logvol;
909}
910
911
912//-------------------------------------------------------------------
913G4VPhysicalVolume*
914G4tgbVolume::ConstructG4PhysVol( const G4tgrPlace* place,
915                                 const G4LogicalVolume* currentLV,
916                                 const G4LogicalVolume* parentLV )
917{
918  G4VPhysicalVolume* physvol = 0;
919  G4int copyNo;
920 
921  //----- Case of placement of top volume
922  if( place == 0 )
923  {
924#ifdef G4VERBOSE
925    if( G4tgrMessenger::GetVerboseLevel() >= 2 )
926    {
927      G4cout << " G4tgbVolume::ConstructG4PhysVol() - World: "
928             << GetName() << G4endl;
929    }
930#endif
931    physvol = new G4PVPlacement(0, G4ThreeVector(),
932                                const_cast<G4LogicalVolume*>(currentLV),
933                                GetName(), 0, false, 0, 
934                                theTgrVolume->GetCheckOverlaps());
935#ifdef G4VERBOSE
936    if( G4tgrMessenger::GetVerboseLevel() >= 1 )
937    {
938      G4cout << " Constructing new : G4PVPlacement " 
939             << physvol->GetName() << G4endl;
940    }
941#endif
942  }
943  else
944  { 
945    copyNo = place->GetCopyNo();
946
947#ifdef G4VERBOSE
948    if( G4tgrMessenger::GetVerboseLevel() >= 2 )
949    {
950      G4cout << " G4tgbVolume::ConstructG4PhysVol() - " << GetName() << G4endl
951             << "   inside " << parentLV->GetName() << " copy No: " << copyNo
952             << " of type= " << theTgrVolume->GetType() << G4endl
953             << "   placement type= " << place->GetType() << G4endl;
954    }
955#endif
956   
957    if( theTgrVolume->GetType() == "VOLSimple" )
958    {
959      //----- Get placement
960#ifdef G4VERBOSE
961      if( G4tgrMessenger::GetVerboseLevel() >= 2 )
962      {
963        G4cout << " G4tgbVolume::ConstructG4PhysVol() - Placement type = "
964               << place->GetType() << G4endl;
965      }
966#endif
967     
968      //--------------- If it is  G4tgrPlaceSimple
969      if( place->GetType() == "PlaceSimple" )
970      {
971        //----- Get rotation matrix
972        G4tgrPlaceSimple* placeSimple = (G4tgrPlaceSimple*)place; 
973        G4String rmName = placeSimple->GetRotMatName();
974
975        G4RotationMatrix* rotmat = G4tgbRotationMatrixMgr::GetInstance()
976                                 ->FindOrBuildG4RotMatrix( rmName );
977        //----- Place volume in mother
978        G4double check = (rotmat->colX().cross(rotmat->colY()))*rotmat->colZ();
979        G4double tol = 1.0e-3;
980        //---- Check that matrix is ortogonal
981        if (1-std::abs(check)>tol)
982        {
983          G4cerr << " Matrix : " << rmName << " " << rotmat->colX()
984                 << " " << rotmat->colY() << " " << rotmat->colZ() << G4endl
985                 << " product x X y * z = " << check << " x X y "
986                 << rotmat->colX().cross(rotmat->colY()) << G4endl;
987          G4String ErrMessage = "Rotation is not ortogonal " + rmName + " !";
988          G4Exception("G4tgbVolume::ConstructG4PhysVol()",
989                      "InvalidSetup", FatalException, ErrMessage);
990          //---- Check if it is reflection
991        }
992        else if (1+check<=tol)
993        {
994          G4Translate3D transl = place->GetPlacement();
995          G4Transform3D trfrm  = transl * G4Rotate3D(*rotmat);
996          physvol = (G4ReflectionFactory::Instance()->Place(trfrm, GetName(),
997                     const_cast<G4LogicalVolume*>(currentLV),
998                     const_cast<G4LogicalVolume*>(parentLV),
999                     false, copyNo, false )).first;
1000        }
1001        else
1002        {
1003#ifdef G4VERBOSE
1004          if( G4tgrMessenger::GetVerboseLevel() >= 1 )
1005          {
1006            G4cout << "Construction new G4VPhysicalVolume"
1007                   << " through G4ReflectionFactory " << GetName() 
1008                   << " in volume " << parentLV->GetName() 
1009                   << " copyNo " << copyNo
1010                   << " position " << place->GetPlacement() 
1011                   << " ROT " << rotmat->colX() 
1012                   << " " << rotmat->colY() 
1013                   << " " << rotmat->colZ() << G4endl;
1014          }
1015#endif
1016          physvol = new G4PVPlacement( rotmat, place->GetPlacement(),
1017                                       const_cast<G4LogicalVolume*>(currentLV),
1018                                       GetName(),
1019                                       const_cast<G4LogicalVolume*>(parentLV),
1020                                       false, copyNo, 
1021                                       theTgrVolume->GetCheckOverlaps());
1022        }
1023       
1024        //--------------- If it is G4tgrPlaceParam
1025      }
1026      else if( place->GetType() == "PlaceParam" )
1027      {
1028        G4tgrPlaceParameterisation* dp = (G4tgrPlaceParameterisation*)(place);
1029
1030        //----- See what parameterisation type
1031#ifdef G4VERBOSE
1032        if( G4tgrMessenger::GetVerboseLevel() >= 2 )
1033        {
1034          G4cout << " G4tgbVolume::ConstructG4PhysVol() -" << G4endl
1035                 << "   param: " << GetName() << " in " <<  parentLV->GetName()
1036                 << " param type= " << dp->GetParamType() << G4endl;
1037        }
1038#endif
1039       
1040        G4tgbPlaceParameterisation * param=0;
1041       
1042        if( (dp->GetParamType() == "CIRCLE")
1043         || (dp->GetParamType() == "CIRCLE_XY")
1044         || (dp->GetParamType() == "CIRCLE_XZ")
1045         || (dp->GetParamType() == "CIRCLE_YZ") )
1046        { 
1047          param = new G4tgbPlaceParamCircle(dp);
1048         
1049        } 
1050        else if( (dp->GetParamType() == "LINEAR")
1051                || (dp->GetParamType() == "LINEAR_X")
1052                || (dp->GetParamType() == "LINEAR_Y")
1053                || (dp->GetParamType() == "LINEAR_Z") )
1054        {   
1055          param = new G4tgbPlaceParamLinear(dp);
1056         
1057        } 
1058        else if( (dp->GetParamType() == "SQUARE")
1059                || (dp->GetParamType() == "SQUARE_XY")
1060                || (dp->GetParamType() == "SQUARE_XZ")
1061                || (dp->GetParamType() == "SQUARE_YZ") )
1062        {
1063          param = new G4tgbPlaceParamSquare(dp);
1064        }
1065        else
1066        {
1067          G4String ErrMessage = "Parameterisation has wrong type, TYPE: "
1068                              + G4String(dp->GetParamType()) + " !";
1069          G4Exception("G4tgbVolume::ConstructG4PhysVol", "WrongArgument",
1070                      FatalException, ErrMessage);
1071        }
1072#ifdef G4VERBOSE
1073        if( G4tgrMessenger::GetVerboseLevel() >= 1 )
1074        {
1075          G4cout << " G4tgbVolume::ConstructG4PhysVol() -" << G4endl
1076                 << "   New G4PVParameterised: " << GetName() << " vol "
1077                 << currentLV->GetName() << " in vol " << parentLV->GetName()
1078                 << " axis " << param->GetAxis() << " nCopies "
1079                 << param->GetNCopies() << G4endl;
1080        }
1081#endif
1082        physvol = new G4PVParameterised(GetName(),
1083                                        const_cast<G4LogicalVolume*>(currentLV),
1084                                        const_cast<G4LogicalVolume*>(parentLV),
1085                                        EAxis(param->GetAxis()),
1086                                        param->GetNCopies(), param);
1087#ifdef G4VERBOSE
1088    if( G4tgrMessenger::GetVerboseLevel() >= 1 )
1089    {
1090      G4cout << " Constructing new G4PVParameterised: " 
1091             << physvol->GetName() << " in volume " << parentLV->GetName() 
1092             << " N copies " << param->GetNCopies() 
1093             << " axis " << param->GetAxis() << G4endl;
1094    }
1095#endif
1096
1097      }
1098      else if( place->GetType() == "PlaceReplica" )
1099      {
1100        //--------------- If it is  PlaceReplica
1101        G4tgrPlaceDivRep* dpr = (G4tgrPlaceDivRep*)place;
1102
1103#ifdef G4VERBOSE
1104        if( G4tgrMessenger::GetVerboseLevel() >= 2 )
1105        {
1106          G4cout << " G4tgbVolume::ConstructG4PhysVol() -" << G4endl
1107                 << "   replica" << " " << currentLV->GetName()
1108                 << " in " <<  parentLV->GetName() 
1109                 << " NDiv " << dpr->GetNDiv() << " Width " << dpr->GetWidth()
1110                 << " offset " << dpr->GetOffset() << G4endl;
1111        }
1112#endif
1113        physvol = new G4PVReplica(GetName(),
1114                                  const_cast<G4LogicalVolume*>(currentLV),
1115                                  const_cast<G4LogicalVolume*>(parentLV),
1116                                  EAxis(dpr->GetAxis()), dpr->GetNDiv(),
1117                                  dpr->GetWidth(), dpr->GetOffset());
1118#ifdef G4VERBOSE
1119    if( G4tgrMessenger::GetVerboseLevel() >= 1 )
1120    {
1121      G4cout << " Constructing new G4PVReplica: " 
1122             << currentLV->GetName()
1123             << " in " <<  parentLV->GetName() 
1124             << " NDiv " << dpr->GetNDiv() << " Width " << dpr->GetWidth()
1125             << " offset " << dpr->GetOffset() << G4endl;
1126    }
1127#endif
1128      }
1129    }
1130    else if( theTgrVolume->GetType() == "VOLDivision" )
1131    {
1132      G4tgrVolumeDivision* volr = (G4tgrVolumeDivision*)theTgrVolume;
1133      G4tgrPlaceDivRep* placeDiv = volr->GetPlaceDivision() ;
1134      G4VSolid* solid = BuildSolidForDivision( parentLV->GetSolid(), placeDiv->GetAxis() );
1135      G4Material* mate = G4tgbMaterialMgr::GetInstance()
1136                  ->FindOrBuildG4Material( theTgrVolume->GetMaterialName() );
1137      G4LogicalVolume* divLV = new G4LogicalVolume(solid,
1138                                                  const_cast<G4Material*>(mate),
1139                                                   GetName() );
1140#ifdef G4VERBOSE
1141      if( G4tgrMessenger::GetVerboseLevel() >= 1 )
1142      {
1143        G4cout << " Constructed new G4LogicalVolume for division: " 
1144               << divLV->GetName() << " mate " << mate->GetName() << G4endl;
1145      }
1146#endif
1147 
1148      G4DivType divType = placeDiv->GetDivType();
1149      switch (divType)
1150      {
1151      case DivByNdiv:
1152        physvol = new G4PVDivision(GetName(), (G4LogicalVolume*)divLV,
1153                                   const_cast<G4LogicalVolume*>(parentLV),
1154                                   placeDiv->GetAxis(), placeDiv->GetNDiv(),
1155                                   placeDiv->GetOffset());
1156#ifdef G4VERBOSE
1157        if( G4tgrMessenger::GetVerboseLevel() >= 1 )
1158        {
1159          G4cout << " Constructing new G4PVDivision by number of divisions: " 
1160                 << GetName() << " in " <<  parentLV->GetName() 
1161                 << " axis " << placeDiv->GetAxis() 
1162                 << " Ndiv " << placeDiv->GetNDiv()
1163                 << " offset " << placeDiv->GetOffset() << G4endl;
1164        }
1165#endif
1166        break;
1167      case DivByWidth:
1168        physvol = new G4PVDivision(GetName(), (G4LogicalVolume*)divLV,
1169                                   const_cast<G4LogicalVolume*>(parentLV),
1170                                   placeDiv->GetAxis(), placeDiv->GetWidth(),
1171                                   placeDiv->GetOffset());
1172#ifdef G4VERBOSE
1173        if( G4tgrMessenger::GetVerboseLevel() >= 1 )
1174        {
1175          G4cout << " Constructing new G4PVDivision by width: " 
1176                 << GetName() << " in " <<  parentLV->GetName() 
1177                 << " axis " << placeDiv->GetAxis() 
1178                 << " width " << placeDiv->GetWidth()
1179                 << " offset " << placeDiv->GetOffset() << G4endl;
1180        }
1181#endif
1182        break;
1183      case DivByNdivAndWidth:
1184        physvol = new G4PVDivision(GetName(), (G4LogicalVolume*)divLV,
1185                                   const_cast<G4LogicalVolume*>(parentLV),
1186                                   placeDiv->GetAxis(), placeDiv->GetNDiv(),
1187                                   placeDiv->GetWidth(),
1188                                   placeDiv->GetOffset());
1189#ifdef G4VERBOSE
1190        if( G4tgrMessenger::GetVerboseLevel() >= 1 )
1191        {
1192          G4cout << " Constructing new G4PVDivision by width"
1193                 << " and number of divisions: " 
1194                 << GetName() << " in " <<  parentLV->GetName() 
1195                 << " axis " << placeDiv->GetAxis() 
1196                 << " Ndiv " << placeDiv->GetNDiv()
1197                 << " width " << placeDiv->GetWidth()
1198                 << " offset " << placeDiv->GetOffset() << G4endl;
1199        }
1200#endif
1201        break;
1202      }
1203    }
1204    else if( theTgrVolume->GetType() == "VOLAssembly" )
1205    {
1206      // Define one layer as one assembly volume
1207      G4tgrVolumeAssembly * tgrAssembly = (G4tgrVolumeAssembly *)theTgrVolume;
1208
1209      if( !theG4AssemblyVolume )
1210      {
1211        theG4AssemblyVolume = new G4AssemblyVolume;
1212#ifdef G4VERBOSE
1213        if( G4tgrMessenger::GetVerboseLevel() >= 1 )
1214          {
1215            G4cout << " Constructing new G4AssemblyVolume: " 
1216                   << " number of assembly components " 
1217                   <<  tgrAssembly->GetNoComponents() << G4endl;
1218          }
1219#endif       
1220        G4tgbVolumeMgr* g4vmgr = G4tgbVolumeMgr::GetInstance();
1221        for( G4int ii = 0; ii < tgrAssembly->GetNoComponents(); ii++ )
1222        {
1223          // Rotation and translation of a plate inside the assembly
1224
1225          G4ThreeVector transl =  tgrAssembly->GetComponentPos(ii);
1226          G4String rmName = tgrAssembly->GetComponentRM(ii);
1227          G4RotationMatrix* rotmat = G4tgbRotationMatrixMgr::GetInstance()
1228                                   ->FindOrBuildG4RotMatrix( rmName );
1229         
1230          //----- Get G4LogicalVolume of component
1231          G4String lvname = tgrAssembly->GetComponentName(ii);
1232          G4LogicalVolume* logvol = g4vmgr->FindG4LogVol( lvname);
1233          if( logvol == 0 )
1234          {
1235            g4vmgr->FindVolume( lvname )->ConstructG4Volumes( 0, 0);
1236            logvol = g4vmgr->FindG4LogVol( lvname, true );
1237          }
1238          // Fill the assembly by the plates
1239          theG4AssemblyVolume->AddPlacedVolume( logvol, transl, rotmat );
1240#ifdef G4VERBOSE
1241          if( G4tgrMessenger::GetVerboseLevel() >= 1 )
1242          {
1243            G4cout << " G4AssemblyVolume->AddPlacedVolume " << ii 
1244                   << " " << logvol->GetName()
1245                   << " translation " << transl 
1246                   << " rotmat " << rotmat->colX() 
1247                   << " " << rotmat->colY() 
1248                   << " " << rotmat->colZ() << G4endl;
1249          }
1250#endif
1251
1252        }
1253      }
1254
1255      // Rotation and Translation of the assembly inside the world
1256
1257      G4tgrPlaceSimple* placeSimple = (G4tgrPlaceSimple*)place; 
1258      G4String rmName = placeSimple->GetRotMatName();
1259      G4RotationMatrix* rotmat = G4tgbRotationMatrixMgr::GetInstance()
1260                               ->FindOrBuildG4RotMatrix( rmName );
1261      G4ThreeVector transl = place->GetPlacement();
1262
1263      G4LogicalVolume* parentLV_nonconst =
1264                       const_cast<G4LogicalVolume*>(parentLV);
1265      theG4AssemblyVolume->MakeImprint( parentLV_nonconst, transl, rotmat );
1266 
1267    }
1268    else   // If it is G4tgrVolumeAssembly
1269    {
1270      G4String ErrMessage = "Volume type not supported: "
1271                          + theTgrVolume->GetType() + ", sorry...";
1272      G4Exception("G4tgbVolume::ConstructG4PhysVol()", "NotImplemented",
1273                  FatalException, ErrMessage);
1274    }   
1275  } 
1276
1277  return physvol;
1278}
1279
1280
1281//-------------------------------------------------------------------
1282G4VSolid* G4tgbVolume::BuildSolidForDivision( G4VSolid* parentSolid, EAxis axis  )
1283{
1284  G4VSolid* solid=0;
1285  G4double redf = (parentSolid->GetExtent().GetXmax()-parentSolid->GetExtent().GetXmin());
1286  redf = std::min(redf,parentSolid->GetExtent().GetYmax()-parentSolid->GetExtent().GetYmin());
1287  redf = std::min(redf,parentSolid->GetExtent().GetZmax()-parentSolid->GetExtent().GetZmin());
1288  redf *= 0.001; //make daugther much smaller, to fit in parent
1289
1290  if( parentSolid->GetEntityType() == "G4Box" )
1291  {
1292    G4Box* psolid = (G4Box*)(parentSolid);
1293    solid = new G4Box(GetName(), psolid->GetXHalfLength()*redf,
1294                                 psolid->GetZHalfLength()*redf,
1295                                 psolid->GetZHalfLength()*redf);
1296  } 
1297  else if ( parentSolid->GetEntityType() == "G4Tubs" )
1298  {
1299    G4Tubs* psolid = (G4Tubs*)(parentSolid);
1300    solid = new G4Tubs( GetName(), psolid->GetInnerRadius()*redf,
1301                                   psolid->GetOuterRadius()*redf,
1302                                   psolid->GetZHalfLength()*redf,
1303                                   psolid->GetSPhi(), psolid->GetDPhi());
1304  } 
1305  else if ( parentSolid->GetEntityType() == "G4Cons" )
1306  {
1307    G4Cons* psolid = (G4Cons*)(parentSolid);
1308    solid = new G4Cons( GetName(), psolid->GetInnerRadiusMinusZ()*redf,
1309                                   psolid->GetOuterRadiusMinusZ()*redf,
1310                                   psolid->GetInnerRadiusPlusZ()*redf,
1311                                   psolid->GetOuterRadiusPlusZ()*redf,
1312                                   psolid->GetZHalfLength()*redf,
1313                                   psolid->GetSPhi(), psolid->GetDPhi());
1314  } 
1315  else if ( parentSolid->GetEntityType() == "G4Trd" )
1316  {
1317    G4Trd* psolid = (G4Trd*)(parentSolid);
1318    G4double mpDx1 = psolid->GetXHalfLength1();
1319    G4double mpDx2 = psolid->GetXHalfLength2();
1320
1321    if( axis == kXAxis && std::fabs(mpDx1 - mpDx2) > G4GeometryTolerance::GetInstance()->GetSurfaceTolerance() )
1322    {
1323      solid = new G4Trap( GetName(), psolid->GetZHalfLength()*redf, 
1324                          psolid->GetYHalfLength1()*redf, 
1325                          psolid->GetXHalfLength2()*redf, 
1326                          psolid->GetXHalfLength1()*redf );
1327    } 
1328    else
1329    {
1330      solid = new G4Trd( GetName(), psolid->GetXHalfLength1()*redf,
1331                         psolid->GetXHalfLength2()*redf,
1332                         psolid->GetYHalfLength1()*redf,
1333                         psolid->GetYHalfLength2()*redf,
1334                         psolid->GetZHalfLength()*redf);
1335    }
1336   
1337  } 
1338  else if ( parentSolid->GetEntityType() == "G4Para" )
1339  {
1340    G4Para* psolid = (G4Para*)(parentSolid);
1341    solid = new G4Para( GetName(), psolid->GetXHalfLength()*redf,
1342                                   psolid->GetYHalfLength()*redf,
1343                                   psolid->GetZHalfLength()*redf,
1344                                   std::atan(psolid->GetTanAlpha()),
1345                                   psolid->GetSymAxis().theta(),
1346                                   psolid->GetSymAxis().phi() ); 
1347  } 
1348  else if ( parentSolid->GetEntityType() == "G4Polycone" )
1349  {
1350    G4Polycone* psolid = (G4Polycone*)(parentSolid);
1351    G4PolyconeHistorical origParam = *(psolid->GetOriginalParameters());
1352    for( G4int ii = 0; ii < origParam.Num_z_planes; ii++ )
1353    {
1354      origParam.Rmin[ii] = origParam.Rmin[ii]*redf;
1355      origParam.Rmax[ii] = origParam.Rmax[ii]*redf;
1356    }
1357    solid = new G4Polycone( GetName(), psolid->GetStartPhi(),
1358                                       psolid->GetEndPhi(),
1359                            origParam.Num_z_planes, origParam.Z_values,
1360                            origParam.Rmin, origParam.Rmax);
1361
1362  } 
1363  else if ( parentSolid->GetEntityType() == "G4Polyhedra" )
1364  {
1365    G4Polyhedra* psolid = (G4Polyhedra*)(parentSolid);
1366    G4PolyhedraHistorical origParam = *(psolid->GetOriginalParameters());
1367    for( G4int ii = 0; ii < origParam.Num_z_planes; ii++ )
1368    {
1369      origParam.Rmin[ii] = origParam.Rmin[ii]*redf;
1370      origParam.Rmax[ii] = origParam.Rmax[ii]*redf;
1371    }
1372    solid = new G4Polyhedra( GetName(), psolid->GetStartPhi(),
1373                                        psolid->GetEndPhi(),
1374                                        psolid->GetNumSide(),
1375                             origParam.Num_z_planes, origParam.Z_values,
1376                             origParam.Rmin, origParam.Rmax);
1377  }
1378  else
1379  { 
1380    G4String ErrMessage = "Solid type not supported. VOLUME= " + GetName()
1381                        + " Solid type= " + parentSolid->GetEntityType() + "\n"
1382                        + "Only supported types are: G4Box, G4Tubs, G4Cons,"
1383                        + " G4Trd, G4Para, G4Polycone, G4Polyhedra.";
1384    G4Exception("G4tgbVolume::BuildSolidForDivision()", "NotImplemented",
1385                FatalException, ErrMessage);
1386  }
1387
1388#ifdef G4VERBOSE
1389    if( G4tgrMessenger::GetVerboseLevel() >= 1 )
1390    {
1391      G4cout << " Constructing new G4Solid for division: " 
1392             << *solid << G4endl;
1393    }
1394#endif
1395  return solid;
1396}
1397 
Note: See TracBrowser for help on using the repository browser.