applications/utilities/mesh/manipulation/transformPoints/transformPoints.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright held by original author
   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 the
  13     Free Software Foundation; either version 2 of the License, or (at your
  14     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, write to the Free Software Foundation,
  23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
  24
  25 Application
  26     transformPoints
  27
  28 Description
  29     Transforms the mesh points in the polyMesh directory according to the
  30     translate, rotate and scale options.
  31
  32 Usage
  33     Options are:
  34
  35     -translate vector
  36         Translates the points by the given vector,
  37
  38     -rotate (vector vector)
  39         Rotates the points from the first vector to the second,
  40
  41      or -yawPitchRoll (yawdegrees pitchdegrees rolldegrees)
  42      or -rollPitchYaw (rolldegrees pitchdegrees yawdegrees)
  43      or -rotateAlongVector (vector and angle)
  44
  45     -scale vector
  46         Scales the points by the given vector.
  47
  48     The any or all of the three options may be specified and are processed
  49     in the above order.
  50
  51     With -rotateFields (in combination with -rotate/yawPitchRoll/rollPitchYaw)
  52     it will also read & transform vector & tensor fields.
  53
  54     Note:
  55     yaw (rotation about z)
  56     pitch (rotation about y)
  57     roll (rotation about x)
  58
  59 \*---------------------------------------------------------------------------*/
  60
  61 #include "argList.H"
  62 #include "objectRegistry.H"
  63 #include "Time.H"
  64 #include "fvMesh.H"
  65 #include "volFields.H"
  66 #include "surfaceFields.H"
  67 #include "ReadFields.H"
  68 #include "pointFields.H"
  69 #include "transformField.H"
  70 #include "transformGeometricField.H"
  71 #include "IStringStream.H"
  72 #include "RodriguesRotation.H"
  73
  74 using namespace Foam;
  75 using namespace Foam::mathematicalConstant;
  76
  77 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  78
  79 template<class GeoField>
  80 void readAndRotateFields
  81 (
  82     PtrList<GeoField>& flds,
  83     const fvMesh& mesh,
  84     const tensor& T,
  85     const IOobjectList& objects
  86 )
  87 {
  88     ReadFields(mesh, objects, flds);
  89     forAll(flds, i)
  90     {
  91         Info<< "Transforming " << flds[i].name() << endl;
  92         dimensionedTensor dimT("t", flds[i].dimensions(), T);
  93         transform(flds[i], dimT, flds[i]);
  94     }
  95 }
  96
  97
  98 void rotateFields(const argList& args, const Time& runTime, const tensor& T)
  99 {
 100 #   include "createNamedMesh.H"
 101
 102     // Read objects in time directory
 103     IOobjectList objects(mesh, runTime.timeName());
 104
 105     // Read vol fields.
 106
 107     PtrList<volScalarField> vsFlds;
 108     readAndRotateFields(vsFlds, mesh, T, objects);
 109
 110     PtrList<volVectorField> vvFlds;
 111     readAndRotateFields(vvFlds, mesh, T, objects);
 112
 113     PtrList<volSphericalTensorField> vstFlds;
 114     readAndRotateFields(vstFlds, mesh, T, objects);
 115
 116     PtrList<volSymmTensorField> vsymtFlds;
 117     readAndRotateFields(vsymtFlds, mesh, T, objects);
 118
 119     PtrList<volTensorField> vtFlds;
 120     readAndRotateFields(vtFlds, mesh, T, objects);
 121
 122     // Read surface fields.
 123
 124     PtrList<surfaceScalarField> ssFlds;
 125     readAndRotateFields(ssFlds, mesh, T, objects);
 126
 127     PtrList<surfaceVectorField> svFlds;
 128     readAndRotateFields(svFlds, mesh, T, objects);
 129
 130     PtrList<surfaceSphericalTensorField> sstFlds;
 131     readAndRotateFields(sstFlds, mesh, T, objects);
 132
 133     PtrList<surfaceSymmTensorField> ssymtFlds;
 134     readAndRotateFields(ssymtFlds, mesh, T, objects);
 135
 136     PtrList<surfaceTensorField> stFlds;
 137     readAndRotateFields(stFlds, mesh, T, objects);
 138
 139     mesh.write();
 140 }
 141
 142
 143 //  Main program:
 144
 145 int main(int argc, char *argv[])
 146 {
 147 #   include "addRegionOption.H"
 148     argList::validOptions.insert("translate", "vector");
 149     argList::validOptions.insert("rotate", "(vector vector)");
 150     argList::validOptions.insert("rotateAlongVector", "(vector angleInDegree)");
 151     argList::validOptions.insert("rollPitchYaw", "(roll pitch yaw)");
 152     argList::validOptions.insert("yawPitchRoll", "(yaw pitch roll)");
 153     argList::validOptions.insert("rotateFields", "");
 154     argList::validOptions.insert("scale", "vector");
 155
 156 #   include "setRootCase.H"
 157 #   include "createTime.H"
 158
 159     word regionName = polyMesh::defaultRegion;
 160     fileName meshDir;
 161
 162     if (args.optionReadIfPresent("region", regionName))
 163     {
 164         meshDir = regionName/polyMesh::meshSubDir;
 165     }
 166     else
 167     {
 168         meshDir = polyMesh::meshSubDir;
 169     }
 170
 171     pointIOField points
 172     (
 173         IOobject
 174         (
 175             "points",
 176             runTime.findInstance(meshDir, "points"),
 177             meshDir,
 178             runTime,
 179             IOobject::MUST_READ,
 180             IOobject::NO_WRITE,
 181             false
 182         )
 183     );
 184
 185
 186     if (args.options().empty())
 187     {
 188         FatalErrorIn(args.executable())
 189             << "No options supplied, please use one or more of "
 190                "-translate, -rotate or -scale options."
 191             << exit(FatalError);
 192     }
 193
 194     if (args.optionFound("translate"))
 195     {
 196         vector transVector(args.optionLookup("translate")());
 197
 198         Info<< "Translating points by " << transVector << endl;
 199
 200         points += transVector;
 201     }
 202
 203     if (args.optionFound("rotate"))
 204     {
 205         Pair<vector> n1n2(args.optionLookup("rotate")());
 206         n1n2[0] /= mag(n1n2[0]);
 207         n1n2[1] /= mag(n1n2[1]);
 208         tensor T = rotationTensor(n1n2[0], n1n2[1]);
 209
 210         Info<< "Rotating points by " << T << endl;
 211
 212         points = transform(T, points);
 213
 214         if (args.optionFound("rotateFields"))
 215         {
 216             rotateFields(args, runTime, T);
 217         }
 218     }
 219     else if (args.optionFound("rollPitchYaw"))
 220     {
 221         vector v(args.optionLookup("rollPitchYaw")());
 222
 223         Info<< "Rotating points by" << nl
 224             << "    roll  " << v.x() << nl
 225             << "    pitch " << v.y() << nl
 226             << "    yaw   " << v.z() << endl;
 227
 228
 229         // Convert to radians
 230         v *= pi/180.0;
 231
 232         quaternion R(v.x(), v.y(), v.z());
 233
 234         Info<< "Rotating points by quaternion " << R << endl;
 235         points = transform(R, points);
 236
 237         if (args.optionFound("rotateFields"))
 238         {
 239             rotateFields(args, runTime, R.R());
 240         }
 241     }
 242     else if (args.optionFound("yawPitchRoll"))
 243     {
 244         vector v(args.optionLookup("yawPitchRoll")());
 245
 246         Info<< "Rotating points by" << nl
 247             << "    yaw   " << v.x() << nl
 248             << "    pitch " << v.y() << nl
 249             << "    roll  " << v.z() << endl;
 250
 251
 252         // Convert to radians
 253         v *= pi/180.0;
 254
 255         scalar yaw = v.x();
 256         scalar pitch = v.y();
 257         scalar roll = v.z();
 258
 259         quaternion R = quaternion(vector(0, 0, 1), yaw);
 260         R *= quaternion(vector(0, 1, 0), pitch);
 261         R *= quaternion(vector(1, 0, 0), roll);
 262
 263         Info<< "Rotating points by quaternion " << R << endl;
 264         points = transform(R, points);
 265
 266         if (args.optionFound("rotateFields"))
 267         {
 268             rotateFields(args, runTime, R.R());
 269         }
 270     }
 271     else if (args.optionFound("rotateAlongVector"))
 272     {
 273         vector rotationAxis;
 274         scalar rotationAngle;
 275
 276         args.optionLookup("rotateAlongVector")()
 277             >> rotationAxis
 278             >> rotationAngle;
 279
 280         tensor T = RodriguesRotation(rotationAxis, rotationAngle);
 281
 282         Info << "Rotating points by " << T << endl;
 283
 284         points = transform(T, points);
 285
 286         if (args.options().found("rotateFields"))
 287         {
 288             rotateFields(args, runTime, T);
 289         }
 290     }
 291
 292
 293
 294     if (args.optionFound("scale"))
 295     {
 296         vector scaleVector(args.optionLookup("scale")());
 297
 298         Info<< "Scaling points by " << scaleVector << endl;
 299
 300         points.replace(vector::X, scaleVector.x()*points.component(vector::X));
 301         points.replace(vector::Y, scaleVector.y()*points.component(vector::Y));
 302         points.replace(vector::Z, scaleVector.z()*points.component(vector::Z));
 303     }
 304
 305     // Set the precision of the points data to 10
 306     IOstream::defaultPrecision(10);
 307
 308     Info << "Writing points into directory " << points.path() << nl << endl;
 309     points.write();
 310
 311     return 0;
 312 }
 313
 314
 315 // ************************************************************************* //