125 ((
TGeoXtru *)
this)->CreateThreadData(1);
134 std::lock_guard<std::mutex> guard(
fMutex);
135 std::vector<ThreadData_t *>::iterator i =
fThreadData.begin();
149 std::lock_guard<std::mutex> guard(
fMutex);
152 for (
Int_t tid = 0; tid < nthreads; tid++) {
161 td.fPoly->SetXY(td.fXc, td.fYc);
162 td.fPoly->FinishPolygon();
163 if (tid == 0 && td.fPoly->IsIllegalCheck()) {
164 Error(
"DefinePolygon",
"Shape %s of type XTRU has an illegal polygon.",
GetName());
224 Error(
"ctor",
"Cannot create TGeoXtru %s with less than 2 Z planes", GetName());
259 SetDimensions(param);
305 area = td.fPoly->Area();
306 for (iz = 0; iz <
fNz - 1; iz++) {
307 dz =
fZ[iz + 1] -
fZ[iz];
312 capacity += (area * dz / 3.) * (sc1 * sc1 + sc1 * sc2 + sc2 * sc2);
324 Error(
"ComputeBBox",
"In shape %s polygon not defined",
GetName());
337 if (td.fXc[j] <
xmin)
339 if (td.fXc[j] >
xmax)
341 if (td.fYc[j] <
ymin)
343 if (td.fYc[j] >
ymax)
349 fOrigin[2] = 0.5 * (zmin + zmax);
352 fDZ = 0.5 * (zmax - zmin);
362 memset(norm, 0, 3 *
sizeof(
Double_t));
363 norm[2] = (dir[2] > 0) ? 1 : -1;
369 Double_t ndotd = norm[0] * dir[0] + norm[1] * dir[1] + norm[2] * dir[2];
385 if (point[2] <
fZ[0])
387 if (point[2] >
fZ[
fNz - 1])
390 if (iz < 0 || iz ==
fNz - 1)
395 if (td.fPoly->Contains(point))
399 return td.fPoly->Contains(point);
402 return td.fPoly->Contains(point);
410 return td.fPoly->Contains(point);
447 snext = (
fZ[iz] - point[2]) / dir[2];
450 pt[0] = point[0] + snext * dir[0];
451 pt[1] = point[1] + snext * dir[1];
452 pt[2] = point[2] + snext * dir[2];
457 if (!td.fPoly->Contains(
pt))
463 Double_t ndotd = norm[0] * dir[0] + norm[1] * dir[1] + norm[2] * dir[2];
467 safe = (vert[0] - point[0]) * norm[0] + (vert[1] - point[1]) * norm[1] + (vert[2] - point[2]) * norm[2];
474 safe = (point[0] - vert[0]) * norm[0] + (point[1] - vert[1]) * norm[1] + (point[2] - vert[2]) * norm[2];
478 snext = safe / ndotd;
481 if (
fZ[iz] <
fZ[iz + 1]) {
482 znew = point[2] + snext * dir[2];
485 if (znew >
fZ[iz + 1])
488 pt[0] = point[0] + snext * dir[0];
489 pt[1] = point[1] + snext * dir[1];
490 pt[2] = point[2] + snext * dir[2];
504 if (iact < 3 && safe) {
508 if (iact == 1 && step < *safe)
536 Bool_t convex = td.fPoly->IsConvex();
542 Int_t iv, ipl, inext;
545 for (iv = 0; iv <
fNvert; iv++) {
561 Int_t incseg = (dir[2] > 0) ? 1 : -1;
564 while (iz >= 0 && iz <
fNz - 1) {
566 ipl = iz + ((incseg + 1) >> 1);
567 inext = ipl + incseg;
568 sz = (
fZ[ipl] - point[2]) / dir[2];
572 pt[0] = point[0] + sz * dir[0];
573 pt[1] = point[1] + sz * dir[1];
575 if (td.fPoly->Contains(
pt)) {
577 if (ipl == 0 || ipl ==
fNz - 1) {
588 if (!td.fPoly->Contains(
pt)) {
603 for (iv = 0; iv <
fNvert; iv++) {
629 if (iact < 3 && safe) {
633 if (iact == 1 && step < *safe)
654 snext = (
fZ[0] - point[2]) / dir[2];
657 for (i = 0; i < 3; i++)
658 pt[i] = point[i] + snext * dir[i];
660 if (td.fPoly->Contains(
pt)) {
671 snext = (
fZ[
fNz - 1] - point[2]) / dir[2];
674 for (i = 0; i < 3; i++)
675 pt[i] = point[i] + snext * dir[i];
677 if (td.fPoly->Contains(
pt)) {
694 for (i = 0; i < 3; i++)
695 pt[i] += dist * dir[i];
699 else if (iz ==
fNz - 1)
707 Bool_t convex = td.fPoly->IsConvex();
712 for (iv = 0; iv <
fNvert; iv++) {
714 if (dist < stepmax) {
717 return (snext + dist);
723 return (snext + stepmax);
727 Int_t incseg = (dir[2] > 0) ? 1 : -1;
728 while (iz >= 0 && iz <
fNz - 1) {
733 for (iv = 0; iv <
fNvert; iv++) {
735 if (dist < stepmax) {
739 return (snext + dist);
745 return (snext + stepmax);
762 Error(
"DefinePolygon",
"In shape %s cannot create polygon with less than 3 vertices",
GetName());
766 for (
Int_t i = 0; i < nvert - 1; i++) {
767 for (
Int_t j = i + 1; j < nvert; j++) {
769 Error(
"DefinePolygon",
"In shape %s 2 vertices cannot be identical",
GetName());
795 if ((snum < 0) || (snum >=
fNz))
802 if (
fZ[snum] <
fZ[snum - 1]) {
804 "In shape: %s, Z position of section "
805 "%i, z=%e, not in increasing order, %i, z=%e",
806 GetName(), snum,
fZ[snum], snum - 1,
fZ[snum - 1]);
810 if (snum == (
fNz - 1)) {
822 if (ipl < 0 || ipl > (
fNz - 1)) {
823 Error(
"GetZ",
"In shape %s, ipl=%i out of range (0,%i)",
GetName(), ipl,
fNz - 1);
836 v1[0] = vert[9] - vert[0];
837 v1[1] = vert[10] - vert[1];
838 v1[2] = vert[11] - vert[2];
839 v2[0] = vert[3] - vert[0];
840 v2[1] = vert[4] - vert[1];
841 v2[2] = vert[5] - vert[2];
842 norm[0] =
v1[1] *
v2[2] -
v1[2] *
v2[1];
843 cross += norm[0] * norm[0];
844 norm[1] =
v1[2] *
v2[0] -
v1[0] *
v2[2];
845 cross += norm[1] * norm[1];
846 norm[2] =
v1[0] *
v2[1] -
v1[1] *
v2[0];
847 cross += norm[2] * norm[2];
851 for (
Int_t i = 0; i < 3; i++)
867 if (td.fPoly->IsClockwise()) {
919 for (
Int_t i = 0; i < 4; i++) {
921 k = 3 * ((i + 1) % 4);
922 v1[0] = point[0] - vert[j];
923 v1[1] = point[1] - vert[j + 1];
924 v1[2] = point[2] - vert[j + 2];
925 v2[0] = vert[k] - vert[j];
926 v2[1] = vert[k + 1] - vert[j + 1];
927 v2[2] = vert[k + 2] - vert[j + 2];
929 (
v1[0] *
v2[1] -
v1[1] *
v2[0]) * norm[2];
941 printf(
"*** Shape %s: TGeoXtru ***\n",
GetName());
942 printf(
" Nz = %i\n",
fNz);
943 printf(
" List of (x,y) of polygon vertices:\n");
945 printf(
" x = %11.5f y = %11.5f\n",
fX[ivert],
fY[ivert]);
946 for (
Int_t ipl = 0; ipl <
fNz; ipl++)
947 printf(
" plane %i: z=%11.5f x0=%11.5f y0=%11.5f scale=%11.5f\n", ipl,
fZ[ipl],
fX0[ipl],
fY0[ipl],
949 printf(
" Bounding box:\n");
961 Int_t nbPnts = nz * nvert;
962 Int_t nbSegs = nvert * (2 * nz - 1);
963 Int_t nbPols = nvert * (nz - 1) + 2;
966 6 * (nbPols - 2) + 2 * (2 + nvert));
985 Int_t indx = 0, indx2, k;
986 for (i = 0; i < nz; i++) {
990 for (j = 0; j < nvert; j++) {
993 buff.
fSegs[indx++] = indx2 + j;
994 buff.
fSegs[indx++] = indx2 + k;
997 for (i = 0; i < nz - 1; i++) {
1001 for (j = 0; j < nvert; j++) {
1004 buff.
fSegs[indx++] = indx2 + j;
1005 buff.
fSegs[indx++] = indx2 + k;
1012 for (i = 0; i < nz - 1; i++) {
1014 for (j = 0; j < nvert; j++) {
1015 k = (j + 1) % nvert;
1016 buff.
fPols[indx++] =
c + j % 3;
1017 buff.
fPols[indx++] = 4;
1018 buff.
fPols[indx++] = indx2 + j;
1019 buff.
fPols[indx++] = nz * nvert + indx2 + k;
1020 buff.
fPols[indx++] = indx2 + nvert + j;
1021 buff.
fPols[indx++] = nz * nvert + indx2 + j;
1024 buff.
fPols[indx++] =
c + 2;
1025 buff.
fPols[indx++] = nvert;
1027 for (j = nvert - 1; j >= 0; --j) {
1028 buff.
fPols[indx++] = indx2 + j;
1032 buff.
fPols[indx++] = nvert;
1033 indx2 = (nz - 1) * nvert;
1035 for (j = 0; j < nvert; j++) {
1036 buff.
fPols[indx++] = indx2 + j;
1055 saf1 = td.fPoly->Safety(point, iseg);
1056 in1 = td.fPoly->Contains(point);
1059 saf2 = td.fPoly->Safety(point, iseg);
1060 in2 = td.fPoly->Contains(point);
1061 if ((in1 & !in2) | (in2 & !in1)) {
1072 safz =
fZ[iz] - point[2];
1076 saf1 = point[2] -
fZ[iz + 1];
1091 for (iseg = 0; iseg <
fNvert; iseg++) {
1094 saf1 = (point[0] - vert[0]) * norm[0] + (point[1] - vert[1]) * norm[1] + (point[2] - vert[2]) * norm[2];
1127 for (iz = 0; iz <
fNz - 1; iz++) {
1140 safz =
fZ[0] - point[2];
1142 if (iz ==
fNz - 1) {
1144 safz = point[2] -
fZ[
fNz - 1];
1149 for (i = iz; i <
fNz - 1; i++) {
1155 for (i = iz - 1; i >= 0; i--) {
1175 out <<
" xvert[" << i <<
"] = " <<
fX[i] <<
"; yvert[" << i <<
"] = " <<
fY[i] <<
";" << std::endl;
1179 out <<
" " <<
GetPointerName() <<
"->DefineSection(" << i <<
", " <<
fZ[i] <<
", " <<
fX0[i] <<
", " <<
fY0[i]
1180 <<
", " <<
fScale[i] <<
");" << std::endl;
1194 a = (
fX0[ind1] *
fZ[ind2] -
fX0[ind2] *
fZ[ind1]) * invdz;
1195 b = (
fX0[ind2] -
fX0[ind1]) * invdz;
1197 a = (
fY0[ind1] *
fZ[ind2] -
fY0[ind2] *
fZ[ind1]) * invdz;
1198 b = (
fY0[ind2] -
fY0[ind1]) * invdz;
1213 td.fXc[i] = scale *
fX[i] + x0;
1214 td.fYc[i] = scale *
fY[i] + y0;
1235 Error(
"SetDimensions",
"Cannot create TGeoXtru %s with less than 2 Z planes",
GetName());
1253 DefineSection(i, param[1 + 4 * i], param[2 + 4 * i], param[3 + 4 * i], param[4 + 4 * i]);
1266 for (i = 0; i <
fNz; i++) {
1268 if (td.fPoly->IsClockwise()) {
1269 for (j = 0; j <
fNvert; j++) {
1270 points[indx++] = td.fXc[j];
1271 points[indx++] = td.fYc[j];
1275 for (j = 0; j <
fNvert; j++) {
1295 for (i = 0; i <
fNz; i++) {
1297 if (td.fPoly->IsClockwise()) {
1298 for (j = 0; j <
fNvert; j++) {
1299 points[indx++] = td.fXc[j];
1300 points[indx++] = td.fYc[j];
1304 for (j = 0; j <
fNvert; j++) {
1322 nsegs = nv * (2 * nz - 1);
1323 npols = nv * (nz - 1) + 2;
1352 Int_t nbPnts = nz * nvert;
1353 Int_t nbSegs = nvert * (2 * nz - 1);
1354 Int_t nbPols = nvert * (nz - 1) + 2;
1355 if (buffer.SetRawSizes(nbPnts, 3 * nbPnts, nbSegs, 3 * nbSegs, nbPols, 6 * (nbPols - 2) + 2 * (2 + nvert))) {
1362 if (!buffer.fLocalFrame) {
1380 for (
Int_t i = 0; i < vecsize; i++)
1391 for (
Int_t i = 0; i < vecsize; i++)
1401 for (
Int_t i = 0; i < vecsize; i++)
1411 for (
Int_t i = 0; i < vecsize; i++)
1422 for (
Int_t i = 0; i < vecsize; i++)
void Error(const char *location, const char *msgfmt,...)
Use this function in case an error occurred.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t points
Generic 3D primitive description class.
void FillBuffer3D(TBuffer3D &buffer, Int_t reqSections, Bool_t localFrame) const override
Fill the supplied buffer, with sections in desired frame See TBuffer3D.h for explanation of sections,...
Double_t DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=nullptr) const override
void InspectShape() const override
Bool_t Contains(const Double_t *point) const override
Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const override
static Int_t ThreadId()
Translates the current thread id to an ordinal number.
An arbitrary polygon defined by vertices.
Int_t GetBasicColor() const
Get the basic color (0-7).
void TransformPoints(Double_t *points, UInt_t NbPoints) const
Tranform a set of points (LocalToMaster)
void SetShapeBit(UInt_t f, Bool_t set)
Equivalent of TObject::SetBit.
static Bool_t IsSameWithinTolerance(Double_t a, Double_t b)
Check if two numbers differ with less than a tolerance.
const char * GetPointerName() const
Provide a pointer name containing uid.
Int_t ShapeDistancetoPrimitive(Int_t numpoints, Int_t px, Int_t py) const
Returns distance to shape primitive mesh.
const char * GetName() const override
Get the shape name.
static Double_t Tolerance()
Bool_t TestShapeBit(UInt_t f) const
ThreadData_t & GetThreadData() const
Int_t GetNmeshVertices() const override
virtual void DefineSection(Int_t snum, Double_t z, Double_t x0=0., Double_t y0=0., Double_t scale=1.)
void DistFromInside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t *step) const override
void ClearThreadData() const override
Double_t Capacity() const override
void ComputeBBox() override
void Contains_v(const Double_t *points, Bool_t *inside, Int_t vecsize) const override
Bool_t Contains(const Double_t *point) const override
Int_t fThreadSize
Navigation data per thread.
void SetCurrentVertices(Double_t x0, Double_t y0, Double_t scale)
void ComputeNormal(const Double_t *point, const Double_t *dir, Double_t *norm) override
void SavePrimitive(std::ostream &out, Option_t *option="") override
Save a primitive as a C++ statement(s) on output stream "out".
void GetPlaneVertices(Int_t iz, Int_t ivert, Double_t *vert) const
void InspectShape() const override
Int_t DistancetoPrimitive(Int_t px, Int_t py) override
Computes distance from point (px,py) to the object.
Bool_t IsPointInsidePlane(const Double_t *point, Double_t *vert, Double_t *norm) const
void Sizeof3D() const override
const TBuffer3D & GetBuffer3D(Int_t reqSections, Bool_t localFrame) const override
Stub implementation to avoid forcing implementation at this stage.
Bool_t DefinePolygon(Int_t nvert, const Double_t *xv, const Double_t *yv)
void SetCurrentZ(Double_t z, Int_t iz)
Double_t Safety(const Double_t *point, Bool_t in=kTRUE) const override
void ComputeNormal_v(const Double_t *points, const Double_t *dirs, Double_t *norms, Int_t vecsize) override
void SetDimensions(Double_t *param) override
void DrawPolygon(Option_t *option="")
void SetPoints(Double_t *points) const override
Double_t SafetyToSector(const Double_t *point, Int_t iz, Double_t safmin, Bool_t in)
Double_t DistFromInside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=nullptr) const override
void DistFromOutside_v(const Double_t *points, const Double_t *dirs, Double_t *dists, Int_t vecsize, Double_t *step) const override
TBuffer3D * MakeBuffer3D() const override
std::vector< ThreadData_t * > fThreadData
Double_t DistToPlane(const Double_t *point, const Double_t *dir, Int_t iz, Int_t ivert, Double_t stepmax, Bool_t in) const
void SetSegsAndPols(TBuffer3D &buff) const override
Double_t DistFromOutside(const Double_t *point, const Double_t *dir, Int_t iact=1, Double_t step=TGeoShape::Big(), Double_t *safe=nullptr) const override
std::mutex fMutex
size of thread-specific array
void GetMeshNumbers(Int_t &nvert, Int_t &nsegs, Int_t &npols) const override
void GetPlaneNormal(const Double_t *vert, Double_t *norm) const
void Safety_v(const Double_t *points, const Bool_t *inside, Double_t *safe, Int_t vecsize) const override
void CreateThreadData(Int_t nthreads) override
R__ALWAYS_INLINE Bool_t TestBit(UInt_t f) const
virtual const char * ClassName() const
Returns name of class to which the object belongs.
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
void SetBit(UInt_t f, Bool_t set)
Set or unset the user status bits as specified in f.
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
double dist(Rotation3D const &r1, Rotation3D const &r2)
Short_t Max(Short_t a, Short_t b)
Returns the largest of a and b.
constexpr Double_t E()
Base of natural log: .
Double_t Sqrt(Double_t x)
Returns the square root of x.
Short_t Min(Short_t a, Short_t b)
Returns the smallest of a and b.
Long64_t BinarySearch(Long64_t n, const T *array, T value)
Binary search in an array of n values to locate value.
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
BVH_ALWAYS_INLINE Vec< T, 3 > cross(const Vec< T, 3 > &a, const Vec< T, 3 > &b)