FreeFOAM The Cross-Platform CFD Toolkit
transformPoints.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  transformPoints
26 
27 Description
28  Transforms the mesh points in the polyMesh directory according to the
29  translate, rotate and scale options.
30 
31 Usage
32  Options are:
33 
34  -translate vector
35  Translates the points by the given vector,
36 
37  -rotate (vector vector)
38  Rotates the points from the first vector to the second,
39 
40  or -yawPitchRoll (yawdegrees pitchdegrees rolldegrees)
41  or -rollPitchYaw (rolldegrees pitchdegrees yawdegrees)
42 
43  -scale vector
44  Scales the points by the given vector.
45 
46  The any or all of the three options may be specified and are processed
47  in the above order.
48 
49  With -rotateFields (in combination with -rotate/yawPitchRoll/rollPitchYaw)
50  it will also read & transform vector & tensor fields.
51 
52 Usage
53 
54  - transformPoints [OPTIONS]
55 
56  @param -translate <vector>\n
57  Translate by given vector.
58 
59  @param -rotateFields \n
60  Also rotate vector/tensor fields.
61 
62  @param -rotate <(vector vector)>\n
63  Rotate from first to second direction.
64 
65  @param -scale <vector (sx sy sz)>\n
66  Scaling factors in x, y and z direction.
67 
68  @param -case <dir>\n
69  Case directory.
70 
71  @param -parallel \n
72  Run in parallel.
73 
74  @param -help \n
75  Display help message.
76 
77  @param -doc \n
78  Display Doxygen API documentation page for this application.
79 
80  @param -srcDoc \n
81  Display Doxygen source documentation page for this application.
82 
83 Note
84  - yaw: rotation about z
85  - pitch: rotation about y
86  - roll: rotation about x
87 \*---------------------------------------------------------------------------*/
88 
89 #include <OpenFOAM/argList.H>
90 #include <OpenFOAM/Time.H>
91 #include <finiteVolume/fvMesh.H>
92 #include <finiteVolume/volFields.H>
94 #include <OpenFOAM/ReadFields.H>
95 #include <OpenFOAM/pointFields.H>
98 #include <OpenFOAM/IStringStream.H>
99 
100 using namespace Foam;
101 using namespace Foam::mathematicalConstant;
102 
103 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
104 
105 template<class GeoField>
106 void readAndRotateFields
107 (
108  PtrList<GeoField>& flds,
109  const fvMesh& mesh,
110  const tensor& T,
111  const IOobjectList& objects
112 )
113 {
114  ReadFields(mesh, objects, flds);
115  forAll(flds, i)
116  {
117  Info<< "Transforming " << flds[i].name() << endl;
118  dimensionedTensor dimT("t", flds[i].dimensions(), T);
119  transform(flds[i], dimT, flds[i]);
120  }
121 }
122 
123 
124 void rotateFields(const argList& args, const Time& runTime, const tensor& T)
125 {
126 # include <OpenFOAM/createNamedMesh.H>
127 
128  // Read objects in time directory
129  IOobjectList objects(mesh, runTime.timeName());
130 
131  // Read vol fields.
132 
134  readAndRotateFields(vsFlds, mesh, T, objects);
135 
137  readAndRotateFields(vvFlds, mesh, T, objects);
138 
140  readAndRotateFields(vstFlds, mesh, T, objects);
141 
142  PtrList<volSymmTensorField> vsymtFlds;
143  readAndRotateFields(vsymtFlds, mesh, T, objects);
144 
146  readAndRotateFields(vtFlds, mesh, T, objects);
147 
148  // Read surface fields.
149 
151  readAndRotateFields(ssFlds, mesh, T, objects);
152 
154  readAndRotateFields(svFlds, mesh, T, objects);
155 
157  readAndRotateFields(sstFlds, mesh, T, objects);
158 
160  readAndRotateFields(ssymtFlds, mesh, T, objects);
161 
163  readAndRotateFields(stFlds, mesh, T, objects);
164 
165  mesh.write();
166 }
167 
168 
169 // Main program:
170 
171 int main(int argc, char *argv[])
172 {
173 # include <OpenFOAM/addRegionOption.H>
174  argList::validOptions.insert("translate", "vector");
175  argList::validOptions.insert("rotate", "(vector vector)");
176  argList::validOptions.insert("rollPitchYaw", "(roll pitch yaw)");
177  argList::validOptions.insert("yawPitchRoll", "(yaw pitch roll)");
178  argList::validOptions.insert("rotateFields", "");
179  argList::validOptions.insert("scale", "vector");
180 
181 # include <OpenFOAM/setRootCase.H>
182 # include <OpenFOAM/createTime.H>
183 
185  fileName meshDir;
186 
187  if (args.optionReadIfPresent("region", regionName))
188  {
189  meshDir = regionName/polyMesh::meshSubDir;
190  }
191  else
192  {
193  meshDir = polyMesh::meshSubDir;
194  }
195 
197  (
198  IOobject
199  (
200  "points",
201  runTime.findInstance(meshDir, "points"),
202  meshDir,
203  runTime,
206  false
207  )
208  );
209 
210 
211  if (args.options().empty())
212  {
213  FatalErrorIn(args.executable())
214  << "No options supplied, please use one or more of "
215  "-translate, -rotate or -scale options."
216  << exit(FatalError);
217  }
218 
219  if (args.optionFound("translate"))
220  {
221  vector transVector(args.optionLookup("translate")());
222 
223  Info<< "Translating points by " << transVector << endl;
224 
225  points += transVector;
226  }
227 
228  if (args.optionFound("rotate"))
229  {
230  Pair<vector> n1n2(args.optionLookup("rotate")());
231  n1n2[0] /= mag(n1n2[0]);
232  n1n2[1] /= mag(n1n2[1]);
233  tensor T = rotationTensor(n1n2[0], n1n2[1]);
234 
235  Info<< "Rotating points by " << T << endl;
236 
237  points = transform(T, points);
238 
239  if (args.optionFound("rotateFields"))
240  {
241  rotateFields(args, runTime, T);
242  }
243  }
244  else if (args.optionFound("rollPitchYaw"))
245  {
246  vector v(args.optionLookup("rollPitchYaw")());
247 
248  Info<< "Rotating points by" << nl
249  << " roll " << v.x() << nl
250  << " pitch " << v.y() << nl
251  << " yaw " << v.z() << endl;
252 
253 
254  // Convert to radians
255  v *= pi/180.0;
256 
257  quaternion R(v.x(), v.y(), v.z());
258 
259  Info<< "Rotating points by quaternion " << R << endl;
260  points = transform(R, points);
261 
262  if (args.optionFound("rotateFields"))
263  {
264  rotateFields(args, runTime, R.R());
265  }
266  }
267  else if (args.optionFound("yawPitchRoll"))
268  {
269  vector v(args.optionLookup("yawPitchRoll")());
270 
271  Info<< "Rotating points by" << nl
272  << " yaw " << v.x() << nl
273  << " pitch " << v.y() << nl
274  << " roll " << v.z() << endl;
275 
276 
277  // Convert to radians
278  v *= pi/180.0;
279 
280  scalar yaw = v.x();
281  scalar pitch = v.y();
282  scalar roll = v.z();
283 
284  quaternion R = quaternion(vector(0, 0, 1), yaw);
285  R *= quaternion(vector(0, 1, 0), pitch);
286  R *= quaternion(vector(1, 0, 0), roll);
287 
288  Info<< "Rotating points by quaternion " << R << endl;
289  points = transform(R, points);
290 
291  if (args.optionFound("rotateFields"))
292  {
293  rotateFields(args, runTime, R.R());
294  }
295  }
296 
297  if (args.optionFound("scale"))
298  {
299  vector scaleVector(args.optionLookup("scale")());
300 
301  Info<< "Scaling points by " << scaleVector << endl;
302 
303  points.replace(vector::X, scaleVector.x()*points.component(vector::X));
304  points.replace(vector::Y, scaleVector.y()*points.component(vector::Y));
305  points.replace(vector::Z, scaleVector.z()*points.component(vector::Z));
306  }
307 
308  // Set the precision of the points data to 10
310 
311  Info << "Writing points into directory " << points.path() << nl << endl;
312  points.write();
313 
314  return 0;
315 }
316 
317 
318 // ************************ vim: set sw=4 sts=4 et: ************************ //