VTK使用轮廓提取技术生成等值面

2022年2月24日 381点热度 1人点赞 0条评论

参考之前记录的标量颜色映射,很多时候为了搞清楚标量在空间或者平面上的分布趋势,通常会通过提取等值面或者等值线的方式将这些数据可视化,这里会记录怎么使用VTK实现这样的功能。

工欲善其事必先利其器,实现这个功能的核心是vtkContourFIlter类,它提供了等值线(等值面)提取算法,只需要准备好数据,然后设置好轮廓参数,就可以实现了。关于过滤器的使用可以参考之前关于VTK渲染流水线的博客。

vtkContourFilter

这个过滤器的功能是从众多的标量值中找出标量值相等的位置,把这些位置连成线或者生成面,当然这个过程中会使用已有的标量值进行插值。它可以将任何数据集作为输入并输出等值面或者等值线(这里可以参考在标量颜色映射中对vtkDataSet类的介绍),输出的数据的形式却决于输入数据的维度,如果是3D单元数据那么输出的就是等值面,2D单元输出等值线,1D单元则输出一些点。

要使用这个过滤器必须为过滤器设置轮廓值(等值面的值),可以使用SetValue()来指定每一个轮廓值,或者使用GenerateValues()来生成一些均匀分布的轮廓值。

生成轮廓值需要耗费大量的时间,可以使用vtkScalarTree来加速这个操作,标量树用于快速定位包好等值表面的单元,尤其是对需要提取多个轮廓的时候效果极佳,如果要使用标量树需要调用UseScalarTreeOn()开启。

等值线提取

将之前标量颜色映射的示例拿过来,增加一个vtkContourFIlter。

//为三角形准备三个顶点
    vtkNew<vtkPoints> points;
    points->InsertNextPoint(1, 0, 0);
    points->InsertNextPoint(0, 0, 0);
    points->InsertNextPoint(0, 1, 0);

    //生成点的标量数据
    vtkNew<vtkIntArray> scalar;
    scalar->InsertNextTuple1(30);
    scalar->InsertNextTuple1(100);
    scalar->InsertNextTuple1(50);

    //使用顶点构建一个三角形 点使用ID传入
    vtkNew<vtkTriangle> triangle;
    triangle->GetPointIds()->SetId(0, 0);
    triangle->GetPointIds()->SetId(1, 1);
    triangle->GetPointIds()->SetId(2, 2);


    //将三角形传入单元集
    vtkNew<vtkCellArray> cellArray;
    cellArray->InsertNextCell(triangle);

    //为多边形数据添加单元集和点集
    vtkNew<vtkPolyData> polydata;
    polydata->SetPolys(cellArray);
    polydata->SetPoints(points);
    //为点数据添加标量
    polydata->GetPointData()->SetScalars(scalar);


    //添加一个过滤器
    vtkNew<vtkContourFilter> contourFilter;
    contourFilter->SetInputData(polydata);
    //设置等值线的值
    contourFilter->SetValue(0, 50);
    contourFilter->SetComputeScalars(true);
    contourFilter->Update();


    vtkNew<vtkPolyDataMapper> mapper;
    mapper->SetInputData(polydata);
    mapper->SetScalarRange(scalar->GetRange());
    mapper->Update();

    vtkNew<vtkPolyDataMapper> contourMapper;
    contourMapper->SetInputData(contourFilter->GetOutput());
    contourMapper->ScalarVisibilityOff();

    //为映射器添加一个颜色映射表
    //不添加颜色映射表会得到一个默认的颜色表
    //如果需要自定义颜色表需要调用SetTableValue()设置
 	vtkNew<vtkLookupTable> lookupTable;
    lookupTable->SetTableRange(scalar->GetRange());
    lookupTable->Build();
 	mapper->SetLookupTable(lookupTable);
    contourMapper->SetLookupTable(lookupTable);

    vtkNew<vtkScalarBarActor> barActor;
    barActor->SetLookupTable(lookupTable);
    barActor->SetNumberOfLabels(6);

    vtkNew<vtkActor> actor;
    actor->SetMapper(mapper);
    actor->GetProperty()->SetColor(255, 0, 0);

    vtkNew<vtkActor> contourActor;
    contourActor->SetMapper(contourMapper);
    contourActor->GetProperty()->SetColor(0,0,0);

    vtkNew<vtkRenderer> renderer;
    vtkNew<vtkRenderWindow> window;
    window->AddRenderer(renderer);

    vtkNew<vtkRenderWindowInteractor> interactor;
    interactor->SetRenderWindow(window);

    renderer->AddActor(actor);
    renderer->AddActor(contourActor);
    renderer->AddActor(barActor);
    window->Render();
    interactor->Start();

 

