Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/geometry/solids/specific/src/G4UPolycone.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  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 // Implementation of G4UPolycone wrapper class
 27 //
 28 // 31.10.13 G.Cosmo, CERN
 29 // --------------------------------------------------------------------
 30 
 31 #include "G4Polycone.hh"
 32 #include "G4UPolycone.hh"
 33 
 34 #if ( defined(G4GEOM_USE_USOLIDS) || defined(G4GEOM_USE_PARTIAL_USOLIDS) )
 35 
 36 #include "G4GeomTools.hh"
 37 #include "G4AffineTransform.hh"
 38 #include "G4VPVParameterisation.hh"
 39 #include "G4BoundingEnvelope.hh"
 40 
 41 using namespace CLHEP;
 42 
 43 ////////////////////////////////////////////////////////////////////////
 44 //
 45 // Constructor (GEANT3 style parameters)
 46 //
 47 G4UPolycone::G4UPolycone( const G4String& name,
 48                               G4double phiStart,
 49                               G4double phiTotal,
 50                               G4int numZPlanes,
 51                         const G4double zPlane[],
 52                         const G4double rInner[],
 53                         const G4double rOuter[]  )
 54   : Base_t(name, phiStart, phiTotal, numZPlanes, zPlane, rInner, rOuter)
 55 {
 56   fGenericPcon = false;
 57   SetOriginalParameters();
 58   wrStart = phiStart;
 59   while (wrStart < 0)
 60   {
 61     wrStart += twopi;
 62   }
 63   wrDelta = phiTotal;
 64   if (wrDelta <= 0 || wrDelta >= twopi*(1-DBL_EPSILON))
 65   {
 66     wrStart = 0;
 67     wrDelta = twopi;
 68   }
 69   rzcorners.resize(0);
 70   for (G4int i=0; i<numZPlanes; ++i)
 71   {
 72     G4double z = zPlane[i];
 73     G4double r = rOuter[i];
 74     rzcorners.emplace_back(r,z);
 75   }
 76   for (G4int i=numZPlanes-1; i>=0; --i)
 77   {
 78     G4double z = zPlane[i];
 79     G4double r = rInner[i];
 80     rzcorners.emplace_back(r,z);
 81   }
 82   std::vector<G4int> iout;
 83   G4GeomTools::RemoveRedundantVertices(rzcorners,iout,2*kCarTolerance);
 84 }
 85 
 86 
 87 ////////////////////////////////////////////////////////////////////////
 88 //
 89 // Constructor (generic parameters)
 90 //
 91 G4UPolycone::G4UPolycone(const G4String& name,
 92                                G4double phiStart,
 93                                G4double phiTotal,
 94                                G4int    numRZ,
 95                          const G4double r[],
 96                          const G4double z[]   )
 97   : Base_t(name, phiStart, phiTotal, numRZ, r, z)
 98 {
 99   fGenericPcon = true;
100   SetOriginalParameters();
101   wrStart = phiStart; while (wrStart < 0) wrStart += twopi;
102   wrDelta = phiTotal;
103   if (wrDelta <= 0 || wrDelta >= twopi*(1-DBL_EPSILON))
104   {
105     wrStart = 0;
106     wrDelta = twopi;
107   }
108   rzcorners.resize(0);
109   for (G4int i=0; i<numRZ; ++i)
110   {
111     rzcorners.emplace_back(r[i],z[i]);
112   }
113   std::vector<G4int> iout;
114   G4GeomTools::RemoveRedundantVertices(rzcorners,iout,2*kCarTolerance);
115 }
116 
117 
118 ////////////////////////////////////////////////////////////////////////
119 //
120 // Fake default constructor - sets only member data and allocates memory
121 //                            for usage restricted to object persistency.
122 //
123 G4UPolycone::G4UPolycone( __void__& a )
124   : Base_t(a)
125 {
126 }
127 
128 
129 ////////////////////////////////////////////////////////////////////////
130 //
131 // Destructor
132 //
133 G4UPolycone::~G4UPolycone() = default;
134 
135 
136 ////////////////////////////////////////////////////////////////////////
137 //
138 // Copy constructor
139 //
140 G4UPolycone::G4UPolycone( const G4UPolycone& source )
141   : Base_t( source )
142 {
143   fGenericPcon = source.fGenericPcon;
144   fOriginalParameters = source.fOriginalParameters;
145   wrStart   = source.wrStart;
146   wrDelta   = source.wrDelta;
147   rzcorners = source.rzcorners;
148 }
149 
150 
151 ////////////////////////////////////////////////////////////////////////
152 //
153 // Assignment operator
154 //
155 G4UPolycone& G4UPolycone::operator=( const G4UPolycone& source )
156 {
157   if (this == &source) return *this;
158 
159   Base_t::operator=( source );
160   fGenericPcon = source.fGenericPcon;
161   fOriginalParameters = source.fOriginalParameters;
162   wrStart   = source.wrStart;
163   wrDelta   = source.wrDelta;
164   rzcorners = source.rzcorners;
165 
166   return *this;
167 }
168 
169 
170 ////////////////////////////////////////////////////////////////////////
171 //
172 // Accessors & modifiers
173 //
174 G4double G4UPolycone::GetStartPhi() const
175 {
176   return wrStart;
177 }
178 G4double G4UPolycone::GetDeltaPhi() const
179 {
180   return wrDelta;
181 }
182 G4double G4UPolycone::GetEndPhi() const
183 {
184   return (wrStart + wrDelta);
185 }
186 G4double G4UPolycone::GetSinStartPhi() const
187 {
188   if (!IsOpen()) return 0.;
189   G4double phi = GetStartPhi();
190   return std::sin(phi);
191 }
192 G4double G4UPolycone::GetCosStartPhi() const
193 {
194   if (!IsOpen()) return 1.;
195   G4double phi = GetStartPhi();
196   return std::cos(phi);
197 }
198 G4double G4UPolycone::GetSinEndPhi() const
199 {
200   if (!IsOpen()) return 0.;
201   G4double phi = GetEndPhi();
202   return std::sin(phi);
203 }
204 G4double G4UPolycone::GetCosEndPhi() const
205 {
206   if (!IsOpen()) return 1.;
207   G4double phi = GetEndPhi();
208   return std::cos(phi);
209 }
210 G4bool G4UPolycone::IsOpen() const
211 {
212   return (wrDelta < twopi);
213 }
214 G4int G4UPolycone::GetNumRZCorner() const
215 {
216   return rzcorners.size();
217 }
218 G4PolyconeSideRZ G4UPolycone::GetCorner(G4int index) const
219 {
220   G4TwoVector rz = rzcorners.at(index);
221   G4PolyconeSideRZ psiderz = { rz.x(), rz.y() };
222 
223   return psiderz;
224 }
225 G4PolyconeHistorical* G4UPolycone::GetOriginalParameters() const
226 {
227   return new G4PolyconeHistorical(fOriginalParameters);
228 }
229 void G4UPolycone::SetOriginalParameters()
230 {
231   vecgeom::PolyconeHistorical* original_parameters = Base_t::GetOriginalParameters();
232 
233   fOriginalParameters.Start_angle   = original_parameters->fHStart_angle;
234   fOriginalParameters.Opening_angle = original_parameters->fHOpening_angle;
235   fOriginalParameters.Num_z_planes  = original_parameters->fHNum_z_planes;
236 
237   delete [] fOriginalParameters.Z_values;
238   delete [] fOriginalParameters.Rmin;
239   delete [] fOriginalParameters.Rmax;
240 
241   G4int numPlanes = fOriginalParameters.Num_z_planes;
242   fOriginalParameters.Z_values = new G4double[numPlanes];
243   fOriginalParameters.Rmin = new G4double[numPlanes];
244   fOriginalParameters.Rmax = new G4double[numPlanes];
245   for (G4int i=0; i<numPlanes; ++i)
246   {
247     fOriginalParameters.Z_values[i] = original_parameters->fHZ_values[i];
248     fOriginalParameters.Rmin[i]     = original_parameters->fHRmin[i];
249     fOriginalParameters.Rmax[i]     = original_parameters->fHRmax[i];
250   }
251 }
252 void G4UPolycone::SetOriginalParameters(G4PolyconeHistorical* pars)
253 {
254   fOriginalParameters = *pars;
255   fRebuildPolyhedron = true;
256   Reset();
257 }
258 
259 G4bool G4UPolycone::Reset()
260 {
261   if (fGenericPcon)
262   {
263     std::ostringstream message;
264     message << "Solid " << GetName() << " built using generic construct."
265             << G4endl << "Not applicable to the generic construct !";
266     G4Exception("G4UPolycone::Reset()", "GeomSolids1001",
267                 JustWarning, message, "Parameters NOT reset.");
268     return true;  // error code set
269   }
270 
271   //
272   // Rebuild polycone based on original parameters
273   //
274   wrStart = fOriginalParameters.Start_angle;
275   while (wrStart < 0.)
276   {
277     wrStart += twopi;
278   }
279   wrDelta = fOriginalParameters.Opening_angle;
280   if (wrDelta <= 0. || wrDelta >= twopi*(1-DBL_EPSILON))
281   {
282     wrStart = 0.;
283     wrDelta = twopi;
284   }
285   rzcorners.resize(0);
286   for (G4int i=0; i<fOriginalParameters.Num_z_planes; ++i)
287     {
288       G4double z = fOriginalParameters.Z_values[i];
289       G4double r = fOriginalParameters.Rmax[i];
290       rzcorners.emplace_back(r,z);
291     }
292   for (G4int i=fOriginalParameters.Num_z_planes-1; i>=0; --i)
293     {
294       G4double z = fOriginalParameters.Z_values[i];
295       G4double r = fOriginalParameters.Rmin[i];
296       rzcorners.emplace_back(r,z);
297     }
298   std::vector<G4int> iout;
299   G4GeomTools::RemoveRedundantVertices(rzcorners,iout,2*kCarTolerance);
300 
301   return false;  // error code unset
302 }
303 
304 ////////////////////////////////////////////////////////////////////////
305 //
306 // Dispatch to parameterisation for replication mechanism dimension
307 // computation & modification.
308 //
309 void G4UPolycone::ComputeDimensions(G4VPVParameterisation* p,
310                                     const G4int n,
311                                     const G4VPhysicalVolume* pRep)
312 {
313   p->ComputeDimensions(*(G4Polycone*)this,n,pRep);
314 }
315 
316 
317 //////////////////////////////////////////////////////////////////////////
318 //
319 // Make a clone of the object
320 
321 G4VSolid* G4UPolycone::Clone() const
322 {
323   return new G4UPolycone(*this);
324 }
325 
326 //////////////////////////////////////////////////////////////////////////
327 //
328 // Get bounding box
329 
330 void G4UPolycone::BoundingLimits(G4ThreeVector& pMin,
331                                  G4ThreeVector& pMax) const
332 {
333   static G4bool checkBBox = true;
334   static G4bool checkPhi  = true;
335 
336   G4double rmin = kInfinity, rmax = -kInfinity;
337   G4double zmin = kInfinity, zmax = -kInfinity;
338 
339   for (G4int i=0; i<GetNumRZCorner(); ++i)
340   {
341     G4PolyconeSideRZ corner = GetCorner(i);
342     if (corner.r < rmin) rmin = corner.r;
343     if (corner.r > rmax) rmax = corner.r;
344     if (corner.z < zmin) zmin = corner.z;
345     if (corner.z > zmax) zmax = corner.z;
346   }
347 
348   if (IsOpen())
349   {
350     G4TwoVector vmin,vmax;
351     G4GeomTools::DiskExtent(rmin,rmax,
352                             GetSinStartPhi(),GetCosStartPhi(),
353                             GetSinEndPhi(),GetCosEndPhi(),
354                             vmin,vmax);
355     pMin.set(vmin.x(),vmin.y(),zmin);
356     pMax.set(vmax.x(),vmax.y(),zmax);
357   }
358   else
359   {
360     pMin.set(-rmax,-rmax, zmin);
361     pMax.set( rmax, rmax, zmax);
362   }
363 
364   // Check correctness of the bounding box
365   //
366   if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
367   {
368     std::ostringstream message;
369     message << "Bad bounding box (min >= max) for solid: "
370             << GetName() << " !"
371             << "\npMin = " << pMin
372             << "\npMax = " << pMax;
373     G4Exception("G4UPolycone::BoundingLimits()", "GeomMgt0001",
374                 JustWarning, message);
375     StreamInfo(G4cout);
376   }
377 
378   // Check consistency of bounding boxes
379   //
380   if (checkBBox)
381   {
382     U3Vector vmin, vmax;
383     Extent(vmin,vmax);
384     if (std::abs(pMin.x()-vmin.x()) > kCarTolerance ||
385         std::abs(pMin.y()-vmin.y()) > kCarTolerance ||
386         std::abs(pMin.z()-vmin.z()) > kCarTolerance ||
387         std::abs(pMax.x()-vmax.x()) > kCarTolerance ||
388         std::abs(pMax.y()-vmax.y()) > kCarTolerance ||
389         std::abs(pMax.z()-vmax.z()) > kCarTolerance)
390     {
391       std::ostringstream message;
392       message << "Inconsistency in bounding boxes for solid: "
393               << GetName() << " !"
394               << "\nBBox min: wrapper = " << pMin << " solid = " << vmin
395               << "\nBBox max: wrapper = " << pMax << " solid = " << vmax;
396       G4Exception("G4UPolycone::BoundingLimits()", "GeomMgt0001",
397                   JustWarning, message);
398       checkBBox = false;
399     }
400   }
401 
402   // Check consistency of angles
403   //
404   if (checkPhi)
405   {
406     if (GetStartPhi() != Base_t::GetStartPhi() ||
407         GetEndPhi()   != Base_t::GetEndPhi()   ||
408         IsOpen()      != (Base_t::GetDeltaPhi() < twopi))
409     {
410       std::ostringstream message;
411       message << "Inconsistency in Phi angles or # of sides for solid: "
412               << GetName() << " !"
413               << "\nPhi start  : wrapper = " << GetStartPhi()
414               << " solid = " <<     Base_t::GetStartPhi()
415               << "\nPhi end    : wrapper = " << GetEndPhi()
416               << " solid = " <<     Base_t::GetEndPhi()
417               << "\nPhi is open: wrapper = " << (IsOpen() ? "true" : "false")
418               << " solid = "
419               << ((Base_t::GetDeltaPhi() < twopi) ? "true" : "false");
420       G4Exception("G4UPolycone::BoundingLimits()", "GeomMgt0001",
421                   JustWarning, message);
422       checkPhi = false;
423     }
424   }
425 }
426 
427 //////////////////////////////////////////////////////////////////////////
428 //
429 // Calculate extent under transform and specified limit
430 
431 G4bool G4UPolycone::CalculateExtent(const EAxis pAxis,
432                                     const G4VoxelLimits& pVoxelLimit,
433                                     const G4AffineTransform& pTransform,
434                                           G4double& pMin, G4double& pMax) const
435 {
436   G4ThreeVector bmin, bmax;
437   G4bool exist;
438 
439   // Check bounding box (bbox)
440   //
441   BoundingLimits(bmin,bmax);
442   G4BoundingEnvelope bbox(bmin,bmax);
443 #ifdef G4BBOX_EXTENT
444   return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
445 #endif
446   if (bbox.BoundingBoxVsVoxelLimits(pAxis,pVoxelLimit,pTransform,pMin,pMax))
447   {
448     return exist = pMin < pMax;
449   }
450 
451   // To find the extent, RZ contour of the polycone is subdivided
452   // in triangles. The extent is calculated as cumulative extent of
453   // all sub-polycones formed by rotation of triangles around Z
454   //
455   G4TwoVectorList contourRZ;
456   G4TwoVectorList triangles;
457   std::vector<G4int> iout;
458   G4double eminlim = pVoxelLimit.GetMinExtent(pAxis);
459   G4double emaxlim = pVoxelLimit.GetMaxExtent(pAxis);
460 
461   // get RZ contour, ensure anticlockwise order of corners
462   for (G4int i=0; i<GetNumRZCorner(); ++i)
463   {
464     G4PolyconeSideRZ corner = GetCorner(i);
465     contourRZ.emplace_back(corner.r,corner.z);
466   }
467   G4GeomTools::RemoveRedundantVertices(contourRZ,iout,2*kCarTolerance);
468   G4double area = G4GeomTools::PolygonArea(contourRZ);
469   if (area < 0.) std::reverse(contourRZ.begin(),contourRZ.end());
470 
471   // triangulate RZ countour
472   if (!G4GeomTools::TriangulatePolygon(contourRZ,triangles))
473   {
474     std::ostringstream message;
475     message << "Triangulation of RZ contour has failed for solid: "
476             << GetName() << " !"
477             << "\nExtent has been calculated using boundary box";
478     G4Exception("G4UPolycone::CalculateExtent()",
479                 "GeomMgt1002", JustWarning, message);
480     return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
481   }
482 
483   // set trigonometric values
484   const G4int NSTEPS = 24;            // number of steps for whole circle
485   G4double astep  = twopi/NSTEPS;     // max angle for one step
486 
487   G4double sphi   = GetStartPhi();
488   G4double ephi   = GetEndPhi();
489   G4double dphi   = IsOpen() ? ephi-sphi : twopi;
490   G4int    ksteps = (dphi <= astep) ? 1 : (G4int)((dphi-deg)/astep) + 1;
491   G4double ang    = dphi/ksteps;
492 
493   G4double sinHalf = std::sin(0.5*ang);
494   G4double cosHalf = std::cos(0.5*ang);
495   G4double sinStep = 2.*sinHalf*cosHalf;
496   G4double cosStep = 1. - 2.*sinHalf*sinHalf;
497 
498   G4double sinStart = GetSinStartPhi();
499   G4double cosStart = GetCosStartPhi();
500   G4double sinEnd   = GetSinEndPhi();
501   G4double cosEnd   = GetCosEndPhi();
502 
503   // define vectors and arrays
504   std::vector<const G4ThreeVectorList *> polygons;
505   polygons.resize(ksteps+2);
506   G4ThreeVectorList pols[NSTEPS+2];
507   for (G4int k=0; k<ksteps+2; ++k) pols[k].resize(6);
508   for (G4int k=0; k<ksteps+2; ++k) polygons[k] = &pols[k];
509   G4double r0[6],z0[6]; // contour with original edges of triangle
510   G4double r1[6];       // shifted radii of external edges of triangle
511 
512   // main loop along triangles
513   pMin = kInfinity;
514   pMax =-kInfinity;
515   G4int ntria = triangles.size()/3;
516   for (G4int i=0; i<ntria; ++i)
517   {
518     G4int i3 = i*3;
519     for (G4int k=0; k<3; ++k)
520     {
521       G4int e0 = i3+k, e1 = (k<2) ? e0+1 : i3;
522       G4int k2 = k*2;
523       // set contour with original edges of triangle
524       r0[k2+0] = triangles[e0].x(); z0[k2+0] = triangles[e0].y();
525       r0[k2+1] = triangles[e1].x(); z0[k2+1] = triangles[e1].y();
526       // set shifted radii
527       r1[k2+0] = r0[k2+0];
528       r1[k2+1] = r0[k2+1];
529       if (z0[k2+1] - z0[k2+0] <= 0) continue;
530       r1[k2+0] /= cosHalf;
531       r1[k2+1] /= cosHalf;
532     }
533 
534     // rotate countour, set sequence of 6-sided polygons
535     G4double sinCur = sinStart*cosHalf + cosStart*sinHalf;
536     G4double cosCur = cosStart*cosHalf - sinStart*sinHalf;
537     for (G4int j=0; j<6; ++j) pols[0][j].set(r0[j]*cosStart,r0[j]*sinStart,z0[j]);
538     for (G4int k=1; k<ksteps+1; ++k)
539     {
540       for (G4int j=0; j<6; ++j) pols[k][j].set(r1[j]*cosCur,r1[j]*sinCur,z0[j]);
541       G4double sinTmp = sinCur;
542       sinCur = sinCur*cosStep + cosCur*sinStep;
543       cosCur = cosCur*cosStep - sinTmp*sinStep;
544     }
545     for (G4int j=0; j<6; ++j) pols[ksteps+1][j].set(r0[j]*cosEnd,r0[j]*sinEnd,z0[j]);
546 
547     // set sub-envelope and adjust extent
548     G4double emin,emax;
549     G4BoundingEnvelope benv(polygons);
550     if (!benv.CalculateExtent(pAxis,pVoxelLimit,pTransform,emin,emax)) continue;
551     if (emin < pMin) pMin = emin;
552     if (emax > pMax) pMax = emax;
553     if (eminlim > pMin && emaxlim < pMax) return true; // max possible extent
554   }
555   return (pMin < pMax);
556 }
557 
558 ////////////////////////////////////////////////////////////////////////
559 //
560 // CreatePolyhedron
561 //
562 G4Polyhedron* G4UPolycone::CreatePolyhedron() const
563 {
564   return new G4PolyhedronPcon(wrStart, wrDelta, rzcorners);
565 }
566 
567 #endif  // G4GEOM_USE_USOLIDS
568