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
37
39
40std::ostream &operator<<(std::ostream &os, TGeoFacet const &facet)
41{
42 os << "{";
43 for (int i = 0; i < facet.GetNvert(); ++i) {
44 os << facet.GetVertex(i);
45 if (i != facet.GetNvert() - 1)
46 os << ", ";
47 }
48 os << "}";
49 return os;
50}
51
52TGeoFacet::TGeoFacet(const TGeoFacet &other) : fVertices(other.fVertices), fNvert(other.fNvert), fShared(other.fShared)
53{
54 memcpy(fIvert, other.fIvert, 4 * sizeof(int));
55 if (!fShared && other.fVertices)
56 fVertices = new VertexVec_t(*other.fVertices);
57}
58
60{
61 if (&other != this) {
62 if (!fShared)
63 delete fVertices;
64 fNvert = other.fNvert;
65 fShared = other.fShared;
66 memcpy(fIvert, other.fIvert, 4 * sizeof(int));
67 if (!fShared && other.fVertices)
68 fVertices = new VertexVec_t(*other.fVertices);
69 else
70 fVertices = other.fVertices;
71 }
72 return *this;
73}
74
75Vertex_t TGeoFacet::ComputeNormal(bool &degenerated) const
76{
77 // Compute normal using non-zero segments
78 constexpr double kTolerance = 1.e-20;
79 degenerated = true;
80 Vertex_t normal;
81 for (int i = 0; i < fNvert - 1; ++i) {
82 Vertex_t e1 = GetVertex(i + 1) - GetVertex(i);
83 if (e1.Mag2() < kTolerance)
84 continue;
85 for (int j = i + 1; j < fNvert; ++j) {
86 Vertex_t e2 = GetVertex((j + 1) % fNvert) - GetVertex(j);
87 if (e2.Mag2() < kTolerance)
88 continue;
89 normal = Vertex_t::Cross(e1, e2);
90 // e1 and e2 may be colinear
91 if (normal.Mag2() < kTolerance)
92 continue;
93 normal.Normalize();
94 degenerated = false;
95 break;
96 }
97 if (!degenerated)
98 break;
99 }
100 return normal;
101}
102
104{
105 constexpr double kTolerance = 1.e-10;
106 bool degenerated = true;
107 ComputeNormal(degenerated);
108 if (degenerated) {
109 std::cout << "Facet: " << *this << " is degenerated\n";
110 return false;
111 }
112
113 // Compute surface area
114 double surfaceArea = 0.;
115 for (int i = 1; i < fNvert - 1; ++i) {
116 Vertex_t e1 = GetVertex(i) - GetVertex(0);
117 Vertex_t e2 = GetVertex(i + 1) - GetVertex(0);
118 surfaceArea += 0.5 * Vertex_t::Cross(e1, e2).Mag();
119 }
120 if (surfaceArea < kTolerance) {
121 std::cout << "Facet: " << *this << " has zero surface area\n";
122 return false;
123 }
124
125 // Center of the tile
126 /*
127 Vertex_t center;
128 for (int i = 0; i < fNvert; ++i)
129 center += GetVertex(i);
130 center /= fNvert;
131 */
132 return true;
133}
134
135int TGeoFacet::CompactFacet(Vertex_t *vert, int nvertices)
136{
137 // Compact the common vertices and return new facet
138 if (nvertices < 2)
139 return nvertices;
140 int nvert = nvertices;
141 int i = 0;
142 while (i < nvert) {
143 if (vert[(i + 1) % nvert] == vert[i]) {
144 // shift last vertices left by one element
145 for (int j = i + 2; j < nvert; ++j)
146 vert[j - 1] = vert[j];
147 nvert--;
148 }
149 i++;
150 }
151 return nvert;
152}
153
154////////////////////////////////////////////////////////////////////////////////
155/// Check if a connected neighbour facet has compatible normal
156
157bool TGeoFacet::IsNeighbour(const TGeoFacet &other, bool &flip) const
158{
159
160 // Find a connecting segment
161 bool neighbour = false;
162 int line1[2], line2[2];
163 int npoints = 0;
164 for (int i = 0; i < fNvert; ++i) {
165 auto ivert = GetVertexIndex(i);
166 // Check if the other facet has the same vertex
167 for (int j = 0; j < other.GetNvert(); ++j) {
168 if (ivert == other.GetVertexIndex(j)) {
169 line1[npoints] = i;
170 line2[npoints] = j;
171 if (++npoints == 2) {
172 neighbour = true;
173 bool order1 = line1[1] == line1[0] + 1;
174 bool order2 = line2[1] == (line2[0] + 1) % other.GetNvert();
175 flip = (order1 == order2);
176 return neighbour;
177 }
178 }
179 }
180 }
181 return neighbour;
182}
183
184////////////////////////////////////////////////////////////////////////////////
185/// Constructor. In case nfacets is zero, it is user's responsibility to
186/// call CloseShape once all faces are defined.
187
188TGeoTessellated::TGeoTessellated(const char *name, int nfacets) : TGeoBBox(name, 0, 0, 0)
189{
190 fNfacets = nfacets;
191 if (nfacets)
192 fFacets.reserve(nfacets);
193}
194
195////////////////////////////////////////////////////////////////////////////////
196/// Constructor providing directly the array of vertices. Facets have to be added
197/// providing vertex indices rather than coordinates.
198
199TGeoTessellated::TGeoTessellated(const char *name, const std::vector<Vertex_t> &vertices) : TGeoBBox(name, 0, 0, 0)
200{
201 fVertices = vertices;
202 fNvert = fVertices.size();
203}
204
205////////////////////////////////////////////////////////////////////////////////
206/// Function to be called after reading tessellated volumes from the geometry file
207
209{
210 // The pointer to the array of vertices is not streamed so update it to facets
211 for (auto facet : fFacets)
212 facet.SetVertices(&fVertices);
213 fDefined = true;
214}
215
216////////////////////////////////////////////////////////////////////////////////
217/// Adding a triangular facet from vertex positions in absolute coordinates
218
219bool TGeoTessellated::AddFacet(const Vertex_t &pt0, const Vertex_t &pt1, const Vertex_t &pt2)
220{
221 if (fDefined) {
222 Error("AddFacet", "Shape %s already fully defined. Not adding", GetName());
223 return false;
224 }
225
226 Vertex_t vert[3];
227 vert[0] = pt0;
228 vert[1] = pt1;
229 vert[2] = pt2;
230 int nvert = TGeoFacet::CompactFacet(vert, 3);
231 if (nvert < 3) {
232 Error("AddFacet", "Triangular facet at index %d degenerated. Not adding.", GetNfacets());
233 return false;
234 }
235 fNvert += 3;
236 fNseg += 3;
237 fFacets.emplace_back(pt0, pt1, pt2);
238
239 if (fNfacets > 0 && GetNfacets() == fNfacets)
240 CloseShape();
241 return true;
242}
243
244////////////////////////////////////////////////////////////////////////////////
245/// Adding a triangular facet from indices of vertices
246
247bool TGeoTessellated::AddFacet(int i0, int i1, int i2)
248{
249 if (fDefined) {
250 Error("AddFacet", "Shape %s already fully defined. Not adding", GetName());
251 return false;
252 }
253 if (fVertices.size() == 0) {
254 Error("AddFacet", "Shape %s Cannot add facets by indices without vertices. Not adding", GetName());
255 return false;
256 }
257
258 fNseg += 3;
259 fFacets.emplace_back(&fVertices, 3, i0, i1, i2);
260 return true;
261}
262
263////////////////////////////////////////////////////////////////////////////////
264/// Adding a quadrilateral facet from vertex positions in absolute coordinates
265
266bool TGeoTessellated::AddFacet(const Vertex_t &pt0, const Vertex_t &pt1, const Vertex_t &pt2, const Vertex_t &pt3)
267{
268 if (fDefined) {
269 Error("AddFacet", "Shape %s already fully defined. Not adding", GetName());
270 return false;
271 }
272 Vertex_t vert[4];
273 vert[0] = pt0;
274 vert[1] = pt1;
275 vert[2] = pt2;
276 vert[3] = pt3;
277 int nvert = TGeoFacet::CompactFacet(vert, 4);
278 if (nvert < 3) {
279 Error("AddFacet", "Quadrilateral facet at index %d degenerated. Not adding.", GetNfacets());
280 return false;
281 }
282
283 fNvert += nvert;
284 fNseg += nvert;
285 if (nvert == 3)
286 fFacets.emplace_back(vert[0], vert[1], vert[2]);
287 else
288 fFacets.emplace_back(vert[0], vert[1], vert[2], vert[3]);
289
290 if (fNfacets > 0 && GetNfacets() == fNfacets)
291 CloseShape(false);
292 return true;
293}
294
295////////////////////////////////////////////////////////////////////////////////
296/// Adding a quadrilateral facet from indices of vertices
297
298bool TGeoTessellated::AddFacet(int i0, int i1, int i2, int i3)
299{
300 if (fDefined) {
301 Error("AddFacet", "Shape %s already fully defined. Not adding", GetName());
302 return false;
303 }
304 if (fVertices.size() == 0) {
305 Error("AddFacet", "Shape %s Cannot add facets by indices without vertices. Not adding", GetName());
306 return false;
307 }
308
309 fNseg += 4;
310 fFacets.emplace_back(&fVertices, 4, i0, i1, i2, i3);
311 return true;
312}
313
314////////////////////////////////////////////////////////////////////////////////
315/// Close the shape: calculate bounding box and compact vertices
316
317void TGeoTessellated::CloseShape(bool check, bool fixFlipped, bool verbose)
318{
319 // Compact the array of vertices
320 constexpr double tolerance = 1.e-10;
321 constexpr int ngrid = 100;
322 ComputeBBox();
323 double minExtent[3];
324 minExtent[0] = fOrigin[0] - fDX - tolerance;
325 minExtent[1] = fOrigin[1] - fDY - tolerance;
326 minExtent[2] = fOrigin[2] - fDZ - tolerance;
327
328 double invExtent[3];
329 invExtent[0] = 0.5 / (fDX + tolerance);
330 invExtent[1] = 0.5 / (fDY + tolerance);
331 invExtent[2] = 0.5 / (fDZ + tolerance);
332
333 if (fVertices.size() > 0) {
334 fDefined = true;
335 if (!check)
336 return;
337
338 // Check facets
339 for (auto &facet : fFacets) {
340 facet.Check();
341 }
342 fClosedBody = CheckClosure(fixFlipped, verbose);
343 return;
344 }
345
346 auto AddVertex = [this](const Vertex_t &vertex) {
347 // Check if vertex exists
348 int ivert = 0;
349 for (const auto &current_vert : fVertices) {
350 if (current_vert == vertex)
351 return ivert;
352 ivert++;
353 }
354 // Vertex new, just add it
355 fVertices.push_back(vertex);
356 return ivert;
357 };
358
359 auto GetHashIndex = [&](const Vertex_t &vertex) {
360 // Get the hash index for a vertex in a 10x10x10 grid in the bounding box
361 int index = 0;
362 for (int i = 0; i < 3; ++i) {
363 int ind = int(ngrid * (vertex[i] - minExtent[i]) * invExtent[i]); // between [0, ngrid-1]
364 assert(ind < (int)ngrid);
365 for (int j = i + 1; j < 3; ++j)
366 ind *= ngrid;
367 index += ind;
368 }
369 return index;
370 };
371
372 // In case the number of vertices is small, just compare with all others
373 int ind[4] = {-1, -1, -1, -1};
374 if (fNvert < 1000) {
375 for (auto &facet : fFacets) {
376 int nvert = facet.GetNvert();
377 for (int i = 0; i < nvert; ++i) {
378 // Check if vertex exists already
379 ind[i] = AddVertex(facet.GetVertex(i));
380 }
381 facet.SetVertices(&fVertices, nvert, ind[0], ind[1], ind[2], ind[3]);
382 }
383 } else {
384 // Use hash index for each vertex
385 using CellVec_t = std::vector<int>;
386 auto grid = new std::array<CellVec_t, ngrid * ngrid * ngrid>;
387 for (auto &facet : fFacets) {
388 int nvert = facet.GetNvert();
389 for (int i = 0; i < nvert; ++i) {
390 // Check if vertex exists already
391 const Vertex_t &vertex = facet.GetVertex(i);
392 int hashind = GetHashIndex(vertex);
393 bool isAdded = false;
394 for (auto ivert : grid->operator[](hashind)) {
395 if (vertex == fVertices[ivert]) {
396 ind[i] = ivert;
397 isAdded = true;
398 break;
399 }
400 }
401 if (!isAdded) {
402 fVertices.push_back(vertex);
403 ind[i] = fVertices.size() - 1;
404 grid->operator[](hashind).push_back(ind[i]);
405 }
406 }
407 facet.SetVertices(&fVertices, nvert, ind[0], ind[1], ind[2], ind[3]);
408 }
409 delete grid;
410 }
411 fNvert = fVertices.size();
412 fNfacets = fFacets.size();
413 fDefined = true;
414 if (!check)
415 return;
416
417 // Check facets
418 for (auto &facet : fFacets) {
419 facet.Check();
420 }
421
422 fClosedBody = CheckClosure(fixFlipped, verbose);
423}
424
425////////////////////////////////////////////////////////////////////////////////
426/// Check closure of the solid and check/fix flipped normals
427
428bool TGeoTessellated::CheckClosure(bool fixFlipped, bool verbose)
429{
430 int *nn = new int[fNfacets];
431 bool *flipped = new bool[fNfacets];
432 bool hasorphans = false;
433 bool hasflipped = false;
434 for (int i = 0; i < fNfacets; ++i) {
435 nn[i] = 0;
436 flipped[i] = false;
437 }
438
439 for (int icrt = 0; icrt < fNfacets; ++icrt) {
440 // all neighbours checked?
441 if (nn[icrt] >= fFacets[icrt].GetNvert())
442 continue;
443 for (int i = icrt + 1; i < fNfacets; ++i) {
444 bool isneighbour = fFacets[icrt].IsNeighbour(fFacets[i], flipped[i]);
445 if (isneighbour) {
446 if (flipped[icrt])
447 flipped[i] = !flipped[i];
448 if (flipped[i])
449 hasflipped = true;
450 nn[icrt]++;
451 nn[i]++;
452 if (nn[icrt] == fFacets[icrt].GetNvert())
453 break;
454 }
455 }
456 if (nn[icrt] < fFacets[icrt].GetNvert())
457 hasorphans = true;
458 }
459
460 if (hasorphans && verbose) {
461 Error("Check", "Tessellated solid %s has following not fully connected facets:", GetName());
462 for (int icrt = 0; icrt < fNfacets; ++icrt) {
463 if (nn[icrt] < fFacets[icrt].GetNvert())
464 std::cout << icrt << " (" << fFacets[icrt].GetNvert() << " edges, " << nn[icrt] << " neighbours)\n";
465 }
466 }
467 fClosedBody = !hasorphans;
468 int nfixed = 0;
469 if (hasflipped) {
470 if (verbose)
471 Warning("Check", "Tessellated solid %s has following facets with flipped normals:", GetName());
472 for (int icrt = 0; icrt < fNfacets; ++icrt) {
473 if (flipped[icrt]) {
474 if (verbose)
475 std::cout << icrt << "\n";
476 if (fixFlipped) {
477 fFacets[icrt].Flip();
478 nfixed++;
479 }
480 }
481 }
482 if (nfixed && verbose)
483 Info("Check", "Automatically flipped %d facets to match first defined facet", nfixed);
484 }
485 delete[] nn;
486 delete[] flipped;
487
488 return !hasorphans;
489}
490
491////////////////////////////////////////////////////////////////////////////////
492/// Compute bounding box
493
495{
496 const double kBig = TGeoShape::Big();
497 double vmin[3] = {kBig, kBig, kBig};
498 double vmax[3] = {-kBig, -kBig, -kBig};
499 for (const auto &facet : fFacets) {
500 for (int i = 0; i < facet.GetNvert(); ++i) {
501 for (int j = 0; j < 3; ++j) {
502 vmin[j] = TMath::Min(vmin[j], facet.GetVertex(i).operator[](j));
503 vmax[j] = TMath::Max(vmax[j], facet.GetVertex(i).operator[](j));
504 }
505 }
506 }
507 fDX = 0.5 * (vmax[0] - vmin[0]);
508 fDY = 0.5 * (vmax[1] - vmin[1]);
509 fDZ = 0.5 * (vmax[2] - vmin[2]);
510 for (int i = 0; i < 3; ++i)
511 fOrigin[i] = 0.5 * (vmax[i] + vmin[i]);
512}
513
514////////////////////////////////////////////////////////////////////////////////
515/// Returns numbers of vertices, segments and polygons composing the shape mesh.
516
517void TGeoTessellated::GetMeshNumbers(int &nvert, int &nsegs, int &npols) const
518{
519 nvert = fNvert;
520 nsegs = fNseg;
521 npols = GetNfacets();
522}
523
524////////////////////////////////////////////////////////////////////////////////
525/// Creates a TBuffer3D describing *this* shape.
526/// Coordinates are in local reference frame.
527
529{
530 const int nvert = fNvert;
531 const int nsegs = fNseg;
532 const int npols = GetNfacets();
533 auto buff = new TBuffer3D(TBuffer3DTypes::kGeneric, nvert, 3 * nvert, nsegs, 3 * nsegs, npols, 6 * npols);
534 if (buff) {
535 SetPoints(buff->fPnts);
536 SetSegsAndPols(*buff);
537 }
538 return buff;
539}
540
541////////////////////////////////////////////////////////////////////////////////
542/// Prints basic info
543
545{
546 std::cout << "=== Tessellated shape " << GetName() << " having " << GetNvertices() << " vertices and "
547 << GetNfacets() << " facets\n";
548}
549
550////////////////////////////////////////////////////////////////////////////////
551/// Fills TBuffer3D structure for segments and polygons.
552
554{
555 const int c = GetBasicColor();
556 int *segs = buff.fSegs;
557 int *pols = buff.fPols;
558
559 int indseg = 0; // segment internal data index
560 int indpol = 0; // polygon internal data index
561 int sind = 0; // segment index
562 for (const auto &facet : fFacets) {
563 auto nvert = facet.GetNvert();
564 pols[indpol++] = c;
565 pols[indpol++] = nvert;
566 for (auto j = 0; j < nvert; ++j) {
567 int k = (j + 1) % nvert;
568 // segment made by next consecutive points
569 segs[indseg++] = c;
570 segs[indseg++] = facet.GetVertexIndex(j);
571 segs[indseg++] = facet.GetVertexIndex(k);
572 // add segment to current polygon and increment segment index
573 pols[indpol + nvert - j - 1] = sind++;
574 }
575 indpol += nvert;
576 }
577}
578
579////////////////////////////////////////////////////////////////////////////////
580/// Fill tessellated points to an array.
581
583{
584 int ind = 0;
585 for (const auto &vertex : fVertices) {
586 vertex.CopyTo(&points[ind]);
587 ind += 3;
588 }
589}
590
591////////////////////////////////////////////////////////////////////////////////
592/// Fill tessellated points in float.
593
595{
596 int ind = 0;
597 for (const auto &vertex : fVertices) {
598 points[ind++] = vertex.x();
599 points[ind++] = vertex.y();
600 points[ind++] = vertex.z();
601 }
602}
603
604////////////////////////////////////////////////////////////////////////////////
605/// Resize the shape by scaling vertices within maxsize and center to origin
606
608{
609 using Vector3_t = Vertex_t;
610
611 if (!fDefined) {
612 Error("ResizeCenter", "Not all faces are defined");
613 return;
614 }
615 Vector3_t origin(fOrigin[0], fOrigin[1], fOrigin[2]);
616 double maxedge = TMath::Max(TMath::Max(fDX, fDY), fDZ);
617 double scale = maxsize / maxedge;
618 for (size_t i = 0; i < fVertices.size(); ++i) {
619 fVertices[i] = scale * (fVertices[i] - origin);
620 }
621 fOrigin[0] = fOrigin[1] = fOrigin[2] = 0;
622 fDX *= scale;
623 fDY *= scale;
624 fDZ *= scale;
625}
626
627////////////////////////////////////////////////////////////////////////////////
628/// Fills a static 3D buffer and returns a reference.
629
630const TBuffer3D &TGeoTessellated::GetBuffer3D(int reqSections, Bool_t localFrame) const
631{
632 static TBuffer3D buffer(TBuffer3DTypes::kGeneric);
633
634 FillBuffer3D(buffer, reqSections, localFrame);
635
636 const int nvert = fNvert;
637 const int nsegs = fNseg;
638 const int npols = GetNfacets();
639
640 if (reqSections & TBuffer3D::kRawSizes) {
641 if (buffer.SetRawSizes(nvert, 3 * nvert, nsegs, 3 * nsegs, npols, 6 * npols)) {
643 }
644 }
645 if ((reqSections & TBuffer3D::kRaw) && buffer.SectionsValid(TBuffer3D::kRawSizes)) {
646 SetPoints(buffer.fPnts);
647 if (!buffer.fLocalFrame) {
648 TransformPoints(buffer.fPnts, buffer.NbPnts());
649 }
650
651 SetSegsAndPols(buffer);
653 }
654
655 return buffer;
656}
657
658////////////////////////////////////////////////////////////////////////////////
659/// Reads a single tessellated solid from an .obj file.
660
661TGeoTessellated *TGeoTessellated::ImportFromObjFormat(const char *objfile, bool check, bool verbose)
662{
663 using namespace std;
664
665 vector<Vertex_t> vertices;
666 vector<string> sfacets;
667
668 struct FacetInd_t {
669 int i0 = -1;
670 int i1 = -1;
671 int i2 = -1;
672 int i3 = -1;
673 int nvert = 0;
674 FacetInd_t(int a, int b, int c)
675 {
676 i0 = a;
677 i1 = b;
678 i2 = c;
679 nvert = 3;
680 };
681 FacetInd_t(int a, int b, int c, int d)
682 {
683 i0 = a;
684 i1 = b;
685 i2 = c;
686 i3 = d;
687 nvert = 4;
688 };
689 };
690
691 vector<FacetInd_t> facets;
692 // List of geometric vertices, with (x, y, z [,w]) coordinates, w is optional and defaults to 1.0.
693 // struct vtx_t { double x = 0; double y = 0; double z = 0; double w = 1; };
694
695 // Texture coordinates in u, [,v ,w]) coordinates, these will vary between 0 and 1. v, w are optional and default to
696 // 0.
697 // struct tex_t { double u; double v; double w; };
698
699 // List of vertex normals in (x,y,z) form; normals might not be unit vectors.
700 // struct vn_t { double x; double y; double z; };
701
702 // Parameter space vertices in ( u [,v] [,w] ) form; free form geometry statement
703 // struct vp_t { double u; double v; double w; };
704
705 // Faces are defined using lists of vertex, texture and normal indices which start at 1.
706 // Polygons such as quadrilaterals can be defined by using more than three vertex/texture/normal indices.
707 // f v1//vn1 v2//vn2 v3//vn3 ...
708
709 // Records starting with the letter "l" specify the order of the vertices which build a polyline.
710 // l v1 v2 v3 v4 v5 v6 ...
711
712 string line;
713 int ind[4] = {0};
714 ifstream file(objfile);
715 if (!file.is_open()) {
716 ::Error("TGeoTessellated::ImportFromObjFormat", "Unable to open %s", objfile);
717 return nullptr;
718 }
719
720 while (getline(file, line)) {
721 stringstream ss(line);
722 string tag;
723
724 // We ignore everything which is not a vertex or a face
725 if (line.rfind("v", 0) == 0 && line.rfind("vt", 0) != 0 && line.rfind("vn", 0) != 0 && line.rfind("vn", 0) != 0) {
726 // Decode the vertex
727 double pos[4] = {0, 0, 0, 1};
728 ss >> tag >> pos[0] >> pos[1] >> pos[2] >> pos[3];
729 vertices.emplace_back(pos[0] * pos[3], pos[1] * pos[3], pos[2] * pos[3]);
730 }
731
732 else if (line.rfind("f", 0) == 0) {
733 // Decode the face
734 ss >> tag;
735 string word;
736 sfacets.clear();
737 while (ss >> word)
738 sfacets.push_back(word);
739 if (sfacets.size() > 4 || sfacets.size() < 3) {
740 ::Error("TGeoTessellated::ImportFromObjFormat", "Detected face having unsupported %zu vertices",
741 sfacets.size());
742 return nullptr;
743 }
744 int nvert = 0;
745 for (auto &sword : sfacets) {
746 stringstream ssword(sword);
747 string token;
748 getline(ssword, token, '/'); // just need the vertex index, which is the first token
749 // Convert string token to integer
750
751 ind[nvert++] = stoi(token) - 1;
752 if (ind[nvert - 1] < 0) {
753 ::Error("TGeoTessellated::ImportFromObjFormat", "Unsupported relative vertex index definition in %s",
754 objfile);
755 return nullptr;
756 }
757 }
758 if (nvert == 3)
759 facets.emplace_back(ind[0], ind[1], ind[2]);
760 else
761 facets.emplace_back(ind[0], ind[1], ind[2], ind[3]);
762 }
763 }
764
765 int nvertices = (int)vertices.size();
766 int nfacets = (int)facets.size();
767 if (nfacets < 3) {
768 ::Error("TGeoTessellated::ImportFromObjFormat", "Not enough faces detected in %s", objfile);
769 return nullptr;
770 }
771
772 string sobjfile(objfile);
773 if (verbose)
774 std::cout << "Read " << nvertices << " vertices and " << nfacets << " facets from " << sobjfile << endl;
775
776 auto tsl = new TGeoTessellated(sobjfile.erase(sobjfile.find_last_of('.')).c_str(), vertices);
777
778 for (int i = 0; i < nfacets; ++i) {
779 auto facet = facets[i];
780 if (facet.nvert == 3)
781 tsl->AddFacet(facet.i0, facet.i1, facet.i2);
782 else
783 tsl->AddFacet(facet.i0, facet.i1, facet.i2, facet.i3);
784 }
785 tsl->CloseShape(check, true, verbose);
786 tsl->Print();
787 return tsl;
788}
#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
Definition RtypesCore.h:57
const char Option_t
Definition RtypesCore.h:66
#define ClassImp(name)
Definition Rtypes.h:364
char name[80]
Definition TGX11.cxx:110
segment * segs
Definition X3DBuffer.c:23
point * points
Definition X3DBuffer.c:22
Generic 3D primitive description class.
Definition TBuffer3D.h:18
Int_t * fPols
Definition TBuffer3D.h:114
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
Int_t * fSegs
Definition TBuffer3D.h:113
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:112
Box class.
Definition TGeoBBox.h:18
Double_t fDX
Definition TGeoBBox.h:21
Double_t fOrigin[3]
Definition TGeoBBox.h:24
Double_t fDY
Definition TGeoBBox.h:22
Double_t fDZ
Definition TGeoBBox.h:23
virtual void FillBuffer3D(TBuffer3D &buffer, Int_t reqSections, Bool_t localFrame) const
Fills the supplied buffer, with sections in desired frame See TBuffer3D.h for explanation of sections...
Tessellated::VertexVec_t VertexVec_t
bool Check() const
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)
int GetVertexIndex(int ivert) const
Vertex_t & GetVertex(int ivert)
VertexVec_t * fVertices
int fNvert
array of vertices
Vertex_t ComputeNormal(bool &degenerated) const
int GetNvert() const
const TGeoFacet & operator=(const TGeoFacet &other)
static Double_t Big()
Definition TGeoShape.h:88
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)
virtual const char * GetName() const
Get the shape name.
Tessellated solid class.
void ResizeCenter(double maxsize)
Resize and center the shape in a box of size maxsize.
virtual void Print(Option_t *option="") const
Prints basic info.
virtual const TBuffer3D & GetBuffer3D(int reqSections, Bool_t localFrame) const
Fills a static 3D buffer and returns a reference.
virtual void SetSegsAndPols(TBuffer3D &buff) const
Fills TBuffer3D structure for segments and polygons.
bool CheckClosure(bool fixFlipped=true, bool verbose=true)
Check closure of the solid and check/fix flipped normals.
int GetNvertices() const
virtual void SetPoints(double *points) const
Fill tessellated points to an array.
virtual TBuffer3D * MakeBuffer3D() const
Creates a TBuffer3D describing this shape.
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
std::vector< TGeoFacet > fFacets
static TGeoTessellated * ImportFromObjFormat(const char *objfile, bool check=false, bool verbose=false)
Reader from .obj format.
void ComputeBBox()
Compute bounding box.
virtual void AfterStreamer()
Function to be called after reading tessellated volumes from the geometry file.
bool AddFacet(const Vertex_t &pt0, const Vertex_t &pt1, const Vertex_t &pt2)
Adding a triangular facet from vertex positions in absolute coordinates.
virtual void GetMeshNumbers(int &nvert, int &nsegs, int &npols) const
Returns numbers of vertices, segments and polygons composing the shape mesh.
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:949
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition TObject.cxx:963
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition TObject.cxx:937
TLine * line
Short_t Max(Short_t a, Short_t b)
Definition TMathBase.h:208
Short_t Min(Short_t a, Short_t b)
Definition TMathBase.h:176
Typedefs used by the geometry group.
Definition file.py:1
static int AddVertex(GLUtesselator *tess, GLdouble coords[3], void *data)
Definition tess.c:355
REAL * vertex
Definition triangle.c:513
void flip(struct mesh *m, struct behavior *b, struct otri *flipedge)
Definition triangle.c:7890