writeMeshObj.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Application
25  writeMeshObj
26 
27 Description
28  For mesh debugging: writes mesh as three separate OBJ files which can
29  be viewed with e.g. javaview.
30 
31  meshPoints_XXX.obj : all points and edges as lines.
32  meshFaceCentres_XXX.obj : all face centres.
33  meshCellCentres_XXX.obj : all cell centres.
34 
35  patch_YYY_XXX.obj : all face centres of patch YYY
36 
37  Optional: - patch faces (as polygons) : patchFaces_YYY_XXX.obj
38  - non-manifold edges : patchEdges_YYY_XXX.obj
39 
40 \*---------------------------------------------------------------------------*/
41 
42 #include "argList.H"
43 #include "timeSelector.H"
44 #include "Time.H"
45 #include "polyMesh.H"
46 #include "OFstream.H"
47 #include "meshTools.H"
48 #include "cellSet.H"
49 #include "faceSet.H"
50 #include "SubField.H"
51 
52 using namespace Foam;
53 
54 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
55 
56 void writeOBJ(const point& pt, Ostream& os)
57 {
58  os << "v " << pt.x() << ' ' << pt.y() << ' ' << pt.z() << nl;
59 }
60 
61 // All edges of mesh
62 void writePoints(const polyMesh& mesh, const fileName& timeName)
63 {
64  label vertI = 0;
65 
66  fileName pointFile(mesh.time().path()/"meshPoints_" + timeName + ".obj");
67 
68  Info<< "Writing mesh points and edges to " << pointFile << endl;
69 
70  OFstream pointStream(pointFile);
71 
72  forAll(mesh.points(), pointi)
73  {
74  writeOBJ(mesh.points()[pointi], pointStream);
75  vertI++;
76  }
77 
78  forAll(mesh.edges(), edgeI)
79  {
80  const edge& e = mesh.edges()[edgeI];
81 
82  pointStream << "l " << e.start() + 1 << ' ' << e.end() + 1 << nl;
83  }
84 }
85 
86 
87 // Edges for subset of cells
88 void writePoints
89 (
90  const polyMesh& mesh,
91  const labelList& cellLabels,
92  const fileName& timeName
93 )
94 {
95  fileName fName(mesh.time().path()/"meshPoints_" + timeName + ".obj");
96 
97  Info<< "Writing mesh points and edges to " << fName << endl;
98 
99  OFstream str(fName);
100 
101  // OBJ file vertex
102  label vertI = 0;
103 
104  // From point to OBJ file vertex
105  Map<label> pointToObj(6*cellLabels.size());
106 
107  forAll(cellLabels, i)
108  {
109  label celli = cellLabels[i];
110 
111  const labelList& cEdges = mesh.cellEdges()[celli];
112 
113  forAll(cEdges, cEdgeI)
114  {
115  const edge& e = mesh.edges()[cEdges[cEdgeI]];
116 
117  label v0;
118 
119  Map<label>::iterator e0Fnd = pointToObj.find(e[0]);
120 
121  if (e0Fnd == pointToObj.end())
122  {
123  meshTools::writeOBJ(str, mesh.points()[e[0]]);
124  v0 = vertI++;
125  pointToObj.insert(e[0], v0);
126  }
127  else
128  {
129  v0 = e0Fnd();
130  }
131 
132  label v1;
133 
134  Map<label>::iterator e1Fnd = pointToObj.find(e[1]);
135 
136  if (e1Fnd == pointToObj.end())
137  {
138  meshTools::writeOBJ(str, mesh.points()[e[1]]);
139  v1 = vertI++;
140  pointToObj.insert(e[1], v1);
141  }
142  else
143  {
144  v1 = e1Fnd();
145  }
146 
147 
148  str << "l " << v0+1 << ' ' << v1+1 << nl;
149  }
150  }
151 }
152 
153 
154 // Edges of single cell
155 void writePoints
156 (
157  const polyMesh& mesh,
158  const label celli,
159  const fileName& timeName
160 )
161 {
162  fileName fName
163  (
164  mesh.time().path()
165  / "meshPoints_" + timeName + '_' + name(celli) + ".obj"
166  );
167 
168  Info<< "Writing mesh points and edges to " << fName << endl;
169 
170  OFstream pointStream(fName);
171 
172  const cell& cFaces = mesh.cells()[celli];
173 
174  meshTools::writeOBJ(pointStream, mesh.faces(), mesh.points(), cFaces);
175 }
176 
177 
178 
179 // All face centres
180 void writeFaceCentres(const polyMesh& mesh,const fileName& timeName)
181 {
182  fileName faceFile
183  (
184  mesh.time().path()
185  / "meshFaceCentres_" + timeName + ".obj"
186  );
187 
188  Info<< "Writing mesh face centres to " << faceFile << endl;
189 
190  OFstream faceStream(faceFile);
191 
192  forAll(mesh.faceCentres(), facei)
193  {
194  writeOBJ(mesh.faceCentres()[facei], faceStream);
195  }
196 }
197 
198 
199 void writeCellCentres(const polyMesh& mesh, const fileName& timeName)
200 {
201  fileName cellFile
202  (
203  mesh.time().path()/"meshCellCentres_" + timeName + ".obj"
204  );
205 
206  Info<< "Writing mesh cell centres to " << cellFile << endl;
207 
208  OFstream cellStream(cellFile);
209 
210  forAll(mesh.cellCentres(), celli)
211  {
212  writeOBJ(mesh.cellCentres()[celli], cellStream);
213  }
214 }
215 
216 
217 void writePatchCentres
218 (
219  const polyMesh& mesh,
220  const fileName& timeName
221 )
222 {
223  const polyBoundaryMesh& patches = mesh.boundaryMesh();
224 
225  forAll(patches, patchi)
226  {
227  const polyPatch& pp = patches[patchi];
228 
229  fileName faceFile
230  (
231  mesh.time().path()/"patch_" + pp.name() + '_' + timeName + ".obj"
232  );
233 
234  Info<< "Writing patch face centres to " << faceFile << endl;
235 
236  OFstream patchFaceStream(faceFile);
237 
238  forAll(pp.faceCentres(), facei)
239  {
240  writeOBJ(pp.faceCentres()[facei], patchFaceStream);
241  }
242  }
243 }
244 
245 
246 void writePatchFaces
247 (
248  const polyMesh& mesh,
249  const fileName& timeName
250 )
251 {
252  const polyBoundaryMesh& patches = mesh.boundaryMesh();
253 
254  forAll(patches, patchi)
255  {
256  const polyPatch& pp = patches[patchi];
257 
258  fileName faceFile
259  (
260  mesh.time().path()
261  / "patchFaces_" + pp.name() + '_' + timeName + ".obj"
262  );
263 
264  Info<< "Writing patch faces to " << faceFile << endl;
265 
266  OFstream patchFaceStream(faceFile);
267 
268  forAll(pp.localPoints(), pointi)
269  {
270  writeOBJ(pp.localPoints()[pointi], patchFaceStream);
271  }
272 
273  forAll(pp.localFaces(), facei)
274  {
275  const face& f = pp.localFaces()[facei];
276 
277  patchFaceStream<< 'f';
278 
279  forAll(f, fp)
280  {
281  patchFaceStream << ' ' << f[fp]+1;
282  }
283  patchFaceStream << nl;
284  }
285  }
286 }
287 
288 
289 void writePatchBoundaryEdges
290 (
291  const polyMesh& mesh,
292  const fileName& timeName
293 )
294 {
295  const polyBoundaryMesh& patches = mesh.boundaryMesh();
296 
297  forAll(patches, patchi)
298  {
299  const polyPatch& pp = patches[patchi];
300 
301  fileName edgeFile
302  (
303  mesh.time().path()
304  / "patchEdges_" + pp.name() + '_' + timeName + ".obj"
305  );
306 
307  Info<< "Writing patch edges to " << edgeFile << endl;
308 
309  OFstream patchEdgeStream(edgeFile);
310 
311  forAll(pp.localPoints(), pointi)
312  {
313  writeOBJ(pp.localPoints()[pointi], patchEdgeStream);
314  }
315 
316  for (label edgeI = pp.nInternalEdges(); edgeI < pp.nEdges(); edgeI++)
317  {
318  if (pp.edgeFaces()[edgeI].size() == 1)
319  {
320  const edge& e = pp.edges()[edgeI];
321 
322  patchEdgeStream<< "l " << e[0]+1 << ' ' << e[1]+1 << nl;
323  }
324  }
325  }
326 }
327 
328 
329 void writePointCells
330 (
331  const polyMesh& mesh,
332  const label pointi,
333  const fileName& timeName
334 )
335 {
336  const labelList& pCells = mesh.pointCells()[pointi];
337 
338  labelHashSet allEdges(6*pCells.size());
339 
340  forAll(pCells, i)
341  {
342  const labelList& cEdges = mesh.cellEdges()[pCells[i]];
343 
344  forAll(cEdges, i)
345  {
346  allEdges.insert(cEdges[i]);
347  }
348  }
349 
350 
351  fileName pFile
352  (
353  mesh.time().path()
354  / "pointEdges_" + timeName + '_' + name(pointi) + ".obj"
355  );
356 
357  Info<< "Writing pointEdges to " << pFile << endl;
358 
359  OFstream pointStream(pFile);
360 
361  label vertI = 0;
362 
363  forAllConstIter(labelHashSet, allEdges, iter)
364  {
365  const edge& e = mesh.edges()[iter.key()];
366 
367  meshTools::writeOBJ(pointStream, mesh.points()[e[0]]); vertI++;
368  meshTools::writeOBJ(pointStream, mesh.points()[e[1]]); vertI++;
369  pointStream<< "l " << vertI-1 << ' ' << vertI << nl;
370  }
371 }
372 
373 
374 
375 int main(int argc, char *argv[])
376 {
378  (
379  "for mesh debugging: write mesh as separate OBJ files"
380  );
381 
384  (
385  "patchFaces",
386  "write patch faces edges"
387  );
389  (
390  "patchEdges",
391  "write patch boundary edges"
392  );
394  (
395  "cell",
396  "int",
397  "write points for the specified cell"
398  );
400  (
401  "face",
402  "int",
403  "write specified face"
404  );
406  (
407  "point",
408  "int",
409  "write specified point"
410  );
412  (
413  "cellSet",
414  "name",
415  "write points for specified cellSet"
416  );
418  (
419  "faceSet",
420  "name",
421  "write points for specified faceSet"
422  );
423  #include "addRegionOption.H"
424 
425  #include "setRootCase.H"
426  #include "createTime.H"
427  runTime.functionObjects().off();
428 
429  const bool patchFaces = args.optionFound("patchFaces");
430  const bool patchEdges = args.optionFound("patchEdges");
431  const bool doCell = args.optionFound("cell");
432  const bool doPoint = args.optionFound("point");
433  const bool doFace = args.optionFound("face");
434  const bool doCellSet = args.optionFound("cellSet");
435  const bool doFaceSet = args.optionFound("faceSet");
436 
437 
438  Info<< "Writing mesh objects as .obj files such that the object"
439  << " numbering" << endl
440  << "(for points, faces, cells) is consistent with"
441  << " Foam numbering (starting from 0)." << endl << endl;
442 
443  instantList timeDirs = timeSelector::select0(runTime, args);
444 
445  #include "createNamedPolyMesh.H"
446 
447  forAll(timeDirs, timeI)
448  {
449  runTime.setTime(timeDirs[timeI], timeI);
450 
451  Info<< "Time = " << runTime.timeName() << endl;
452 
453  polyMesh::readUpdateState state = mesh.readUpdate();
454 
455  if (!timeI || state != polyMesh::UNCHANGED)
456  {
457  if (patchFaces)
458  {
459  writePatchFaces(mesh, runTime.timeName());
460  }
461  if (patchEdges)
462  {
463  writePatchBoundaryEdges(mesh, runTime.timeName());
464  }
465  if (doCell)
466  {
467  label celli = args.optionRead<label>("cell");
468 
469  writePoints(mesh, celli, runTime.timeName());
470  }
471  if (doPoint)
472  {
473  label pointi = args.optionRead<label>("point");
474 
475  writePointCells(mesh, pointi, runTime.timeName());
476  }
477  if (doFace)
478  {
479  label facei = args.optionRead<label>("face");
480 
481  fileName fName
482  (
483  mesh.time().path()
484  / "meshPoints_"
485  + runTime.timeName()
486  + '_'
487  + name(facei)
488  + ".obj"
489  );
490 
491  Info<< "Writing mesh points and edges to " << fName << endl;
492 
493  OFstream str(fName);
494 
495  const face& f = mesh.faces()[facei];
496 
497  meshTools::writeOBJ(str, faceList(1, f), mesh.points());
498  }
499  if (doCellSet)
500  {
501  const word setName = args["cellSet"];
502 
503  cellSet cells(mesh, setName);
504 
505  Info<< "Read " << cells.size() << " cells from set " << setName
506  << endl;
507 
508  writePoints(mesh, cells.toc(), runTime.timeName());
509  }
510  if (doFaceSet)
511  {
512  const word setName = args["faceSet"];
513 
514  faceSet faces(mesh, setName);
515 
516  Info<< "Read " << faces.size() << " faces from set " << setName
517  << endl;
518 
519  fileName fName
520  (
521  mesh.time().path()
522  / "meshPoints_"
523  + runTime.timeName()
524  + '_'
525  + setName
526  + ".obj"
527  );
528 
529  Info<< "Writing mesh points and edges to " << fName << endl;
530 
531  OFstream str(fName);
532 
534  (
535  str,
536  mesh.faces(),
537  mesh.points(),
538  faces.toc()
539  );
540  }
541  else if
542  (
543  !patchFaces
544  && !patchEdges
545  && !doCell
546  && !doPoint
547  && !doFace
548  && !doCellSet
549  && !doFaceSet
550  )
551  {
552  // points & edges
553  writePoints(mesh, runTime.timeName());
554 
555  // face centres
556  writeFaceCentres(mesh, runTime.timeName());
557 
558  // cell centres
559  writeCellCentres(mesh, runTime.timeName());
560 
561  // Patch face centres
562  writePatchCentres(mesh, runTime.timeName());
563  }
564  }
565  else
566  {
567  Info<< "No mesh." << endl;
568  }
569 
570  Info<< nl << endl;
571  }
572 
573 
574  Info<< "End\n" << endl;
575 
576  return 0;
577 }
578 
579 
580 // ************************************************************************* //
label end() const
Return end vertex label.
Definition: edgeI.H:92
const Time & time() const
Return time.
const Cmpt & z() const
Definition: VectorI.H:87
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
A class for handling file names.
Definition: fileName.H:69
T optionRead(const word &opt) const
Read a value from the named option.
Definition: argListI.H:187
A list of face labels.
Definition: faceSet.H:48
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
static iteratorEnd end()
iteratorEnd set to beyond the end of any HashTable
Definition: HashTable.H:106
const Cmpt & x() const
Definition: VectorI.H:75
const List< Face > & localFaces() const
Return patch faces addressing into local point list.
Output to file stream.
Definition: OFstream.H:81
const edgeList & edges() const
Return mesh edges. Uses calcEdges.
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
List< face > faceList
Definition: faceListFwd.H:43
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
void writeOBJ(Ostream &os, const point &pt)
Write obj representation of point.
Definition: meshTools.C:203
const vectorField & faceCentres() const
const cellList & cells() const
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:979
iterator find(const label &)
Find and return an iterator set at the hashedEntry.
const labelListList & pointCells() const
const cellShapeList & cells
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:108
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
A class for handling words, derived from string.
Definition: word.H:59
const edgeList & edges() const
Return list of edges, address into LOCAL point list.
static void addOption(const word &opt, const string &param="", const string &usage="")
Add to an option to validOptions with usage information.
Definition: argList.C:93
const Cmpt & y() const
Definition: VectorI.H:81
const vectorField & cellCentres() const
const word & name() const
Return name.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
Foam::polyBoundaryMesh.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
label start() const
Return start vertex label.
Definition: edgeI.H:81
static const char nl
Definition: Ostream.H:262
const labelListList & cellEdges() const
static instantList select0(Time &runTime, const argList &args)
Return the set of times selected based on the argList options.
Definition: timeSelector.C:253
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
label nEdges() const
Return number of edges in patch.
label patchi
A cell is defined as a list of faces with extra functionality.
Definition: cell.H:56
A collection of cell labels.
Definition: cellSet.H:48
label nInternalEdges() const
Number of internal edges.
messageStream Info
Inserts points at locations specified in a pointFile into the surfaces to be conformed to of the conf...
Definition: pointFile.H:52
const labelListList & edgeFaces() const
Return edge-face addressing.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
static void addBoolOption(const word &opt, const string &usage="")
Add to a bool option to validOptions with usage information.
Definition: argList.C:83
const vectorField::subField faceCentres() const
Return face centres.
Definition: polyPatch.C:284
static void addNote(const string &)
Add extra notes for the usage information.
Definition: argList.C:124
readUpdateState
Enumeration defining the state of the mesh after a read update.
Definition: polyMesh.H:88
virtual readUpdateState readUpdate()
Update the mesh based on the mesh files saved in.
Definition: polyMeshIO.C:65
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1004
fileName path() const
Return path.
Definition: Time.H:269
Foam::argList args(argc, argv)
const Field< PointType > & localPoints() const
Return pointField of points in patch.
static void addOptions(const bool constant=true, const bool withZero=false)
Add the options handled by timeSelector to argList::validOptions.
Definition: timeSelector.C:114
Namespace for OpenFOAM.