Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
HepPolyhedron.cc
Go to the documentation of this file.
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// $Id$
28//
29//
30//
31// G4 Polyhedron library
32//
33// History:
34// 23.07.96 E.Chernyaev <[email protected]> - initial version
35//
36// 30.09.96 E.Chernyaev
37// - added GetNextVertexIndex, GetVertex by Yasuhide Sawada
38// - added GetNextUnitNormal, GetNextEdgeIndeces, GetNextEdge
39//
40// 15.12.96 E.Chernyaev
41// - added GetNumberOfRotationSteps, RotateEdge, RotateAroundZ, SetReferences
42// - rewritten G4PolyhedronCons;
43// - added G4PolyhedronPara, ...Trap, ...Pgon, ...Pcon, ...Sphere, ...Torus
44//
45// 01.06.97 E.Chernyaev
46// - modified RotateAroundZ, added SetSideFacets
47//
48// 19.03.00 E.Chernyaev
49// - implemented boolean operations (add, subtract, intersect) on polyhedra;
50//
51// 25.05.01 E.Chernyaev
52// - added GetSurfaceArea() and GetVolume();
53//
54// 05.11.02 E.Chernyaev
55// - added createTwistedTrap() and createPolyhedron();
56//
57// 20.06.05 G.Cosmo
58// - added HepPolyhedronEllipsoid;
59//
60// 18.07.07 T.Nikitin
61// - added HepParaboloid;
62
63#include "HepPolyhedron.h"
64#include "G4PhysicalConstants.hh"
65#include "G4Vector3D.hh"
66
67#include <cstdlib> // Required on some compilers for std::abs(int) ...
68#include <cmath>
69
70using CLHEP::perMillion;
71using CLHEP::deg;
72using CLHEP::pi;
73using CLHEP::twopi;
74
75/***********************************************************************
76 * *
77 * Name: HepPolyhedron operator << Date: 09.05.96 *
78 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
79 * *
80 * Function: Print contents of G4 polyhedron *
81 * *
82 ***********************************************************************/
83std::ostream & operator<<(std::ostream & ostr, const G4Facet & facet) {
84 for (G4int k=0; k<4; k++) {
85 ostr << " " << facet.edge[k].v << "/" << facet.edge[k].f;
86 }
87 return ostr;
88}
89
90std::ostream & operator<<(std::ostream & ostr, const HepPolyhedron & ph) {
91 ostr << std::endl;
92 ostr << "Nverteces=" << ph.nvert << ", Nfacets=" << ph.nface << std::endl;
93 G4int i;
94 for (i=1; i<=ph.nvert; i++) {
95 ostr << "xyz(" << i << ")="
96 << ph.pV[i].x() << ' ' << ph.pV[i].y() << ' ' << ph.pV[i].z()
97 << std::endl;
98 }
99 for (i=1; i<=ph.nface; i++) {
100 ostr << "face(" << i << ")=" << ph.pF[i] << std::endl;
101 }
102 return ostr;
103}
104
105HepPolyhedron::HepPolyhedron(const HepPolyhedron &from)
106/***********************************************************************
107 * *
108 * Name: HepPolyhedron copy constructor Date: 23.07.96 *
109 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
110 * *
111 ***********************************************************************/
112: nvert(0), nface(0), pV(0), pF(0)
113{
114 AllocateMemory(from.nvert, from.nface);
115 for (G4int i=1; i<=nvert; i++) pV[i] = from.pV[i];
116 for (G4int k=1; k<=nface; k++) pF[k] = from.pF[k];
117}
118
119HepPolyhedron & HepPolyhedron::operator=(const HepPolyhedron &from)
120/***********************************************************************
121 * *
122 * Name: HepPolyhedron operator = Date: 23.07.96 *
123 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
124 * *
125 * Function: Copy contents of one polyhedron to another *
126 * *
127 ***********************************************************************/
128{
129 if (this != &from) {
130 AllocateMemory(from.nvert, from.nface);
131 for (G4int i=1; i<=nvert; i++) pV[i] = from.pV[i];
132 for (G4int k=1; k<=nface; k++) pF[k] = from.pF[k];
133 }
134 return *this;
135}
136
137G4int
138HepPolyhedron::FindNeighbour(G4int iFace, G4int iNode, G4int iOrder) const
139/***********************************************************************
140 * *
141 * Name: HepPolyhedron::FindNeighbour Date: 22.11.99 *
142 * Author: E.Chernyaev Revised: *
143 * *
144 * Function: Find neighbouring face *
145 * *
146 ***********************************************************************/
147{
148 G4int i;
149 for (i=0; i<4; i++) {
150 if (iNode == std::abs(pF[iFace].edge[i].v)) break;
151 }
152 if (i == 4) {
153 std::cerr
154 << "HepPolyhedron::FindNeighbour: face " << iFace
155 << " has no node " << iNode
156 << std::endl;
157 return 0;
158 }
159 if (iOrder < 0) {
160 if ( --i < 0) i = 3;
161 if (pF[iFace].edge[i].v == 0) i = 2;
162 }
163 return (pF[iFace].edge[i].v > 0) ? 0 : pF[iFace].edge[i].f;
164}
165
166G4Normal3D HepPolyhedron::FindNodeNormal(G4int iFace, G4int iNode) const
167/***********************************************************************
168 * *
169 * Name: HepPolyhedron::FindNodeNormal Date: 22.11.99 *
170 * Author: E.Chernyaev Revised: *
171 * *
172 * Function: Find normal at given node *
173 * *
174 ***********************************************************************/
175{
176 G4Normal3D normal = GetUnitNormal(iFace);
177 G4int k = iFace, iOrder = 1, n = 1;
178
179 for(;;) {
180 k = FindNeighbour(k, iNode, iOrder);
181 if (k == iFace) break;
182 if (k > 0) {
183 n++;
184 normal += GetUnitNormal(k);
185 }else{
186 if (iOrder < 0) break;
187 k = iFace;
188 iOrder = -iOrder;
189 }
190 }
191 return normal.unit();
192}
193
194G4int HepPolyhedron::GetNumberOfRotationSteps()
195/***********************************************************************
196 * *
197 * Name: HepPolyhedron::GetNumberOfRotationSteps Date: 24.06.97 *
198 * Author: J.Allison (Manchester University) Revised: *
199 * *
200 * Function: Get number of steps for whole circle *
201 * *
202 ***********************************************************************/
203{
204 return fNumberOfRotationSteps;
205}
206
207void HepPolyhedron::SetNumberOfRotationSteps(G4int n)
208/***********************************************************************
209 * *
210 * Name: HepPolyhedron::SetNumberOfRotationSteps Date: 24.06.97 *
211 * Author: J.Allison (Manchester University) Revised: *
212 * *
213 * Function: Set number of steps for whole circle *
214 * *
215 ***********************************************************************/
216{
217 const G4int nMin = 3;
218 if (n < nMin) {
219 std::cerr
220 << "HepPolyhedron::SetNumberOfRotationSteps: attempt to set the\n"
221 << "number of steps per circle < " << nMin << "; forced to " << nMin
222 << std::endl;
223 fNumberOfRotationSteps = nMin;
224 }else{
225 fNumberOfRotationSteps = n;
226 }
227}
228
229void HepPolyhedron::ResetNumberOfRotationSteps()
230/***********************************************************************
231 * *
232 * Name: HepPolyhedron::GetNumberOfRotationSteps Date: 24.06.97 *
233 * Author: J.Allison (Manchester University) Revised: *
234 * *
235 * Function: Reset number of steps for whole circle to default value *
236 * *
237 ***********************************************************************/
238{
239 fNumberOfRotationSteps = DEFAULT_NUMBER_OF_STEPS;
240}
241
242void HepPolyhedron::AllocateMemory(G4int Nvert, G4int Nface)
243/***********************************************************************
244 * *
245 * Name: HepPolyhedron::AllocateMemory Date: 19.06.96 *
246 * Author: E.Chernyaev (IHEP/Protvino) Revised: 05.11.02 *
247 * *
248 * Function: Allocate memory for GEANT4 polyhedron *
249 * *
250 * Input: Nvert - number of nodes *
251 * Nface - number of faces *
252 * *
253 ***********************************************************************/
254{
255 if (nvert == Nvert && nface == Nface) return;
256 if (pV != 0) delete [] pV;
257 if (pF != 0) delete [] pF;
258 if (Nvert > 0 && Nface > 0) {
259 nvert = Nvert;
260 nface = Nface;
261 pV = new G4Point3D[nvert+1];
262 pF = new G4Facet[nface+1];
263 }else{
264 nvert = 0; nface = 0; pV = 0; pF = 0;
265 }
266}
267
268void HepPolyhedron::CreatePrism()
269/***********************************************************************
270 * *
271 * Name: HepPolyhedron::CreatePrism Date: 15.07.96 *
272 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
273 * *
274 * Function: Set facets for a prism *
275 * *
276 ***********************************************************************/
277{
278 enum {DUMMY, BOTTOM, LEFT, BACK, RIGHT, FRONT, TOP};
279
280 pF[1] = G4Facet(1,LEFT, 4,BACK, 3,RIGHT, 2,FRONT);
281 pF[2] = G4Facet(5,TOP, 8,BACK, 4,BOTTOM, 1,FRONT);
282 pF[3] = G4Facet(8,TOP, 7,RIGHT, 3,BOTTOM, 4,LEFT);
283 pF[4] = G4Facet(7,TOP, 6,FRONT, 2,BOTTOM, 3,BACK);
284 pF[5] = G4Facet(6,TOP, 5,LEFT, 1,BOTTOM, 2,RIGHT);
285 pF[6] = G4Facet(5,FRONT, 6,RIGHT, 7,BACK, 8,LEFT);
286}
287
288void HepPolyhedron::RotateEdge(G4int k1, G4int k2, G4double r1, G4double r2,
289 G4int v1, G4int v2, G4int vEdge,
290 G4bool ifWholeCircle, G4int nds, G4int &kface)
291/***********************************************************************
292 * *
293 * Name: HepPolyhedron::RotateEdge Date: 05.12.96 *
294 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
295 * *
296 * Function: Create set of facets by rotation of an edge around Z-axis *
297 * *
298 * Input: k1, k2 - end vertices of the edge *
299 * r1, r2 - radiuses of the end vertices *
300 * v1, v2 - visibility of edges produced by rotation of the end *
301 * vertices *
302 * vEdge - visibility of the edge *
303 * ifWholeCircle - is true in case of whole circle rotation *
304 * nds - number of discrete steps *
305 * r[] - r-coordinates *
306 * kface - current free cell in the pF array *
307 * *
308 ***********************************************************************/
309{
310 if (r1 == 0. && r2 == 0) return;
311
312 G4int i;
313 G4int i1 = k1;
314 G4int i2 = k2;
315 G4int ii1 = ifWholeCircle ? i1 : i1+nds;
316 G4int ii2 = ifWholeCircle ? i2 : i2+nds;
317 G4int vv = ifWholeCircle ? vEdge : 1;
318
319 if (nds == 1) {
320 if (r1 == 0.) {
321 pF[kface++] = G4Facet(i1,0, v2*i2,0, (i2+1),0);
322 }else if (r2 == 0.) {
323 pF[kface++] = G4Facet(i1,0, i2,0, v1*(i1+1),0);
324 }else{
325 pF[kface++] = G4Facet(i1,0, v2*i2,0, (i2+1),0, v1*(i1+1),0);
326 }
327 }else{
328 if (r1 == 0.) {
329 pF[kface++] = G4Facet(vv*i1,0, v2*i2,0, vEdge*(i2+1),0);
330 for (i2++,i=1; i<nds-1; i2++,i++) {
331 pF[kface++] = G4Facet(vEdge*i1,0, v2*i2,0, vEdge*(i2+1),0);
332 }
333 pF[kface++] = G4Facet(vEdge*i1,0, v2*i2,0, vv*ii2,0);
334 }else if (r2 == 0.) {
335 pF[kface++] = G4Facet(vv*i1,0, vEdge*i2,0, v1*(i1+1),0);
336 for (i1++,i=1; i<nds-1; i1++,i++) {
337 pF[kface++] = G4Facet(vEdge*i1,0, vEdge*i2,0, v1*(i1+1),0);
338 }
339 pF[kface++] = G4Facet(vEdge*i1,0, vv*i2,0, v1*ii1,0);
340 }else{
341 pF[kface++] = G4Facet(vv*i1,0, v2*i2,0, vEdge*(i2+1),0,v1*(i1+1),0);
342 for (i1++,i2++,i=1; i<nds-1; i1++,i2++,i++) {
343 pF[kface++] = G4Facet(vEdge*i1,0, v2*i2,0, vEdge*(i2+1),0,v1*(i1+1),0);
344 }
345 pF[kface++] = G4Facet(vEdge*i1,0, v2*i2,0, vv*ii2,0, v1*ii1,0);
346 }
347 }
348}
349
350void HepPolyhedron::SetSideFacets(G4int ii[4], G4int vv[4],
351 G4int *kk, G4double *r,
352 G4double dphi, G4int nds, G4int &kface)
353/***********************************************************************
354 * *
355 * Name: HepPolyhedron::SetSideFacets Date: 20.05.97 *
356 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
357 * *
358 * Function: Set side facets for the case of incomplete rotation *
359 * *
360 * Input: ii[4] - indeces of original verteces *
361 * vv[4] - visibility of edges *
362 * kk[] - indeces of nodes *
363 * r[] - radiuses *
364 * dphi - delta phi *
365 * nds - number of discrete steps *
366 * kface - current free cell in the pF array *
367 * *
368 ***********************************************************************/
369{
370 G4int k1, k2, k3, k4;
371
372 if (std::abs((G4double)(dphi-pi)) < perMillion) { // half a circle
373 for (G4int i=0; i<4; i++) {
374 k1 = ii[i];
375 k2 = (i == 3) ? ii[0] : ii[i+1];
376 if (r[k1] == 0. && r[k2] == 0.) vv[i] = -1;
377 }
378 }
379
380 if (ii[1] == ii[2]) {
381 k1 = kk[ii[0]];
382 k2 = kk[ii[2]];
383 k3 = kk[ii[3]];
384 pF[kface++] = G4Facet(vv[0]*k1,0, vv[2]*k2,0, vv[3]*k3,0);
385 if (r[ii[0]] != 0.) k1 += nds;
386 if (r[ii[2]] != 0.) k2 += nds;
387 if (r[ii[3]] != 0.) k3 += nds;
388 pF[kface++] = G4Facet(vv[2]*k3,0, vv[0]*k2,0, vv[3]*k1,0);
389 }else if (kk[ii[0]] == kk[ii[1]]) {
390 k1 = kk[ii[0]];
391 k2 = kk[ii[2]];
392 k3 = kk[ii[3]];
393 pF[kface++] = G4Facet(vv[1]*k1,0, vv[2]*k2,0, vv[3]*k3,0);
394 if (r[ii[0]] != 0.) k1 += nds;
395 if (r[ii[2]] != 0.) k2 += nds;
396 if (r[ii[3]] != 0.) k3 += nds;
397 pF[kface++] = G4Facet(vv[2]*k3,0, vv[1]*k2,0, vv[3]*k1,0);
398 }else if (kk[ii[2]] == kk[ii[3]]) {
399 k1 = kk[ii[0]];
400 k2 = kk[ii[1]];
401 k3 = kk[ii[2]];
402 pF[kface++] = G4Facet(vv[0]*k1,0, vv[1]*k2,0, vv[3]*k3,0);
403 if (r[ii[0]] != 0.) k1 += nds;
404 if (r[ii[1]] != 0.) k2 += nds;
405 if (r[ii[2]] != 0.) k3 += nds;
406 pF[kface++] = G4Facet(vv[1]*k3,0, vv[0]*k2,0, vv[3]*k1,0);
407 }else{
408 k1 = kk[ii[0]];
409 k2 = kk[ii[1]];
410 k3 = kk[ii[2]];
411 k4 = kk[ii[3]];
412 pF[kface++] = G4Facet(vv[0]*k1,0, vv[1]*k2,0, vv[2]*k3,0, vv[3]*k4,0);
413 if (r[ii[0]] != 0.) k1 += nds;
414 if (r[ii[1]] != 0.) k2 += nds;
415 if (r[ii[2]] != 0.) k3 += nds;
416 if (r[ii[3]] != 0.) k4 += nds;
417 pF[kface++] = G4Facet(vv[2]*k4,0, vv[1]*k3,0, vv[0]*k2,0, vv[3]*k1,0);
418 }
419}
420
421void HepPolyhedron::RotateAroundZ(G4int nstep, G4double phi, G4double dphi,
422 G4int np1, G4int np2,
423 const G4double *z, G4double *r,
424 G4int nodeVis, G4int edgeVis)
425/***********************************************************************
426 * *
427 * Name: HepPolyhedron::RotateAroundZ Date: 27.11.96 *
428 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
429 * *
430 * Function: Create HepPolyhedron for a solid produced by rotation of *
431 * two polylines around Z-axis *
432 * *
433 * Input: nstep - number of discrete steps, if 0 then default *
434 * phi - starting phi angle *
435 * dphi - delta phi *
436 * np1 - number of points in external polyline *
437 * (must be negative in case of closed polyline) *
438 * np2 - number of points in internal polyline (may be 1) *
439 * z[] - z-coordinates (+z >>> -z for both polylines) *
440 * r[] - r-coordinates *
441 * nodeVis - how to Draw edges joing consecutive positions of *
442 * node during rotation *
443 * edgeVis - how to Draw edges *
444 * *
445 ***********************************************************************/
446{
447 static G4double wholeCircle = twopi;
448
449 // S E T R O T A T I O N P A R A M E T E R S
450
451 G4bool ifWholeCircle = (std::abs(dphi-wholeCircle) < perMillion) ? true : false;
452 G4double delPhi = ifWholeCircle ? wholeCircle : dphi;
453 G4int nSphi = (nstep > 0) ?
454 nstep : G4int(delPhi*GetNumberOfRotationSteps()/wholeCircle+.5);
455 if (nSphi == 0) nSphi = 1;
456 G4int nVphi = ifWholeCircle ? nSphi : nSphi+1;
457 G4bool ifClosed = np1 > 0 ? false : true;
458
459 // C O U N T V E R T E C E S
460
461 G4int absNp1 = std::abs(np1);
462 G4int absNp2 = std::abs(np2);
463 G4int i1beg = 0;
464 G4int i1end = absNp1-1;
465 G4int i2beg = absNp1;
466 G4int i2end = absNp1+absNp2-1;
467 G4int i, j, k;
468
469 for(i=i1beg; i<=i2end; i++) {
470 if (std::abs(r[i]) < perMillion) r[i] = 0.;
471 }
472
473 j = 0; // external nodes
474 for (i=i1beg; i<=i1end; i++) {
475 j += (r[i] == 0.) ? 1 : nVphi;
476 }
477
478 G4bool ifSide1 = false; // internal nodes
479 G4bool ifSide2 = false;
480
481 if (r[i2beg] != r[i1beg] || z[i2beg] != z[i1beg]) {
482 j += (r[i2beg] == 0.) ? 1 : nVphi;
483 ifSide1 = true;
484 }
485
486 for(i=i2beg+1; i<i2end; i++) {
487 j += (r[i] == 0.) ? 1 : nVphi;
488 }
489
490 if (r[i2end] != r[i1end] || z[i2end] != z[i1end]) {
491 if (absNp2 > 1) j += (r[i2end] == 0.) ? 1 : nVphi;
492 ifSide2 = true;
493 }
494
495 // C O U N T F A C E S
496
497 k = ifClosed ? absNp1*nSphi : (absNp1-1)*nSphi; // external faces
498
499 if (absNp2 > 1) { // internal faces
500 for(i=i2beg; i<i2end; i++) {
501 if (r[i] > 0. || r[i+1] > 0.) k += nSphi;
502 }
503
504 if (ifClosed) {
505 if (r[i2end] > 0. || r[i2beg] > 0.) k += nSphi;
506 }
507 }
508
509 if (!ifClosed) { // side faces
510 if (ifSide1 && (r[i1beg] > 0. || r[i2beg] > 0.)) k += nSphi;
511 if (ifSide2 && (r[i1end] > 0. || r[i2end] > 0.)) k += nSphi;
512 }
513
514 if (!ifWholeCircle) { // phi_side faces
515 k += ifClosed ? 2*absNp1 : 2*(absNp1-1);
516 }
517
518 // A L L O C A T E M E M O R Y
519
520 AllocateMemory(j, k);
521
522 // G E N E R A T E V E R T E C E S
523
524 G4int *kk;
525 kk = new G4int[absNp1+absNp2];
526
527 k = 1;
528 for(i=i1beg; i<=i1end; i++) {
529 kk[i] = k;
530 if (r[i] == 0.)
531 { pV[k++] = G4Point3D(0, 0, z[i]); } else { k += nVphi; }
532 }
533
534 i = i2beg;
535 if (ifSide1) {
536 kk[i] = k;
537 if (r[i] == 0.)
538 { pV[k++] = G4Point3D(0, 0, z[i]); } else { k += nVphi; }
539 }else{
540 kk[i] = kk[i1beg];
541 }
542
543 for(i=i2beg+1; i<i2end; i++) {
544 kk[i] = k;
545 if (r[i] == 0.)
546 { pV[k++] = G4Point3D(0, 0, z[i]); } else { k += nVphi; }
547 }
548
549 if (absNp2 > 1) {
550 i = i2end;
551 if (ifSide2) {
552 kk[i] = k;
553 if (r[i] == 0.) pV[k] = G4Point3D(0, 0, z[i]);
554 }else{
555 kk[i] = kk[i1end];
556 }
557 }
558
559 G4double cosPhi, sinPhi;
560
561 for(j=0; j<nVphi; j++) {
562 cosPhi = std::cos(phi+j*delPhi/nSphi);
563 sinPhi = std::sin(phi+j*delPhi/nSphi);
564 for(i=i1beg; i<=i2end; i++) {
565 if (r[i] != 0.)
566 pV[kk[i]+j] = G4Point3D(r[i]*cosPhi,r[i]*sinPhi,z[i]);
567 }
568 }
569
570 // G E N E R A T E E X T E R N A L F A C E S
571
572 G4int v1,v2;
573
574 k = 1;
575 v2 = ifClosed ? nodeVis : 1;
576 for(i=i1beg; i<i1end; i++) {
577 v1 = v2;
578 if (!ifClosed && i == i1end-1) {
579 v2 = 1;
580 }else{
581 v2 = (r[i] == r[i+1] && r[i+1] == r[i+2]) ? -1 : nodeVis;
582 }
583 RotateEdge(kk[i], kk[i+1], r[i], r[i+1], v1, v2,
584 edgeVis, ifWholeCircle, nSphi, k);
585 }
586 if (ifClosed) {
587 RotateEdge(kk[i1end], kk[i1beg], r[i1end],r[i1beg], nodeVis, nodeVis,
588 edgeVis, ifWholeCircle, nSphi, k);
589 }
590
591 // G E N E R A T E I N T E R N A L F A C E S
592
593 if (absNp2 > 1) {
594 v2 = ifClosed ? nodeVis : 1;
595 for(i=i2beg; i<i2end; i++) {
596 v1 = v2;
597 if (!ifClosed && i==i2end-1) {
598 v2 = 1;
599 }else{
600 v2 = (r[i] == r[i+1] && r[i+1] == r[i+2]) ? -1 : nodeVis;
601 }
602 RotateEdge(kk[i+1], kk[i], r[i+1], r[i], v2, v1,
603 edgeVis, ifWholeCircle, nSphi, k);
604 }
605 if (ifClosed) {
606 RotateEdge(kk[i2beg], kk[i2end], r[i2beg], r[i2end], nodeVis, nodeVis,
607 edgeVis, ifWholeCircle, nSphi, k);
608 }
609 }
610
611 // G E N E R A T E S I D E F A C E S
612
613 if (!ifClosed) {
614 if (ifSide1) {
615 RotateEdge(kk[i2beg], kk[i1beg], r[i2beg], r[i1beg], 1, 1,
616 -1, ifWholeCircle, nSphi, k);
617 }
618 if (ifSide2) {
619 RotateEdge(kk[i1end], kk[i2end], r[i1end], r[i2end], 1, 1,
620 -1, ifWholeCircle, nSphi, k);
621 }
622 }
623
624 // G E N E R A T E S I D E F A C E S for the case of incomplete circle
625
626 if (!ifWholeCircle) {
627
628 G4int ii[4], vv[4];
629
630 if (ifClosed) {
631 for (i=i1beg; i<=i1end; i++) {
632 ii[0] = i;
633 ii[3] = (i == i1end) ? i1beg : i+1;
634 ii[1] = (absNp2 == 1) ? i2beg : ii[0]+absNp1;
635 ii[2] = (absNp2 == 1) ? i2beg : ii[3]+absNp1;
636 vv[0] = -1;
637 vv[1] = 1;
638 vv[2] = -1;
639 vv[3] = 1;
640 SetSideFacets(ii, vv, kk, r, dphi, nSphi, k);
641 }
642 }else{
643 for (i=i1beg; i<i1end; i++) {
644 ii[0] = i;
645 ii[3] = i+1;
646 ii[1] = (absNp2 == 1) ? i2beg : ii[0]+absNp1;
647 ii[2] = (absNp2 == 1) ? i2beg : ii[3]+absNp1;
648 vv[0] = (i == i1beg) ? 1 : -1;
649 vv[1] = 1;
650 vv[2] = (i == i1end-1) ? 1 : -1;
651 vv[3] = 1;
652 SetSideFacets(ii, vv, kk, r, dphi, nSphi, k);
653 }
654 }
655 }
656
657 delete [] kk;
658
659 if (k-1 != nface) {
660 std::cerr
661 << "Polyhedron::RotateAroundZ: number of generated faces ("
662 << k-1 << ") is not equal to the number of allocated faces ("
663 << nface << ")"
664 << std::endl;
665 }
666}
667
668void HepPolyhedron::SetReferences()
669/***********************************************************************
670 * *
671 * Name: HepPolyhedron::SetReferences Date: 04.12.96 *
672 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
673 * *
674 * Function: For each edge set reference to neighbouring facet *
675 * *
676 ***********************************************************************/
677{
678 if (nface <= 0) return;
679
680 struct edgeListMember {
681 edgeListMember *next;
682 G4int v2;
683 G4int iface;
684 G4int iedge;
685 } *edgeList, *freeList, **headList;
686
687
688 // A L L O C A T E A N D I N I T I A T E L I S T S
689
690 edgeList = new edgeListMember[2*nface];
691 headList = new edgeListMember*[nvert];
692
693 G4int i;
694 for (i=0; i<nvert; i++) {
695 headList[i] = 0;
696 }
697 freeList = edgeList;
698 for (i=0; i<2*nface-1; i++) {
699 edgeList[i].next = &edgeList[i+1];
700 }
701 edgeList[2*nface-1].next = 0;
702
703 // L O O P A L O N G E D G E S
704
705 G4int iface, iedge, nedge, i1, i2, k1, k2;
706 edgeListMember *prev, *cur;
707
708 for(iface=1; iface<=nface; iface++) {
709 nedge = (pF[iface].edge[3].v == 0) ? 3 : 4;
710 for (iedge=0; iedge<nedge; iedge++) {
711 i1 = iedge;
712 i2 = (iedge < nedge-1) ? iedge+1 : 0;
713 i1 = std::abs(pF[iface].edge[i1].v);
714 i2 = std::abs(pF[iface].edge[i2].v);
715 k1 = (i1 < i2) ? i1 : i2; // k1 = ::min(i1,i2);
716 k2 = (i1 > i2) ? i1 : i2; // k2 = ::max(i1,i2);
717
718 // check head of the List corresponding to k1
719 cur = headList[k1];
720 if (cur == 0) {
721 headList[k1] = freeList;
722 freeList = freeList->next;
723 cur = headList[k1];
724 cur->next = 0;
725 cur->v2 = k2;
726 cur->iface = iface;
727 cur->iedge = iedge;
728 continue;
729 }
730
731 if (cur->v2 == k2) {
732 headList[k1] = cur->next;
733 cur->next = freeList;
734 freeList = cur;
735 pF[iface].edge[iedge].f = cur->iface;
736 pF[cur->iface].edge[cur->iedge].f = iface;
737 i1 = (pF[iface].edge[iedge].v < 0) ? -1 : 1;
738 i2 = (pF[cur->iface].edge[cur->iedge].v < 0) ? -1 : 1;
739 if (i1 != i2) {
740 std::cerr
741 << "Polyhedron::SetReferences: different edge visibility "
742 << iface << "/" << iedge << "/"
743 << pF[iface].edge[iedge].v << " and "
744 << cur->iface << "/" << cur->iedge << "/"
745 << pF[cur->iface].edge[cur->iedge].v
746 << std::endl;
747 }
748 continue;
749 }
750
751 // check List itself
752 for (;;) {
753 prev = cur;
754 cur = prev->next;
755 if (cur == 0) {
756 prev->next = freeList;
757 freeList = freeList->next;
758 cur = prev->next;
759 cur->next = 0;
760 cur->v2 = k2;
761 cur->iface = iface;
762 cur->iedge = iedge;
763 break;
764 }
765
766 if (cur->v2 == k2) {
767 prev->next = cur->next;
768 cur->next = freeList;
769 freeList = cur;
770 pF[iface].edge[iedge].f = cur->iface;
771 pF[cur->iface].edge[cur->iedge].f = iface;
772 i1 = (pF[iface].edge[iedge].v < 0) ? -1 : 1;
773 i2 = (pF[cur->iface].edge[cur->iedge].v < 0) ? -1 : 1;
774 if (i1 != i2) {
775 std::cerr
776 << "Polyhedron::SetReferences: different edge visibility "
777 << iface << "/" << iedge << "/"
778 << pF[iface].edge[iedge].v << " and "
779 << cur->iface << "/" << cur->iedge << "/"
780 << pF[cur->iface].edge[cur->iedge].v
781 << std::endl;
782 }
783 break;
784 }
785 }
786 }
787 }
788
789 // C H E C K T H A T A L L L I S T S A R E E M P T Y
790
791 for (i=0; i<nvert; i++) {
792 if (headList[i] != 0) {
793 std::cerr
794 << "Polyhedron::SetReferences: List " << i << " is not empty"
795 << std::endl;
796 }
797 }
798
799 // F R E E M E M O R Y
800
801 delete [] edgeList;
802 delete [] headList;
803}
804
805void HepPolyhedron::InvertFacets()
806/***********************************************************************
807 * *
808 * Name: HepPolyhedron::InvertFacets Date: 01.12.99 *
809 * Author: E.Chernyaev Revised: *
810 * *
811 * Function: Invert the order of the nodes in the facets *
812 * *
813 ***********************************************************************/
814{
815 if (nface <= 0) return;
816 G4int i, k, nnode, v[4],f[4];
817 for (i=1; i<=nface; i++) {
818 nnode = (pF[i].edge[3].v == 0) ? 3 : 4;
819 for (k=0; k<nnode; k++) {
820 v[k] = (k+1 == nnode) ? pF[i].edge[0].v : pF[i].edge[k+1].v;
821 if (v[k] * pF[i].edge[k].v < 0) v[k] = -v[k];
822 f[k] = pF[i].edge[k].f;
823 }
824 for (k=0; k<nnode; k++) {
825 pF[i].edge[nnode-1-k].v = v[k];
826 pF[i].edge[nnode-1-k].f = f[k];
827 }
828 }
829}
830
831HepPolyhedron & HepPolyhedron::Transform(const G4Transform3D &t)
832/***********************************************************************
833 * *
834 * Name: HepPolyhedron::Transform Date: 01.12.99 *
835 * Author: E.Chernyaev Revised: *
836 * *
837 * Function: Make transformation of the polyhedron *
838 * *
839 ***********************************************************************/
840{
841 if (nvert > 0) {
842 for (G4int i=1; i<=nvert; i++) { pV[i] = t * pV[i]; }
843
844 // C H E C K D E T E R M I N A N T A N D
845 // I N V E R T F A C E T S I F I T I S N E G A T I V E
846
847 G4Vector3D d = t * G4Vector3D(0,0,0);
848 G4Vector3D x = t * G4Vector3D(1,0,0) - d;
849 G4Vector3D y = t * G4Vector3D(0,1,0) - d;
850 G4Vector3D z = t * G4Vector3D(0,0,1) - d;
851 if ((x.cross(y))*z < 0) InvertFacets();
852 }
853 return *this;
854}
855
856G4bool HepPolyhedron::GetNextVertexIndex(G4int &index, G4int &edgeFlag) const
857/***********************************************************************
858 * *
859 * Name: HepPolyhedron::GetNextVertexIndex Date: 03.09.96 *
860 * Author: Yasuhide Sawada Revised: *
861 * *
862 * Function: *
863 * *
864 ***********************************************************************/
865{
866 static G4int iFace = 1;
867 static G4int iQVertex = 0;
868 G4int vIndex = pF[iFace].edge[iQVertex].v;
869
870 edgeFlag = (vIndex > 0) ? 1 : 0;
871 index = std::abs(vIndex);
872
873 if (iQVertex >= 3 || pF[iFace].edge[iQVertex+1].v == 0) {
874 iQVertex = 0;
875 if (++iFace > nface) iFace = 1;
876 return false; // Last Edge
877 }else{
878 ++iQVertex;
879 return true; // not Last Edge
880 }
881}
882
883G4Point3D HepPolyhedron::GetVertex(G4int index) const
884/***********************************************************************
885 * *
886 * Name: HepPolyhedron::GetVertex Date: 03.09.96 *
887 * Author: Yasuhide Sawada Revised: 17.11.99 *
888 * *
889 * Function: Get vertex of the index. *
890 * *
891 ***********************************************************************/
892{
893 if (index <= 0 || index > nvert) {
894 std::cerr
895 << "HepPolyhedron::GetVertex: irrelevant index " << index
896 << std::endl;
897 return G4Point3D();
898 }
899 return pV[index];
900}
901
902G4bool
903HepPolyhedron::GetNextVertex(G4Point3D &vertex, G4int &edgeFlag) const
904/***********************************************************************
905 * *
906 * Name: HepPolyhedron::GetNextVertex Date: 22.07.96 *
907 * Author: John Allison Revised: *
908 * *
909 * Function: Get vertices of the quadrilaterals in order for each *
910 * face in face order. Returns false when finished each *
911 * face. *
912 * *
913 ***********************************************************************/
914{
915 G4int index;
916 G4bool rep = GetNextVertexIndex(index, edgeFlag);
917 vertex = pV[index];
918 return rep;
919}
920
921G4bool HepPolyhedron::GetNextVertex(G4Point3D &vertex, G4int &edgeFlag,
922 G4Normal3D &normal) const
923/***********************************************************************
924 * *
925 * Name: HepPolyhedron::GetNextVertex Date: 26.11.99 *
926 * Author: E.Chernyaev Revised: *
927 * *
928 * Function: Get vertices with normals of the quadrilaterals in order *
929 * for each face in face order. *
930 * Returns false when finished each face. *
931 * *
932 ***********************************************************************/
933{
934 static G4int iFace = 1;
935 static G4int iNode = 0;
936
937 if (nface == 0) return false; // empty polyhedron
938
939 G4int k = pF[iFace].edge[iNode].v;
940 if (k > 0) { edgeFlag = 1; } else { edgeFlag = -1; k = -k; }
941 vertex = pV[k];
942 normal = FindNodeNormal(iFace,k);
943 if (iNode >= 3 || pF[iFace].edge[iNode+1].v == 0) {
944 iNode = 0;
945 if (++iFace > nface) iFace = 1;
946 return false; // last node
947 }else{
948 ++iNode;
949 return true; // not last node
950 }
951}
952
953G4bool HepPolyhedron::GetNextEdgeIndeces(G4int &i1, G4int &i2, G4int &edgeFlag,
954 G4int &iface1, G4int &iface2) const
955/***********************************************************************
956 * *
957 * Name: HepPolyhedron::GetNextEdgeIndeces Date: 30.09.96 *
958 * Author: E.Chernyaev Revised: 17.11.99 *
959 * *
960 * Function: Get indeces of the next edge together with indeces of *
961 * of the faces which share the edge. *
962 * Returns false when the last edge. *
963 * *
964 ***********************************************************************/
965{
966 static G4int iFace = 1;
967 static G4int iQVertex = 0;
968 static G4int iOrder = 1;
969 G4int k1, k2, kflag, kface1, kface2;
970
971 if (iFace == 1 && iQVertex == 0) {
972 k2 = pF[nface].edge[0].v;
973 k1 = pF[nface].edge[3].v;
974 if (k1 == 0) k1 = pF[nface].edge[2].v;
975 if (std::abs(k1) > std::abs(k2)) iOrder = -1;
976 }
977
978 do {
979 k1 = pF[iFace].edge[iQVertex].v;
980 kflag = k1;
981 k1 = std::abs(k1);
982 kface1 = iFace;
983 kface2 = pF[iFace].edge[iQVertex].f;
984 if (iQVertex >= 3 || pF[iFace].edge[iQVertex+1].v == 0) {
985 iQVertex = 0;
986 k2 = std::abs(pF[iFace].edge[iQVertex].v);
987 iFace++;
988 }else{
989 iQVertex++;
990 k2 = std::abs(pF[iFace].edge[iQVertex].v);
991 }
992 } while (iOrder*k1 > iOrder*k2);
993
994 i1 = k1; i2 = k2; edgeFlag = (kflag > 0) ? 1 : 0;
995 iface1 = kface1; iface2 = kface2;
996
997 if (iFace > nface) {
998 iFace = 1; iOrder = 1;
999 return false;
1000 }else{
1001 return true;
1002 }
1003}
1004
1005G4bool
1006HepPolyhedron::GetNextEdgeIndeces(G4int &i1, G4int &i2, G4int &edgeFlag) const
1007/***********************************************************************
1008 * *
1009 * Name: HepPolyhedron::GetNextEdgeIndeces Date: 17.11.99 *
1010 * Author: E.Chernyaev Revised: *
1011 * *
1012 * Function: Get indeces of the next edge. *
1013 * Returns false when the last edge. *
1014 * *
1015 ***********************************************************************/
1016{
1017 G4int kface1, kface2;
1018 return GetNextEdgeIndeces(i1, i2, edgeFlag, kface1, kface2);
1019}
1020
1021G4bool
1022HepPolyhedron::GetNextEdge(G4Point3D &p1,
1023 G4Point3D &p2,
1024 G4int &edgeFlag) const
1025/***********************************************************************
1026 * *
1027 * Name: HepPolyhedron::GetNextEdge Date: 30.09.96 *
1028 * Author: E.Chernyaev Revised: *
1029 * *
1030 * Function: Get next edge. *
1031 * Returns false when the last edge. *
1032 * *
1033 ***********************************************************************/
1034{
1035 G4int i1,i2;
1036 G4bool rep = GetNextEdgeIndeces(i1,i2,edgeFlag);
1037 p1 = pV[i1];
1038 p2 = pV[i2];
1039 return rep;
1040}
1041
1042G4bool
1043HepPolyhedron::GetNextEdge(G4Point3D &p1, G4Point3D &p2,
1044 G4int &edgeFlag, G4int &iface1, G4int &iface2) const
1045/***********************************************************************
1046 * *
1047 * Name: HepPolyhedron::GetNextEdge Date: 17.11.99 *
1048 * Author: E.Chernyaev Revised: *
1049 * *
1050 * Function: Get next edge with indeces of the faces which share *
1051 * the edge. *
1052 * Returns false when the last edge. *
1053 * *
1054 ***********************************************************************/
1055{
1056 G4int i1,i2;
1057 G4bool rep = GetNextEdgeIndeces(i1,i2,edgeFlag,iface1,iface2);
1058 p1 = pV[i1];
1059 p2 = pV[i2];
1060 return rep;
1061}
1062
1063void HepPolyhedron::GetFacet(G4int iFace, G4int &n, G4int *iNodes,
1064 G4int *edgeFlags, G4int *iFaces) const
1065/***********************************************************************
1066 * *
1067 * Name: HepPolyhedron::GetFacet Date: 15.12.99 *
1068 * Author: E.Chernyaev Revised: *
1069 * *
1070 * Function: Get face by index *
1071 * *
1072 ***********************************************************************/
1073{
1074 if (iFace < 1 || iFace > nface) {
1075 std::cerr
1076 << "HepPolyhedron::GetFacet: irrelevant index " << iFace
1077 << std::endl;
1078 n = 0;
1079 }else{
1080 G4int i, k;
1081 for (i=0; i<4; i++) {
1082 k = pF[iFace].edge[i].v;
1083 if (k == 0) break;
1084 if (iFaces != 0) iFaces[i] = pF[iFace].edge[i].f;
1085 if (k > 0) {
1086 iNodes[i] = k;
1087 if (edgeFlags != 0) edgeFlags[i] = 1;
1088 }else{
1089 iNodes[i] = -k;
1090 if (edgeFlags != 0) edgeFlags[i] = -1;
1091 }
1092 }
1093 n = i;
1094 }
1095}
1096
1097void HepPolyhedron::GetFacet(G4int index, G4int &n, G4Point3D *nodes,
1098 G4int *edgeFlags, G4Normal3D *normals) const
1099/***********************************************************************
1100 * *
1101 * Name: HepPolyhedron::GetFacet Date: 17.11.99 *
1102 * Author: E.Chernyaev Revised: *
1103 * *
1104 * Function: Get face by index *
1105 * *
1106 ***********************************************************************/
1107{
1108 G4int iNodes[4];
1109 GetFacet(index, n, iNodes, edgeFlags);
1110 if (n != 0) {
1111 for (G4int i=0; i<n; i++) {
1112 nodes[i] = pV[iNodes[i]];
1113 if (normals != 0) normals[i] = FindNodeNormal(index,iNodes[i]);
1114 }
1115 }
1116}
1117
1118G4bool
1119HepPolyhedron::GetNextFacet(G4int &n, G4Point3D *nodes,
1120 G4int *edgeFlags, G4Normal3D *normals) const
1121/***********************************************************************
1122 * *
1123 * Name: HepPolyhedron::GetNextFacet Date: 19.11.99 *
1124 * Author: E.Chernyaev Revised: *
1125 * *
1126 * Function: Get next face with normals of unit length at the nodes. *
1127 * Returns false when finished all faces. *
1128 * *
1129 ***********************************************************************/
1130{
1131 static G4int iFace = 1;
1132
1133 if (edgeFlags == 0) {
1134 GetFacet(iFace, n, nodes);
1135 }else if (normals == 0) {
1136 GetFacet(iFace, n, nodes, edgeFlags);
1137 }else{
1138 GetFacet(iFace, n, nodes, edgeFlags, normals);
1139 }
1140
1141 if (++iFace > nface) {
1142 iFace = 1;
1143 return false;
1144 }else{
1145 return true;
1146 }
1147}
1148
1149G4Normal3D HepPolyhedron::GetNormal(G4int iFace) const
1150/***********************************************************************
1151 * *
1152 * Name: HepPolyhedron::GetNormal Date: 19.11.99 *
1153 * Author: E.Chernyaev Revised: *
1154 * *
1155 * Function: Get normal of the face given by index *
1156 * *
1157 ***********************************************************************/
1158{
1159 if (iFace < 1 || iFace > nface) {
1160 std::cerr
1161 << "HepPolyhedron::GetNormal: irrelevant index " << iFace
1162 << std::endl;
1163 return G4Normal3D();
1164 }
1165
1166 G4int i0 = std::abs(pF[iFace].edge[0].v);
1167 G4int i1 = std::abs(pF[iFace].edge[1].v);
1168 G4int i2 = std::abs(pF[iFace].edge[2].v);
1169 G4int i3 = std::abs(pF[iFace].edge[3].v);
1170 if (i3 == 0) i3 = i0;
1171 return (pV[i2] - pV[i0]).cross(pV[i3] - pV[i1]);
1172}
1173
1174G4Normal3D HepPolyhedron::GetUnitNormal(G4int iFace) const
1175/***********************************************************************
1176 * *
1177 * Name: HepPolyhedron::GetNormal Date: 19.11.99 *
1178 * Author: E.Chernyaev Revised: *
1179 * *
1180 * Function: Get unit normal of the face given by index *
1181 * *
1182 ***********************************************************************/
1183{
1184 if (iFace < 1 || iFace > nface) {
1185 std::cerr
1186 << "HepPolyhedron::GetUnitNormal: irrelevant index " << iFace
1187 << std::endl;
1188 return G4Normal3D();
1189 }
1190
1191 G4int i0 = std::abs(pF[iFace].edge[0].v);
1192 G4int i1 = std::abs(pF[iFace].edge[1].v);
1193 G4int i2 = std::abs(pF[iFace].edge[2].v);
1194 G4int i3 = std::abs(pF[iFace].edge[3].v);
1195 if (i3 == 0) i3 = i0;
1196 return ((pV[i2] - pV[i0]).cross(pV[i3] - pV[i1])).unit();
1197}
1198
1199G4bool HepPolyhedron::GetNextNormal(G4Normal3D &normal) const
1200/***********************************************************************
1201 * *
1202 * Name: HepPolyhedron::GetNextNormal Date: 22.07.96 *
1203 * Author: John Allison Revised: 19.11.99 *
1204 * *
1205 * Function: Get normals of each face in face order. Returns false *
1206 * when finished all faces. *
1207 * *
1208 ***********************************************************************/
1209{
1210 static G4int iFace = 1;
1211 normal = GetNormal(iFace);
1212 if (++iFace > nface) {
1213 iFace = 1;
1214 return false;
1215 }else{
1216 return true;
1217 }
1218}
1219
1220G4bool HepPolyhedron::GetNextUnitNormal(G4Normal3D &normal) const
1221/***********************************************************************
1222 * *
1223 * Name: HepPolyhedron::GetNextUnitNormal Date: 16.09.96 *
1224 * Author: E.Chernyaev Revised: *
1225 * *
1226 * Function: Get normals of unit length of each face in face order. *
1227 * Returns false when finished all faces. *
1228 * *
1229 ***********************************************************************/
1230{
1231 G4bool rep = GetNextNormal(normal);
1232 normal = normal.unit();
1233 return rep;
1234}
1235
1236G4double HepPolyhedron::GetSurfaceArea() const
1237/***********************************************************************
1238 * *
1239 * Name: HepPolyhedron::GetSurfaceArea Date: 25.05.01 *
1240 * Author: E.Chernyaev Revised: *
1241 * *
1242 * Function: Returns area of the surface of the polyhedron. *
1243 * *
1244 ***********************************************************************/
1245{
1246 G4double srf = 0.;
1247 for (G4int iFace=1; iFace<=nface; iFace++) {
1248 G4int i0 = std::abs(pF[iFace].edge[0].v);
1249 G4int i1 = std::abs(pF[iFace].edge[1].v);
1250 G4int i2 = std::abs(pF[iFace].edge[2].v);
1251 G4int i3 = std::abs(pF[iFace].edge[3].v);
1252 if (i3 == 0) i3 = i0;
1253 srf += ((pV[i2] - pV[i0]).cross(pV[i3] - pV[i1])).mag();
1254 }
1255 return srf/2.;
1256}
1257
1258G4double HepPolyhedron::GetVolume() const
1259/***********************************************************************
1260 * *
1261 * Name: HepPolyhedron::GetVolume Date: 25.05.01 *
1262 * Author: E.Chernyaev Revised: *
1263 * *
1264 * Function: Returns volume of the polyhedron. *
1265 * *
1266 ***********************************************************************/
1267{
1268 G4double v = 0.;
1269 for (G4int iFace=1; iFace<=nface; iFace++) {
1270 G4int i0 = std::abs(pF[iFace].edge[0].v);
1271 G4int i1 = std::abs(pF[iFace].edge[1].v);
1272 G4int i2 = std::abs(pF[iFace].edge[2].v);
1273 G4int i3 = std::abs(pF[iFace].edge[3].v);
1274 G4Point3D pt;
1275 if (i3 == 0) {
1276 i3 = i0;
1277 pt = (pV[i0]+pV[i1]+pV[i2]) * (1./3.);
1278 }else{
1279 pt = (pV[i0]+pV[i1]+pV[i2]+pV[i3]) * 0.25;
1280 }
1281 v += ((pV[i2] - pV[i0]).cross(pV[i3] - pV[i1])).dot(pt);
1282 }
1283 return v/6.;
1284}
1285
1286G4int
1287HepPolyhedron::createTwistedTrap(G4double Dz,
1288 const G4double xy1[][2],
1289 const G4double xy2[][2])
1290/***********************************************************************
1291 * *
1292 * Name: createTwistedTrap Date: 05.11.02 *
1293 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
1294 * *
1295 * Function: Creates polyhedron for twisted trapezoid *
1296 * *
1297 * Input: Dz - half-length along Z 8----7 *
1298 * xy1[2,4] - quadrilateral at Z=-Dz 5----6 ! *
1299 * xy2[2,4] - quadrilateral at Z=+Dz ! 4-!--3 *
1300 * 1----2 *
1301 * *
1302 ***********************************************************************/
1303{
1304 AllocateMemory(12,18);
1305
1306 pV[ 1] = G4Point3D(xy1[0][0],xy1[0][1],-Dz);
1307 pV[ 2] = G4Point3D(xy1[1][0],xy1[1][1],-Dz);
1308 pV[ 3] = G4Point3D(xy1[2][0],xy1[2][1],-Dz);
1309 pV[ 4] = G4Point3D(xy1[3][0],xy1[3][1],-Dz);
1310
1311 pV[ 5] = G4Point3D(xy2[0][0],xy2[0][1], Dz);
1312 pV[ 6] = G4Point3D(xy2[1][0],xy2[1][1], Dz);
1313 pV[ 7] = G4Point3D(xy2[2][0],xy2[2][1], Dz);
1314 pV[ 8] = G4Point3D(xy2[3][0],xy2[3][1], Dz);
1315
1316 pV[ 9] = (pV[1]+pV[2]+pV[5]+pV[6])/4.;
1317 pV[10] = (pV[2]+pV[3]+pV[6]+pV[7])/4.;
1318 pV[11] = (pV[3]+pV[4]+pV[7]+pV[8])/4.;
1319 pV[12] = (pV[4]+pV[1]+pV[8]+pV[5])/4.;
1320
1321 enum {DUMMY, BOTTOM,
1322 LEFT_BOTTOM, LEFT_FRONT, LEFT_TOP, LEFT_BACK,
1323 BACK_BOTTOM, BACK_LEFT, BACK_TOP, BACK_RIGHT,
1324 RIGHT_BOTTOM, RIGHT_BACK, RIGHT_TOP, RIGHT_FRONT,
1325 FRONT_BOTTOM, FRONT_RIGHT, FRONT_TOP, FRONT_LEFT,
1326 TOP};
1327
1328 pF[ 1]=G4Facet(1,LEFT_BOTTOM, 4,BACK_BOTTOM, 3,RIGHT_BOTTOM, 2,FRONT_BOTTOM);
1329
1330 pF[ 2]=G4Facet(4,BOTTOM, -1,LEFT_FRONT, -12,LEFT_BACK, 0,0);
1331 pF[ 3]=G4Facet(1,FRONT_LEFT, -5,LEFT_TOP, -12,LEFT_BOTTOM, 0,0);
1332 pF[ 4]=G4Facet(5,TOP, -8,LEFT_BACK, -12,LEFT_FRONT, 0,0);
1333 pF[ 5]=G4Facet(8,BACK_LEFT, -4,LEFT_BOTTOM, -12,LEFT_TOP, 0,0);
1334
1335 pF[ 6]=G4Facet(3,BOTTOM, -4,BACK_LEFT, -11,BACK_RIGHT, 0,0);
1336 pF[ 7]=G4Facet(4,LEFT_BACK, -8,BACK_TOP, -11,BACK_BOTTOM, 0,0);
1337 pF[ 8]=G4Facet(8,TOP, -7,BACK_RIGHT, -11,BACK_LEFT, 0,0);
1338 pF[ 9]=G4Facet(7,RIGHT_BACK, -3,BACK_BOTTOM, -11,BACK_TOP, 0,0);
1339
1340 pF[10]=G4Facet(2,BOTTOM, -3,RIGHT_BACK, -10,RIGHT_FRONT, 0,0);
1341 pF[11]=G4Facet(3,BACK_RIGHT, -7,RIGHT_TOP, -10,RIGHT_BOTTOM, 0,0);
1342 pF[12]=G4Facet(7,TOP, -6,RIGHT_FRONT, -10,RIGHT_BACK, 0,0);
1343 pF[13]=G4Facet(6,FRONT_RIGHT,-2,RIGHT_BOTTOM,-10,RIGHT_TOP, 0,0);
1344
1345 pF[14]=G4Facet(1,BOTTOM, -2,FRONT_RIGHT, -9,FRONT_LEFT, 0,0);
1346 pF[15]=G4Facet(2,RIGHT_FRONT,-6,FRONT_TOP, -9,FRONT_BOTTOM, 0,0);
1347 pF[16]=G4Facet(6,TOP, -5,FRONT_LEFT, -9,FRONT_RIGHT, 0,0);
1348 pF[17]=G4Facet(5,LEFT_FRONT, -1,FRONT_BOTTOM, -9,FRONT_TOP, 0,0);
1349
1350 pF[18]=G4Facet(5,FRONT_TOP, 6,RIGHT_TOP, 7,BACK_TOP, 8,LEFT_TOP);
1351
1352 return 0;
1353}
1354
1355G4int
1356HepPolyhedron::createPolyhedron(G4int Nnodes, G4int Nfaces,
1357 const G4double xyz[][3],
1358 const G4int faces[][4])
1359/***********************************************************************
1360 * *
1361 * Name: createPolyhedron Date: 05.11.02 *
1362 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
1363 * *
1364 * Function: Creates user defined polyhedron *
1365 * *
1366 * Input: Nnodes - number of nodes *
1367 * Nfaces - number of faces *
1368 * nodes[][3] - node coordinates *
1369 * faces[][4] - faces *
1370 * *
1371 ***********************************************************************/
1372{
1373 AllocateMemory(Nnodes, Nfaces);
1374 if (nvert == 0) return 1;
1375
1376 for (G4int i=0; i<Nnodes; i++) {
1377 pV[i+1] = G4Point3D(xyz[i][0], xyz[i][1], xyz[i][2]);
1378 }
1379 for (G4int k=0; k<Nfaces; k++) {
1380 pF[k+1] = G4Facet(faces[k][0],0,faces[k][1],0,faces[k][2],0,faces[k][3],0);
1381 }
1382 SetReferences();
1383 return 0;
1384}
1385
1386HepPolyhedronTrd2::HepPolyhedronTrd2(G4double Dx1, G4double Dx2,
1387 G4double Dy1, G4double Dy2,
1388 G4double Dz)
1389/***********************************************************************
1390 * *
1391 * Name: HepPolyhedronTrd2 Date: 22.07.96 *
1392 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
1393 * *
1394 * Function: Create GEANT4 TRD2-trapezoid *
1395 * *
1396 * Input: Dx1 - half-length along X at -Dz 8----7 *
1397 * Dx2 - half-length along X ay +Dz 5----6 ! *
1398 * Dy1 - half-length along Y ay -Dz ! 4-!--3 *
1399 * Dy2 - half-length along Y ay +Dz 1----2 *
1400 * Dz - half-length along Z *
1401 * *
1402 ***********************************************************************/
1403{
1404 AllocateMemory(8,6);
1405
1406 pV[1] = G4Point3D(-Dx1,-Dy1,-Dz);
1407 pV[2] = G4Point3D( Dx1,-Dy1,-Dz);
1408 pV[3] = G4Point3D( Dx1, Dy1,-Dz);
1409 pV[4] = G4Point3D(-Dx1, Dy1,-Dz);
1410 pV[5] = G4Point3D(-Dx2,-Dy2, Dz);
1411 pV[6] = G4Point3D( Dx2,-Dy2, Dz);
1412 pV[7] = G4Point3D( Dx2, Dy2, Dz);
1413 pV[8] = G4Point3D(-Dx2, Dy2, Dz);
1414
1415 CreatePrism();
1416}
1417
1418HepPolyhedronTrd2::~HepPolyhedronTrd2() {}
1419
1420HepPolyhedronTrd1::HepPolyhedronTrd1(G4double Dx1, G4double Dx2,
1421 G4double Dy, G4double Dz)
1422 : HepPolyhedronTrd2(Dx1, Dx2, Dy, Dy, Dz) {}
1423
1424HepPolyhedronTrd1::~HepPolyhedronTrd1() {}
1425
1426HepPolyhedronBox::HepPolyhedronBox(G4double Dx, G4double Dy, G4double Dz)
1427 : HepPolyhedronTrd2(Dx, Dx, Dy, Dy, Dz) {}
1428
1429HepPolyhedronBox::~HepPolyhedronBox() {}
1430
1431HepPolyhedronTrap::HepPolyhedronTrap(G4double Dz,
1432 G4double Theta,
1433 G4double Phi,
1434 G4double Dy1,
1435 G4double Dx1,
1436 G4double Dx2,
1437 G4double Alp1,
1438 G4double Dy2,
1439 G4double Dx3,
1440 G4double Dx4,
1441 G4double Alp2)
1442/***********************************************************************
1443 * *
1444 * Name: HepPolyhedronTrap Date: 20.11.96 *
1445 * Author: E.Chernyaev Revised: *
1446 * *
1447 * Function: Create GEANT4 TRAP-trapezoid *
1448 * *
1449 * Input: DZ - half-length in Z *
1450 * Theta,Phi - polar angles of the line joining centres of the *
1451 * faces at Z=-Dz and Z=+Dz *
1452 * Dy1 - half-length in Y of the face at Z=-Dz *
1453 * Dx1 - half-length in X of low edge of the face at Z=-Dz *
1454 * Dx2 - half-length in X of top edge of the face at Z=-Dz *
1455 * Alp1 - angle between Y-axis and the median joining top and *
1456 * low edges of the face at Z=-Dz *
1457 * Dy2 - half-length in Y of the face at Z=+Dz *
1458 * Dx3 - half-length in X of low edge of the face at Z=+Dz *
1459 * Dx4 - half-length in X of top edge of the face at Z=+Dz *
1460 * Alp2 - angle between Y-axis and the median joining top and *
1461 * low edges of the face at Z=+Dz *
1462 * *
1463 ***********************************************************************/
1464{
1465 G4double DzTthetaCphi = Dz*std::tan(Theta)*std::cos(Phi);
1466 G4double DzTthetaSphi = Dz*std::tan(Theta)*std::sin(Phi);
1467 G4double Dy1Talp1 = Dy1*std::tan(Alp1);
1468 G4double Dy2Talp2 = Dy2*std::tan(Alp2);
1469
1470 AllocateMemory(8,6);
1471
1472 pV[1] = G4Point3D(-DzTthetaCphi-Dy1Talp1-Dx1,-DzTthetaSphi-Dy1,-Dz);
1473 pV[2] = G4Point3D(-DzTthetaCphi-Dy1Talp1+Dx1,-DzTthetaSphi-Dy1,-Dz);
1474 pV[3] = G4Point3D(-DzTthetaCphi+Dy1Talp1+Dx2,-DzTthetaSphi+Dy1,-Dz);
1475 pV[4] = G4Point3D(-DzTthetaCphi+Dy1Talp1-Dx2,-DzTthetaSphi+Dy1,-Dz);
1476 pV[5] = G4Point3D( DzTthetaCphi-Dy2Talp2-Dx3, DzTthetaSphi-Dy2, Dz);
1477 pV[6] = G4Point3D( DzTthetaCphi-Dy2Talp2+Dx3, DzTthetaSphi-Dy2, Dz);
1478 pV[7] = G4Point3D( DzTthetaCphi+Dy2Talp2+Dx4, DzTthetaSphi+Dy2, Dz);
1479 pV[8] = G4Point3D( DzTthetaCphi+Dy2Talp2-Dx4, DzTthetaSphi+Dy2, Dz);
1480
1481 CreatePrism();
1482}
1483
1484HepPolyhedronTrap::~HepPolyhedronTrap() {}
1485
1486HepPolyhedronPara::HepPolyhedronPara(G4double Dx, G4double Dy, G4double Dz,
1487 G4double Alpha, G4double Theta,
1488 G4double Phi)
1489 : HepPolyhedronTrap(Dz, Theta, Phi, Dy, Dx, Dx, Alpha, Dy, Dx, Dx, Alpha) {}
1490
1491HepPolyhedronPara::~HepPolyhedronPara() {}
1492
1493HepPolyhedronParaboloid::HepPolyhedronParaboloid(G4double r1,
1494 G4double r2,
1495 G4double dz,
1496 G4double sPhi,
1497 G4double dPhi)
1498/***********************************************************************
1499 * *
1500 * Name: HepPolyhedronParaboloid Date: 28.06.07 *
1501 * Author: L.Lindroos, T.Nikitina (CERN), July 2007 Revised: 28.06.07 *
1502 * *
1503 * Function: Constructor for paraboloid *
1504 * *
1505 * Input: r1 - inside and outside radiuses at -Dz *
1506 * r2 - inside and outside radiuses at +Dz *
1507 * dz - half length in Z *
1508 * sPhi - starting angle of the segment *
1509 * dPhi - segment range *
1510 * *
1511 ***********************************************************************/
1512{
1513 static G4double wholeCircle=twopi;
1514
1515 // C H E C K I N P U T P A R A M E T E R S
1516
1517 G4int k = 0;
1518 if (r1 < 0. || r2 <= 0.) k = 1;
1519
1520 if (dz <= 0.) k += 2;
1521
1522 G4double phi1, phi2, dphi;
1523
1524 if(dPhi < 0.)
1525 {
1526 phi2 = sPhi; phi1 = phi2 + dPhi;
1527 }
1528 else if(dPhi == 0.)
1529 {
1530 phi1 = sPhi; phi2 = phi1 + wholeCircle;
1531 }
1532 else
1533 {
1534 phi1 = sPhi; phi2 = phi1 + dPhi;
1535 }
1536 dphi = phi2 - phi1;
1537
1538 if (std::abs(dphi-wholeCircle) < perMillion) dphi = wholeCircle;
1539 if (dphi > wholeCircle) k += 4;
1540
1541 if (k != 0) {
1542 std::cerr << "HepPolyhedronParaboloid: error in input parameters";
1543 if ((k & 1) != 0) std::cerr << " (radiuses)";
1544 if ((k & 2) != 0) std::cerr << " (half-length)";
1545 if ((k & 4) != 0) std::cerr << " (angles)";
1546 std::cerr << std::endl;
1547 std::cerr << " r1=" << r1;
1548 std::cerr << " r2=" << r2;
1549 std::cerr << " dz=" << dz << " sPhi=" << sPhi << " dPhi=" << dPhi
1550 << std::endl;
1551 return;
1552 }
1553
1554 // P R E P A R E T W O P O L Y L I N E S
1555
1556 G4int n = GetNumberOfRotationSteps();
1557 G4double dl = (r2 - r1) / n;
1558 G4double k1 = (r2*r2 - r1*r1) / 2 / dz;
1559 G4double k2 = (r2*r2 + r1*r1) / 2;
1560
1561 G4double *zz = new G4double[n + 2], *rr = new G4double[n + 2];
1562
1563 zz[0] = dz;
1564 rr[0] = r2;
1565
1566 for(G4int i = 1; i < n - 1; i++)
1567 {
1568 rr[i] = rr[i-1] - dl;
1569 zz[i] = (rr[i]*rr[i] - k2) / k1;
1570 if(rr[i] < 0)
1571 {
1572 rr[i] = 0;
1573 zz[i] = 0;
1574 }
1575 }
1576
1577 zz[n-1] = -dz;
1578 rr[n-1] = r1;
1579
1580 zz[n] = dz;
1581 rr[n] = 0;
1582
1583 zz[n+1] = -dz;
1584 rr[n+1] = 0;
1585
1586 // R O T A T E P O L Y L I N E S
1587
1588 RotateAroundZ(0, phi1, dphi, n, 2, zz, rr, -1, -1);
1589 SetReferences();
1590
1591 delete [] zz;
1592 delete [] rr;
1593}
1594
1595HepPolyhedronParaboloid::~HepPolyhedronParaboloid() {}
1596
1597HepPolyhedronHype::HepPolyhedronHype(G4double r1,
1598 G4double r2,
1599 G4double sqrtan1,
1600 G4double sqrtan2,
1601 G4double halfZ)
1602/***********************************************************************
1603 * *
1604 * Name: HepPolyhedronHype Date: 14.04.08 *
1605 * Author: Tatiana Nikitina (CERN) Revised: 14.04.08 *
1606 * *
1607 * Function: Constructor for Hype *
1608 * *
1609 * Input: r1 - inside radius at z=0 *
1610 * r2 - outside radiuses at z=0 *
1611 * sqrtan1 - sqr of tan of Inner Stereo Angle *
1612 * sqrtan2 - sqr of tan of Outer Stereo Angle *
1613 * halfZ - half length in Z *
1614 * *
1615 ***********************************************************************/
1616{
1617 static G4double wholeCircle=twopi;
1618
1619 // C H E C K I N P U T P A R A M E T E R S
1620
1621 G4int k = 0;
1622 if (r2 < 0. || r1 < 0. ) k = 1;
1623 if (r1 > r2 ) k = 1;
1624 if (r1 == r2) k = 1;
1625
1626 if (halfZ <= 0.) k += 2;
1627
1628 if (sqrtan1<0.||sqrtan2<0.) k += 4;
1629
1630 if (k != 0)
1631 {
1632 std::cerr << "HepPolyhedronHype: error in input parameters";
1633 if ((k & 1) != 0) std::cerr << " (radiuses)";
1634 if ((k & 2) != 0) std::cerr << " (half-length)";
1635 if ((k & 4) != 0) std::cerr << " (angles)";
1636 std::cerr << std::endl;
1637 std::cerr << " r1=" << r1 << " r2=" << r2;
1638 std::cerr << " halfZ=" << halfZ << " sqrTan1=" << sqrtan1
1639 << " sqrTan2=" << sqrtan2
1640 << std::endl;
1641 return;
1642 }
1643
1644 // P R E P A R E T W O P O L Y L I N E S
1645
1646 G4int n = GetNumberOfRotationSteps();
1647 G4double dz = 2.*halfZ / n;
1648 G4double k1 = r1*r1;
1649 G4double k2 = r2*r2;
1650
1651 G4double *zz = new G4double[n+n+1], *rr = new G4double[n+n+1];
1652
1653 zz[0] = halfZ;
1654 rr[0] = std::sqrt(sqrtan2*halfZ*halfZ+k2);
1655
1656 for(G4int i = 1; i < n-1; i++)
1657 {
1658 zz[i] = zz[i-1] - dz;
1659 rr[i] =std::sqrt(sqrtan2*zz[i]*zz[i]+k2);
1660 }
1661
1662 zz[n-1] = -halfZ;
1663 rr[n-1] = rr[0];
1664
1665 zz[n] = halfZ;
1666 rr[n] = std::sqrt(sqrtan1*halfZ*halfZ+k1);
1667
1668 for(G4int i = n+1; i < n+n; i++)
1669 {
1670 zz[i] = zz[i-1] - dz;
1671 rr[i] =std::sqrt(sqrtan1*zz[i]*zz[i]+k1);
1672 }
1673 zz[n+n] = -halfZ;
1674 rr[n+n] = rr[n];
1675
1676 // R O T A T E P O L Y L I N E S
1677
1678 RotateAroundZ(0, 0., wholeCircle, n, n, zz, rr, -1, -1);
1679 SetReferences();
1680
1681 delete [] zz;
1682 delete [] rr;
1683}
1684
1685HepPolyhedronHype::~HepPolyhedronHype() {}
1686
1687HepPolyhedronCons::HepPolyhedronCons(G4double Rmn1,
1688 G4double Rmx1,
1689 G4double Rmn2,
1690 G4double Rmx2,
1691 G4double Dz,
1692 G4double Phi1,
1693 G4double Dphi)
1694/***********************************************************************
1695 * *
1696 * Name: HepPolyhedronCons::HepPolyhedronCons Date: 15.12.96 *
1697 * Author: E.Chernyaev (IHEP/Protvino) Revised: 15.12.96 *
1698 * *
1699 * Function: Constructor for CONS, TUBS, CONE, TUBE *
1700 * *
1701 * Input: Rmn1, Rmx1 - inside and outside radiuses at -Dz *
1702 * Rmn2, Rmx2 - inside and outside radiuses at +Dz *
1703 * Dz - half length in Z *
1704 * Phi1 - starting angle of the segment *
1705 * Dphi - segment range *
1706 * *
1707 ***********************************************************************/
1708{
1709 static G4double wholeCircle=twopi;
1710
1711 // C H E C K I N P U T P A R A M E T E R S
1712
1713 G4int k = 0;
1714 if (Rmn1 < 0. || Rmx1 < 0. || Rmn2 < 0. || Rmx2 < 0.) k = 1;
1715 if (Rmn1 > Rmx1 || Rmn2 > Rmx2) k = 1;
1716 if (Rmn1 == Rmx1 && Rmn2 == Rmx2) k = 1;
1717
1718 if (Dz <= 0.) k += 2;
1719
1720 G4double phi1, phi2, dphi;
1721 if (Dphi < 0.) {
1722 phi2 = Phi1; phi1 = phi2 - Dphi;
1723 }else if (Dphi == 0.) {
1724 phi1 = Phi1; phi2 = phi1 + wholeCircle;
1725 }else{
1726 phi1 = Phi1; phi2 = phi1 + Dphi;
1727 }
1728 dphi = phi2 - phi1;
1729 if (std::abs(dphi-wholeCircle) < perMillion) dphi = wholeCircle;
1730 if (dphi > wholeCircle) k += 4;
1731
1732 if (k != 0) {
1733 std::cerr << "HepPolyhedronCone(s)/Tube(s): error in input parameters";
1734 if ((k & 1) != 0) std::cerr << " (radiuses)";
1735 if ((k & 2) != 0) std::cerr << " (half-length)";
1736 if ((k & 4) != 0) std::cerr << " (angles)";
1737 std::cerr << std::endl;
1738 std::cerr << " Rmn1=" << Rmn1 << " Rmx1=" << Rmx1;
1739 std::cerr << " Rmn2=" << Rmn2 << " Rmx2=" << Rmx2;
1740 std::cerr << " Dz=" << Dz << " Phi1=" << Phi1 << " Dphi=" << Dphi
1741 << std::endl;
1742 return;
1743 }
1744
1745 // P R E P A R E T W O P O L Y L I N E S
1746
1747 G4double zz[4], rr[4];
1748 zz[0] = Dz;
1749 zz[1] = -Dz;
1750 zz[2] = Dz;
1751 zz[3] = -Dz;
1752 rr[0] = Rmx2;
1753 rr[1] = Rmx1;
1754 rr[2] = Rmn2;
1755 rr[3] = Rmn1;
1756
1757 // R O T A T E P O L Y L I N E S
1758
1759 RotateAroundZ(0, phi1, dphi, 2, 2, zz, rr, -1, -1);
1760 SetReferences();
1761}
1762
1763HepPolyhedronCons::~HepPolyhedronCons() {}
1764
1765HepPolyhedronCone::HepPolyhedronCone(G4double Rmn1, G4double Rmx1,
1766 G4double Rmn2, G4double Rmx2,
1767 G4double Dz) :
1768 HepPolyhedronCons(Rmn1, Rmx1, Rmn2, Rmx2, Dz, 0*deg, 360*deg) {}
1769
1770HepPolyhedronCone::~HepPolyhedronCone() {}
1771
1772HepPolyhedronTubs::HepPolyhedronTubs(G4double Rmin, G4double Rmax,
1773 G4double Dz,
1774 G4double Phi1, G4double Dphi)
1775 : HepPolyhedronCons(Rmin, Rmax, Rmin, Rmax, Dz, Phi1, Dphi) {}
1776
1777HepPolyhedronTubs::~HepPolyhedronTubs() {}
1778
1779HepPolyhedronTube::HepPolyhedronTube (G4double Rmin, G4double Rmax,
1780 G4double Dz)
1781 : HepPolyhedronCons(Rmin, Rmax, Rmin, Rmax, Dz, 0*deg, 360*deg) {}
1782
1783HepPolyhedronTube::~HepPolyhedronTube () {}
1784
1785HepPolyhedronPgon::HepPolyhedronPgon(G4double phi,
1786 G4double dphi,
1787 G4int npdv,
1788 G4int nz,
1789 const G4double *z,
1790 const G4double *rmin,
1791 const G4double *rmax)
1792/***********************************************************************
1793 * *
1794 * Name: HepPolyhedronPgon Date: 09.12.96 *
1795 * Author: E.Chernyaev Revised: *
1796 * *
1797 * Function: Constructor of polyhedron for PGON, PCON *
1798 * *
1799 * Input: phi - initial phi *
1800 * dphi - delta phi *
1801 * npdv - number of steps along phi *
1802 * nz - number of z-planes (at least two) *
1803 * z[] - z coordinates of the slices *
1804 * rmin[] - smaller r at the slices *
1805 * rmax[] - bigger r at the slices *
1806 * *
1807 ***********************************************************************/
1808{
1809 // C H E C K I N P U T P A R A M E T E R S
1810
1811 if (dphi <= 0. || dphi > twopi) {
1812 std::cerr
1813 << "HepPolyhedronPgon/Pcon: wrong delta phi = " << dphi
1814 << std::endl;
1815 return;
1816 }
1817
1818 if (nz < 2) {
1819 std::cerr
1820 << "HepPolyhedronPgon/Pcon: number of z-planes less than two = " << nz
1821 << std::endl;
1822 return;
1823 }
1824
1825 if (npdv < 0) {
1826 std::cerr
1827 << "HepPolyhedronPgon/Pcon: error in number of phi-steps =" << npdv
1828 << std::endl;
1829 return;
1830 }
1831
1832 G4int i;
1833 for (i=0; i<nz; i++) {
1834 if (rmin[i] < 0. || rmax[i] < 0. || rmin[i] > rmax[i]) {
1835 std::cerr
1836 << "HepPolyhedronPgon: error in radiuses rmin[" << i << "]="
1837 << rmin[i] << " rmax[" << i << "]=" << rmax[i]
1838 << std::endl;
1839 return;
1840 }
1841 }
1842
1843 // P R E P A R E T W O P O L Y L I N E S
1844
1845 G4double *zz, *rr;
1846 zz = new G4double[2*nz];
1847 rr = new G4double[2*nz];
1848
1849 if (z[0] > z[nz-1]) {
1850 for (i=0; i<nz; i++) {
1851 zz[i] = z[i];
1852 rr[i] = rmax[i];
1853 zz[i+nz] = z[i];
1854 rr[i+nz] = rmin[i];
1855 }
1856 }else{
1857 for (i=0; i<nz; i++) {
1858 zz[i] = z[nz-i-1];
1859 rr[i] = rmax[nz-i-1];
1860 zz[i+nz] = z[nz-i-1];
1861 rr[i+nz] = rmin[nz-i-1];
1862 }
1863 }
1864
1865 // R O T A T E P O L Y L I N E S
1866
1867 RotateAroundZ(npdv, phi, dphi, nz, nz, zz, rr, -1, (npdv == 0) ? -1 : 1);
1868 SetReferences();
1869
1870 delete [] zz;
1871 delete [] rr;
1872}
1873
1874HepPolyhedronPgon::~HepPolyhedronPgon() {}
1875
1876HepPolyhedronPcon::HepPolyhedronPcon(G4double phi, G4double dphi, G4int nz,
1877 const G4double *z,
1878 const G4double *rmin,
1879 const G4double *rmax)
1880 : HepPolyhedronPgon(phi, dphi, 0, nz, z, rmin, rmax) {}
1881
1882HepPolyhedronPcon::~HepPolyhedronPcon() {}
1883
1884HepPolyhedronSphere::HepPolyhedronSphere(G4double rmin, G4double rmax,
1885 G4double phi, G4double dphi,
1886 G4double the, G4double dthe)
1887/***********************************************************************
1888 * *
1889 * Name: HepPolyhedronSphere Date: 11.12.96 *
1890 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
1891 * *
1892 * Function: Constructor of polyhedron for SPHERE *
1893 * *
1894 * Input: rmin - internal radius *
1895 * rmax - external radius *
1896 * phi - initial phi *
1897 * dphi - delta phi *
1898 * the - initial theta *
1899 * dthe - delta theta *
1900 * *
1901 ***********************************************************************/
1902{
1903 // C H E C K I N P U T P A R A M E T E R S
1904
1905 if (dphi <= 0. || dphi > twopi) {
1906 std::cerr
1907 << "HepPolyhedronSphere: wrong delta phi = " << dphi
1908 << std::endl;
1909 return;
1910 }
1911
1912 if (the < 0. || the > pi) {
1913 std::cerr
1914 << "HepPolyhedronSphere: wrong theta = " << the
1915 << std::endl;
1916 return;
1917 }
1918
1919 if (dthe <= 0. || dthe > pi) {
1920 std::cerr
1921 << "HepPolyhedronSphere: wrong delta theta = " << dthe
1922 << std::endl;
1923 return;
1924 }
1925
1926 if (the+dthe > pi) {
1927 std::cerr
1928 << "HepPolyhedronSphere: wrong theta + delta theta = "
1929 << the << " " << dthe
1930 << std::endl;
1931 return;
1932 }
1933
1934 if (rmin < 0. || rmin >= rmax) {
1935 std::cerr
1936 << "HepPolyhedronSphere: error in radiuses"
1937 << " rmin=" << rmin << " rmax=" << rmax
1938 << std::endl;
1939 return;
1940 }
1941
1942 // P R E P A R E T W O P O L Y L I N E S
1943
1944 G4int nds = (GetNumberOfRotationSteps() + 1) / 2;
1945 G4int np1 = G4int(dthe*nds/pi+.5) + 1;
1946 if (np1 <= 1) np1 = 2;
1947 G4int np2 = rmin < perMillion ? 1 : np1;
1948
1949 G4double *zz, *rr;
1950 zz = new G4double[np1+np2];
1951 rr = new G4double[np1+np2];
1952
1953 G4double a = dthe/(np1-1);
1954 G4double cosa, sina;
1955 for (G4int i=0; i<np1; i++) {
1956 cosa = std::cos(the+i*a);
1957 sina = std::sin(the+i*a);
1958 zz[i] = rmax*cosa;
1959 rr[i] = rmax*sina;
1960 if (np2 > 1) {
1961 zz[i+np1] = rmin*cosa;
1962 rr[i+np1] = rmin*sina;
1963 }
1964 }
1965 if (np2 == 1) {
1966 zz[np1] = 0.;
1967 rr[np1] = 0.;
1968 }
1969
1970 // R O T A T E P O L Y L I N E S
1971
1972 RotateAroundZ(0, phi, dphi, np1, np2, zz, rr, -1, -1);
1973 SetReferences();
1974
1975 delete [] zz;
1976 delete [] rr;
1977}
1978
1979HepPolyhedronSphere::~HepPolyhedronSphere() {}
1980
1981HepPolyhedronTorus::HepPolyhedronTorus(G4double rmin,
1982 G4double rmax,
1983 G4double rtor,
1984 G4double phi,
1985 G4double dphi)
1986/***********************************************************************
1987 * *
1988 * Name: HepPolyhedronTorus Date: 11.12.96 *
1989 * Author: E.Chernyaev (IHEP/Protvino) Revised: *
1990 * *
1991 * Function: Constructor of polyhedron for TORUS *
1992 * *
1993 * Input: rmin - internal radius *
1994 * rmax - external radius *
1995 * rtor - radius of torus *
1996 * phi - initial phi *
1997 * dphi - delta phi *
1998 * *
1999 ***********************************************************************/
2000{
2001 // C H E C K I N P U T P A R A M E T E R S
2002
2003 if (dphi <= 0. || dphi > twopi) {
2004 std::cerr
2005 << "HepPolyhedronTorus: wrong delta phi = " << dphi
2006 << std::endl;
2007 return;
2008 }
2009
2010 if (rmin < 0. || rmin >= rmax || rmax >= rtor) {
2011 std::cerr
2012 << "HepPolyhedronTorus: error in radiuses"
2013 << " rmin=" << rmin << " rmax=" << rmax << " rtorus=" << rtor
2014 << std::endl;
2015 return;
2016 }
2017
2018 // P R E P A R E T W O P O L Y L I N E S
2019
2020 G4int np1 = GetNumberOfRotationSteps();
2021 G4int np2 = rmin < perMillion ? 1 : np1;
2022
2023 G4double *zz, *rr;
2024 zz = new G4double[np1+np2];
2025 rr = new G4double[np1+np2];
2026
2027 G4double a = twopi/np1;
2028 G4double cosa, sina;
2029 for (G4int i=0; i<np1; i++) {
2030 cosa = std::cos(i*a);
2031 sina = std::sin(i*a);
2032 zz[i] = rmax*cosa;
2033 rr[i] = rtor+rmax*sina;
2034 if (np2 > 1) {
2035 zz[i+np1] = rmin*cosa;
2036 rr[i+np1] = rtor+rmin*sina;
2037 }
2038 }
2039 if (np2 == 1) {
2040 zz[np1] = 0.;
2041 rr[np1] = rtor;
2042 np2 = -1;
2043 }
2044
2045 // R O T A T E P O L Y L I N E S
2046
2047 RotateAroundZ(0, phi, dphi, -np1, -np2, zz, rr, -1,-1);
2048 SetReferences();
2049
2050 delete [] zz;
2051 delete [] rr;
2052}
2053
2054HepPolyhedronTorus::~HepPolyhedronTorus() {}
2055
2056HepPolyhedronEllipsoid::HepPolyhedronEllipsoid(G4double ax, G4double by,
2057 G4double cz, G4double zCut1,
2058 G4double zCut2)
2059/***********************************************************************
2060 * *
2061 * Name: HepPolyhedronEllipsoid Date: 25.02.05 *
2062 * Author: G.Guerrieri Revised: *
2063 * *
2064 * Function: Constructor of polyhedron for ELLIPSOID *
2065 * *
2066 * Input: ax - semiaxis x *
2067 * by - semiaxis y *
2068 * cz - semiaxis z *
2069 * zCut1 - lower cut plane level (solid lies above this plane) *
2070 * zCut2 - upper cut plane level (solid lies below this plane) *
2071 * *
2072 ***********************************************************************/
2073{
2074 // C H E C K I N P U T P A R A M E T E R S
2075
2076 if (zCut1 >= cz || zCut2 <= -cz || zCut1 > zCut2) {
2077 std::cerr << "HepPolyhedronEllipsoid: wrong zCut1 = " << zCut1
2078 << " zCut2 = " << zCut2
2079 << " for given cz = " << cz << std::endl;
2080 return;
2081 }
2082 if (cz <= 0.0) {
2083 std::cerr << "HepPolyhedronEllipsoid: bad z semi-axis: cz = " << cz
2084 << std::endl;
2085 return;
2086 }
2087
2088 G4double dthe;
2089 G4double sthe;
2090 G4int cutflag;
2091 cutflag= 0;
2092 if (zCut2 >= cz)
2093 {
2094 sthe= 0.0;
2095 }
2096 else
2097 {
2098 sthe= std::acos(zCut2/cz);
2099 cutflag++;
2100 }
2101 if (zCut1 <= -cz)
2102 {
2103 dthe= pi - sthe;
2104 }
2105 else
2106 {
2107 dthe= std::acos(zCut1/cz)-sthe;
2108 cutflag++;
2109 }
2110
2111 // P R E P A R E T W O P O L Y L I N E S
2112 // generate sphere of radius cz first, then rescale x and y later
2113
2114 G4int nds = (GetNumberOfRotationSteps() + 1) / 2;
2115 G4int np1 = G4int(dthe*nds/pi) + 2 + cutflag;
2116
2117 G4double *zz, *rr;
2118 zz = new G4double[np1+1];
2119 rr = new G4double[np1+1];
2120 if (!zz || !rr)
2121 {
2122 G4Exception("HepPolyhedronEllipsoid::HepPolyhedronEllipsoid",
2123 "greps1002", FatalException, "Out of memory");
2124 }
2125
2126 G4double a = dthe/(np1-cutflag-1);
2127 G4double cosa, sina;
2128 G4int j=0;
2129 if (sthe > 0.0)
2130 {
2131 zz[j]= zCut2;
2132 rr[j]= 0.;
2133 j++;
2134 }
2135 for (G4int i=0; i<np1-cutflag; i++) {
2136 cosa = std::cos(sthe+i*a);
2137 sina = std::sin(sthe+i*a);
2138 zz[j] = cz*cosa;
2139 rr[j] = cz*sina;
2140 j++;
2141 }
2142 if (j < np1)
2143 {
2144 zz[j]= zCut1;
2145 rr[j]= 0.;
2146 j++;
2147 }
2148 if (j > np1)
2149 {
2150 std::cerr << "Logic error in HepPolyhedronEllipsoid, memory corrupted!"
2151 << std::endl;
2152 }
2153 if (j < np1)
2154 {
2155 std::cerr << "Warning: logic error in HepPolyhedronEllipsoid."
2156 << std::endl;
2157 np1= j;
2158 }
2159 zz[j] = 0.;
2160 rr[j] = 0.;
2161
2162
2163 // R O T A T E P O L Y L I N E S
2164
2165 RotateAroundZ(0, 0.0, twopi, np1, 1, zz, rr, -1, 1);
2166 SetReferences();
2167
2168 delete [] zz;
2169 delete [] rr;
2170
2171 // rescale x and y vertex coordinates
2172 {
2173 G4Point3D * p= pV;
2174 for (G4int i=0; i<nvert; i++, p++) {
2175 p->setX( p->x() * ax/cz );
2176 p->setY( p->y() * by/cz );
2177 }
2178 }
2179}
2180
2181HepPolyhedronEllipsoid::~HepPolyhedronEllipsoid() {}
2182
2183HepPolyhedronEllipticalCone::HepPolyhedronEllipticalCone(G4double ax,
2184 G4double ay,
2185 G4double h,
2186 G4double zTopCut)
2187/***********************************************************************
2188 * *
2189 * Name: HepPolyhedronEllipticalCone Date: 8.9.2005 *
2190 * Author: D.Anninos Revised: 9.9.2005 *
2191 * *
2192 * Function: Constructor for EllipticalCone *
2193 * *
2194 * Input: ax, ay - X & Y semi axes at z = 0 *
2195 * h - height of full cone *
2196 * zTopCut - Top Cut in Z Axis *
2197 * *
2198 ***********************************************************************/
2199{
2200 // C H E C K I N P U T P A R A M E T E R S
2201
2202 G4int k = 0;
2203 if ( (ax <= 0.) || (ay <= 0.) || (h <= 0.) || (zTopCut <= 0.) ) { k = 1; }
2204
2205 if (k != 0) {
2206 std::cerr << "HepPolyhedronCone: error in input parameters";
2207 std::cerr << std::endl;
2208 return;
2209 }
2210
2211 // P R E P A R E T W O P O L Y L I N E S
2212
2213 zTopCut = (h >= zTopCut ? zTopCut : h);
2214
2215 G4double *zz, *rr;
2216 zz = new G4double[4];
2217 rr = new G4double[4];
2218 zz[0] = zTopCut;
2219 zz[1] = -zTopCut;
2220 zz[2] = zTopCut;
2221 zz[3] = -zTopCut;
2222 rr[0] = (h-zTopCut);
2223 rr[1] = (h+zTopCut);
2224 rr[2] = 0.;
2225 rr[3] = 0.;
2226
2227 // R O T A T E P O L Y L I N E S
2228
2229 RotateAroundZ(0, 0., twopi, 2, 2, zz, rr, -1, -1);
2230 SetReferences();
2231
2232 delete [] zz;
2233 delete [] rr;
2234
2235 // rescale x and y vertex coordinates
2236 {
2237 G4Point3D * p= pV;
2238 for (G4int i=0; i<nvert; i++, p++) {
2239 p->setX( p->x() * ax );
2240 p->setY( p->y() * ay );
2241 }
2242 }
2243}
2244
2245HepPolyhedronEllipticalCone::~HepPolyhedronEllipticalCone() {}
2246
2247G4int HepPolyhedron::fNumberOfRotationSteps = DEFAULT_NUMBER_OF_STEPS;
2248/***********************************************************************
2249 * *
2250 * Name: HepPolyhedron::fNumberOfRotationSteps Date: 24.06.97 *
2251 * Author: J.Allison (Manchester University) Revised: *
2252 * *
2253 * Function: Number of steps for whole circle *
2254 * *
2255 ***********************************************************************/
2256
2257#include "BooleanProcessor.src"
2258
2259HepPolyhedron HepPolyhedron::add(const HepPolyhedron & p) const
2260/***********************************************************************
2261 * *
2262 * Name: HepPolyhedron::add Date: 19.03.00 *
2263 * Author: E.Chernyaev Revised: *
2264 * *
2265 * Function: Boolean "union" of two polyhedra *
2266 * *
2267 ***********************************************************************/
2268{
2269 G4int ierr;
2270 BooleanProcessor processor;
2271 return processor.execute(OP_UNION, *this, p,ierr);
2272}
2273
2274HepPolyhedron HepPolyhedron::intersect(const HepPolyhedron & p) const
2275/***********************************************************************
2276 * *
2277 * Name: HepPolyhedron::intersect Date: 19.03.00 *
2278 * Author: E.Chernyaev Revised: *
2279 * *
2280 * Function: Boolean "intersection" of two polyhedra *
2281 * *
2282 ***********************************************************************/
2283{
2284 G4int ierr;
2285 BooleanProcessor processor;
2286 return processor.execute(OP_INTERSECTION, *this, p,ierr);
2287}
2288
2289HepPolyhedron HepPolyhedron::subtract(const HepPolyhedron & p) const
2290/***********************************************************************
2291 * *
2292 * Name: HepPolyhedron::add Date: 19.03.00 *
2293 * Author: E.Chernyaev Revised: *
2294 * *
2295 * Function: Boolean "subtraction" of "p" from "this" *
2296 * *
2297 ***********************************************************************/
2298{
2299 G4int ierr;
2300 BooleanProcessor processor;
2301 return processor.execute(OP_SUBTRACTION, *this, p,ierr);
2302}
2303
2304//NOTE : include the code of HepPolyhedronProcessor here
2305// since there is no BooleanProcessor.h
2306
2307#undef INTERSECTION
2308
2309#include "HepPolyhedronProcessor.src"
2310
FatalException
@ FatalException
Definition: G4ExceptionSeverity.hh:60
G4Normal3D
HepGeom::Normal3D< G4double > G4Normal3D
Definition: G4Normal3D.hh:35
G4Point3D
HepGeom::Point3D< G4double > G4Point3D
Definition: G4Point3D.hh:35
G4double
double G4double
Definition: G4Types.hh:64
G4int
int G4int
Definition: G4Types.hh:66
G4bool
bool G4bool
Definition: G4Types.hh:67
G4Vector3D
HepGeom::Vector3D< G4double > G4Vector3D
Definition: G4Vector3D.hh:35
operator<<
std::ostream & operator<<(std::ostream &ostr, const G4Facet &facet)
Definition: HepPolyhedron.cc:83
DEFAULT_NUMBER_OF_STEPS
#define DEFAULT_NUMBER_OF_STEPS
Definition: HepPolyhedron.h:177
HepGeom::BasicVector3D::cross
BasicVector3D< T > cross(const BasicVector3D< T > &v) const
Definition: BasicVector3D.h:275
HepGeom::BasicVector3D::unit
BasicVector3D< T > unit() const
Definition: BasicVector3D.h:304
HepPolyhedronBox::~HepPolyhedronBox
virtual ~HepPolyhedronBox()
Definition: HepPolyhedron.cc:1429
HepPolyhedronBox::HepPolyhedronBox
HepPolyhedronBox(G4double Dx, G4double Dy, G4double Dz)
Definition: HepPolyhedron.cc:1426
HepPolyhedronCone::~HepPolyhedronCone
virtual ~HepPolyhedronCone()
Definition: HepPolyhedron.cc:1770
HepPolyhedronCone::HepPolyhedronCone
HepPolyhedronCone(G4double Rmn1, G4double Rmx1, G4double Rmn2, G4double Rmx2, G4double Dz)
Definition: HepPolyhedron.cc:1765
HepPolyhedronCons::~HepPolyhedronCons
virtual ~HepPolyhedronCons()
Definition: HepPolyhedron.cc:1763
HepPolyhedronCons::HepPolyhedronCons
HepPolyhedronCons(G4double Rmn1, G4double Rmx1, G4double Rmn2, G4double Rmx2, G4double Dz, G4double Phi1, G4double Dphi)
Definition: HepPolyhedron.cc:1687
HepPolyhedronEllipsoid::HepPolyhedronEllipsoid
HepPolyhedronEllipsoid(G4double dx, G4double dy, G4double dz, G4double zcut1, G4double zcut2)
Definition: HepPolyhedron.cc:2056
HepPolyhedronEllipticalCone::HepPolyhedronEllipticalCone
HepPolyhedronEllipticalCone(G4double dx, G4double dy, G4double z, G4double zcut1)
Definition: HepPolyhedron.cc:2183
HepPolyhedronHype::HepPolyhedronHype
HepPolyhedronHype(G4double r1, G4double r2, G4double tan1, G4double tan2, G4double halfZ)
Definition: HepPolyhedron.cc:1597
HepPolyhedronHype::~HepPolyhedronHype
virtual ~HepPolyhedronHype()
Definition: HepPolyhedron.cc:1685
HepPolyhedronPara::~HepPolyhedronPara
virtual ~HepPolyhedronPara()
Definition: HepPolyhedron.cc:1491
HepPolyhedronPara::HepPolyhedronPara
HepPolyhedronPara(G4double Dx, G4double Dy, G4double Dz, G4double Alpha, G4double Theta, G4double Phi)
Definition: HepPolyhedron.cc:1486
HepPolyhedronParaboloid::HepPolyhedronParaboloid
HepPolyhedronParaboloid(G4double r1, G4double r2, G4double dz, G4double Phi1, G4double Dphi)
Definition: HepPolyhedron.cc:1493
HepPolyhedronPcon::~HepPolyhedronPcon
virtual ~HepPolyhedronPcon()
Definition: HepPolyhedron.cc:1882
HepPolyhedronPcon::HepPolyhedronPcon
HepPolyhedronPcon(G4double phi, G4double dphi, G4int nz, const G4double *z, const G4double *rmin, const G4double *rmax)
Definition: HepPolyhedron.cc:1876
HepPolyhedronPgon::~HepPolyhedronPgon
virtual ~HepPolyhedronPgon()
Definition: HepPolyhedron.cc:1874
HepPolyhedronPgon::HepPolyhedronPgon
HepPolyhedronPgon(G4double phi, G4double dphi, G4int npdv, G4int nz, const G4double *z, const G4double *rmin, const G4double *rmax)
Definition: HepPolyhedron.cc:1785
HepPolyhedronSphere::~HepPolyhedronSphere
virtual ~HepPolyhedronSphere()
Definition: HepPolyhedron.cc:1979
HepPolyhedronSphere::HepPolyhedronSphere
HepPolyhedronSphere(G4double rmin, G4double rmax, G4double phi, G4double dphi, G4double the, G4double dthe)
Definition: HepPolyhedron.cc:1884
HepPolyhedronTorus::HepPolyhedronTorus
HepPolyhedronTorus(G4double rmin, G4double rmax, G4double rtor, G4double phi, G4double dphi)
Definition: HepPolyhedron.cc:1981
HepPolyhedronTorus::~HepPolyhedronTorus
virtual ~HepPolyhedronTorus()
Definition: HepPolyhedron.cc:2054
HepPolyhedronTrap::HepPolyhedronTrap
HepPolyhedronTrap(G4double Dz, G4double Theta, G4double Phi, G4double Dy1, G4double Dx1, G4double Dx2, G4double Alp1, G4double Dy2, G4double Dx3, G4double Dx4, G4double Alp2)
Definition: HepPolyhedron.cc:1431
HepPolyhedronTrap::~HepPolyhedronTrap
virtual ~HepPolyhedronTrap()
Definition: HepPolyhedron.cc:1484
HepPolyhedronTrd1::HepPolyhedronTrd1
HepPolyhedronTrd1(G4double Dx1, G4double Dx2, G4double Dy, G4double Dz)
Definition: HepPolyhedron.cc:1420
HepPolyhedronTrd1::~HepPolyhedronTrd1
virtual ~HepPolyhedronTrd1()
Definition: HepPolyhedron.cc:1424
HepPolyhedronTrd2::~HepPolyhedronTrd2
virtual ~HepPolyhedronTrd2()
Definition: HepPolyhedron.cc:1418
HepPolyhedronTrd2::HepPolyhedronTrd2
HepPolyhedronTrd2(G4double Dx1, G4double Dx2, G4double Dy1, G4double Dy2, G4double Dz)
Definition: HepPolyhedron.cc:1386
HepPolyhedronTube::~HepPolyhedronTube
virtual ~HepPolyhedronTube()
Definition: HepPolyhedron.cc:1783
HepPolyhedronTube::HepPolyhedronTube
HepPolyhedronTube(G4double Rmin, G4double Rmax, G4double Dz)
Definition: HepPolyhedron.cc:1779
HepPolyhedronTubs::HepPolyhedronTubs
HepPolyhedronTubs(G4double Rmin, G4double Rmax, G4double Dz, G4double Phi1, G4double Dphi)
Definition: HepPolyhedron.cc:1772
HepPolyhedronTubs::~HepPolyhedronTubs
virtual ~HepPolyhedronTubs()
Definition: HepPolyhedron.cc:1777
HepPolyhedron::GetNextEdgeIndeces
G4bool GetNextEdgeIndeces(G4int &i1, G4int &i2, G4int &edgeFlag, G4int &iface1, G4int &iface2) const
Definition: HepPolyhedron.cc:953
HepPolyhedron::SetNumberOfRotationSteps
static void SetNumberOfRotationSteps(G4int n)
Definition: HepPolyhedron.cc:207
HepPolyhedron::RotateAroundZ
void RotateAroundZ(G4int nstep, G4double phi, G4double dphi, G4int np1, G4int np2, const G4double *z, G4double *r, G4int nodeVis, G4int edgeVis)
Definition: HepPolyhedron.cc:421
HepPolyhedron::GetUnitNormal
G4Normal3D GetUnitNormal(G4int iFace) const
Definition: HepPolyhedron.cc:1174
HepPolyhedron::GetFacet
void GetFacet(G4int iFace, G4int &n, G4int *iNodes, G4int *edgeFlags=0, G4int *iFaces=0) const
Definition: HepPolyhedron.cc:1063
HepPolyhedron::createTwistedTrap
G4int createTwistedTrap(G4double Dz, const G4double xy1[][2], const G4double xy2[][2])
Definition: HepPolyhedron.cc:1287
HepPolyhedron::pV
G4Point3D * pV
Definition: HepPolyhedron.h:201
HepPolyhedron::SetReferences
void SetReferences()
Definition: HepPolyhedron.cc:668
HepPolyhedron::fNumberOfRotationSteps
static G4int fNumberOfRotationSteps
Definition: HepPolyhedron.h:199
HepPolyhedron::GetNextNormal
G4bool GetNextNormal(G4Normal3D &normal) const
Definition: HepPolyhedron.cc:1199
HepPolyhedron::operator=
HepPolyhedron & operator=(const HepPolyhedron &from)
Definition: HepPolyhedron.cc:119
HepPolyhedron::FindNodeNormal
G4Normal3D FindNodeNormal(G4int iFace, G4int iNode) const
Definition: HepPolyhedron.cc:166
HepPolyhedron::GetNumberOfRotationSteps
static G4int GetNumberOfRotationSteps()
Definition: HepPolyhedron.cc:194
HepPolyhedron::GetNormal
G4Normal3D GetNormal(G4int iFace) const
Definition: HepPolyhedron.cc:1149
HepPolyhedron::FindNeighbour
G4int FindNeighbour(G4int iFace, G4int iNode, G4int iOrder) const
Definition: HepPolyhedron.cc:138
HepPolyhedron::GetVertex
G4Point3D GetVertex(G4int index) const
Definition: HepPolyhedron.cc:883
HepPolyhedron::GetSurfaceArea
G4double GetSurfaceArea() const
Definition: HepPolyhedron.cc:1236
HepPolyhedron::subtract
HepPolyhedron subtract(const HepPolyhedron &p) const
Definition: HepPolyhedron.cc:2289
HepPolyhedron::GetNextVertex
G4bool GetNextVertex(G4Point3D &vertex, G4int &edgeFlag) const
Definition: HepPolyhedron.cc:903
HepPolyhedron::GetVolume
G4double GetVolume() const
Definition: HepPolyhedron.cc:1258
HepPolyhedron::intersect
HepPolyhedron intersect(const HepPolyhedron &p) const
Definition: HepPolyhedron.cc:2274
HepPolyhedron::Transform
HepPolyhedron & Transform(const G4Transform3D &t)
Definition: HepPolyhedron.cc:831
HepPolyhedron::createPolyhedron
G4int createPolyhedron(G4int Nnodes, G4int Nfaces, const G4double xyz[][3], const G4int faces[][4])
Definition: HepPolyhedron.cc:1356
HepPolyhedron::AllocateMemory
void AllocateMemory(G4int Nvert, G4int Nface)
Definition: HepPolyhedron.cc:242
HepPolyhedron::GetNextUnitNormal
G4bool GetNextUnitNormal(G4Normal3D &normal) const
Definition: HepPolyhedron.cc:1220
HepPolyhedron::ResetNumberOfRotationSteps
static void ResetNumberOfRotationSteps()
Definition: HepPolyhedron.cc:229
HepPolyhedron::RotateEdge
void RotateEdge(G4int k1, G4int k2, G4double r1, G4double r2, G4int v1, G4int v2, G4int vEdge, G4bool ifWholeCircle, G4int ns, G4int &kface)
Definition: HepPolyhedron.cc:288
HepPolyhedron::GetNextFacet
G4bool GetNextFacet(G4int &n, G4Point3D *nodes, G4int *edgeFlags=0, G4Normal3D *normals=0) const
Definition: HepPolyhedron.cc:1119
HepPolyhedron::GetNextVertexIndex
G4bool GetNextVertexIndex(G4int &index, G4int &edgeFlag) const
Definition: HepPolyhedron.cc:856
HepPolyhedron::pF
G4Facet * pF
Definition: HepPolyhedron.h:202
HepPolyhedron::add
HepPolyhedron add(const HepPolyhedron &p) const
Definition: HepPolyhedron.cc:2259
HepPolyhedron::GetNextEdge
G4bool GetNextEdge(G4Point3D &p1, G4Point3D &p2, G4int &edgeFlag) const
Definition: HepPolyhedron.cc:1022
HepPolyhedron::SetSideFacets
void SetSideFacets(G4int ii[4], G4int vv[4], G4int *kk, G4double *r, G4double dphi, G4int ns, G4int &kface)
Definition: HepPolyhedron.cc:350
G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
processor
#define processor
Definition: xmlparse.cc:600
const
#define const
Definition: zconf.h:118