58 if (f[fp] < 0 || f[fp] >= surf.
points().
size())
61 <<
"triangle " << facei <<
" vertices " << f
62 <<
" uses point indices outside point range 0.." 68 if ((f[0] == f[1]) || (f[0] == f[2]) || (f[1] == f[2]))
71 <<
"triangle " << facei
72 <<
" uses non-unique vertices " << f
86 label nbrFacei = fFaces[i];
88 if (nbrFacei <= facei)
98 ((f[0] == nbrF[0]) || (f[0] == nbrF[1]) || (f[0] == nbrF[2]))
99 && ((f[1] == nbrF[0]) || (f[1] == nbrF[1]) || (f[1] == nbrF[2]))
100 && ((f[2] == nbrF[0]) || (f[2] == nbrF[1]) || (f[2] == nbrF[2]))
104 <<
"triangle " << facei <<
" vertices " << f
105 <<
" has the same vertices as triangle " << nbrFacei
106 <<
" vertices " << nbrF
125 scalar dist = nBins/(max -
min);
131 scalar val = vals[i];
139 else if (val >= max - SMALL)
145 index =
label((val - min)*dist);
147 if ((index < 0) || (index >= nBins))
150 <<
"value " << val <<
" at index " << i
151 <<
" outside range " << min <<
" .. " << max <<
endl;
172 int main(
int argc,
char *argv[])
178 "checkSelfIntersection",
179 "also check for self-intersection" 184 "split surface along non-manifold edges" 185 " (default split is fully disconnected)" 195 "write vertices/blocks for blockMeshDict" 201 const bool checkSelfIntersect =
args.
optionFound(
"checkSelfIntersection");
205 Info<<
"Reading surface from " << surfFileName <<
" ..." << nl <<
endl;
223 Info<<
"// blockMeshDict info" << nl <<
nl;
225 Info<<
"vertices\n(" <<
nl;
228 Info<<
" " << cornerPts[ptI] <<
nl;
235 <<
" hex (0 1 2 3 4 5 6 7) (10 10 10) simpleGrading (1 1 1)\n" 238 Info<<
"edges\n();" << nl
239 <<
"patches\n();" <<
endl;
241 Info<< nl <<
"// end blockMeshDict info" << nl <<
endl;
253 label region = surf[facei].region();
255 if (region < 0 || region >= regionSize.size())
258 <<
"Triangle " << facei <<
" vertices " << surf[facei]
259 <<
" has region " << region <<
" which is outside the range" 265 regionSize[region]++;
269 Info<<
"Region\tSize" << nl
270 <<
"------\t----" <<
nl;
288 if (!validTri(verbose, surf, facei))
290 illegalFaces.append(facei);
294 if (illegalFaces.size())
296 Info<<
"Surface has " << illegalFaces.size()
297 <<
" illegal triangles." <<
endl;
300 Info<<
"Dumping conflicting face labels to " << str.name() << endl
301 <<
"Paste this into the input for surfaceSubset" <<
endl;
306 Info<<
"Surface has no illegal triangles." <<
endl;
322 if (f[0] == f[1] || f[0] == f[2] || f[1] == f[2])
339 labelList binCount = countBins(0, 1, 20, triQ);
341 Info<<
"Triangle quality (equilateral=1, collapsed=0):" 348 scalar dist = (1.0 - 0.0)/20.0;
352 Info<<
" " << min <<
" .. " << min+dist <<
" : " 353 << 1.0/surf.
size() * binCount[binI]
362 Info<<
" min " << triQ[minIndex] <<
" for triangle " << minIndex
364 <<
" max " << triQ[maxIndex] <<
" for triangle " << maxIndex
369 if (triQ[minIndex] < SMALL)
372 << triQ[minIndex] <<
". This might give problems in" 373 <<
" self-intersection testing later on." <<
endl;
382 if (triQ[facei] < 1e-11)
384 problemFaces.append(facei);
388 if (!problemFaces.empty())
392 Info<<
"Dumping bad quality faces to " << str.name() << endl
393 <<
"Paste this into the input for surfaceSubset" << nl
413 edgeMag[edgeI] = edges[edgeI].mag(localPoints);
419 const edge& minE = edges[minEdgeI];
420 const edge& maxE = edges[maxEdgeI];
423 Info<<
"Edges:" << nl
424 <<
" min " << edgeMag[minEdgeI] <<
" for edge " << minEdgeI
425 <<
" points " << localPoints[minE[0]] << localPoints[minE[1]]
427 <<
" max " << edgeMag[maxEdgeI] <<
" for edge " << maxEdgeI
428 <<
" points " << localPoints[maxE[0]] << localPoints[maxE[1]]
442 scalar smallDim = 1e-6 * bb.
mag();
444 Info<<
"Checking for points less than 1e-6 of bounding box (" 445 << bb.span() <<
" metre) apart." 453 for (
label i = 1; i < sortedMag.size(); i++)
455 label ptI = sortedMag.indices()[i];
457 label prevPtI = sortedMag.indices()[i-1];
459 if (
mag(localPoints[ptI] - localPoints[prevPtI]) < smallDim)
468 const edge& e = edges[pEdges[i]];
470 if (e[0] == prevPtI || e[1] == prevPtI)
483 Info<<
" close unconnected points " 484 << ptI <<
' ' << localPoints[ptI]
485 <<
" and " << prevPtI <<
' ' 486 << localPoints[prevPtI]
488 <<
mag(localPoints[ptI] - localPoints[prevPtI])
493 Info<<
" small edge between points " 494 << ptI <<
' ' << localPoints[ptI]
495 <<
" and " << prevPtI <<
' ' 496 << localPoints[prevPtI]
498 <<
mag(localPoints[ptI] - localPoints[prevPtI])
504 Info<<
"Found " << nClose <<
" nearby points." << nl
517 label nSingleEdges = 0;
520 const labelList& myFaces = eFaces[edgeI];
522 if (myFaces.
size() == 1)
524 problemFaces.
append(myFaces[0]);
530 label nMultEdges = 0;
533 const labelList& myFaces = eFaces[edgeI];
535 if (myFaces.
size() > 2)
539 problemFaces.append(myFaces[myFacei]);
545 problemFaces.shrink();
547 if ((nSingleEdges != 0) || (nMultEdges != 0))
549 Info<<
"Surface is not closed since not all edges connected to " 550 <<
"two faces:" << endl
551 <<
" connected to one face : " << nSingleEdges << endl
552 <<
" connected to >2 faces : " << nMultEdges <<
endl;
554 Info<<
"Conflicting face labels:" << problemFaces.size() <<
endl;
558 Info<<
"Dumping conflicting face labels to " << str.name() << endl
559 <<
"Paste this into the input for surfaceSubset" <<
endl;
565 Info<<
"Surface is closed. All edges connected to two faces." <<
endl;
576 if (splitNonManifold)
581 if (eFaces[edgeI].size() > 2)
583 borderEdge[edgeI] =
true;
591 Info<<
"Number of unconnected parts : " << numZones <<
endl;
595 Info<<
"Splitting surface into parts ..." << endl <<
endl;
598 const word fileType = surfFileNameBase.ext();
600 surfFileNameBase = surfFileNameBase.lessExt();
602 if (fileType ==
"gz")
604 surfFileNameBase = surfFileNameBase.lessExt();
609 Info<<
"Writing zoning to " 615 + vtkSurfaceWriter::typeName
617 <<
"..." << endl <<
endl;
623 scalarFaceZone[i] = faceZone[i];
628 faces[i] = surf[i].triFaceFace();
650 if (faceZone[facei] ==
zone)
652 includeMap[facei] =
true;
671 Info<<
"writing part " <<
zone <<
" size " << subSurf.
size()
672 <<
" to " << subName <<
endl;
674 subSurf.
write(subName);
694 <<
"Number of zones (connected area with consistent normal) : " 695 << numNormalZones <<
endl;
697 if (numNormalZones > 1)
699 Info<<
"More than one normal orientation." <<
endl;
708 if (checkSelfIntersect)
710 Info<<
"Checking self-intersection." <<
endl;
716 OBJstream intStream(
"selfInterPoints.obj");
740 intStream.write(hitInfo.hitPoint());
747 Info<<
"Surface is not self-intersecting" <<
endl;
751 Info<<
"Surface is self-intersecting at " << nInt
752 <<
" locations." <<
endl;
753 Info<<
"Writing intersection points to " << intStream.name()
const labelListList & pointEdges() const
Return point-edge addressing.
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
A class for handling file names.
A list that is sorted upon construction or when explicitly requested with the sort() method...
void size(const label)
Override size to be inconsistent with allocated storage.
Ostream & endl(Ostream &os)
Add newline and flush stream.
bool optionFound(const word &opt) const
Return true if the named option is found.
static void noParallel()
Remove the parallel options.
static SLList< string > validArgs
A list of valid (mandatory) arguments.
label findMin(const ListType &, const label start=0)
Find index of min element (and less than given element).
A bounding box defined in terms of the points at its extremities.
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
const labelList & meshPoints() const
Return labelList of mesh points in patch. They are constructed.
Pair< int > faceMap(const label facePi, const face &faceP, const label faceNi, const face &faceN)
Helper class to search on triSurface.
const labelListList & faceFaces() const
Return face-face addressing.
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
void writeStats(Ostream &) const
Write some statistics.
A class for handling words, derived from string.
pointField points(const pointField &) const
Return the points corresponding to this face.
Extract command arguments and options from the supplied argc and argv parameters. ...
word name() const
Return file name (part beyond last /)
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
const Field< PointType > & points() const
Return reference to global points.
const labelListList & edgeFaces() const
Return edge-face addressing.
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
scalar mag(const pointField &) const
Return scalar magnitude.
Triangle with additional region number.
label findMax(const ListType &, const label start=0)
Find index of max element (and larger than given element).
const Field< PointType > & localPoints() const
Return pointField of points in patch.
label nEdges() const
Return number of edges in patch.
OFstream which keeps track of vertices.
virtual void write(const fileName &outputDir, const fileName &surfaceName, const pointField &points, const faceList &faces, const bool verbose=false) const
Write single surface geometry to file.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
word name(const complex &)
Return a string representation of a complex.
const geometricSurfacePatchList & patches() const
virtual int width() const
Get width of output field.
Non-pointer based hierarchical recursive searching.
#define WarningInFunction
Report a warning using Foam::Warning.
A surfaceWriter for VTK legacy format.
triangle< point, const point & > triPointRef
dimensioned< scalar > mag(const dimensioned< Type > &)
fileName path() const
Return directory path name (part before last /)
virtual Ostream & write(const token &)=0
Write next token to stream.
static void addBoolOption(const word &opt, const string &usage="")
Add to a bool option to validOptions with usage information.
label markZones(const boolList &borderEdge, labelList &faceZone) const
(size and) fills faceZone with zone of face. Zone is area
triSurface subsetMesh(const boolList &include, labelList &pointMap, labelList &faceMap) const
Return new surface. Returns pointMap, faceMap from.
Triangulated surface description with patch information.
Foam::argList args(argc, argv)