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) 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  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  Note:
53  yaw (rotation about z)
54  pitch (rotation about y)
55  roll (rotation about x)
56 
57 \*---------------------------------------------------------------------------*/
58 
59 #include "argList.H"
60 #include "Time.H"
61 #include "fvMesh.H"
62 #include "volFields.H"
63 #include "surfaceFields.H"
64 #include "ReadFields.H"
65 #include "pointFields.H"
66 #include "transformField.H"
68 #include "IStringStream.H"
69 #include "mathematicalConstants.H"
70 
71 using namespace Foam;
72 using namespace Foam::constant::mathematical;
73 
74 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
75 
76 template<class GeoField>
77 void readAndRotateFields
78 (
79  PtrList<GeoField>& flds,
80  const fvMesh& mesh,
81  const tensor& T,
82  const IOobjectList& objects
83 )
84 {
85  ReadFields(mesh, objects, flds);
86  forAll(flds, i)
87  {
88  Info<< "Transforming " << flds[i].name() << endl;
89  dimensionedTensor dimT("t", flds[i].dimensions(), T);
90  transform(flds[i], dimT, flds[i]);
91  }
92 }
93 
94 
95 void rotateFields(const argList& args, const Time& runTime, const tensor& T)
96 {
97  #include "createNamedMesh.H"
98 
99  // Read objects in time directory
100  IOobjectList objects(mesh, runTime.timeName());
101 
102  // Read vol fields.
103 
105  readAndRotateFields(vsFlds, mesh, T, objects);
106 
108  readAndRotateFields(vvFlds, mesh, T, objects);
109 
111  readAndRotateFields(vstFlds, mesh, T, objects);
112 
113  PtrList<volSymmTensorField> vsymtFlds;
114  readAndRotateFields(vsymtFlds, mesh, T, objects);
115 
117  readAndRotateFields(vtFlds, mesh, T, objects);
118 
119  // Read surface fields.
120 
122  readAndRotateFields(ssFlds, mesh, T, objects);
123 
125  readAndRotateFields(svFlds, mesh, T, objects);
126 
128  readAndRotateFields(sstFlds, mesh, T, objects);
129 
131  readAndRotateFields(ssymtFlds, mesh, T, objects);
132 
134  readAndRotateFields(stFlds, mesh, T, objects);
135 
136  mesh.write();
137 }
138 
139 
140 
141 int main(int argc, char *argv[])
142 {
144  (
145  "translate",
146  "vector",
147  "translate by the specified <vector> - eg, '(1 0 0)'"
148  );
150  (
151  "rotate",
152  "(vectorA vectorB)",
153  "transform in terms of a rotation between <vectorA> and <vectorB> "
154  "- eg, '( (1 0 0) (0 0 1) )'"
155  );
157  (
158  "rollPitchYaw",
159  "vector",
160  "transform in terms of '(roll pitch yaw)' in degrees"
161  );
163  (
164  "yawPitchRoll",
165  "vector",
166  "transform in terms of '(yaw pitch roll)' in degrees"
167  );
169  (
170  "rotateFields",
171  "read and transform vector and tensor fields too"
172  );
174  (
175  "scale",
176  "vector",
177  "scale by the specified amount - eg, '(0.001 0.001 0.001)' for a "
178  "uniform [mm] to [m] scaling"
179  );
180 
181  #include "addRegionOption.H"
182  #include "setRootCase.H"
183  #include "createTime.H"
184 
185  word regionName = polyMesh::defaultRegion;
186  fileName meshDir;
187 
188  if (args.optionReadIfPresent("region", regionName))
189  {
190  meshDir = regionName/polyMesh::meshSubDir;
191  }
192  else
193  {
194  meshDir = polyMesh::meshSubDir;
195  }
196 
198  (
199  IOobject
200  (
201  "points",
202  runTime.findInstance(meshDir, "points"),
203  meshDir,
204  runTime,
207  false
208  )
209  );
210 
211  const bool doRotateFields = args.optionFound("rotateFields");
212 
213  // this is not actually stringent enough:
214  if (args.options().empty())
215  {
217  << "No options supplied, please use one or more of "
218  "-translate, -rotate or -scale options."
219  << exit(FatalError);
220  }
221 
222  vector v;
223  if (args.optionReadIfPresent("translate", v))
224  {
225  Info<< "Translating points by " << v << endl;
226 
227  points += v;
228  }
229 
230  if (args.optionFound("rotate"))
231  {
232  Pair<vector> n1n2
233  (
234  args.optionLookup("rotate")()
235  );
236  n1n2[0] /= mag(n1n2[0]);
237  n1n2[1] /= mag(n1n2[1]);
238  tensor T = rotationTensor(n1n2[0], n1n2[1]);
239 
240  Info<< "Rotating points by " << T << endl;
241 
242  points = transform(T, points);
243 
244  if (doRotateFields)
245  {
246  rotateFields(args, runTime, T);
247  }
248  }
249  else if (args.optionReadIfPresent("rollPitchYaw", v))
250  {
251  Info<< "Rotating points by" << nl
252  << " roll " << v.x() << nl
253  << " pitch " << v.y() << nl
254  << " yaw " << v.z() << nl;
255 
256  // Convert to radians
257  v *= pi/180.0;
258 
259  quaternion R(quaternion::rotationSequence::XYZ, v);
260 
261  Info<< "Rotating points by quaternion " << R << endl;
262  points = transform(R, points);
263 
264  if (doRotateFields)
265  {
266  rotateFields(args, runTime, R.R());
267  }
268  }
269  else if (args.optionReadIfPresent("yawPitchRoll", v))
270  {
271  Info<< "Rotating points by" << nl
272  << " yaw " << v.x() << nl
273  << " pitch " << v.y() << nl
274  << " roll " << v.z() << nl;
275 
276  // Convert to radians
277  v *= pi/180.0;
278 
279  scalar yaw = v.x();
280  scalar pitch = v.y();
281  scalar roll = v.z();
282 
283  quaternion R = quaternion(vector(0, 0, 1), yaw);
284  R *= quaternion(vector(0, 1, 0), pitch);
285  R *= quaternion(vector(1, 0, 0), roll);
286 
287  Info<< "Rotating points by quaternion " << R << endl;
288  points = transform(R, points);
289 
290  if (doRotateFields)
291  {
292  rotateFields(args, runTime, R.R());
293  }
294  }
295 
296  if (args.optionReadIfPresent("scale", v))
297  {
298  Info<< "Scaling points by " << v << endl;
299 
303  }
304 
305  // Set the precision of the points data to 10
307 
308  Info<< "Writing points into directory " << points.path() << nl << endl;
309  points.write();
310 
311  Info<< "End\n" << endl;
312 
313  return 0;
314 }
315 
316 
317 // ************************************************************************* //
Foam::surfaceFields.
const Cmpt & z() const
Definition: VectorI.H:87
wordList ReadFields(const Mesh &mesh, const IOobjectList &objects, PtrList< GeoField > &fields, const bool syncPar=true)
Read all fields of the specified type.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
A class for handling file names.
Definition: fileName.H:69
List of IOobjects with searching and retrieving facilities.
Definition: IOobjectList.H:50
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
const Cmpt & x() const
Definition: VectorI.H:75
error FatalError
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
Spatial transformation functions for FieldFields.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:309
static unsigned int defaultPrecision()
Return the default precision.
Definition: IOstream.H:461
bool optionReadIfPresent(const word &opt, T &) const
Read a value from the named option if present.
Definition: argListI.H:198
tensor R() const
The rotation tensor corresponding the quaternion.
Definition: quaternionI.H:323
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:306
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
Field reading functions for post-processing utilities.
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
Generic dimensioned Type class.
tensor rotationTensor(const vector &n1, const vector &n2)
Rotational transformation tensor from vector n1 to n2.
Definition: transform.H:47
void replace(const direction, const UList< cmptType > &)
Replace a component field of the field.
Definition: Field.C:662
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
const Foam::HashTable< string > & options() const
Return options.
Definition: argListI.H:90
Spatial transformation functions for primitive fields.
word findInstance(const fileName &dir, const word &name=word::null, const IOobject::readOption rOpt=IOobject::MUST_READ, const word &stopInstance=word::null) const
Return the location of "dir" containing the file "name".
Definition: findInstance.C:38
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
const pointField & points
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:108
mathematical constants.
A class for handling words, derived from string.
Definition: word.H:59
Extract command arguments and options from the supplied argc and argv parameters. ...
Definition: argList.H:102
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
Quaternion class used to perform rotations in 3D space.
Definition: quaternion.H:60
static const char nl
Definition: Ostream.H:262
tmp< Field< cmptType > > component(const direction) const
Return a component field of the field.
Definition: Field.C:650
#define R(A, B, C, D, E, F, K, M)
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:62
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
IStringStream optionLookup(const word &opt) const
Return an IStringStream from the named option.
Definition: argListI.H:114
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
static void addBoolOption(const word &opt, const string &usage="")
Add to a bool option to validOptions with usage information.
Definition: argList.C:83
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
virtual bool write() const
Write mesh using IO settings from time.
Definition: fvMesh.C:870
Namespace for OpenFOAM.
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:465