Creates user defined polyhedron. This function allows to the user to define arbitrary polyhedron. The faces of the polyhedron should be either triangles or planar quadrilateral. Nodes of a face are defined by indexes pointing to the elements in the xyz array. Numeration of the elements in the array starts from 1 (like in fortran). The indexes can be positive or negative. Negative sign means that the corresponding edge is invisible. The normal of the face should be directed to exterior of the polyhedron.
691 const G4int numSide =
698 typedef G4int int4[4];
705 std::vector<G4bool> chopped(
numCorner,
false);
706 std::vector<G4int*> triQuads;
709 while (remaining >= 3)
714 G4int iStepper = iStarter;
717 if (A < 0) { A = iStepper; }
718 else if (
B < 0) {
B = iStepper; }
719 else if (
C < 0) {
C = iStepper; }
722 if (++iStepper >=
numCorner) { iStepper = 0; }
724 while (chopped[iStepper]);
726 while (
C < 0 && iStepper != iStarter);
741 triQuads.push_back(tq);
749 if (++iStarter >=
numCorner) { iStarter = 0; }
751 while (chopped[iStarter]);
757 nFaces = numSide *
numCorner + 2 * triQuads.size();
758 faces_vec =
new int4[nFaces];
762 for (
G4int iEnd = 0; iEnd < 2; ++iEnd)
764 for (
size_t i = 0; i < triQuads.size(); ++i)
777 a = triQuads[i][0] + addition;
778 b = triQuads[i][2] + addition;
779 c = triQuads[i][1] + addition;
782 G4int bc = std::abs(c - b);
783 G4int ca = std::abs(a - c);
784 faces_vec[iface][0] = (ab == 1 || ab == d)? a: -a;
785 faces_vec[iface][1] = (bc == 1 || bc == d)? b: -b;
786 faces_vec[iface][2] = (ca == 1 || ca == d)? c: -c;
787 faces_vec[iface][3] = 0;
794 xyz =
new double3[nNodes];
798 for (
G4int iSide = 0; iSide < numSide; ++iSide)
802 xyz[ixyz][0] =
corners[iCorner].
r * std::cos(phi);
803 xyz[ixyz][1] =
corners[iCorner].
r * std::sin(phi);
807 if (iCorner < numCorner - 1)
809 faces_vec[iface][0] = ixyz + 1;
810 faces_vec[iface][1] = -(ixyz + numCorner + 1);
811 faces_vec[iface][2] = ixyz + numCorner + 2;
812 faces_vec[iface][3] = ixyz + 2;
816 faces_vec[iface][0] = ixyz + 1;
817 faces_vec[iface][1] = -(ixyz + numCorner + 1);
818 faces_vec[iface][2] = ixyz + 2;
819 faces_vec[iface][3] = ixyz - numCorner + 2;
822 else if (iSide == numSide - 1)
824 if (iCorner < numCorner - 1)
826 faces_vec[iface][0] = ixyz + 1;
827 faces_vec[iface][1] = ixyz + numCorner + 1;
828 faces_vec[iface][2] = ixyz + numCorner + 2;
829 faces_vec[iface][3] = -(ixyz + 2);
833 faces_vec[iface][0] = ixyz + 1;
834 faces_vec[iface][1] = ixyz + numCorner + 1;
835 faces_vec[iface][2] = ixyz + 2;
836 faces_vec[iface][3] = -(ixyz - numCorner + 2);
841 if (iCorner < numCorner - 1)
843 faces_vec[iface][0] = ixyz + 1;
844 faces_vec[iface][1] = -(ixyz + numCorner + 1);
845 faces_vec[iface][2] = ixyz + numCorner + 2;
846 faces_vec[iface][3] = -(ixyz + 2);
850 faces_vec[iface][0] = ixyz + 1;
851 faces_vec[iface][1] = -(ixyz + numCorner + 1);
852 faces_vec[iface][2] = ixyz + 2;
853 faces_vec[iface][3] = -(ixyz - numCorner + 2);
866 xyz[ixyz][0] =
corners[iCorner].
r * std::cos(phi);
867 xyz[ixyz][1] =
corners[iCorner].
r * std::sin(phi);
876 xyz =
new double3[nNodes];
877 faces_vec =
new int4[nFaces];
880 G4int ixyz = 0, iface = 0;
881 for (
G4int iSide = 0; iSide < numSide; ++iSide)
885 xyz[ixyz][0] =
corners[iCorner].
r * std::cos(phi);
886 xyz[ixyz][1] =
corners[iCorner].
r * std::sin(phi);
889 if (iSide < numSide - 1)
891 if (iCorner < numCorner - 1)
893 faces_vec[iface][0] = ixyz + 1;
894 faces_vec[iface][1] = -(ixyz + numCorner + 1);
895 faces_vec[iface][2] = ixyz + numCorner + 2;
896 faces_vec[iface][3] = -(ixyz + 2);
900 faces_vec[iface][0] = ixyz + 1;
901 faces_vec[iface][1] = -(ixyz + numCorner + 1);
902 faces_vec[iface][2] = ixyz + 2;
903 faces_vec[iface][3] = -(ixyz - numCorner + 2);
908 if (iCorner < numCorner - 1)
910 faces_vec[iface][0] = ixyz + 1;
911 faces_vec[iface][1] = -(ixyz + numCorner - nFaces + 1);
912 faces_vec[iface][2] = ixyz + numCorner - nFaces + 2;
913 faces_vec[iface][3] = -(ixyz + 2);
917 faces_vec[iface][0] = ixyz + 1;
918 faces_vec[iface][1] = -(ixyz - nFaces + numCorner + 1);
919 faces_vec[iface][2] = ixyz - nFaces + 2;
920 faces_vec[iface][3] = -(ixyz - numCorner + 2);
935 std::ostringstream message;
936 message <<
"Problem creating G4Polyhedron for: " <<
GetName();
937 G4Exception(
"G4GenericPolycone::CreatePolyhedron()",
"GeomSolids1002",
G4PolyconeSideRZ * corners
G4int createPolyhedron(G4int Nnodes, G4int Nfaces, const G4double xyz[][3], const G4int faces[][4])
double B(double temperature)
static constexpr double twopi
double A(double temperature)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
static G4int GetNumberOfRotationSteps()