メインコンテンツに移動

非構造格子とポリデータ

vtu_1は, 非構造格子を読み取り, ポリデータに変換し, 領域の境界を判定する例です.

たとえば元の VTK UnstructuredGrid データが

a1

こうであったとしましょう.  これをまあ普通に読み込んで,

#include <iostream>
#include <boost/filesystem.hpp>
//Load VTU file
#include <vtkSmartPointer.h>
#include <vtkNew.h>
#include <vtkUnstructuredGrid.h>
#include <vtkXMLUnstructuredGridReader.h>
...
    //Load VTU file
    boost::filesystem::path cpwd=getenv("HOME");
    cpwd/="CGAL_test";
    boost::filesystem::path filename = cpwd/"2dmesh_test.vtu";
    vtkNew<vtkXMLUnstructuredGridReader> reader;
    reader->SetFileName(filename.c_str());

PolyDataに変換するには,

//Convert to PolyData
#include <vtkGeometryFilter.h>
#include <vtkPolyData.h>
...
    //Convert to vtkPolyData
    //    https://kitware.github.io/vtk-examples/site/Cxx/PolyData/GeometryFilter/
    auto geometryFilter = vtkSmartPointer<vtkGeometryFilter>::New();
    geometryFilter->SetInputConnection(reader->GetOutputPort());

PolyDataでは, FeatureEdgeフィルターが使えるので, いろいろ取得できます.たとえば領域の縁をゲットしたい場合

//Extract Boundary
#include <vtkFeatureEdges.h>    //in vtkFiltersCore, vtkCommonExecutionModel
...
    //Get Boundary of the data
    //    https://kitware.github.io/vtk-examples/site/Cxx/Meshes/BoundaryEdges/
    vtkNew<vtkFeatureEdges> featureEdges;
    featureEdges->SetInputConnection(geometryFilter->GetOutputPort());
    featureEdges->BoundaryEdgesOn();
    featureEdges->FeatureEdgesOff();
    featureEdges->ManifoldEdgesOff();
    featureEdges->NonManifoldEdgesOff();
    featureEdges->ColoringOn();
    featureEdges->Update();
    vtkPolyData& polydata = *geometryFilter->GetOutput();
    vtkPolyData& polydata_B = *featureEdges->GetOutput();

物理量を追加することも可能:

//Add function
#include <vtkDoubleArray.h>
#include <vtkCellData.h>
#include <vtkCell.h>
...
    //Add function to polydata
    vtkDoubleArray& pressure = *vtkDoubleArray::New();        //Create array
    pressure.SetNumberOfComponents(0);                        //It is scalar
    pressure.SetName("p");                                    //Its name
    pressure.SetNumberOfTuples(polydata.GetNumberOfCells());  //Its number
    polydata.GetCellData()->AddArray(&pressure);              //Add it !
    //Add sample values
    for (int cid = 0; cid < polydata.GetNumberOfCells(); cid++)
    {
        double p=cid;
        pressure.SetValue(cid,p);
    }

もちろんVTPファイルに保存できる:

    //Save VTP file
    {
        vtkNew<vtkXMLPolyDataWriter> writer;
        boost::filesystem::path p_name = cpwd/"2dmesh_p.vtp";
        writer->SetFileName(p_name.c_str());
        writer->SetInputData(&polydata);
        writer->Write();
    }
    {
        vtkNew<vtkXMLPolyDataWriter> writer;
        boost::filesystem::path p_name = cpwd/"2dmesh_b.vtp";
        writer->SetFileName(p_name.c_str());
        writer->SetInputData(&polydata_B);
        writer->Write();
    }

で,こうなる:

f2