OBJstream.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2012-2016 OpenFOAM Foundation
9  Copyright (C) 2020-2021 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "OBJstream.H"
30 #include "primitivePatch.H"
31 #include "treeBoundBox.H"
32 
33 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 defineTypeNameAndDebug(OBJstream, 0);
38 }
39 
40 
41 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
42 
43 void Foam::OBJstream::writeAndCheck(const char c)
44 {
45  if (c == '\n')
46  {
47  startOfLine_ = true;
48  }
49  else if (startOfLine_)
50  {
51  startOfLine_ = false;
52  if (c == 'v')
53  {
54  ++nVertices_;
55  }
56  }
57 
59 }
60 
61 
62 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
63 
65 (
66  const fileName& pathname,
67  IOstreamOption streamOpt
68 )
69 :
70  OFstream(pathname, streamOpt),
71  startOfLine_(true),
72  nVertices_(0)
73 {}
74 
75 
76 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
77 
79 {
80  writeAndCheck(c);
81  return *this;
82 }
83 
84 
86 {
87  for (const char* iter = str; *iter; ++iter)
88  {
89  writeAndCheck(*iter);
90  }
91  return *this;
92 }
93 
94 
96 {
97  return writeQuoted(str, false);
98 }
99 
100 
102 {
103  return writeQuoted(str, true);
104 }
105 
106 
108 (
109  const std::string& str,
110  const bool quoted
111 )
112 {
113  if (!quoted)
114  {
115  // Output unquoted, only advance line number on newline
116  for (auto iter = str.cbegin(); iter != str.cend(); ++iter)
117  {
118  writeAndCheck(*iter);
119  }
120  return *this;
121  }
122 
123  // Output with surrounding quotes and backslash escaping
124  OFstream::write(static_cast<char>(token::DQUOTE));
125 
126  unsigned backslash = 0;
127  for (auto iter = str.cbegin(); iter != str.cend(); ++iter)
128  {
129  const char c = *iter;
130 
131  if (c == '\\')
132  {
133  ++backslash;
134  continue; // only output after escaped character is known
135  }
136  else if (c == token::NL)
137  {
138  ++lineNumber_;
139  ++backslash; // backslash escape for newline
140  }
141  else if (c == token::DQUOTE)
142  {
143  ++backslash; // backslash escape for quote
144  }
145 
146  // output all pending backslashes
147  while (backslash)
148  {
149  OFstream::write('\\');
150  --backslash;
151  }
152 
153  writeAndCheck(c);
154  }
155 
156  // silently drop any trailing backslashes
157  // they would otherwise appear like an escaped end-quote
158  OFstream::write(static_cast<char>(token::DQUOTE));
159 
160  return *this;
161 }
162 
163 
165 {
166  write("v ") << pt.x() << ' ' << pt.y() << ' ' << pt.z() << nl;
167  return *this;
168 }
169 
170 
172 {
173  write(pt);
174  OFstream::write("vn ") << n.x() << ' ' << n.y() << ' ' << n.z() << nl;
175  return *this;
176 }
177 
178 
180 {
181  write(points[e[0]]);
182  write(points[e[1]]);
183  write("l ") << nVertices_-1 << ' ' << nVertices_ << nl;
184  return *this;
185 }
186 
187 
189 {
190  write(ln.start());
191  write(ln.end());
192  write("l ") << nVertices_-1 << ' ' << nVertices_ << nl;
193  return *this;
194 }
195 
196 
198 (
199  const linePointRef& ln,
200  const vector& n0,
201  const vector& n1
202 )
203 {
204  write(ln.start(), n0);
205  write(ln.end(), n1);
206  write("l ") << nVertices_-1 << ' ' << nVertices_ << nl;
207  return *this;
208 }
209 
210 
212 (
213  const triPointRef& f,
214  const bool lines
215 )
216 {
217  const label start = nVertices_+1; // 1-offset for obj included here
218  write(f.a());
219  write(f.b());
220  write(f.c());
221  if (lines)
222  {
223  write('l');
224  for (int i = 0; i < 3; i++)
225  {
226  write(' ') << i+start;
227  }
228  write(' ') << start << '\n';
229  }
230  else
231  {
232  write('f');
233  for (int i = 0; i < 3; i++)
234  {
235  write(' ') << i+start;
236  }
237  write('\n');
238  }
239  return *this;
240 }
241 
242 
244 (
245  const face& f,
246  const UList<point>& points,
247  const bool lines
248 )
249 {
250  const label start = nVertices_+1; // 1-offset for obj included here
251  forAll(f, i)
252  {
253  write(points[f[i]]);
254  }
255  if (lines)
256  {
257  write('l');
258  forAll(f, i)
259  {
260  write(' ') << i+start;
261  }
262  write(' ') << start << '\n';
263  }
264  else
265  {
266  write('f');
267  forAll(f, i)
268  {
269  write(' ') << i+start;
270  }
271  write('\n');
272  }
273  return *this;
274 }
275 
276 
278 (
279  const UList<face>& faces,
280  const pointField& points,
281  const bool lines
282 )
283 {
284  SubList<face> allFcs(faces, faces.size());
285 
286  primitivePatch pp(allFcs, points);
287 
288  const pointField& localPoints = pp.localPoints();
289  const faceList& localFaces = pp.localFaces();
290 
291  const label start = nVertices_+1; // 1-offset for obj included here
292 
293  forAll(localPoints, i)
294  {
295  write(localPoints[i]);
296  }
297 
298  if (lines)
299  {
300  const edgeList& edges = pp.edges();
301  forAll(edges, edgeI)
302  {
303  const edge& e = edges[edgeI];
304 
305  write("l ") << e[0]+start << ' ' << e[1]+start << nl;
306  }
307  }
308  else
309  {
310  forAll(localFaces, facei)
311  {
312  const face& f = localFaces[facei];
313  write('f');
314  forAll(f, i)
315  {
316  write(' ') << f[i]+start;
317  }
318  write('\n');
319  }
320  }
321  return *this;
322 }
323 
324 
326 (
327  const UList<edge>& edges,
328  const UList<point>& points,
329  const bool compact
330 )
331 {
332  if (compact)
333  {
334  // Code similar to PrimitivePatch::calcMeshData()
335  // Unsorted version
336 
337  label objPointId = nVertices_+1; // 1-offset for obj included here
338 
339  Map<label> markedPoints(2*edges.size());
340  forAll(edges, edgei)
341  {
342  const edge& e = edges[edgei];
343 
344  if (markedPoints.insert(e[0], objPointId))
345  {
346  write(points[e[0]]);
347  ++objPointId;
348  }
349  if (markedPoints.insert(e[1], objPointId))
350  {
351  write(points[e[1]]);
352  ++objPointId;
353  }
354  }
355 
356  forAll(edges, edgei)
357  {
358  const edge& e = edges[edgei];
359 
360  write("l ")
361  << markedPoints[e[0]] << ' '
362  << markedPoints[e[1]] << nl;
363  }
364  }
365  else
366  {
367  const label start = nVertices_+1; // 1-offset for obj included here
368 
369  forAll(points, i)
370  {
371  write(points[i]);
372  }
373 
374  forAll(edges, edgei)
375  {
376  const edge& e = edges[edgei];
377 
378  write("l ")
379  << e[0]+start << ' ' << e[1]+start << nl;
380  }
381  }
382 
383  return *this;
384 }
385 
386 
388 (
389  const treeBoundBox& bb,
390  const bool lines
391 )
392 {
393  const label start = nVertices_+1; // 1-offset for obj included here
394 
395  pointField points(bb.points());
396  forAll(points, i)
397  {
398  write(points[i]);
399  }
400 
401  if (lines)
402  {
404  {
405  const edge& e = treeBoundBox::edges[edgei];
406 
407  write("l ") << e[0]+start << ' ' << e[1]+start << nl;
408  }
409  }
410  else
411  {
413  {
414  const face& f = treeBoundBox::faces[facei];
415 
416  write('f');
417  forAll(f, i)
418  {
419  write(' ') << f[i]+start;
420  }
421  write('\n');
422  }
423  }
424 
425  return *this;
426 }
427 
428 
429 // ************************************************************************* //
Foam::Vector::x
const Cmpt & x() const
Access to the vector x component.
Definition: VectorI.H:73
Foam::OSstream::write
virtual bool write(const token &tok)
Write token to stream or otherwise handle it.
Definition: OSstream.C:36
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::fileName
A class for handling file names.
Definition: fileName.H:73
Foam::treeBoundBox::points
tmp< pointField > points() const
Vertex coordinates. In octant coding.
Definition: treeBoundBox.C:92
Foam::SubList
A List obtained as a section of another List.
Definition: SubList.H:54
Foam::OBJstream::writeQuoted
virtual Ostream & writeQuoted(const std::string &str, const bool quoted=true)
Write std::string surrounded by quotes.
Definition: OBJstream.C:108
Foam::treeBoundBox
Standard boundBox with extra functionality for use in octree.
Definition: treeBoundBox.H:86
Foam::edge
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:63
Foam::Map< label >
Foam::treeBoundBox::edges
static const edgeList edges
Edge to point addressing.
Definition: treeBoundBox.H:154
Foam::OBJstream::write
virtual Ostream & write(const char c)
Write character.
Definition: OBJstream.C:78
Foam::Vector::z
const Cmpt & z() const
Access to the vector z component.
Definition: VectorI.H:85
writeQuoted
os writeQuoted(("# "+outputName+"\n"), false)
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::triangle
A triangle primitive used to calculate face normals and swept volumes.
Definition: triangle.H:59
n
label n
Definition: TABSMDCalcMethod2.H:31
Foam::token::DQUOTE
Double quote.
Definition: token.H:135
Foam::Field< vector >
treeBoundBox.H
Foam::treeBoundBox::faces
static const faceList faces
Face to point addressing.
Definition: treeBoundBox.H:151
Foam::IOstreamOption
The IOstreamOption is a simple container for options an IOstream can normally have.
Definition: IOstreamOption.H:63
Foam::token::NL
Newline [isspace].
Definition: token.H:124
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::OFstream
Output to file stream, using an OSstream.
Definition: OFstream.H:53
Foam::Vector::y
const Cmpt & y() const
Access to the vector y component.
Definition: VectorI.H:79
Foam::nl
constexpr char nl
Definition: Ostream.H:404
f
labelList f(nPoints)
Foam::Vector< scalar >
Foam::List< face >
Foam::UList
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:103
points
const pointField & points
Definition: gmvOutputHeader.H:1
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
Foam::line
A line primitive.
Definition: line.H:53
Foam::vtk::write
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
Definition: foamVtkOutputTemplates.C:36
Foam::constant::universal::c
const dimensionedScalar c
Speed of light in a vacuum.
Foam::ln
bool ln(const fileName &src, const fileName &dst)
Create a softlink. dst should not exist. Returns true if successful.
Definition: MSwindows.C:925
Foam::face
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:72
Foam::UList::size
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::defineTypeNameAndDebug
defineTypeNameAndDebug(combustionModel, 0)
primitivePatch.H
OBJstream.H
Foam::PrimitivePatch
A list of faces which address into the list of points.
Definition: PrimitivePatch.H:79
Foam::OBJstream::OBJstream
OBJstream(const fileName &pathname, IOstreamOption streamOpt=IOstreamOption())
Construct from pathname.
Definition: OBJstream.C:65