Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
TGeoTessellated.cxx
Go to the documentation of this file.
1// @(#)root/geom:$Id$// Author: Andrei Gheata 24/10/01
2
3// Contains() and DistFromOutside/Out() implemented by Mihaela Gheata
4
5/*************************************************************************
6 * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
7 * All rights reserved. *
8 * *
9 * For the licensing terms see $ROOTSYS/LICENSE. *
10 * For the list of contributors see $ROOTSYS/README/CREDITS. *
11 *************************************************************************/
12
13/** \class TGeoTessellated
14\ingroup Geometry_classes
15
16Tessellated solid class. It is composed by a set of planar faces having triangular or
17quadrilateral shape. The class does not provide navigation functionality, it just wraps the data
18for the composing faces.
19*/
20
21#include <iostream>
22#include <sstream>
23
24#include "TGeoManager.h"
25#include "TGeoMatrix.h"
26#include "TGeoVolume.h"
27#include "TVirtualGeoPainter.h"
28#include "TGeoTessellated.h"
29#include "TBuffer3D.h"
30#include "TBuffer3DTypes.h"
31#include "TMath.h"
32
33#include <array>
34#include <vector>
35
36
38
39////////////////////////////////////////////////////////////////////////////////
40/// Compact consecutive equal vertices
41
43{
44 // Compact the common vertices and return new facet
45 if (nvertices < 2)
46 return nvertices;
47 int nvert = nvertices;
48 int i = 0;
49 while (i < nvert) {
50 if (vert[(i + 1) % nvert] == vert[i]) {
51 // shift last vertices left by one element
52 for (int j = i + 2; j < nvert; ++j)
53 vert[j - 1] = vert[j];
54 nvert--;
55 }
56 i++;
57 }
58 return nvert;
59}
60
61////////////////////////////////////////////////////////////////////////////////
62/// Check if a connected neighbour facet has compatible normal
63
64bool TGeoFacet::IsNeighbour(const TGeoFacet &other, bool &flip) const
65{
66
67 // Find a connecting segment
68 bool neighbour = false;
69 int line1[2], line2[2];
70 int npoints = 0;
71 for (int i = 0; i < fNvert; ++i) {
72 auto ivert = fIvert[i];
73 // Check if the other facet has the same vertex
74 for (int j = 0; j < other.GetNvert(); ++j) {
75 if (ivert == other[j]) {
76 line1[npoints] = i;
77 line2[npoints] = j;
78 if (++npoints == 2) {
79 neighbour = true;
80 bool order1 = line1[1] == line1[0] + 1;
81 bool order2 = line2[1] == (line2[0] + 1) % other.GetNvert();
82 flip = (order1 == order2);
83 return neighbour;
84 }
85 }
86 }
87 }
88 return neighbour;
89}
90
91////////////////////////////////////////////////////////////////////////////////
92/// Constructor. In case nfacets is zero, it is user's responsibility to
93/// call CloseShape once all faces are defined.
94
96{
98 if (nfacets)
99 fFacets.reserve(nfacets);
100}
101
102////////////////////////////////////////////////////////////////////////////////
103/// Constructor providing directly the array of vertices. Facets have to be added
104/// providing vertex indices rather than coordinates.
105
106TGeoTessellated::TGeoTessellated(const char *name, const std::vector<Vertex_t> &vertices) : TGeoBBox(name, 0, 0, 0)
107{
108 fVertices = vertices;
109 fNvert = fVertices.size();
110}
111
112////////////////////////////////////////////////////////////////////////////////
113/// Add a vertex checking for duplicates, returning the vertex index
114
116{
117 constexpr double tolerance = 1.e-10;
118 auto vertexHash = [&](Vertex_t const &vertex) {
119 // Compute hash for the vertex
120 long hash = 0;
121 // helper function to generate hash from integer numbers
122 auto hash_combine = [](long seed, const long value) {
123 return seed ^ (std::hash<long>{}(value) + 0x9e3779b9 + (seed << 6) + (seed >> 2));
124 };
125 for (int i = 0; i < 3; i++) {
126 // use tolerance to generate int with the desired precision from a real number for hashing
127 hash = hash_combine(hash, std::roundl(vertex[i] / tolerance));
128 }
129 return hash;
130 };
131
132 auto hash = vertexHash(vert);
133 bool isAdded = false;
134 int ivert = -1;
135 // Get the compatible vertices
136 auto range = fVerticesMap.equal_range(hash);
137 for (auto it = range.first; it != range.second; ++it) {
138 ivert = it->second;
139 if (fVertices[ivert] == vert) {
140 isAdded = true;
141 break;
142 }
143 }
144 if (!isAdded) {
145 ivert = fVertices.size();
146 fVertices.push_back(vert);
147 fVerticesMap.insert(std::make_pair(hash, ivert));
148 }
149 return ivert;
150}
151
152////////////////////////////////////////////////////////////////////////////////
153/// Adding a triangular facet from vertex positions in absolute coordinates
154
156{
157 if (fDefined) {
158 Error("AddFacet", "Shape %s already fully defined. Not adding", GetName());
159 return false;
160 }
161
162 Vertex_t vert[3];
163 vert[0] = pt0;
164 vert[1] = pt1;
165 vert[2] = pt2;
167 if (nvert < 3) {
168 Error("AddFacet", "Triangular facet at index %d degenerated. Not adding.", GetNfacets());
169 return false;
170 }
171 int ind[3];
172 for (auto i = 0; i < 3; ++i)
173 ind[i] = AddVertex(vert[i]);
174 fNseg += 3;
175 fFacets.emplace_back(ind[0], ind[1], ind[2]);
176
177 if (fNfacets > 0 && GetNfacets() == fNfacets)
178 CloseShape();
179 return true;
180}
181
182////////////////////////////////////////////////////////////////////////////////
183/// Adding a triangular facet from indices of vertices
184
186{
187 if (fDefined) {
188 Error("AddFacet", "Shape %s already fully defined. Not adding", GetName());
189 return false;
190 }
191 if (fVertices.empty()) {
192 Error("AddFacet", "Shape %s Cannot add facets by indices without vertices. Not adding", GetName());
193 return false;
194 }
195
196 fNseg += 3;
197 fFacets.emplace_back(i0, i1, i2);
198 return true;
199}
200
201////////////////////////////////////////////////////////////////////////////////
202/// Adding a quadrilateral facet from vertex positions in absolute coordinates
203
205{
206 if (fDefined) {
207 Error("AddFacet", "Shape %s already fully defined. Not adding", GetName());
208 return false;
209 }
210 Vertex_t vert[4];
211 vert[0] = pt0;
212 vert[1] = pt1;
213 vert[2] = pt2;
214 vert[3] = pt3;
216 if (nvert < 3) {
217 Error("AddFacet", "Quadrilateral facet at index %d degenerated. Not adding.", GetNfacets());
218 return false;
219 }
220
221 int ind[4];
222 for (auto i = 0; i < nvert; ++i)
223 ind[i] = AddVertex(vert[i]);
224 fNseg += nvert;
225 if (nvert == 3)
226 fFacets.emplace_back(ind[0], ind[1], ind[2]);
227 else
228 fFacets.emplace_back(ind[0], ind[1], ind[2], ind[3]);
229
230 if (fNfacets > 0 && GetNfacets() == fNfacets)
231 CloseShape(false);
232 return true;
233}
234
235////////////////////////////////////////////////////////////////////////////////
236/// Adding a quadrilateral facet from indices of vertices
237
238bool TGeoTessellated::AddFacet(int i0, int i1, int i2, int i3)
239{
240 if (fDefined) {
241 Error("AddFacet", "Shape %s already fully defined. Not adding", GetName());
242 return false;
243 }
244 if (fVertices.empty()) {
245 Error("AddFacet", "Shape %s Cannot add facets by indices without vertices. Not adding", GetName());
246 return false;
247 }
248
249 fNseg += 4;
250 fFacets.emplace_back(i0, i1, i2, i3);
251 return true;
252}
253
254////////////////////////////////////////////////////////////////////////////////
255/// Compute normal for a given facet
256
258{
259 // Compute normal using non-zero segments
260 constexpr double kTolerance = 1.e-20;
261 auto const &facet = fFacets[ifacet];
262 int nvert = facet.GetNvert();
263 degenerated = true;
265 for (int i = 0; i < nvert - 1; ++i) {
266 Vertex_t e1 = fVertices[facet[i + 1]] - fVertices[facet[i]];
267 if (e1.Mag2() < kTolerance)
268 continue;
269 for (int j = i + 1; j < nvert; ++j) {
270 Vertex_t e2 = fVertices[facet[(j + 1) % nvert]] - fVertices[facet[j]];
271 if (e2.Mag2() < kTolerance)
272 continue;
274 // e1 and e2 may be colinear
275 if (normal.Mag2() < kTolerance)
276 continue;
277 normal.Normalize();
278 degenerated = false;
279 break;
280 }
281 if (!degenerated)
282 break;
283 }
284 return normal;
285}
286
287////////////////////////////////////////////////////////////////////////////////
288/// Check validity of facet
289
291{
292 constexpr double kTolerance = 1.e-10;
293 auto const &facet = fFacets[ifacet];
294 int nvert = facet.GetNvert();
295 bool degenerated = true;
297 if (degenerated) {
298 std::cout << "Facet: " << ifacet << " is degenerated\n";
299 return false;
300 }
301
302 // Compute surface area
303 double surfaceArea = 0.;
304 for (int i = 1; i < nvert - 1; ++i) {
306 Vertex_t e2 = fVertices[facet[i + 1]] - fVertices[facet[0]];
307 surfaceArea += 0.5 * Vertex_t::Cross(e1, e2).Mag();
308 }
309 if (surfaceArea < kTolerance) {
310 std::cout << "Facet: " << ifacet << " has zero surface area\n";
311 return false;
312 }
313
314 return true;
315}
316
317////////////////////////////////////////////////////////////////////////////////
318/// Close the shape: calculate bounding box and compact vertices
319
320void TGeoTessellated::CloseShape(bool check, bool fixFlipped, bool verbose)
321{
322 // Compute bounding box
323 fDefined = true;
324 fNvert = fVertices.size();
325 fNfacets = fFacets.size();
326 ComputeBBox();
327 // Cleanup the vertex map
328 std::multimap<long, int>().swap(fVerticesMap);
329
330 if (fVertices.size() > 0) {
331 if (!check)
332 return;
333
334 // Check facets
335 for (auto i = 0; i < fNfacets; ++i)
336 FacetCheck(i);
337
339 }
340}
341
342////////////////////////////////////////////////////////////////////////////////
343/// Check closure of the solid and check/fix flipped normals
344
346{
347 int *nn = new int[fNfacets];
348 bool *flipped = new bool[fNfacets];
349 bool hasorphans = false;
350 bool hasflipped = false;
351 for (int i = 0; i < fNfacets; ++i) {
352 nn[i] = 0;
353 flipped[i] = false;
354 }
355
356 for (int icrt = 0; icrt < fNfacets; ++icrt) {
357 // all neighbours checked?
358 if (nn[icrt] >= fFacets[icrt].GetNvert())
359 continue;
360 for (int i = icrt + 1; i < fNfacets; ++i) {
361 bool isneighbour = fFacets[icrt].IsNeighbour(fFacets[i], flipped[i]);
362 if (isneighbour) {
363 if (flipped[icrt])
364 flipped[i] = !flipped[i];
365 if (flipped[i])
366 hasflipped = true;
367 nn[icrt]++;
368 nn[i]++;
369 if (nn[icrt] == fFacets[icrt].GetNvert())
370 break;
371 }
372 }
373 if (nn[icrt] < fFacets[icrt].GetNvert())
374 hasorphans = true;
375 }
376
377 if (hasorphans && verbose) {
378 Error("Check", "Tessellated solid %s has following not fully connected facets:", GetName());
379 for (int icrt = 0; icrt < fNfacets; ++icrt) {
380 if (nn[icrt] < fFacets[icrt].GetNvert())
381 std::cout << icrt << " (" << fFacets[icrt].GetNvert() << " edges, " << nn[icrt] << " neighbours)\n";
382 }
383 }
385 int nfixed = 0;
386 if (hasflipped) {
387 if (verbose)
388 Warning("Check", "Tessellated solid %s has following facets with flipped normals:", GetName());
389 for (int icrt = 0; icrt < fNfacets; ++icrt) {
390 if (flipped[icrt]) {
391 if (verbose)
392 std::cout << icrt << "\n";
393 if (fixFlipped) {
394 fFacets[icrt].Flip();
395 nfixed++;
396 }
397 }
398 }
399 if (nfixed && verbose)
400 Info("Check", "Automatically flipped %d facets to match first defined facet", nfixed);
401 }
402 delete[] nn;
403 delete[] flipped;
404
405 return !hasorphans;
406}
407
408////////////////////////////////////////////////////////////////////////////////
409/// Compute bounding box
410
412{
413 const double kBig = TGeoShape::Big();
414 double vmin[3] = {kBig, kBig, kBig};
415 double vmax[3] = {-kBig, -kBig, -kBig};
416 for (const auto &facet : fFacets) {
417 for (int i = 0; i < facet.GetNvert(); ++i) {
418 for (int j = 0; j < 3; ++j) {
419 vmin[j] = TMath::Min(vmin[j], fVertices[facet[i]].operator[](j));
420 vmax[j] = TMath::Max(vmax[j], fVertices[facet[i]].operator[](j));
421 }
422 }
423 }
424 fDX = 0.5 * (vmax[0] - vmin[0]);
425 fDY = 0.5 * (vmax[1] - vmin[1]);
426 fDZ = 0.5 * (vmax[2] - vmin[2]);
427 for (int i = 0; i < 3; ++i)
428 fOrigin[i] = 0.5 * (vmax[i] + vmin[i]);
429}
430
431////////////////////////////////////////////////////////////////////////////////
432/// Returns numbers of vertices, segments and polygons composing the shape mesh.
433
435{
436 nvert = fNvert;
437 nsegs = fNseg;
438 npols = GetNfacets();
439}
440
441////////////////////////////////////////////////////////////////////////////////
442/// Creates a TBuffer3D describing *this* shape.
443/// Coordinates are in local reference frame.
444
446{
447 const int nvert = fNvert;
448 const int nsegs = fNseg;
449 const int npols = GetNfacets();
451 if (buff) {
452 SetPoints(buff->fPnts);
454 }
455 return buff;
456}
457
458////////////////////////////////////////////////////////////////////////////////
459/// Prints basic info
460
462{
463 std::cout << "=== Tessellated shape " << GetName() << " having " << GetNvertices() << " vertices and "
464 << GetNfacets() << " facets\n";
465}
466
467////////////////////////////////////////////////////////////////////////////////
468/// Fills TBuffer3D structure for segments and polygons.
469
471{
472 const int c = GetBasicColor();
473 int *segs = buff.fSegs;
474 int *pols = buff.fPols;
475
476 int indseg = 0; // segment internal data index
477 int indpol = 0; // polygon internal data index
478 int sind = 0; // segment index
479 for (const auto &facet : fFacets) {
480 auto nvert = facet.GetNvert();
481 pols[indpol++] = c;
482 pols[indpol++] = nvert;
483 for (auto j = 0; j < nvert; ++j) {
484 int k = (j + 1) % nvert;
485 // segment made by next consecutive points
486 segs[indseg++] = c;
487 segs[indseg++] = facet[j];
488 segs[indseg++] = facet[k];
489 // add segment to current polygon and increment segment index
490 pols[indpol + nvert - j - 1] = sind++;
491 }
492 indpol += nvert;
493 }
494}
495
496////////////////////////////////////////////////////////////////////////////////
497/// Fill tessellated points to an array.
498
500{
501 int ind = 0;
502 for (const auto &vertex : fVertices) {
503 vertex.CopyTo(&points[ind]);
504 ind += 3;
505 }
506}
507
508////////////////////////////////////////////////////////////////////////////////
509/// Fill tessellated points in float.
510
512{
513 int ind = 0;
514 for (const auto &vertex : fVertices) {
515 points[ind++] = vertex.x();
516 points[ind++] = vertex.y();
517 points[ind++] = vertex.z();
518 }
519}
520
521////////////////////////////////////////////////////////////////////////////////
522/// Resize the shape by scaling vertices within maxsize and center to origin
523
525{
526 using Vector3_t = Vertex_t;
527
528 if (!fDefined) {
529 Error("ResizeCenter", "Not all faces are defined");
530 return;
531 }
533 double maxedge = TMath::Max(TMath::Max(fDX, fDY), fDZ);
534 double scale = maxsize / maxedge;
535 for (size_t i = 0; i < fVertices.size(); ++i) {
536 fVertices[i] = scale * (fVertices[i] - origin);
537 }
538 fOrigin[0] = fOrigin[1] = fOrigin[2] = 0;
539 fDX *= scale;
540 fDY *= scale;
541 fDZ *= scale;
542}
543
544////////////////////////////////////////////////////////////////////////////////
545/// Fills a static 3D buffer and returns a reference.
546
548{
549 static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
550
552
553 const int nvert = fNvert;
554 const int nsegs = fNseg;
555 const int npols = GetNfacets();
556
558 if (buffer.SetRawSizes(nvert, 3 * nvert, nsegs, 3 * nsegs, npols, 6 * npols)) {
560 }
561 }
563 SetPoints(buffer.fPnts);
564 if (!buffer.fLocalFrame) {
565 TransformPoints(buffer.fPnts, buffer.NbPnts());
566 }
567
568 SetSegsAndPols(buffer);
570 }
571
572 return buffer;
573}
574
575////////////////////////////////////////////////////////////////////////////////
576/// Reads a single tessellated solid from an .obj file.
577
579{
580 using std::vector, std::string, std::ifstream, std::stringstream, std::endl;
581
582 vector<Vertex_t> vertices;
584
585 struct FacetInd_t {
586 int i0 = -1;
587 int i1 = -1;
588 int i2 = -1;
589 int i3 = -1;
590 int nvert = 0;
591 FacetInd_t(int a, int b, int c)
592 {
593 i0 = a;
594 i1 = b;
595 i2 = c;
596 nvert = 3;
597 };
598 FacetInd_t(int a, int b, int c, int d)
599 {
600 i0 = a;
601 i1 = b;
602 i2 = c;
603 i3 = d;
604 nvert = 4;
605 };
606 };
607
609 // List of geometric vertices, with (x, y, z [,w]) coordinates, w is optional and defaults to 1.0.
610 // struct vtx_t { double x = 0; double y = 0; double z = 0; double w = 1; };
611
612 // Texture coordinates in u, [,v ,w]) coordinates, these will vary between 0 and 1. v, w are optional and default to
613 // 0.
614 // struct tex_t { double u; double v; double w; };
615
616 // List of vertex normals in (x,y,z) form; normals might not be unit vectors.
617 // struct vn_t { double x; double y; double z; };
618
619 // Parameter space vertices in ( u [,v] [,w] ) form; free form geometry statement
620 // struct vp_t { double u; double v; double w; };
621
622 // Faces are defined using lists of vertex, texture and normal indices which start at 1.
623 // Polygons such as quadrilaterals can be defined by using more than three vertex/texture/normal indices.
624 // f v1//vn1 v2//vn2 v3//vn3 ...
625
626 // Records starting with the letter "l" specify the order of the vertices which build a polyline.
627 // l v1 v2 v3 v4 v5 v6 ...
628
629 string line;
630 int ind[4] = {0};
631 ifstream file(objfile);
632 if (!file.is_open()) {
633 ::Error("TGeoTessellated::ImportFromObjFormat", "Unable to open %s", objfile);
634 return nullptr;
635 }
636
637 while (getline(file, line)) {
638 stringstream ss(line);
639 string tag;
640
641 // We ignore everything which is not a vertex or a face
642 if (line.rfind('v', 0) == 0 && line.rfind("vt", 0) != 0 && line.rfind("vn", 0) != 0 && line.rfind("vn", 0) != 0) {
643 // Decode the vertex
644 double pos[4] = {0, 0, 0, 1};
645 ss >> tag >> pos[0] >> pos[1] >> pos[2] >> pos[3];
646 vertices.emplace_back(pos[0] * pos[3], pos[1] * pos[3], pos[2] * pos[3]);
647 }
648
649 else if (line.rfind('f', 0) == 0) {
650 // Decode the face
651 ss >> tag;
652 string word;
653 sfacets.clear();
654 while (ss >> word)
655 sfacets.push_back(word);
656 if (sfacets.size() > 4 || sfacets.size() < 3) {
657 ::Error("TGeoTessellated::ImportFromObjFormat", "Detected face having unsupported %zu vertices",
658 sfacets.size());
659 return nullptr;
660 }
661 int nvert = 0;
662 for (auto &sword : sfacets) {
663 stringstream ssword(sword);
664 string token;
665 getline(ssword, token, '/'); // just need the vertex index, which is the first token
666 // Convert string token to integer
667
668 ind[nvert++] = stoi(token) - 1;
669 if (ind[nvert - 1] < 0) {
670 ::Error("TGeoTessellated::ImportFromObjFormat", "Unsupported relative vertex index definition in %s",
671 objfile);
672 return nullptr;
673 }
674 }
675 if (nvert == 3)
676 facets.emplace_back(ind[0], ind[1], ind[2]);
677 else
678 facets.emplace_back(ind[0], ind[1], ind[2], ind[3]);
679 }
680 }
681
682 int nvertices = (int)vertices.size();
683 int nfacets = (int)facets.size();
684 if (nfacets < 3) {
685 ::Error("TGeoTessellated::ImportFromObjFormat", "Not enough faces detected in %s", objfile);
686 return nullptr;
687 }
688
689 string sobjfile(objfile);
690 if (verbose)
691 std::cout << "Read " << nvertices << " vertices and " << nfacets << " facets from " << sobjfile << endl;
692
693 auto tsl = new TGeoTessellated(sobjfile.erase(sobjfile.find_last_of('.')).c_str(), vertices);
694
695 for (int i = 0; i < nfacets; ++i) {
696 auto facet = facets[i];
697 if (facet.nvert == 3)
698 tsl->AddFacet(facet.i0, facet.i1, facet.i2);
699 else
700 tsl->AddFacet(facet.i0, facet.i1, facet.i2, facet.i3);
701 }
702 tsl->CloseShape(check, true, verbose);
703 tsl->Print();
704 return tsl;
705}
#define d(i)
Definition RSha256.hxx:102
#define b(i)
Definition RSha256.hxx:100
#define c(i)
Definition RSha256.hxx:101
#define a(i)
Definition RSha256.hxx:99
float Float_t
Float 4 bytes (float)
Definition RtypesCore.h:71
const char Option_t
Option string (const char)
Definition RtypesCore.h:80
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void value
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t points
char name[80]
Definition TGX11.cxx:110
segment * segs
Definition X3DBuffer.c:23
Generic 3D primitive description class.
Definition TBuffer3D.h:18
UInt_t NbPnts() const
Definition TBuffer3D.h:80
Bool_t SectionsValid(UInt_t mask) const
Definition TBuffer3D.h:67
void SetSectionsValid(UInt_t mask)
Definition TBuffer3D.h:65
Bool_t fLocalFrame
Definition TBuffer3D.h:90
Bool_t SetRawSizes(UInt_t reqPnts, UInt_t reqPntsCapacity, UInt_t reqSegs, UInt_t reqSegsCapacity, UInt_t reqPols, UInt_t reqPolsCapacity)
Set kRaw tessellation section of buffer with supplied sizes.
Double_t * fPnts
Definition TBuffer3D.h:113
Box class.
Definition TGeoBBox.h:17
void FillBuffer3D(TBuffer3D &buffer, Int_t reqSections, Bool_t localFrame) const override
Fills the supplied buffer, with sections in desired frame See TBuffer3D.h for explanation of sections...
Double_t fDX
Definition TGeoBBox.h:20
Double_t fOrigin[3]
Definition TGeoBBox.h:23
Double_t fDY
Definition TGeoBBox.h:21
Double_t fDZ
Definition TGeoBBox.h:22
bool IsNeighbour(const TGeoFacet &other, bool &flip) const
Check if a connected neighbour facet has compatible normal.
static int CompactFacet(Vertex_t *vert, int nvertices)
Compact consecutive equal vertices.
static Double_t Big()
Definition TGeoShape.h:94
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)
const char * GetName() const override
Get the shape name.
Tessellated solid class.
void ResizeCenter(double maxsize)
Resize and center the shape in a box of size maxsize.
int AddVertex(const Vertex_t &vert)
Add a vertex checking for duplicates, returning the vertex index.
bool FacetCheck(int ifacet) const
Check validity of facet.
void Print(Option_t *option="") const override
Prints basic info.
void SetSegsAndPols(TBuffer3D &buff) const override
Fills TBuffer3D structure for segments and polygons.
const TBuffer3D & GetBuffer3D(int reqSections, Bool_t localFrame) const override
Fills a static 3D buffer and returns a reference.
void SetPoints(double *points) const override
Fill tessellated points to an array.
bool CheckClosure(bool fixFlipped=true, bool verbose=true)
Check closure of the solid and check/fix flipped normals.
int GetNvertices() const
Vertex_t FacetComputeNormal(int ifacet, bool &degenerated) const
Compute normal for a given facet.
void CloseShape(bool check=true, bool fixFlipped=true, bool verbose=true)
Close the shape: calculate bounding box and compact vertices.
Tessellated::Vertex_t Vertex_t
void GetMeshNumbers(int &nvert, int &nsegs, int &npols) const override
Returns numbers of vertices, segments and polygons composing the shape mesh.
std::vector< TGeoFacet > fFacets
static TGeoTessellated * ImportFromObjFormat(const char *objfile, bool check=false, bool verbose=false)
Reader from .obj format.
TBuffer3D * MakeBuffer3D() const override
Creates a TBuffer3D describing this shape.
void ComputeBBox() override
Compute bounding box.
std::multimap< long, int > fVerticesMap
bool AddFacet(const Vertex_t &pt0, const Vertex_t &pt1, const Vertex_t &pt2)
Adding a triangular facet from vertex positions in absolute coordinates.
std::vector< Vertex_t > fVertices
int GetNfacets() const
bool fClosedBody
Shape fully defined.
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition TObject.cxx:1057
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:1071
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition TObject.cxx:1045
TLine * line
Short_t Max(Short_t a, Short_t b)
Returns the largest of a and b.
Definition TMathBase.h:251
Short_t Min(Short_t a, Short_t b)
Returns the smallest of a and b.
Definition TMathBase.h:199
ROOT::Geom::Vertex_t Vertex_t
static Vertex_t Cross(Vertex_t const &left, Vertex_t const &right)
The cross (vector) product of two Vector3D<T> objects.