FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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 Usage
41 
42  - writeMeshObj [OPTIONS]
43 
44  @param -face <face label>\n
45  Export face with specified label.
46 
47  @param -patchFaces \n
48  Also export patch faces.
49 
50  @param -cell <cell label>\n
51  Export cell with specified label.
52 
53  @param -point <point label>\n
54  Export point with specified label.
55 
56  @param -cellSet <cellSet name>\n
57  Only apply to named cell set.
58 
59  @param -faceSet <faceSet name>\n
60  Only apply to named face set.
61 
62  @param -noZero \n
63  Ignore timestep 0.
64 
65  @param -constant \n
66  Include the constant directory.
67 
68  @param -time <time>\n
69  Apply only to specific time.
70 
71  @param -latestTime \n
72  Only apply to latest time step.
73 
74  @param -case <dir>\n
75  Case directory.
76 
77  @param -parallel \n
78  Run in parallel.
79 
80  @param -help \n
81  Display help message.
82 
83  @param -doc \n
84  Display Doxygen API documentation page for this application.
85 
86  @param -srcDoc \n
87  Display Doxygen source documentation page for this application.
88 
89 \*---------------------------------------------------------------------------*/
90 
91 #include <OpenFOAM/argList.H>
92 #include <OpenFOAM/timeSelector.H>
93 #include <OpenFOAM/Time.H>
94 #include <OpenFOAM/polyMesh.H>
95 #include <OpenFOAM/OFstream.H>
96 #include <meshTools/meshTools.H>
97 #include <meshTools/cellSet.H>
98 #include <meshTools/faceSet.H>
99 
100 using namespace Foam;
101 
102 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
103 
104 void writeOBJ(const point& pt, Ostream& os)
105 {
106  os << "v " << pt.x() << ' ' << pt.y() << ' ' << pt.z() << nl;
107 }
108 
109 // All edges of mesh
110 void writePoints(const polyMesh& mesh, const fileName& timeName)
111 {
112  label vertI = 0;
113 
114  fileName pointFile(mesh.time().path()/"meshPoints_" + timeName + ".obj");
115 
116  Info << "Writing mesh points and edges to " << pointFile << endl;
117 
118  OFstream pointStream(pointFile);
119 
120  forAll(mesh.points(), pointI)
121  {
122  writeOBJ(mesh.points()[pointI], pointStream);
123  vertI++;
124  }
125 
126  forAll(mesh.edges(), edgeI)
127  {
128  const edge& e = mesh.edges()[edgeI];
129 
130  pointStream << "l " << e.start() + 1 << ' ' << e.end() + 1 << nl;
131  }
132 }
133 
134 
135 // Edges for subset of cells
136 void writePoints
137 (
138  const polyMesh& mesh,
139  const labelList& cellLabels,
140  const fileName& timeName
141 )
142 {
143  fileName fName(mesh.time().path()/"meshPoints_" + timeName + ".obj");
144 
145  Info << "Writing mesh points and edges to " << fName << endl;
146 
147  OFstream str(fName);
148 
149  // OBJ file vertex
150  label vertI = 0;
151 
152  // From point to OBJ file vertex
153  Map<label> pointToObj(6*cellLabels.size());
154 
155  forAll(cellLabels, i)
156  {
157  label cellI = cellLabels[i];
158 
159  const labelList& cEdges = mesh.cellEdges()[cellI];
160 
161  forAll(cEdges, cEdgeI)
162  {
163  const edge& e = mesh.edges()[cEdges[cEdgeI]];
164 
165  label v0;
166 
167  Map<label>::iterator e0Fnd = pointToObj.find(e[0]);
168 
169  if (e0Fnd == pointToObj.end())
170  {
171  meshTools::writeOBJ(str, mesh.points()[e[0]]);
172  v0 = vertI++;
173  pointToObj.insert(e[0], v0);
174  }
175  else
176  {
177  v0 = e0Fnd();
178  }
179 
180  label v1;
181 
182  Map<label>::iterator e1Fnd = pointToObj.find(e[1]);
183 
184  if (e1Fnd == pointToObj.end())
185  {
186  meshTools::writeOBJ(str, mesh.points()[e[1]]);
187  v1 = vertI++;
188  pointToObj.insert(e[1], v1);
189  }
190  else
191  {
192  v1 = e1Fnd();
193  }
194 
195 
196  str << "l " << v0+1 << ' ' << v1+1 << nl;
197  }
198  }
199 }
200 
201 
202 // Edges of single cell
203 void writePoints
204 (
205  const polyMesh& mesh,
206  const label cellI,
207  const fileName& timeName
208 )
209 {
210  fileName fName
211  (
212  mesh.time().path()
213  / "meshPoints_" + timeName + '_' + name(cellI) + ".obj"
214  );
215 
216  Info << "Writing mesh points and edges to " << fName << endl;
217 
218  OFstream pointStream(fName);
219 
220  const cell& cFaces = mesh.cells()[cellI];
221 
222  meshTools::writeOBJ(pointStream, mesh.faces(), mesh.points(), cFaces);
223 }
224 
225 
226 
227 // All face centres
228 void writeFaceCentres(const polyMesh& mesh,const fileName& timeName)
229 {
230  fileName faceFile
231  (
232  mesh.time().path()
233  / "meshFaceCentres_" + timeName + ".obj"
234  );
235 
236  Info << "Writing mesh face centres to " << faceFile << endl;
237 
238  OFstream faceStream(faceFile);
239 
240  forAll(mesh.faceCentres(), faceI)
241  {
242  writeOBJ(mesh.faceCentres()[faceI], faceStream);
243  }
244 }
245 
246 
247 void writeCellCentres(const polyMesh& mesh, const fileName& timeName)
248 {
249  fileName cellFile
250  (
251  mesh.time().path()/"meshCellCentres_" + timeName + ".obj"
252  );
253 
254  Info << "Writing mesh cell centres to " << cellFile << endl;
255 
256  OFstream cellStream(cellFile);
257 
258  forAll(mesh.cellCentres(), cellI)
259  {
260  writeOBJ(mesh.cellCentres()[cellI], cellStream);
261  }
262 }
263 
264 
265 void writePatchCentres
266 (
267  const polyMesh& mesh,
268  const fileName& timeName
269 )
270 {
271  const polyBoundaryMesh& patches = mesh.boundaryMesh();
272 
273  forAll(patches, patchI)
274  {
275  const polyPatch& pp = patches[patchI];
276 
277  fileName faceFile
278  (
279  mesh.time().path()/"patch_" + pp.name() + '_' + timeName + ".obj"
280  );
281 
282  Info << "Writing patch face centres to " << faceFile << endl;
283 
284  OFstream patchFaceStream(faceFile);
285 
286  forAll(pp.faceCentres(), faceI)
287  {
288  writeOBJ(pp.faceCentres()[faceI], patchFaceStream);
289  }
290  }
291 }
292 
293 
294 void writePatchFaces
295 (
296  const polyMesh& mesh,
297  const fileName& timeName
298 )
299 {
300  const polyBoundaryMesh& patches = mesh.boundaryMesh();
301 
302  forAll(patches, patchI)
303  {
304  const polyPatch& pp = patches[patchI];
305 
306  fileName faceFile
307  (
308  mesh.time().path()
309  / "patchFaces_" + pp.name() + '_' + timeName + ".obj"
310  );
311 
312  Info << "Writing patch faces to " << faceFile << endl;
313 
314  OFstream patchFaceStream(faceFile);
315 
316  forAll(pp.localPoints(), pointI)
317  {
318  writeOBJ(pp.localPoints()[pointI], patchFaceStream);
319  }
320 
321  forAll(pp.localFaces(), faceI)
322  {
323  const face& f = pp.localFaces()[faceI];
324 
325  patchFaceStream<< 'f';
326 
327  forAll(f, fp)
328  {
329  patchFaceStream << ' ' << f[fp]+1;
330  }
331  patchFaceStream << nl;
332  }
333  }
334 }
335 
336 
337 void writePatchBoundaryEdges
338 (
339  const polyMesh& mesh,
340  const fileName& timeName
341 )
342 {
343  const polyBoundaryMesh& patches = mesh.boundaryMesh();
344 
345  forAll(patches, patchI)
346  {
347  const polyPatch& pp = patches[patchI];
348 
349  fileName edgeFile
350  (
351  mesh.time().path()
352  / "patchEdges_" + pp.name() + '_' + timeName + ".obj"
353  );
354 
355  Info << "Writing patch edges to " << edgeFile << endl;
356 
357  OFstream patchEdgeStream(edgeFile);
358 
359  forAll(pp.localPoints(), pointI)
360  {
361  writeOBJ(pp.localPoints()[pointI], patchEdgeStream);
362  }
363 
364  for (label edgeI = pp.nInternalEdges(); edgeI < pp.nEdges(); edgeI++)
365  {
366  if (pp.edgeFaces()[edgeI].size() == 1)
367  {
368  const edge& e = pp.edges()[edgeI];
369 
370  patchEdgeStream<< "l " << e[0]+1 << ' ' << e[1]+1 << nl;
371  }
372  }
373  }
374 }
375 
376 
377 void writePointCells
378 (
379  const polyMesh& mesh,
380  const label pointI,
381  const fileName& timeName
382 )
383 {
384  const labelList& pCells = mesh.pointCells()[pointI];
385 
386  labelHashSet allEdges(6*pCells.size());
387 
388  forAll(pCells, i)
389  {
390  const labelList& cEdges = mesh.cellEdges()[pCells[i]];
391 
392  forAll(cEdges, i)
393  {
394  allEdges.insert(cEdges[i]);
395  }
396  }
397 
398 
399  fileName pFile
400  (
401  mesh.time().path()
402  / "pointEdges_" + timeName + '_' + name(pointI) + ".obj"
403  );
404 
405  Info << "Writing pointEdges to " << pFile << endl;
406 
407  OFstream pointStream(pFile);
408 
409  label vertI = 0;
410 
411  for
412  (
413  labelHashSet::const_iterator iter = allEdges.begin();
414  iter != allEdges.end();
415  ++iter
416  )
417  {
418  const edge& e = mesh.edges()[iter.key()];
419 
420  meshTools::writeOBJ(pointStream, mesh.points()[e[0]]); vertI++;
421  meshTools::writeOBJ(pointStream, mesh.points()[e[1]]); vertI++;
422  pointStream<< "l " << vertI-1 << ' ' << vertI << nl;
423  }
424 }
425 
426 
427 // Main program:
428 
429 int main(int argc, char *argv[])
430 {
432  argList::validOptions.insert("patchFaces", "");
433  argList::validOptions.insert("patchEdges", "");
434  argList::validOptions.insert("cell", "cellI");
435  argList::validOptions.insert("face", "faceI");
436  argList::validOptions.insert("point", "pointI");
437  argList::validOptions.insert("cellSet", "setName");
438  argList::validOptions.insert("faceSet", "setName");
439 # include <OpenFOAM/addRegionOption.H>
440 
441 # include <OpenFOAM/setRootCase.H>
442 # include <OpenFOAM/createTime.H>
443  runTime.functionObjects().off();
444 
445  bool patchFaces = args.optionFound("patchFaces");
446  bool patchEdges = args.optionFound("patchEdges");
447  bool doCell = args.optionFound("cell");
448  bool doPoint = args.optionFound("point");
449  bool doFace = args.optionFound("face");
450  bool doCellSet = args.optionFound("cellSet");
451  bool doFaceSet = args.optionFound("faceSet");
452 
453 
454  Info<< "Writing mesh objects as .obj files such that the object"
455  << " numbering" << endl
456  << "(for points, faces, cells) is consistent with"
457  << " Foam numbering (starting from 0)." << endl << endl;
458 
460 
462 
463  forAll(timeDirs, timeI)
464  {
465  runTime.setTime(timeDirs[timeI], timeI);
466 
467  Info<< "Time = " << runTime.timeName() << endl;
468 
469  polyMesh::readUpdateState state = mesh.readUpdate();
470 
471  if (!timeI || state != polyMesh::UNCHANGED)
472  {
473  if (patchFaces)
474  {
475  writePatchFaces(mesh, runTime.timeName());
476  }
477  if (patchEdges)
478  {
479  writePatchBoundaryEdges(mesh, runTime.timeName());
480  }
481  if (doCell)
482  {
483  label cellI = args.optionRead<label>("cell");
484 
485  writePoints(mesh, cellI, runTime.timeName());
486  }
487  if (doPoint)
488  {
489  label pointI = args.optionRead<label>("point");
490 
491  writePointCells(mesh, pointI, runTime.timeName());
492  }
493  if (doFace)
494  {
495  label faceI = args.optionRead<label>("face");
496 
497  fileName fName
498  (
499  mesh.time().path()
500  / "meshPoints_"
501  + runTime.timeName()
502  + '_'
503  + name(faceI)
504  + ".obj"
505  );
506 
507  Info<< "Writing mesh points and edges to " << fName << endl;
508 
509  OFstream str(fName);
510 
511  const face& f = mesh.faces()[faceI];
512 
513  meshTools::writeOBJ(str, faceList(1, f), mesh.points());
514  }
515  if (doCellSet)
516  {
517  word setName(args.option("cellSet"));
518 
519  cellSet cells(mesh, setName);
520 
521  Info<< "Read " << cells.size() << " cells from set " << setName
522  << endl;
523 
524  writePoints(mesh, cells.toc(), runTime.timeName());
525 
526  }
527  if (doFaceSet)
528  {
529  word setName(args.option("faceSet"));
530 
531  faceSet faces(mesh, setName);
532 
533  Info<< "Read " << faces.size() << " faces from set " << setName
534  << endl;
535 
536  fileName fName
537  (
538  mesh.time().path()
539  / "meshPoints_"
540  + runTime.timeName()
541  + '_'
542  + setName
543  + ".obj"
544  );
545 
546  Info << "Writing mesh points and edges to " << fName << endl;
547 
548  OFstream str(fName);
549 
551  (
552  str,
553  mesh.faces(),
554  mesh.points(),
555  faces.toc()
556  );
557  }
558  else if
559  (
560  !patchFaces
561  && !patchEdges
562  && !doCell
563  && !doPoint
564  && !doFace
565  && !doCellSet
566  && !doFaceSet
567  )
568  {
569  // points & edges
570  writePoints(mesh, runTime.timeName());
571 
572  // face centres
573  writeFaceCentres(mesh, runTime.timeName());
574 
575  // cell centres
576  writeCellCentres(mesh, runTime.timeName());
577 
578  // Patch face centres
579  writePatchCentres(mesh, runTime.timeName());
580  }
581  }
582  else
583  {
584  Info << "No mesh." << endl;
585  }
586 
587  Info << nl << endl;
588  }
589 
590 
591  Info << "End\n" << endl;
592 
593  return 0;
594 }
595 
596 
597 // ************************ vim: set sw=4 sts=4 et: ************************ //