VTK的网格模型切割算法

2020-03-20 21:39:02|?次阅读|上传:wustguangh【已有?条评论】发表评论

关键词:C/C++, VTK|来源:唯设编程网

VTK提供了vtkClip相关的类来实现网格模型分割,其中vtkClipDataSet使用用户定义的隐函数或者point的标量数据来切割任意类型的网格模型,vtkClipVolume使用隐函数或者point的标量数据来切割volume类型数据,vtkClipPolyData使用隐函数或者point的标量数据来切割vtkPolyData类型的网格数据。

隐函数通过vtkImplicitFunction的子类来定义,vtkImplicitFunction类提供了两个重要的虚函数,一个是计算给定坐标位置函数值的EvaluateFunction函数,另一个是计算给定坐标位置梯度向量的EvaluateGradient函数,这两个虚函数需要在子类中准确定义。VTK提供了许多可以直接使用的隐函数类,经常使用到的包括包括平面函数vtkPlane(函数值是对应外置坐标到平面的距离,外侧为正值,内侧为负值),球面函数vtkSphere,圆锥面函数vtkCone,圆柱面函数vtkCylinder,长方体函数vtkBox,通用二次曲面函数vtkQuadric。此外,用户还可以使用vtkImplicitBoolean通过集合运算(求交、并、取反)将简单函数进行组合,从而得到各种更加复杂的函数。

如何使用隐函数来执行切割运算,vtkClip切割算法的首先根据隐函数求得每个节点处隐函数的函数值,并将其作为point的scalar值(如何没有设定隐函数则直接使用该point的scalar值),如果一个cell每个point的scalar值均小于vtkClip算法设定的value值(默认为0),则该cell完全在应该被切割的一侧,反之在应该保留的一侧。

如果一个cell的某条edge的两个端点scalar值一个大于value,而另一个小于value,则在该edge上通过插值算法找出等于value值的位置,使用这些点所在的平面形成针对该cell的切割平面,调用对应类型vtkCell3D类的clip函数完成该cell的切割。

vtkClip算法默认只保留scalar值大于value的一侧,也可以调用GenerateClippedOutputOn函数使其保留被切割的一侧,被切割一侧的结果可以使用GetClippedOutput函数得到。

