OBJedgeFormat.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-2013 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 \*---------------------------------------------------------------------------*/
25 
26 #include "OBJedgeFormat.H"
27 #include "clock.H"
28 #include "IFstream.H"
29 #include "IStringStream.H"
30 #include "Ostream.H"
31 #include "OFstream.H"
32 #include "ListOps.H"
33 
34 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
35 
36 Foam::fileFormats::OBJedgeFormat::OBJedgeFormat
37 (
38  const fileName& filename
39 )
40 {
41  read(filename);
42 }
43 
44 
45 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
46 
47 void Foam::fileFormats::OBJedgeFormat::readVertices
48 (
49  const string& line,
50  string::size_type& endNum,
51  DynamicList<label>& dynVertices
52 )
53 {
54  dynVertices.clear();
55  while (true)
56  {
57  string::size_type startNum =
58  line.find_first_not_of(' ', endNum);
59 
60  if (startNum == string::npos)
61  {
62  break;
63  }
64 
65  endNum = line.find(' ', startNum);
66 
67  string vertexSpec;
68  if (endNum != string::npos)
69  {
70  vertexSpec = line.substr(startNum, endNum-startNum);
71  }
72  else
73  {
74  vertexSpec = line.substr(startNum, line.size() - startNum);
75  }
76 
77  string::size_type slashPos = vertexSpec.find('/');
78 
79  label vertI = 0;
80  if (slashPos != string::npos)
81  {
82  IStringStream intStream(vertexSpec.substr(0, slashPos));
83 
84  intStream >> vertI;
85  }
86  else
87  {
88  IStringStream intStream(vertexSpec);
89 
90  intStream >> vertI;
91  }
92  dynVertices.append(vertI - 1);
93  }
94 }
95 
96 
98 {
99  clear();
100 
101  IFstream is(filename);
102  if (!is.good())
103  {
105  (
106  "fileFormats::OBJedgeFormat::read(const fileName&)"
107  )
108  << "Cannot read file " << filename
109  << exit(FatalError);
110  }
111 
112  DynamicList<point> dynPoints;
113  DynamicList<edge> dynEdges;
114  DynamicList<label> dynUsedPoints;
115 
116  DynamicList<label> dynVertices;
117 
118  while (is.good())
119  {
120  string line = this->getLineNoComment(is);
121 
122  // handle continuations
123  if (line[line.size()-1] == '\\')
124  {
125  line.substr(0, line.size()-1);
126  line += this->getLineNoComment(is);
127  }
128 
129  // Read first word
130  IStringStream lineStream(line);
131  word cmd;
132  lineStream >> cmd;
133 
134  if (cmd == "v")
135  {
136  scalar x, y, z;
137  lineStream >> x >> y >> z;
138 
139  dynPoints.append(point(x, y, z));
140  dynUsedPoints.append(-1);
141  }
142  else if (cmd == "l")
143  {
144  // Assume 'l' is followed by space.
145  string::size_type endNum = 1;
146 
147  readVertices
148  (
149  line,
150  endNum,
151  dynVertices
152  );
153 
154 
155  for (label i = 1; i < dynVertices.size(); i++)
156  {
157  edge edgeRead(dynVertices[i-1], dynVertices[i]);
158 
159  dynUsedPoints[edgeRead[0]] = edgeRead[0];
160  dynUsedPoints[edgeRead[1]] = edgeRead[1];
161 
162  dynEdges.append(edgeRead);
163  }
164  }
165  else if (cmd == "f")
166  {
167  // Support for faces with 2 vertices
168 
169  // Assume 'f' is followed by space.
170  string::size_type endNum = 1;
171 
172  readVertices
173  (
174  line,
175  endNum,
176  dynVertices
177  );
178 
179  if (dynVertices.size() == 2)
180  {
181  for (label i = 1; i < dynVertices.size(); i++)
182  {
183  edge edgeRead(dynVertices[i-1], dynVertices[i]);
184 
185  dynUsedPoints[edgeRead[0]] = edgeRead[0];
186  dynUsedPoints[edgeRead[1]] = edgeRead[1];
187 
188  dynEdges.append(edgeRead);
189  }
190  }
191  }
192  }
193 
194  // cull unused points
195  label nUsed = 0;
196 
197  forAll(dynPoints, pointI)
198  {
199  if (dynUsedPoints[pointI] >= 0)
200  {
201  if (nUsed != pointI)
202  {
203  dynPoints[nUsed] = dynPoints[pointI];
204  dynUsedPoints[pointI] = nUsed; // new position
205  }
206  ++nUsed;
207  }
208  }
209 
210  dynPoints.setSize(nUsed);
211 
212  // transfer to normal lists
213  storedPoints().transfer(dynPoints);
214 
215  // renumber edge vertices
216  if (nUsed != dynUsedPoints.size())
217  {
218  forAll(dynEdges, edgeI)
219  {
220  edge& e = dynEdges[edgeI];
221 
222  e[0] = dynUsedPoints[e[0]];
223  e[1] = dynUsedPoints[e[1]];
224  }
225  }
226  storedEdges().transfer(dynEdges);
227 
228  return true;
229 }
230 
231 
233 (
234  const fileName& filename,
235  const edgeMesh& mesh
236 )
237 {
238  const pointField& pointLst = mesh.points();
239  const edgeList& edgeLst = mesh.edges();
240 
241  OFstream os(filename);
242  if (!os.good())
243  {
245  (
246  "fileFormats::OBJedgeFormat::write"
247  "(const fileName&, const edgeMesh&)"
248  )
249  << "Cannot open file for writing " << filename
250  << exit(FatalError);
251  }
252 
253 
254  os << "# Wavefront OBJ file written " << clock::dateTime().c_str() << nl
255  << "o " << os.name().lessExt().name() << nl
256  << nl
257  << "# points : " << pointLst.size() << nl
258  << "# lines : " << edgeLst.size() << nl;
259 
260  os << nl
261  << "# <points count=\"" << pointLst.size() << "\">" << nl;
262 
263  // Write vertex coords
264  forAll(pointLst, ptI)
265  {
266  const point& p = pointLst[ptI];
267 
268  os << "v " << p.x() << ' ' << p.y() << ' ' << p.z() << nl;
269  }
270 
271  os << "# </points>" << nl
272  << nl
273  << "# <edges count=\"" << edgeLst.size() << "\">" << endl;
274 
275  // Write line connectivity
276  forAll(edgeLst, edgeI)
277  {
278  const edge& e = edgeLst[edgeI];
279 
280  os << "l " << (e[0] + 1) << " " << (e[1] + 1) << nl;
281  }
282  os << "# </edges>" << endl;
283 }
284 
285 
286 // ************************************************************************* //
Input from memory buffer stream.
Definition: IStringStream.H:49
Output to file stream.
Definition: OFstream.H:81
word name() const
Return file name (part beyond last /)
Definition: fileName.C:206
edgeList & storedEdges()
Non-const access to the edges.
Definition: edgeMeshI.H:67
vector point
Point is a vector.
Definition: point.H:41
void transfer(List< T > &)
Transfer the contents of the argument List into this list.
Definition: List.C:390
void clear()
Clear the addressed list, i.e. set the size to zero.
Definition: DynamicListI.H:242
A class for handling words, derived from string.
Definition: word.H:59
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
const pointField & points() const
Return points.
Definition: edgeMeshI.H:39
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:76
DynamicList< T, SizeInc, SizeMult, SizeDiv > & append(const T &)
Append an element at the end of the list.
Definition: DynamicListI.H:310
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
static string getLineNoComment(IFstream &)
Read non-comment line.
Various functions to operate on Lists.
Input from file stream.
Definition: IFstream.H:81
An edge is a list of two point labels. The functionality it provides supports the discretisation on a...
Definition: edge.H:58
static string dateTime()
Return the current wall-clock date/time as a string.
Definition: clock.C:57
graph_traits< Graph >::vertices_size_type size_type
Definition: SloanRenumber.C:73
virtual void clear()
Clear all storage.
Definition: edgeMesh.C:169
static const char nl
Definition: Ostream.H:260
const double e
Elementary charge.
Definition: doubleFloat.H:78
const Cmpt & y() const
Definition: VectorI.H:71
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
volScalarField & p
Definition: createFields.H:51
const fileName & name() const
Return the name of the stream.
Definition: OFstream.H:118
#define forAll(list, i)
Definition: UList.H:421
bool good() const
Return true if next operation might succeed.
Definition: IOstream.H:333
const Cmpt & x() const
Definition: VectorI.H:65
void setSize(const label)
Alter the addressed list size.
Definition: DynamicListI.H:167
const Cmpt & z() const
Definition: VectorI.H:77
static void write(const fileName &, const edgeMesh &)
Write surface mesh components by proxy.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
Points connected by edges.
Definition: edgeMesh.H:69
error FatalError
scalar y
virtual bool read(const fileName &)
Read from file.
Definition: OBJedgeFormat.C:97
A class for handling file names.
Definition: fileName.H:69
fileName lessExt() const
Return file name without extension (part before last .)
Definition: fileName.C:307
bool read(const char *, int32_t &)
Definition: int32IO.C:87
pointField & storedPoints()
Non-const access to global points.
Definition: edgeMeshI.H:61
const edgeList & edges() const
Return edges.
Definition: edgeMeshI.H:45