Geant4 Cross Reference |
1 // 1 // 2 // ******************************************* 2 // ******************************************************************** 3 // * License and Disclaimer 3 // * License and Disclaimer * 4 // * 4 // * * 5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. 9 // * include a list of copyright holders. * 10 // * 10 // * * 11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. * 17 // * 17 // * * 18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************* 24 // ******************************************************************** >> 25 // >> 26 // >> 27 // $Id:$ >> 28 // 25 // 29 // 26 // Implementation of G4UExtrudedSolid wrapper 30 // Implementation of G4UExtrudedSolid wrapper class 27 // << 28 // 17.11.17 G.Cosmo, CERN << 29 // ------------------------------------------- 31 // -------------------------------------------------------------------- 30 32 31 #include "G4ExtrudedSolid.hh" 33 #include "G4ExtrudedSolid.hh" 32 #include "G4UExtrudedSolid.hh" 34 #include "G4UExtrudedSolid.hh" 33 35 34 #if ( defined(G4GEOM_USE_USOLIDS) || defined(G 36 #if ( defined(G4GEOM_USE_USOLIDS) || defined(G4GEOM_USE_PARTIAL_USOLIDS) ) 35 37 36 #include "G4GeomTools.hh" 38 #include "G4GeomTools.hh" 37 #include "G4AffineTransform.hh" 39 #include "G4AffineTransform.hh" 38 #include "G4BoundingEnvelope.hh" 40 #include "G4BoundingEnvelope.hh" 39 41 >> 42 #include "G4PolyhedronArbitrary.hh" >> 43 40 ////////////////////////////////////////////// 44 //////////////////////////////////////////////////////////////////////// 41 // 45 // 42 // Constructors 46 // Constructors 43 // 47 // 44 G4UExtrudedSolid::G4UExtrudedSolid(const G4Str 48 G4UExtrudedSolid::G4UExtrudedSolid(const G4String& name, 45 const std:: << 49 std::vector<G4TwoVector> polygon, 46 const std:: << 50 std::vector<ZSection> zsections) 47 : Base_t(name) // General constructor << 51 : G4USolid(name, new UExtrudedSolid()) 48 { 52 { 49 unsigned int nVertices = polygon.size(); << 53 GetShape()->SetName(name); 50 unsigned int nSections = zsections.size(); << 54 std::vector<UVector2> pvec; 51 << 55 for (unsigned int i=0; i<polygon.size(); ++i) 52 auto vertices = new vecgeom::XtruVertex2[nVe << 53 auto sections = new vecgeom::XtruSection[nSe << 54 << 55 for (unsigned int i = 0; i < nVertices; ++i) << 56 { 56 { 57 vertices[i].x = polygon[i].x(); << 57 pvec.push_back(UVector2(polygon[i].x(), polygon[i].y())); 58 vertices[i].y = polygon[i].y(); << 59 } 58 } 60 for (unsigned int i = 0; i < nSections; ++i) << 59 std::vector<UExtrudedSolid::ZSection> svec; >> 60 for (unsigned int i=0; i<zsections.size(); ++i) 61 { 61 { 62 sections[i].fOrigin.Set(zsections[i].fOffs << 62 ZSection sec = zsections[i]; 63 zsections[i].fOffs << 63 svec.push_back(UExtrudedSolid::ZSection(sec.fZ, 64 zsections[i].fZ); << 64 UVector2(sec.fOffset.x(), sec.fOffset.y()), sec.fScale)); 65 sections[i].fScale = zsections[i].fScale; << 66 } 65 } 67 Base_t::Initialize(nVertices, vertices, nSec << 66 GetShape()->Initialise(pvec, svec); 68 delete[] vertices; << 69 delete[] sections; << 70 } 67 } 71 68 72 69 73 G4UExtrudedSolid::G4UExtrudedSolid(const G4Str 70 G4UExtrudedSolid::G4UExtrudedSolid(const G4String& name, 74 const std:: << 71 std::vector<G4TwoVector> polygon, 75 G4double 72 G4double halfZ, 76 const G4Two << 73 G4TwoVector off1, G4double scale1, 77 const G4Two << 74 G4TwoVector off2, G4double scale2) 78 : Base_t(name) // Special constructor for 2 << 75 : G4USolid(name, new UExtrudedSolid()) 79 { << 76 { 80 unsigned int nVertices = polygon.size(); << 77 GetShape()->SetName(name); 81 unsigned int nSections = 2; << 78 std::vector<UVector2> pvec; 82 << 79 for (unsigned int i=0; i<polygon.size(); ++i) 83 auto vertices = new vecgeom::XtruVertex2[nVe << 80 { 84 auto sections = new vecgeom::XtruSection[nSe << 81 pvec.push_back(UVector2(polygon[i].x(), polygon[i].y())); 85 << 82 } 86 for (unsigned int i = 0; i < nVertices; ++i) << 83 GetShape()->Initialise(pvec, halfZ, UVector2(off1.x(), off1.y()), scale1, 87 { << 84 UVector2(off2.x(), off2.y()), scale2); 88 vertices[i].x = polygon[i].x(); << 89 vertices[i].y = polygon[i].y(); << 90 } << 91 sections[0].fOrigin.Set(off1.x(), off1.y(), << 92 sections[0].fScale = scale1; << 93 sections[1].fOrigin.Set(off2.x(), off2.y(), << 94 sections[1].fScale = scale2; << 95 Base_t::Initialize(nVertices, vertices, nSec << 96 delete[] vertices; << 97 delete[] sections; << 98 } 85 } 99 86 100 ////////////////////////////////////////////// 87 //////////////////////////////////////////////////////////////////////// 101 // 88 // 102 // Fake default constructor - sets only member 89 // Fake default constructor - sets only member data and allocates memory 103 // for usage restri 90 // for usage restricted to object persistency. 104 // 91 // 105 G4UExtrudedSolid::G4UExtrudedSolid(__void__& a 92 G4UExtrudedSolid::G4UExtrudedSolid(__void__& a) 106 : Base_t(a) << 93 : G4USolid(a) 107 { 94 { 108 } 95 } 109 96 110 97 111 ////////////////////////////////////////////// 98 ////////////////////////////////////////////////////////////////////////// 112 // 99 // 113 // Destructor 100 // Destructor 114 // 101 // 115 G4UExtrudedSolid::~G4UExtrudedSolid() = defaul << 102 G4UExtrudedSolid::~G4UExtrudedSolid() >> 103 { >> 104 } 116 105 117 106 118 ////////////////////////////////////////////// 107 ////////////////////////////////////////////////////////////////////////// 119 // 108 // 120 // Copy constructor 109 // Copy constructor 121 // 110 // 122 G4UExtrudedSolid::G4UExtrudedSolid(const G4UEx 111 G4UExtrudedSolid::G4UExtrudedSolid(const G4UExtrudedSolid &source) 123 : Base_t(source) << 112 : G4USolid(source) 124 { 113 { 125 } 114 } 126 115 127 116 128 ////////////////////////////////////////////// 117 ////////////////////////////////////////////////////////////////////////// 129 // 118 // 130 // Assignment operator 119 // Assignment operator 131 // 120 // 132 G4UExtrudedSolid& 121 G4UExtrudedSolid& 133 G4UExtrudedSolid::operator=(const G4UExtrudedS 122 G4UExtrudedSolid::operator=(const G4UExtrudedSolid &source) 134 { 123 { 135 if (this == &source) return *this; 124 if (this == &source) return *this; 136 125 137 Base_t::operator=( source ); << 126 G4USolid::operator=( source ); 138 127 139 return *this; 128 return *this; 140 } 129 } 141 130 142 131 143 ////////////////////////////////////////////// 132 ////////////////////////////////////////////////////////////////////////// 144 // 133 // 145 // Accessors 134 // Accessors 146 135 147 G4int G4UExtrudedSolid::GetNofVertices() const 136 G4int G4UExtrudedSolid::GetNofVertices() const 148 { 137 { 149 return Base_t::GetNVertices(); << 138 return GetShape()->GetNofVertices(); 150 } 139 } 151 << 152 G4TwoVector G4UExtrudedSolid::GetVertex(G4int 140 G4TwoVector G4UExtrudedSolid::GetVertex(G4int i) const 153 { 141 { 154 G4double xx, yy; << 142 UVector2 v = GetShape()->GetVertex(i); 155 Base_t::GetVertex(i, xx, yy); << 143 return G4TwoVector(v.x, v.y); 156 return { xx, yy }; << 157 } 144 } 158 << 159 std::vector<G4TwoVector> G4UExtrudedSolid::Get 145 std::vector<G4TwoVector> G4UExtrudedSolid::GetPolygon() const 160 { 146 { 161 std::vector<G4TwoVector> pol; << 147 std::vector<UVector2> pol = GetShape()->GetPolygon(); 162 for (unsigned int i = 0; i < Base_t::GetNVer << 148 std::vector<G4TwoVector> v; >> 149 for (unsigned int i=0; i<pol.size(); ++i) 163 { 150 { 164 pol.push_back(GetVertex(i)); << 151 v.push_back(G4TwoVector(pol[i].x, pol[i].y)); 165 } 152 } 166 return pol; << 153 return v; 167 } 154 } 168 << 169 G4int G4UExtrudedSolid::GetNofZSections() cons 155 G4int G4UExtrudedSolid::GetNofZSections() const 170 { 156 { 171 return Base_t::GetNSections(); << 157 return GetShape()->GetNofZSections(); 172 } 158 } 173 << 174 G4UExtrudedSolid::ZSection G4UExtrudedSolid::G 159 G4UExtrudedSolid::ZSection G4UExtrudedSolid::GetZSection(G4int i) const 175 { 160 { 176 vecgeom::XtruSection sect = Base_t::GetSecti << 161 return ZSection(GetShape()->GetZSection(i)); 177 return { sect.fOrigin[2], << 178 G4TwoVector(sect.fOrigin[0], sect.f << 179 sect.fScale }; << 180 } 162 } 181 << 182 std::vector<G4UExtrudedSolid::ZSection> G4UExt 163 std::vector<G4UExtrudedSolid::ZSection> G4UExtrudedSolid::GetZSections() const 183 { 164 { 184 std::vector<G4UExtrudedSolid::ZSection> sect << 165 std::vector<UExtrudedSolid::ZSection> sv = GetShape()->GetZSections(); 185 for (unsigned int i = 0; i < Base_t::GetNSec << 166 std::vector<G4UExtrudedSolid::ZSection> vec; >> 167 for (unsigned int i=0; i<sv.size(); ++i) 186 { 168 { 187 vecgeom::XtruSection sect = Base_t::GetSec << 169 vec.push_back(ZSection(sv[i])); 188 sections.emplace_back(sect.fOrigin[2], << 189 G4TwoVector(sect.fOr << 190 sect.fScale); << 191 } 170 } 192 return sections; << 171 return vec; 193 } 172 } 194 173 195 174 196 ////////////////////////////////////////////// 175 /////////////////////////////////////////////////////////////////////////////// 197 // 176 // 198 // Get bounding box 177 // Get bounding box 199 178 200 void G4UExtrudedSolid::BoundingLimits(G4ThreeV << 179 void G4UExtrudedSolid::Extent(G4ThreeVector& pMin, G4ThreeVector& pMax) const 201 G4ThreeV << 202 { 180 { 203 static G4bool checkBBox = true; 181 static G4bool checkBBox = true; 204 182 205 G4double xmin0 = kInfinity, xmax0 = -kInfini 183 G4double xmin0 = kInfinity, xmax0 = -kInfinity; 206 G4double ymin0 = kInfinity, ymax0 = -kInfini 184 G4double ymin0 = kInfinity, ymax0 = -kInfinity; 207 185 208 for (G4int i=0; i<GetNofVertices(); ++i) 186 for (G4int i=0; i<GetNofVertices(); ++i) 209 { 187 { 210 G4TwoVector vertex = GetVertex(i); 188 G4TwoVector vertex = GetVertex(i); 211 G4double x = vertex.x(); 189 G4double x = vertex.x(); 212 if (x < xmin0) xmin0 = x; 190 if (x < xmin0) xmin0 = x; 213 if (x > xmax0) xmax0 = x; 191 if (x > xmax0) xmax0 = x; 214 G4double y = vertex.y(); 192 G4double y = vertex.y(); 215 if (y < ymin0) ymin0 = y; 193 if (y < ymin0) ymin0 = y; 216 if (y > ymax0) ymax0 = y; 194 if (y > ymax0) ymax0 = y; 217 } 195 } 218 196 219 G4double xmin = kInfinity, xmax = -kInfinity 197 G4double xmin = kInfinity, xmax = -kInfinity; 220 G4double ymin = kInfinity, ymax = -kInfinity 198 G4double ymin = kInfinity, ymax = -kInfinity; 221 199 222 G4int nsect = GetNofZSections(); 200 G4int nsect = GetNofZSections(); 223 for (G4int i=0; i<nsect; ++i) 201 for (G4int i=0; i<nsect; ++i) 224 { 202 { 225 ZSection zsect = GetZSection(i); 203 ZSection zsect = GetZSection(i); 226 G4double dx = zsect.fOffset.x(); 204 G4double dx = zsect.fOffset.x(); 227 G4double dy = zsect.fOffset.y(); 205 G4double dy = zsect.fOffset.y(); 228 G4double scale = zsect.fScale; 206 G4double scale = zsect.fScale; 229 xmin = std::min(xmin,xmin0*scale+dx); 207 xmin = std::min(xmin,xmin0*scale+dx); 230 xmax = std::max(xmax,xmax0*scale+dx); 208 xmax = std::max(xmax,xmax0*scale+dx); 231 ymin = std::min(ymin,ymin0*scale+dy); 209 ymin = std::min(ymin,ymin0*scale+dy); 232 ymax = std::max(ymax,ymax0*scale+dy); 210 ymax = std::max(ymax,ymax0*scale+dy); 233 } 211 } 234 212 235 G4double zmin = GetZSection(0).fZ; 213 G4double zmin = GetZSection(0).fZ; 236 G4double zmax = GetZSection(nsect-1).fZ; 214 G4double zmax = GetZSection(nsect-1).fZ; 237 215 238 pMin.set(xmin,ymin,zmin); 216 pMin.set(xmin,ymin,zmin); 239 pMax.set(xmax,ymax,zmax); 217 pMax.set(xmax,ymax,zmax); 240 218 241 // Check correctness of the bounding box 219 // Check correctness of the bounding box 242 // 220 // 243 if (pMin.x() >= pMax.x() || pMin.y() >= pMax 221 if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z()) 244 { 222 { 245 std::ostringstream message; 223 std::ostringstream message; 246 message << "Bad bounding box (min >= max) 224 message << "Bad bounding box (min >= max) for solid: " 247 << GetName() << " !" 225 << GetName() << " !" 248 << "\npMin = " << pMin 226 << "\npMin = " << pMin 249 << "\npMax = " << pMax; 227 << "\npMax = " << pMax; 250 G4Exception("G4UExtrudedSolid::BoundingLim << 228 G4Exception("G4UExtrudedSolid::Extent()", "GeomMgt0001", 251 JustWarning, message); 229 JustWarning, message); 252 StreamInfo(G4cout); 230 StreamInfo(G4cout); 253 } 231 } 254 232 255 // Check consistency of bounding boxes 233 // Check consistency of bounding boxes 256 // 234 // 257 if (checkBBox) 235 if (checkBBox) 258 { 236 { 259 U3Vector vmin, vmax; << 237 UVector3 vmin, vmax; 260 Base_t::Extent(vmin,vmax); << 238 GetShape()->Extent(vmin,vmax); 261 if (std::abs(pMin.x()-vmin.x()) > kCarTole 239 if (std::abs(pMin.x()-vmin.x()) > kCarTolerance || 262 std::abs(pMin.y()-vmin.y()) > kCarTole 240 std::abs(pMin.y()-vmin.y()) > kCarTolerance || 263 std::abs(pMin.z()-vmin.z()) > kCarTole 241 std::abs(pMin.z()-vmin.z()) > kCarTolerance || 264 std::abs(pMax.x()-vmax.x()) > kCarTole 242 std::abs(pMax.x()-vmax.x()) > kCarTolerance || 265 std::abs(pMax.y()-vmax.y()) > kCarTole 243 std::abs(pMax.y()-vmax.y()) > kCarTolerance || 266 std::abs(pMax.z()-vmax.z()) > kCarTole 244 std::abs(pMax.z()-vmax.z()) > kCarTolerance) 267 { 245 { 268 std::ostringstream message; 246 std::ostringstream message; 269 message << "Inconsistency in bounding bo 247 message << "Inconsistency in bounding boxes for solid: " 270 << GetName() << " !" 248 << GetName() << " !" 271 << "\nBBox min: wrapper = " << p 249 << "\nBBox min: wrapper = " << pMin << " solid = " << vmin 272 << "\nBBox max: wrapper = " << p 250 << "\nBBox max: wrapper = " << pMax << " solid = " << vmax; 273 G4Exception("G4UExtrudedSolid::BoundingL << 251 G4Exception("G4UExtrudedSolid::Extent()", "GeomMgt0001", JustWarning, message); 274 JustWarning, message); << 275 checkBBox = false; 252 checkBBox = false; 276 } 253 } 277 } 254 } 278 } 255 } 279 256 280 257 281 ////////////////////////////////////////////// 258 ////////////////////////////////////////////////////////////////////////////// 282 // 259 // 283 // Calculate extent under transform and specif 260 // Calculate extent under transform and specified limit 284 261 285 G4bool 262 G4bool 286 G4UExtrudedSolid::CalculateExtent(const EAxis 263 G4UExtrudedSolid::CalculateExtent(const EAxis pAxis, 287 const G4Voxe 264 const G4VoxelLimits& pVoxelLimit, 288 const G4Affi 265 const G4AffineTransform& pTransform, 289 G4doub 266 G4double& pMin, G4double& pMax) const 290 { 267 { 291 G4ThreeVector bmin, bmax; 268 G4ThreeVector bmin, bmax; 292 G4bool exist; 269 G4bool exist; 293 270 294 // Check bounding box (bbox) 271 // Check bounding box (bbox) 295 // 272 // 296 BoundingLimits(bmin,bmax); << 273 Extent(bmin,bmax); 297 G4BoundingEnvelope bbox(bmin,bmax); 274 G4BoundingEnvelope bbox(bmin,bmax); 298 #ifdef G4BBOX_EXTENT 275 #ifdef G4BBOX_EXTENT 299 return bbox.CalculateExtent(pAxis,pVoxelLimi << 276 if (true) return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax); 300 #endif 277 #endif 301 if (bbox.BoundingBoxVsVoxelLimits(pAxis,pVox 278 if (bbox.BoundingBoxVsVoxelLimits(pAxis,pVoxelLimit,pTransform,pMin,pMax)) 302 { 279 { 303 return exist = pMin < pMax; << 280 return exist = (pMin < pMax) ? true : false; 304 } 281 } 305 282 306 // To find the extent, the base polygon is s 283 // To find the extent, the base polygon is subdivided in triangles. 307 // The extent is calculated as cumulative ex 284 // The extent is calculated as cumulative extent of the parts 308 // formed by extrusion of the triangles 285 // formed by extrusion of the triangles 309 // 286 // 310 G4TwoVectorList basePolygon = GetPolygon(); 287 G4TwoVectorList basePolygon = GetPolygon(); 311 G4TwoVectorList triangles; 288 G4TwoVectorList triangles; 312 G4double eminlim = pVoxelLimit.GetMinExtent( 289 G4double eminlim = pVoxelLimit.GetMinExtent(pAxis); 313 G4double emaxlim = pVoxelLimit.GetMaxExtent( 290 G4double emaxlim = pVoxelLimit.GetMaxExtent(pAxis); 314 291 315 // triangulate the base polygon 292 // triangulate the base polygon 316 if (!G4GeomTools::TriangulatePolygon(basePol 293 if (!G4GeomTools::TriangulatePolygon(basePolygon,triangles)) 317 { 294 { 318 std::ostringstream message; 295 std::ostringstream message; 319 message << "Triangulation of the base poly 296 message << "Triangulation of the base polygon has failed for solid: " 320 << GetName() << " !" 297 << GetName() << " !" 321 << "\nExtent has been calculated u 298 << "\nExtent has been calculated using boundary box"; 322 G4Exception("G4UExtrudedSolid::CalculateEx 299 G4Exception("G4UExtrudedSolid::CalculateExtent()", 323 "GeomMgt1002",JustWarning,mess 300 "GeomMgt1002",JustWarning,message); 324 return bbox.CalculateExtent(pAxis,pVoxelLi 301 return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax); 325 } 302 } 326 303 327 // allocate vector lists 304 // allocate vector lists 328 G4int nsect = GetNofZSections(); 305 G4int nsect = GetNofZSections(); 329 std::vector<const G4ThreeVectorList *> polyg 306 std::vector<const G4ThreeVectorList *> polygons; 330 polygons.resize(nsect); 307 polygons.resize(nsect); 331 for (G4int k=0; k<nsect; ++k) { polygons[k] 308 for (G4int k=0; k<nsect; ++k) { polygons[k] = new G4ThreeVectorList(3); } 332 309 333 // main loop along triangles 310 // main loop along triangles 334 pMin = kInfinity; 311 pMin = kInfinity; 335 pMax = -kInfinity; 312 pMax = -kInfinity; 336 G4int ntria = triangles.size()/3; 313 G4int ntria = triangles.size()/3; 337 for (G4int i=0; i<ntria; ++i) 314 for (G4int i=0; i<ntria; ++i) 338 { 315 { 339 G4int i3 = i*3; 316 G4int i3 = i*3; 340 for (G4int k=0; k<nsect; ++k) // extrude t 317 for (G4int k=0; k<nsect; ++k) // extrude triangle 341 { 318 { 342 ZSection zsect = GetZSection(k); 319 ZSection zsect = GetZSection(k); 343 G4double z = zsect.fZ; 320 G4double z = zsect.fZ; 344 G4double dx = zsect.fOffset.x(); 321 G4double dx = zsect.fOffset.x(); 345 G4double dy = zsect.fOffset.y(); 322 G4double dy = zsect.fOffset.y(); 346 G4double scale = zsect.fScale; 323 G4double scale = zsect.fScale; 347 324 348 auto ptr = const_cast<G4ThreeVectorList* << 325 G4ThreeVectorList* ptr = const_cast<G4ThreeVectorList*>(polygons[k]); 349 auto iter = ptr->begin(); << 326 G4ThreeVectorList::iterator iter = ptr->begin(); 350 G4double x0 = triangles[i3+0].x()*scale+ 327 G4double x0 = triangles[i3+0].x()*scale+dx; 351 G4double y0 = triangles[i3+0].y()*scale+ 328 G4double y0 = triangles[i3+0].y()*scale+dy; 352 iter->set(x0,y0,z); 329 iter->set(x0,y0,z); 353 iter++; 330 iter++; 354 G4double x1 = triangles[i3+1].x()*scale+ 331 G4double x1 = triangles[i3+1].x()*scale+dx; 355 G4double y1 = triangles[i3+1].y()*scale+ 332 G4double y1 = triangles[i3+1].y()*scale+dy; 356 iter->set(x1,y1,z); 333 iter->set(x1,y1,z); 357 iter++; 334 iter++; 358 G4double x2 = triangles[i3+2].x()*scale+ 335 G4double x2 = triangles[i3+2].x()*scale+dx; 359 G4double y2 = triangles[i3+2].y()*scale+ 336 G4double y2 = triangles[i3+2].y()*scale+dy; 360 iter->set(x2,y2,z); 337 iter->set(x2,y2,z); 361 } 338 } 362 339 363 // set sub-envelope and adjust extent 340 // set sub-envelope and adjust extent 364 G4double emin,emax; 341 G4double emin,emax; 365 G4BoundingEnvelope benv(polygons); 342 G4BoundingEnvelope benv(polygons); 366 if (!benv.CalculateExtent(pAxis,pVoxelLimi 343 if (!benv.CalculateExtent(pAxis,pVoxelLimit,pTransform,emin,emax)) continue; 367 if (emin < pMin) pMin = emin; 344 if (emin < pMin) pMin = emin; 368 if (emax > pMax) pMax = emax; 345 if (emax > pMax) pMax = emax; 369 if (eminlim > pMin && emaxlim < pMax) brea 346 if (eminlim > pMin && emaxlim < pMax) break; // max possible extent 370 } 347 } 371 // free memory 348 // free memory 372 for (G4int k=0; k<nsect; ++k) { delete polyg << 349 for (G4int k=0; k<nsect; ++k) { delete polygons[k]; polygons[k]=0;} 373 return (pMin < pMax); 350 return (pMin < pMax); 374 } 351 } 375 352 376 353 377 ////////////////////////////////////////////// 354 /////////////////////////////////////////////////////////////////////////////// 378 // 355 // 379 // CreatePolyhedron() 356 // CreatePolyhedron() 380 // 357 // 381 G4Polyhedron* G4UExtrudedSolid::CreatePolyhedr 358 G4Polyhedron* G4UExtrudedSolid::CreatePolyhedron () const 382 { 359 { 383 unsigned int nFacets = Base_t::GetStruct().f << 360 G4int nFacets = GetShape()->GetNumberOfFacets(); 384 unsigned int nVertices = Base_t::GetStruct() << 361 G4int nVertices = 0; 385 << 362 for (G4int l = 0; l<nFacets; ++l) // compute total number of vertices first 386 auto polyhedron = new G4Polyhedron(nVertices << 387 << 388 // Copy vertices << 389 for (unsigned int i = 0; i < nVertices; ++i) << 390 { 363 { 391 U3Vector v = Base_t::GetStruct().fTslHelpe << 364 VUFacet* facet = GetShape()->GetFacet(l); 392 polyhedron->SetVertex(i+1, G4ThreeVector(v << 365 G4int n = facet->GetNumberOfVertices(); >> 366 nVertices += n; 393 } 367 } 394 368 395 // Copy facets << 369 G4PolyhedronArbitrary *polyhedron = 396 for (unsigned int i = 0; i < nFacets; ++i) << 370 new G4PolyhedronArbitrary (nVertices,nFacets); >> 371 >> 372 for (G4int i = 0; i<nFacets; ++i) 397 { 373 { 398 // Facets are only triangular in VecGeom << 374 VUFacet* facet = GetShape()->GetFacet(i); 399 G4int i1 = Base_t::GetStruct().fTslHelper. << 375 G4int v[4]; 400 G4int i2 = Base_t::GetStruct().fTslHelper. << 376 G4int n = facet->GetNumberOfVertices(); 401 G4int i3 = Base_t::GetStruct().fTslHelper. << 377 for (G4int m = 0; m<n; ++m) 402 polyhedron->SetFacet(i+1, i1, i2, i3); << 378 { >> 379 UVector3 vtx = facet->GetVertex(m); >> 380 polyhedron->AddVertex(G4ThreeVector(vtx.x(), vtx.y(), vtx.z())); >> 381 } >> 382 if (n > 4) n = 4; >> 383 else if (n == 3) v[3] = 0; >> 384 for (G4int j=0; j<n; ++j) >> 385 { >> 386 G4int k = facet->GetVertexIndex(j); >> 387 v[j] = k+1; >> 388 } >> 389 polyhedron->AddFacet(v[0],v[1],v[2],v[3]); 403 } 390 } 404 polyhedron->SetReferences(); << 391 polyhedron->SetReferences(); 405 392 406 return polyhedron; << 393 return (G4Polyhedron*) polyhedron; 407 } 394 } 408 395 409 #endif // G4GEOM_USE_USOLIDS 396 #endif // G4GEOM_USE_USOLIDS 410 397