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 (C) 2011 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"
  67 #include "transformGeometricField.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
 104     PtrList<volScalarField> vsFlds;
 105     readAndRotateFields(vsFlds, mesh, T, objects);
 106
 107     PtrList<volVectorField> vvFlds;
 108     readAndRotateFields(vvFlds, mesh, T, objects);
 109
 110     PtrList<volSphericalTensorField> vstFlds;
 111     readAndRotateFields(vstFlds, mesh, T, objects);
 112
 113     PtrList<volSymmTensorField> vsymtFlds;
 114     readAndRotateFields(vsymtFlds, mesh, T, objects);
 115
 116     PtrList<volTensorField> vtFlds;
 117     readAndRotateFields(vtFlds, mesh, T, objects);
 118
 119     // Read surface fields.
 120
 121     PtrList<surfaceScalarField> ssFlds;
 122     readAndRotateFields(ssFlds, mesh, T, objects);
 123
 124     PtrList<surfaceVectorField> svFlds;
 125     readAndRotateFields(svFlds, mesh, T, objects);
 126
 127     PtrList<surfaceSphericalTensorField> sstFlds;
 128     readAndRotateFields(sstFlds, mesh, T, objects);
 129
 130     PtrList<surfaceSymmTensorField> ssymtFlds;
 131     readAndRotateFields(ssymtFlds, mesh, T, objects);
 132
 133     PtrList<surfaceTensorField> stFlds;
 134     readAndRotateFields(stFlds, mesh, T, objects);
 135
 136     mesh.write();
 137 }
 138
 139
 140 //  Main program:
 141
 142 int main(int argc, char *argv[])
 143 {
 144     argList::addOption
 145     (
 146         "translate",
 147         "vector",
 148         "translate by the specified <vector> - eg, '(1 0 0)'"
 149     );
 150     argList::addOption
 151     (
 152         "rotate",
 153         "(vectorA vectorB)",
 154         "transform in terms of a rotation between <vectorA> and <vectorB> "
 155         "- eg, '( (1 0 0) (0 0 1) )'"
 156     );
 157     argList::addOption
 158     (
 159         "rollPitchYaw",
 160         "vector",
 161         "transform in terms of '(roll pitch yaw)' in degrees"
 162     );
 163     argList::addOption
 164     (
 165         "yawPitchRoll",
 166         "vector",
 167         "transform in terms of '(yaw pitch roll)' in degrees"
 168     );
 169     argList::addBoolOption
 170     (
 171         "rotateFields",
 172         "read and transform vector and tensor fields too"
 173     );
 174     argList::addOption
 175     (
 176         "scale",
 177         "vector",
 178         "scale by the specified amount - eg, '(0.001 0.001 0.001)' for a "
 179         "uniform [mm] to [m] scaling"
 180     );
 181
 182 #   include "addRegionOption.H"
 183 #   include "setRootCase.H"
 184 #   include "createTime.H"
 185
 186     word regionName = polyMesh::defaultRegion;
 187     fileName meshDir;
 188
 189     if (args.optionReadIfPresent("region", regionName))
 190     {
 191         meshDir = regionName/polyMesh::meshSubDir;
 192     }
 193     else
 194     {
 195         meshDir = polyMesh::meshSubDir;
 196     }
 197
 198     pointIOField points
 199     (
 200         IOobject
 201         (
 202             "points",
 203             runTime.findInstance(meshDir, "points"),
 204             meshDir,
 205             runTime,
 206             IOobject::MUST_READ,
 207             IOobject::NO_WRITE,
 208             false
 209         )
 210     );
 211
 212     const bool doRotateFields = args.optionFound("rotateFields");
 213
 214     // this is not actually stringent enough:
 215     if (args.options().empty())
 216     {
 217         FatalErrorIn(args.executable())
 218             << "No options supplied, please use one or more of "
 219                "-translate, -rotate or -scale options."
 220             << exit(FatalError);
 221     }
 222
 223     vector v;
 224     if (args.optionReadIfPresent("translate", v))
 225     {
 226         Info<< "Translating points by " << v << endl;
 227
 228         points += v;
 229     }
 230
 231     if (args.optionFound("rotate"))
 232     {
 233         Pair<vector> n1n2
 234         (
 235             args.optionLookup("rotate")()
 236         );
 237         n1n2[0] /= mag(n1n2[0]);
 238         n1n2[1] /= mag(n1n2[1]);
 239         tensor T = rotationTensor(n1n2[0], n1n2[1]);
 240
 241         Info<< "Rotating points by " << T << endl;
 242
 243         points = transform(T, points);
 244
 245         if (doRotateFields)
 246         {
 247             rotateFields(args, runTime, T);
 248         }
 249     }
 250     else if (args.optionReadIfPresent("rollPitchYaw", v))
 251     {
 252         Info<< "Rotating points by" << nl
 253             << "    roll  " << v.x() << nl
 254             << "    pitch " << v.y() << nl
 255             << "    yaw   " << v.z() << nl;
 256
 257         // Convert to radians
 258         v *= pi/180.0;
 259
 260         quaternion R(v.x(), v.y(), v.z());
 261
 262         Info<< "Rotating points by quaternion " << R << endl;
 263         points = transform(R, points);
 264
 265         if (doRotateFields)
 266         {
 267             rotateFields(args, runTime, R.R());
 268         }
 269     }
 270     else if (args.optionReadIfPresent("yawPitchRoll", v))
 271     {
 272         Info<< "Rotating points by" << nl
 273             << "    yaw   " << v.x() << nl
 274             << "    pitch " << v.y() << nl
 275             << "    roll  " << v.z() << nl;
 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 (doRotateFields)
 292         {
 293             rotateFields(args, runTime, R.R());
 294         }
 295     }
 296
 297     if (args.optionReadIfPresent("scale", v))
 298     {
 299         Info<< "Scaling points by " << v << endl;
 300
 301         points.replace(vector::X, v.x()*points.component(vector::X));
 302         points.replace(vector::Y, v.y()*points.component(vector::Y));
 303         points.replace(vector::Z, v.z()*points.component(vector::Z));
 304     }
 305
 306     // Set the precision of the points data to 10
 307     IOstream::defaultPrecision(10);
 308
 309     Info<< "Writing points into directory " << points.path() << nl << endl;
 310     points.write();
 311
 312     return 0;
 313 }
 314
 315
 316 // ************************************************************************* //