等值面轮廓提取

提取等值面与提取等值线的过程如出一辙,首先得提供一个三维点数据与标量数据。

(x^2 + (9/4) y^2 + z^2 - 1)^3 - x^2*z^3 - (9/80) y^2*z^3 =F, {x, -2, 2}, {y, -2, 2}, {z, -2, 2}

我在网上找到了这么一个函数,据说这个F=0的时候,xyz的解会形成一个小心心。于是我生成了一个40*40*40的三维网格数据,然后每个网格的标量值使用这个函数进行赋值,最后再求F=0时的等值面。

vtkNew<vtkStructuredGrid> structuredGrid;
    vtkNew<vtkPoints> points;
    vtkNew<vtkFloatArray> scalars;
    unsigned int gridSize = 40;
    //构造一个点云,使用函数为点云构造标量值
    for (unsigned int i = 0; i < gridSize; i++)
    {
        for (unsigned int j = 0; j < gridSize; j++)
        {
            for (unsigned int k = 0; k < gridSize; k++)
            {
                points->InsertNextPoint(i, j, k);
                //实现函数 (x^2 + (9/4) y^2 + z^2 - 1)^3 - x^2*z^3 - (9/80) y^2*z^3 =F, {x, -2, 2}, {y, -2, 2}, {z, -2, 2}
                double x = (i - (gridSize / 2.0)) / 10.0;
                double y = (j - (gridSize / 2.0)) / 10.0;
                double z = (k - (gridSize / 2.0)) / 10.0;
                scalars->InsertNextTuple1( pow(x*x + (9.0 / 4.0)*y*y + z*z- 1.0,3) - x*x*pow(z,3) - (9.0 / 80.0)*y*y * pow(z,3) );
            }
        }
    }



    //设置网格大小
    structuredGrid->SetDimensions(gridSize, gridSize, gridSize);
    //附加数据
    structuredGrid->SetPoints(points);
    structuredGrid->GetPointData()->SetScalars(scalars);


    //设置映射器
    vtkNew<vtkDataSetMapper> gridMapper;
    gridMapper->SetInputData(structuredGrid.Get());
    gridMapper->ScalarVisibilityOff();


    vtkNew<vtkActor> gridActor;
    gridActor->SetMapper(gridMapper.Get());
    gridActor->GetProperty()->SetOpacity(0.3);


    //设置等值线过滤器
    vtkNew<vtkContourFilter> contourFilter;
    contourFilter->SetInputData(structuredGrid.Get());
    //为标量为0时取一个等值面
    contourFilter->SetValue(0, 0);
    contourFilter->Update();

    //创建颜色表 这里只输入颜色值的数量,然后自动生成表
    vtkNew<vtkLookupTable> lut;
    lut->SetTableRange(scalars->GetRange());
    lut->Build();

    //为生成的等值面添加法向
    vtkNew<vtkPolyDataNormals> normals;
    normals->SetInputConnection(contourFilter->GetOutputPort());
    normals->SetFeatureAngle(45);
    normals->Update();

    vtkNew<vtkDataSetMapper> contourMaper;
    contourMaper->SetInputData(normals->GetOutput());
    contourMaper->SetLookupTable(lut);
    contourMaper->SetScalarRange(scalars->GetRange());
    vtkNew<vtkActor> contourActor;
    contourActor->SetMapper(contourMaper);


    vtkNew<vtkRenderer> renderer;
    vtkNew<vtkRenderWindow> renderWindow;
    renderWindow->AddRenderer(renderer);


    vtkNew<vtkScalarBarActor> barActor;
    barActor->SetLookupTable(lut);


    vtkNew<vtkRenderWindowInteractor> renderWindowInteractor;
    renderWindowInteractor->SetRenderWindow(renderWindow);


    
    renderer->AddActor(contourActor);
    renderer->AddActor(barActor);
    renderer->AddActor(gridActor);

    renderWindow->Render();

    renderWindowInteractor->Start();

哇!真的有小心心欸。

赶紧送给对象,什么?你没有对象?那感觉new一个吧。

小结

以上,即是记录实现等值线面的提取的全过程,过程十分简单,不过有些细节需要十分注意。

  1. 映射器一定要设置好标量范围,否则无法正常完成标量值的颜色映射
  2. 三维等值面提取时一定要重新计算多边形的法向

第二点需要特别注意,我在写示例时就遇到了这个问题,轮廓提取之后的多边形是没有法向信息的,这就会导致它不能正常的反射光照,即生成的等值面颜色不对,或者直接和背景黑色融为一体,所以一定要重新计算法向!!

大脸猫

这个人虽然很勤快,但什么也没有留下!

文章评论