applications/utilities/mesh/generation/cfMesh/FMSToVTK/FMSToVTK.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | foam-extend: Open Source CFD
   4    \\    /   O peration     | Version:     3.2
   5     \\  /    A nd           | Web:         http://www.foam-extend.org
   6      \\/     M anipulation  | For copyright notice see file Copyright
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of foam-extend.
  10
  11     foam-extend 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 3 of the License, or (at your
  14     option) any later version.
  15
  16     foam-extend is distributed in the hope that it will be useful, but
  17     WITHOUT ANY WARRANTY; without even the implied warranty of
  18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19     General Public License for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 Description
  25     cfMesh utility to convert a surface file to VTK multiblock dataset
  26     format, including the patches, feature edges and surface features.
  27
  28 Author
  29     Ivor Clifford <ivor.clifford@psi.ch>
  30
  31 \*---------------------------------------------------------------------------*/
  32
  33 #include "argList.H"
  34 #include "autoPtr.H"
  35 #include "OFstream.H"
  36 #include "triSurf.H"
  37 #include "triSurfModifier.H"
  38 #include "xmlTag.H"
  39
  40 using namespace Foam;
  41
  42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  43
  44 //- Write the supplied pointList in serial vtkPolyData format
  45 void writePointsToVTK
  46 (
  47     const fileName& fn,
  48     const string& title,
  49     const UList<point>& points
  50 )
  51 {
  52     xmlTag xmlRoot("VTKFile");
  53     xmlRoot.addAttribute("type", "PolyData");
  54
  55     xmlTag& xmlPolyData = xmlRoot.addChild("PolyData");
  56
  57     xmlTag& xmlPiece = xmlPolyData.addChild("Piece");
  58     xmlPiece.addAttribute("NumberOfPoints", points.size());
  59
  60     xmlTag& xmlPoints = xmlPiece.addChild("Points");
  61
  62     xmlTag& xmlPointData = xmlPoints.addChild("DataArray");
  63     xmlPointData.addAttribute("type", "Float32");
  64     xmlPointData.addAttribute("NumberOfComponents", 3);
  65     xmlPointData.addAttribute("format", "ascii");
  66     xmlPointData << points;
  67
  68     OFstream os(fn);
  69     os << xmlRoot << endl;
  70
  71     Info << "Created " << fn << endl;
  72 }
  73
  74
  75 //- Write the supplied addressed pointList in serial vtkPolyData format
  76 void writePointsToVTK
  77 (
  78     const fileName& fn,
  79     const string& title,
  80     const UList<point>& points,
  81     unallocLabelList& addr
  82 )
  83 {
  84     // Create subaddressed points
  85     pointField newPoints(addr.size());
  86
  87     forAll(addr, i)
  88     {
  89         newPoints[i] = points[addr[i]];
  90     }
  91
  92     writePointsToVTK
  93     (
  94         fn,
  95         title,
  96         newPoints
  97     );
  98 }
  99
 100
 101 //- Write the supplied edgeList in serial vtkPolyData format
 102 void writeEdgesToVTK
 103 (
 104     const fileName& fn,
 105     const string& title,
 106     const UList<point>& points,
 107     const LongList<edge>& edges
 108 )
 109 {
 110     labelList connectivity(edges.size());
 111
 112     forAll(edges, edgeI)
 113     {
 114         connectivity[edgeI] = 2*(edgeI+1);
 115     }
 116
 117     xmlTag xmlRoot("VTKFile");
 118     xmlRoot.addAttribute("type", "PolyData");
 119
 120     xmlTag& xmlPolyData = xmlRoot.addChild("PolyData");
 121
 122     xmlTag& xmlPiece = xmlPolyData.addChild("Piece");
 123     xmlPiece.addAttribute("NumberOfPoints", points.size());
 124     xmlPiece.addAttribute("NumberOfLines", edges.size());
 125
 126     xmlTag& xmlPoints = xmlPiece.addChild("Points");
 127
 128     xmlTag& xmlPointData = xmlPoints.addChild("DataArray");
 129     xmlPointData.addAttribute("type", "Float32");
 130     xmlPointData.addAttribute("NumberOfComponents", 3);
 131     xmlPointData.addAttribute("format", "ascii");
 132     xmlPointData << points;
 133
 134     xmlTag& xmlEdges = xmlPiece.addChild("Lines");
 135
 136     xmlTag& xmlEdgeData = xmlEdges.addChild("DataArray");
 137     xmlEdgeData.addAttribute("type", "Int32");
 138     xmlEdgeData.addAttribute("Name", "connectivity");
 139     xmlEdgeData.addAttribute("format", "ascii");
 140     xmlEdgeData << edges;
 141
 142     xmlTag& xmlConnectData = xmlEdges.addChild("DataArray");
 143     xmlConnectData.addAttribute("type", "Int32");
 144     xmlConnectData.addAttribute("Name", "offsets");
 145     xmlConnectData.addAttribute("format", "ascii");
 146     xmlConnectData << connectivity;
 147
 148     OFstream os(fn);
 149     os << xmlRoot << endl;
 150
 151     Info << "Created " << fn << endl;
 152 }
 153
 154 //- Write the supplied edgeList subset in serial vtkPolyData format
 155 void writeEdgesToVTK
 156 (
 157     const fileName& fn,
 158     const string& title,
 159     const UList<point>& points,
 160     const LongList<edge>& edges,
 161     const unallocLabelList& addr
 162 )
 163 {
 164     // Remove unused points and create subaddressed edges
 165     DynamicList<point> newPoints;
 166     labelList newPointAddr(points.size(), -1);
 167     LongList<edge> newEdges(addr.size());
 168
 169     forAll(addr, addrI)
 170     {
 171         label edgeI = addr[addrI];
 172
 173         const edge& curEdge = edges[edgeI];
 174         edge& newEdge = newEdges[addrI];
 175
 176         forAll(curEdge, i)
 177         {
 178             label pointId = curEdge[i];
 179
 180             if (newPointAddr[pointId] == -1)
 181             {
 182                 newPoints.append(points[pointId]);
 183                 newPointAddr[pointId] = newPoints.size()-1;
 184             }
 185
 186             newEdge[i] = newPointAddr[pointId];
 187         }
 188     }
 189
 190     writeEdgesToVTK
 191     (
 192         fn,
 193         title,
 194         newPoints,
 195         newEdges
 196     );
 197 }
 198
 199
 200 //- Write the supplied facet list in serial vtkPolyData format
 201 void writeFacetsToVTK
 202 (
 203     const fileName& fn,
 204     const string& title,
 205     const UList<point>& points,
 206     const LongList<labelledTri>& facets
 207 )
 208 {
 209     labelList connectivity(facets.size());
 210
 211     forAll(facets, faceI)
 212     {
 213         connectivity[faceI] = 3*(faceI+1);
 214     }
 215
 216     labelList regionData(facets.size());
 217
 218     forAll(facets, faceI)
 219     {
 220         regionData[faceI] =  facets[faceI].region();
 221     }
 222
 223     xmlTag xmlRoot("VTKFile");
 224     xmlRoot.addAttribute("type", "PolyData");
 225
 226     xmlTag& xmlPolyData = xmlRoot.addChild("PolyData");
 227
 228     xmlTag& xmlPiece = xmlPolyData.addChild("Piece");
 229     xmlPiece.addAttribute("NumberOfPoints", points.size());
 230     xmlPiece.addAttribute("NumberOfPolys", facets.size());
 231
 232     xmlTag& xmlPoints = xmlPiece.addChild("Points");
 233
 234     xmlTag& xmlPointData = xmlPoints.addChild("DataArray");
 235     xmlPointData.addAttribute("type", "Float32");
 236     xmlPointData.addAttribute("NumberOfComponents", 3);
 237     xmlPointData.addAttribute("format", "ascii");
 238     xmlPointData << points;
 239
 240     xmlTag& xmlPolys = xmlPiece.addChild("Polys");
 241
 242     xmlTag& xmlPolyDataArray = xmlPolys.addChild("DataArray");
 243     xmlPolyDataArray.addAttribute("type", "Int32");
 244     xmlPolyDataArray.addAttribute("Name", "connectivity");
 245     xmlPolyDataArray.addAttribute("format", "ascii");
 246     xmlPolyDataArray << facets;
 247
 248     xmlTag& xmlConnectData = xmlPolys.addChild("DataArray");
 249     xmlConnectData.addAttribute("type", "Int32");
 250     xmlConnectData.addAttribute("Name", "offsets");
 251     xmlConnectData.addAttribute("format", "ascii");
 252     xmlConnectData << connectivity;
 253
 254     xmlTag& xmlCellData = xmlPiece.addChild("CellData");
 255
 256     xmlTag& xmlCellDataArray = xmlCellData.addChild("DataArray");
 257     xmlCellDataArray.addAttribute("type", "Int32");
 258     xmlCellDataArray.addAttribute("Name", "region");
 259     xmlCellDataArray.addAttribute("format", "ascii");
 260     xmlCellDataArray << regionData;
 261
 262     OFstream os(fn);
 263     os << xmlRoot << endl;
 264
 265     Info << "Created " << fn << endl;
 266 }
 267
 268
 269 //- Write an addressed subset of the supplied facet list
 270 //-  in serial vtkPolyData format
 271 void writeFacetsToVTK
 272 (
 273     const fileName& fn,
 274     const string& title,
 275     const pointField& points,
 276     const LongList<labelledTri>& facets,
 277     const unallocLabelList& addr
 278 )
 279 {
 280     // Remove unused points and create subaddressed facets
 281     DynamicList<point> newPoints;
 282     labelList newPointAddr(points.size(), -1);
 283     LongList<labelledTri> newFacets(addr.size());
 284
 285     forAll(addr, addrI)
 286     {
 287         label faceI = addr[addrI];
 288
 289         const labelledTri& facet = facets[faceI];
 290         const FixedList<label, 3>& pointIds = facet;
 291         FixedList<label, 3> newPointIds;
 292
 293         forAll(pointIds, i)
 294         {
 295             label pointId = pointIds[i];
 296
 297             if (newPointAddr[pointId] == -1)
 298             {
 299                 newPoints.append(points[pointId]);
 300                 newPointAddr[pointId] = newPoints.size()-1;
 301             }
 302
 303             newPointIds[i] = newPointAddr[pointId];
 304         }
 305
 306         newFacets[addrI] = labelledTri
 307         (
 308             newPointIds[0],
 309             newPointIds[1],
 310             newPointIds[2],
 311             facet.region()
 312         );
 313     }
 314
 315     writeFacetsToVTK
 316     (
 317         fn,
 318         title,
 319         newPoints,
 320         newFacets
 321     );
 322 }
 323
 324
 325 int main(int argc, char *argv[])
 326 {
 327     argList::noParallel();
 328     argList::validArgs.clear();
 329
 330     argList::validArgs.append("input surface file");
 331     argList::validArgs.append("output prefix");
 332     argList args(argc, argv);
 333
 334     // Process commandline arguments
 335     fileName inFileName(args.args()[1]);
 336     fileName outPrefix(args.args()[2]);
 337
 338     // Read original surface
 339     triSurf origSurf(inFileName);
 340
 341     const pointField& points = origSurf.points();
 342     const LongList<labelledTri>& facets = origSurf.facets();
 343     const LongList<edge>& edges = origSurf.featureEdges();
 344     const geometricSurfacePatchList& patches = origSurf.patches();
 345
 346     label index = 0;
 347
 348     // Create file structure for multiblock dataset
 349     mkDir(outPrefix);
 350     mkDir(outPrefix + "/patches");
 351     mkDir(outPrefix + "/pointSubsets");
 352     mkDir(outPrefix + "/edgeSubsets");
 353     mkDir(outPrefix + "/faceSubsets");
 354
 355     // Create VTK multiblock dataset file
 356     xmlTag xmlRoot("VTKFile");
 357     xmlRoot.addAttribute("type", "vtkMultiBlockDataSet");
 358     xmlRoot.addAttribute("version", "1.0");
 359     xmlRoot.addAttribute("byte_order", "LittleEndian");
 360
 361     xmlTag& xmlDataSet = xmlRoot.addChild("vtkMultiBlockDataSet");
 362
 363     // Write faces and feature edges
 364     {
 365         fileName fn = outPrefix / "facets.vtp";
 366
 367         writeFacetsToVTK
 368         (
 369             outPrefix / "facets.vtp",
 370             outPrefix,
 371             points,
 372             facets
 373         );
 374
 375         xmlTag& tag = xmlDataSet.addChild("DataSet");
 376         tag.addAttribute("index", Foam::name(index++));
 377         tag.addAttribute("name", "facets");
 378         tag.addAttribute("file", fn);
 379     }
 380
 381     {
 382         fileName fn = outPrefix / "featureEdges.vtp";
 383
 384         writeEdgesToVTK
 385         (
 386             outPrefix / "featureEdges.vtp",
 387             "featureEdges",
 388             points,
 389             edges
 390         );
 391
 392         xmlTag& tag = xmlDataSet.addChild("DataSet");
 393         tag.addAttribute("index", Foam::name(index++));
 394         tag.addAttribute("name", "featureEdges");
 395         tag.addAttribute("file", fn);
 396     }
 397
 398     // Write patches
 399     // Create patch addressing
 400     List<DynamicList<label> > patchAddr(patches.size());
 401
 402     forAll(facets, faceI)
 403     {
 404         patchAddr[facets[faceI].region()].append(faceI);
 405     }
 406
 407     {
 408         xmlTag& xmlBlock = xmlDataSet.addChild("Block");
 409         xmlBlock.addAttribute("index", Foam::name(index++));
 410         xmlBlock.addAttribute("name", "patches");
 411
 412         forAll(patches, patchI)
 413         {
 414             word patchName = patches[patchI].name();
 415
 416             fileName fn = outPrefix / "patches" / patchName + ".vtp";
 417
 418             writeFacetsToVTK
 419             (
 420                 fn,
 421                 patchName,
 422                 points,
 423                 facets,
 424                 patchAddr[patchI]
 425             );
 426
 427             xmlTag& tag = xmlBlock.addChild("DataSet");
 428             tag.addAttribute("index", Foam::name(patchI));
 429             tag.addAttribute("name", patchName);
 430             tag.addAttribute("file", fn);
 431         }
 432     }
 433
 434     // Write point subsets
 435     {
 436         xmlTag& xmlBlock = xmlDataSet.addChild("Block");
 437         xmlBlock.addAttribute("index", Foam::name(index++));
 438         xmlBlock.addAttribute("name", "pointSubsets");
 439
 440         DynList<label> subsetIndices;
 441         labelList subsetAddr;
 442
 443         origSurf.pointSubsetIndices(subsetIndices);
 444
 445         forAll(subsetIndices, id)
 446         {
 447             word subsetName = origSurf.pointSubsetName(id);
 448             origSurf.pointsInSubset(id, subsetAddr);
 449
 450             fileName fn = outPrefix / "pointSubsets" / subsetName + ".vtp";
 451
 452             writePointsToVTK
 453             (
 454                 fn,
 455                 subsetName,
 456                 points,
 457                 subsetAddr
 458             );
 459
 460             xmlTag& tag = xmlBlock.addChild("DataSet");
 461             tag.addAttribute("index", Foam::name(id));
 462             tag.addAttribute("name", subsetName);
 463             tag.addAttribute("file", fn);
 464         }
 465     }
 466
 467     // Write edge subsets
 468     {
 469         xmlTag& xmlBlock = xmlDataSet.addChild("Block");
 470         xmlBlock.addAttribute("index", Foam::name(index++));
 471         xmlBlock.addAttribute("name", "edgeSubsets");
 472
 473         DynList<label> subsetIndices;
 474         labelList subsetAddr;
 475
 476         origSurf.edgeSubsetIndices(subsetIndices);
 477
 478         forAll(subsetIndices, id)
 479         {
 480             word subsetName = origSurf.edgeSubsetName(id);
 481             origSurf.edgesInSubset(id, subsetAddr);
 482
 483             fileName fn = outPrefix / "edgeSubsets" / subsetName + ".vtp";
 484
 485             writeEdgesToVTK
 486             (
 487                 fn,
 488                 subsetName,
 489                 points,
 490                 edges,
 491                 subsetAddr
 492             );
 493
 494             xmlTag& tag = xmlBlock.addChild("DataSet");
 495             tag.addAttribute("index", Foam::name(id));
 496             tag.addAttribute("name", subsetName);
 497             tag.addAttribute("file", fn);
 498         }
 499     }
 500
 501     // Write facet subsets
 502     {
 503         xmlTag& xmlBlock = xmlDataSet.addChild("Block");
 504         xmlBlock.addAttribute("index", Foam::name(index++));
 505         xmlBlock.addAttribute("name", "faceSubsets");
 506
 507         DynList<label> subsetIndices;
 508         labelList subsetAddr;
 509
 510         origSurf.facetSubsetIndices(subsetIndices);
 511
 512         forAll(subsetIndices, id)
 513         {
 514             word subsetName = origSurf.facetSubsetName(id);
 515             origSurf.facetsInSubset(id, subsetAddr);
 516
 517             fileName fn = outPrefix / "faceSubsets"  / subsetName + ".vtp";
 518
 519             writeFacetsToVTK
 520             (
 521                 fn,
 522                 subsetName,
 523                 points,
 524                 facets,
 525                 subsetAddr
 526             );
 527
 528             xmlTag& tag = xmlBlock.addChild("DataSet");
 529             tag.addAttribute("index", Foam::name(id));
 530             tag.addAttribute("name", subsetName);
 531             tag.addAttribute("file", fn);
 532         }
 533     }
 534
 535     OFstream os(outPrefix + ".vtm");
 536     os << xmlRoot << endl;
 537
 538     Info << "Created " << outPrefix + ".vtm" << endl;
 539
 540     Info << "End\n" << endl;
 541
 542     return 0;
 543 }
 544
 545 // ************************************************************************* //