下面是使用vtkClipDataSet实现网格模型平面和长方体组合切割的示例代码:

    vtkSmartPointer<vtkImageData> imgData = vtkSmartPointer<vtkImageData>::New();
    imgData->SetDimensions(15, 15, 15);
    imgData->SetSpacing(1.0, 1.0, 1.0);
    imgData->SetOrigin(0.0, 0.0, 0.0);

    vtkSmartPointer<vtkFloatArray> xCoordes = vtkSmartPointer<vtkFloatArray>::New();
    xCoordes->SetNumberOfComponents(1);
    for(int i=0; i<15; ++i)
        xCoordes->InsertNextTuple1(i);

    vtkSmartPointer<vtkFloatArray> yCoordes = vtkSmartPointer<vtkFloatArray>::New();
    yCoordes->SetNumberOfComponents(1);
    for(int i=0; i<15; ++i)
        yCoordes->InsertNextTuple1(i);

    vtkSmartPointer<vtkFloatArray> zCoordes = vtkSmartPointer<vtkFloatArray>::New();
    zCoordes->SetNumberOfComponents(1);
    for(int i=0; i<15; ++i){
        if(i%2)
            zCoordes->InsertNextTuple1(i);
        else
            zCoordes->InsertNextTuple1(i+0.25);
    }

    //使用线性网格
    vtkSmartPointer<vtkRectilinearGrid> recGrid = vtkSmartPointer<vtkRectilinearGrid>::New();
    recGrid->SetDimensions(xCoordes->GetNumberOfTuples(),
                           yCoordes->GetNumberOfTuples(),
                           zCoordes->GetNumberOfTuples());
    recGrid->SetXCoordinates(xCoordes);
    recGrid->SetYCoordinates(yCoordes);
    recGrid->SetZCoordinates(zCoordes);

    //CELL数据
    vtkSmartPointer<vtkIntArray> cellDataArray = vtkSmartPointer<vtkIntArray>::New();
    cellDataArray->SetName("cell_scalar");
    for(int i=0; i<imgData->GetNumberOfCells(); ++i){
        cellDataArray->InsertNextTuple1(i);
    }
    imgData->GetCellData()->SetScalars(cellDataArray);
    recGrid->GetCellData()->SetScalars(cellDataArray);

    //Point数据
    vtkSmartPointer<vtkIntArray> pointDataArray = vtkSmartPointer<vtkIntArray>::New();
    pointDataArray->SetName("point_scalar");
    for(int i=0; i<imgData->GetNumberOfPoints(); ++i){
        pointDataArray->InsertNextTuple1(i);
    }
    imgData->GetPointData()->SetScalars(pointDataArray);
    recGrid->GetPointData()->SetScalars(pointDataArray);

    //切割平面
    vtkSmartPointer<vtkPlane> planeYOZ = vtkSmartPointer<vtkPlane>::New();
    planeYOZ->SetOrigin(7.5, 7.5, 7.5);
    planeYOZ->SetNormal(1.0, 0.0, 0.0);

    vtkSmartPointer<vtkPlane> planeXOZ = vtkSmartPointer<vtkPlane>::New();
    planeXOZ->SetOrigin(7.5, 7.5, 7.5);
    planeXOZ->SetNormal(0.0, 1.0, 0.0);

    vtkSmartPointer<vtkGeometryFilter> imgGeoFilter = vtkSmartPointer<vtkGeometryFilter>::New();
    imgGeoFilter->SetInputData(imgData);
    imgGeoFilter->Update();

    vtkSmartPointer<vtkClipDataSet> clipYOZ = vtkSmartPointer<vtkClipDataSet>::New();
    clipYOZ->SetInputData(recGrid);
    clipYOZ->SetClipFunction(planeYOZ);
    clipYOZ->GenerateClippedOutputOn();
    clipYOZ->Update();

    vtkSmartPointer<vtkAppendFilter> appendFilter1 = vtkSmartPointer<vtkAppendFilter>::New();
    appendFilter1->AddInputConnection(clipYOZ->GetOutputPort());
    appendFilter1->AddInputData(clipYOZ->GetClippedOutput());
    appendFilter1->Update();

    vtkSmartPointer<vtkClipDataSet> clipXOZ = vtkSmartPointer<vtkClipDataSet>::New();
    clipXOZ->SetInputData(appendFilter1->GetOutput());
    clipXOZ->SetClipFunction(planeXOZ);
    clipXOZ->GenerateClippedOutputOn();
    clipXOZ->Update();

    vtkSmartPointer<vtkAppendFilter> appendFilter2 = vtkSmartPointer<vtkAppendFilter>::New();
    appendFilter2->AddInputConnection(clipXOZ->GetOutputPort());
    appendFilter2->AddInputData(clipXOZ->GetClippedOutput());
    appendFilter2->Update();

    //开始造型切割
    vtkSmartPointer<vtkBox> box1 = vtkSmartPointer<vtkBox>::New();
    box1->SetBounds(-5.0, 7.5, -5.0, 7.5, -5.0, 16.0);
    //切割
    vtkSmartPointer<vtkClipDataSet> clipper1 = vtkSmartPointer<vtkClipDataSet>::New();
    clipper1->SetInputConnection(appendFilter2->GetOutputPort());
    clipper1->SetClipFunction(box1);

    //第二轮切割
    vtkSmartPointer<vtkBox> box2 = vtkSmartPointer<vtkBox>::New();
    box2->SetBounds(-5.0, 7.5, -5.0, 7.5, 5.0, 16.0);

    vtkSmartPointer<vtkTransform> trans = vtkSmartPointer<vtkTransform>::New();
    trans->Translate(7.5, 7.5, 0.0);
    trans->RotateZ(-45);
    trans->Translate(-7.5, -7.5, 0.0);

    box2->SetTransform(trans);

    //切割
    vtkSmartPointer<vtkClipDataSet> clipper2 = vtkSmartPointer<vtkClipDataSet>::New();
    clipper2->SetInputConnection(clipper1->GetOutputPort());
    clipper2->SetClipFunction(box2);

    //第三轮切割
    vtkSmartPointer<vtkBox> box3 = vtkSmartPointer<vtkBox>::New();
    box3->SetBounds(7.5, 15.0, -5.0, 7.5, 10.0, 16.0);

    //切割
    vtkSmartPointer<vtkClipDataSet> clipper3 = vtkSmartPointer<vtkClipDataSet>::New();
    clipper3->SetInputConnection(clipper2->GetOutputPort());
    clipper3->SetClipFunction(box3);

    //保存到文件中
    vtkSmartPointer<vtkXMLUnstructuredGridWriter> writer = vtkSmartPointer<vtkXMLUnstructuredGridWriter>::New();
    writer->SetFileName("D:/VTKTestData/clipper.vtu");
    writer->SetInputConnection(clipper1->GetOutputPort());
    writer->Update();

    vtkSmartPointer<vtkDataSetMapper> mapper = vtkSmartPointer<vtkDataSetMapper>::New();
    mapper->AddInputConnection(clipper1->GetOutputPort());

    vtkSmartPointer<vtkActor> actor = vtkSmartPointer<vtkActor>::New();
    actor->SetMapper(mapper);
    actor->GetProperty()->EdgeVisibilityOn();

最终实现的效果如下:

VTK切割算法

发表评论0条 】
网友评论(共?条评论)..
VTK的网格模型切割算法