1 VTK(Visualization Toolkit)是一个开源的跨平台软件系统,用于三维计算机图形学、图像处理和可视化。学习VTK的主要目的有:
- 3D可视化: VTK提供了丰富的工具和算法,可以用来可视化各种科学数据,包括医学图像、有限元分析结果、流体模拟数据等。
- 图像处理: VTK支持多种图像处理操作,如滤波、分割、配准等,可以对医学图像进行预处理和分析。
- 科学计算可视化: VTK可以与其他科学计算软件(如MATLAB、Python)结合,实现科学计算结果的可视化。
- 自定义可视化: VTK提供了灵活的接口,可以自定义各种可视化组件和交互方式,满足个性化的需求。
对于有科研需求的小伙伴,可以利用VTK将你的数据和实验结果进行可视化,对数据有一个感官的理解。比如下面的图像就是一个细胞分裂增殖的模拟图像,使用VTK对数据进行了可视化。
3图1 细胞克隆增殖的模拟
VTK就是一个工具,不用刻意的去学习它。假如你不做数据的显示,也不需要做3D模型的渲染,你就不用去学习它,因为学习任何事物都是需要学习成本的。
2 学习方法
学习VTK的方法有很多种,以下是一些建议:
- 官方文档: VTK官方文档是学习VTK最权威的资料,提供了详细的类、函数和示例。
- 教程和书籍: 网上有很多关于VTK的教程和书籍,可以帮助你快速入门。
- 示例代码: VTK提供了大量的示例代码,可以帮助你理解各种概念和用法。
- 社区: VTK有一个活跃的社区,你可以通过论坛、邮件列表等方式向其他用户提问和交流。
- 实践: 最好的学习方式是实践。通过编写自己的代码,尝试实现各种功能,可以加深对VTK的理解。
3 一个比较有意思的例子
我们模拟一下太阳光线照射到地球上,然后有反射的情况,其实就是一个光线追踪算法。
图2 射线追踪算法
先看代码:
#include <iostream>
#include <vtkSmartPointer.h>
#include <vtkSphereSource.h>
#include <vtkActor.h>
#include <vtkConeSource.h>
#include <vtkRenderer.h>
#include <vtkRenderWindow.h>
#include <vtkPolyDataMapper.h>
#include <vtkProperty.h>
#include <vtkRenderWindowInteractor.h>
#include <vtkOBBTree.h>
#include <vtkMath.h>
#include <vtkTransform.h>
#include <vtkTransformFilter.h>
#include <vtkMatrix4x4.h>
#include <vtkPolyDataNormals.h>
#include <vtkPlaneSource.h>
#include <vtkPoints.h>
#include <vtkInteractorStyleTrackballCamera.h>
#include <vtkCamera.h>
#include <vtkLight.h>
#include <vtkJPEGReader.h>
#include <vtkTexture.h>
#include <vtkTextureMapToSphere.h>
#include <vtkCellCenters.h>
#include <vtkPolyDataNormals.h>
#include <vtkVertexGlyphFilter.h>
#include <vtkGlyph3D.h>
#include <vtkArrowSource.h>
#include <vtkCellData.h>
#include <vtkPointData.h>
#include <vtkIdList.h>
#include <vtkDoubleArray.h>
#include <vtkLineSource.h>
using namespace std;
vtkSmartPointer<vtkActor> point2Actor(vtkSmartPointer<vtkPoints> points){
vtkNew<vtkPolyData> polyData;
polyData->SetPoints(points);
vtkNew<vtkVertexGlyphFilter> filter;
filter->SetInputData(polyData);
filter->Update();
vtkNew<vtkPolyDataMapper> mapper;
mapper->SetInputConnection(filter->GetOutputPort());
vtkSmartPointer<vtkActor> actor=vtkSmartPointer<vtkActor>::New();
actor->SetMapper(mapper);
return actor;
}
int main()
{
// Sun Options
// Radius of the sun half-sphere
double RadiusSun=10.0;
// Distance of the sun's center from (0,0,0)
double DistanceSun=50.0;
// Phi & Theta Resolution of sun
int ResolutionSun=10;
// Earth Options
// Radius of the earth sphere
double RadiusEarth=150.0;
// Phi & Theta Resolution of earth
int ResolutionEarth=120;
// Ray Options
double RayCastLength=500.0;
// Color Options
double ColorSun[3]={1.0,1.0,0.0};
double ColorSunEdge[3]={0.0,0.0,0.0};
double ColorEarthEdge[3]={1.0,1.0,1.0};
double ColorBackground[]={0.0,0.0,0.0};
double ColorLight[]={1.0,1.0,0.0};
double ColorSunPoints[]={1.0,1.0,0.0};
double ColorSunGlyphs[]={1.0,1.0,0.0};
double ColorRayHit[]={1.0,1.0,0.0};
double ColorRayMiss[]={1.0,1.0,1.0};
double OpacityRayMiss=0.5;
double ColorEarthGlyphs[]={0.0,0.0,1.0};
double ColorRayReflected[]={1.0,1.0,0.0};
// Load a JPEG file with an 'earth' texture
vtkNew<vtkJPEGReader> reader;
reader->SetFileName("../data/R-C.jpeg");
reader->Update();
vtkNew<vtkSphereSource> sun;
sun->SetCenter(0.0,DistanceSun,0.0);
sun->SetRadius(RadiusSun);
sun->SetThetaResolution(ResolutionSun);
sun->SetPhiResolution(ResolutionSun);
sun->SetStartTheta(180); // create a half-sphere
// Create mapper
vtkNew<vtkPolyDataMapper> sunMapper;
sunMapper->SetInputConnection(sun->GetOutputPort());
// Create actor
vtkNew<vtkActor> sunActor;
sunActor->SetMapper(sunMapper);
sunActor->GetProperty()->SetColor(ColorSun);
sunActor->GetProperty()->EdgeVisibilityOn(); // show edges/wireframe
sunActor->GetProperty()->SetEdgeColor(ColorSunEdge);
// Create and configure the earth spere
vtkNew<vtkSphereSource> earth;
earth->SetCenter(0.0,-RadiusEarth,0.0);
earth->SetThetaResolution(ResolutionEarth);
earth->SetPhiResolution(ResolutionEarth);
earth->SetRadius(RadiusEarth);
// texture to earth
vtkNew<vtkTexture> texture;
texture->SetInputConnection(reader->GetOutputPort());
texture->Update();
vtkNew<vtkTextureMapToSphere> map_to_sphere;
map_to_sphere->SetInputConnection(earth->GetOutputPort());
map_to_sphere->PreventSeamOn();
map_to_sphere->Update();
vtkNew<vtkPolyDataMapper> earthMapper;
earthMapper->SetInputConnection(map_to_sphere->GetOutputPort());
vtkNew<vtkActor> earthActor;
earthActor->SetMapper(earthMapper);
earthActor->SetTexture(texture);
earthActor->GetProperty()->EdgeVisibilityOn();
earthActor->GetProperty()->SetEdgeColor(ColorEarthEdge);
// Renderer to display
vtkNew<vtkRenderer> ren;
ren->AddActor(sunActor);
ren->AddActor(earthActor);
ren->SetBackground(ColorBackground);
vtkNew<vtkCellCenters> cellCenterCalcSun;
cellCenterCalcSun->SetInputConnection(sun->GetOutputPort());
cellCenterCalcSun->Update();
//cellCenterCalcSun->Print(cout);
//cout<<cellCenterCalcSun->GetOutput()->GetNumberOfPoints()<<endl;
vtkPoints* pointsInCellCenterCalcSun=cellCenterCalcSun->GetOutput()->GetPoints();
/*
for(vtkIdType i=0;i<pointsInCellCenterCalcSun->GetNumberOfPoints();++i){
double point[3];
pointsInCellCenterCalcSun->GetPoint(i,point);
cout<<point[0]<<","<<point[1]<<","<<point[2]<<endl;
}
*/
for(vtkIdType i=0;i<cellCenterCalcSun->GetOutput()->GetNumberOfPoints();++i){
/*
cout<<cellCenterCalcSun->GetOutput()->GetPoint(i)[0]<<","
<<cellCenterCalcSun->GetOutput()->GetPoint(i)[1]<<","
<<cellCenterCalcSun->GetOutput()->GetPoint(i)[2]<<endl;
*/
vtkNew<vtkPoints> points;
points->InsertNextPoint(cellCenterCalcSun->GetOutput()->GetPoint(i)[0],
cellCenterCalcSun->GetOutput()->GetPoint(i)[1],
cellCenterCalcSun->GetOutput()->GetPoint(i)[2]);
vtkSmartPointer<vtkActor> actor=point2Actor(points);
actor->GetProperty()->SetColor(ColorSunPoints);
actor->GetProperty()->SetPointSize(5);
ren->AddActor(actor);
}
vtkNew<vtkPolyDataNormals> normalsCalcSun;
normalsCalcSun->SetInputConnection(sun->GetOutputPort());
// Disable normal calculation at cell vertices
normalsCalcSun->ComputePointNormalsOff();
// Enable normal calculation at cell centers
normalsCalcSun->ComputeCellNormalsOn();
// Disable splitting of sharp edges
normalsCalcSun->SplittingOff();
// Disable global flipping of normal orientation
normalsCalcSun->FlipNormalsOff();
// Enable automatic determination of correct normal orientation
normalsCalcSun->AutoOrientNormalsOn();
normalsCalcSun->Update();
// Create a 'dummy' vtkCellCenters to force the glyphs to the cell-centers
vtkNew<vtkCellCenters> dummy_cellCenterCalcSun;
dummy_cellCenterCalcSun->VertexCellsOn();
dummy_cellCenterCalcSun->SetInputConnection(normalsCalcSun->GetOutputPort());
// Create a new 'default' arrow to use as a glyph
vtkNew<vtkArrowSource> arrow;
vtkNew<vtkGlyph3D> glyphSun;
glyphSun->SetInputConnection(dummy_cellCenterCalcSun->GetOutputPort());
glyphSun->SetSourceConnection(arrow->GetOutputPort());
glyphSun->SetVectorModeToUseNormal();
glyphSun->SetScaleFactor(5);
vtkNew<vtkPolyDataMapper> glyphMapperSun;
glyphMapperSun->SetInputConnection(glyphSun->GetOutputPort());
vtkNew<vtkActor> glyphActorSun;
glyphActorSun->SetMapper(glyphMapperSun);
glyphActorSun->GetProperty()->SetColor(ColorSunGlyphs);
ren->AddActor(glyphActorSun);
// Prepare for ray-tracing
vtkNew<vtkOBBTree> obbEarth;
obbEarth->SetDataSet(earth->GetOutput());
obbEarth->BuildLocator();
double tol=1.e-8;
obbEarth->SetTolerance(tol);
vtkNew<vtkPolyDataNormals> normalsCalcEarth;
normalsCalcEarth->SetInputConnection(earth->GetOutputPort());
normalsCalcEarth->ComputePointNormalsOff();
normalsCalcEarth->ComputeCellNormalsOn();
normalsCalcEarth->SplittingOff();
normalsCalcEarth->FlipNormalsOff();
normalsCalcEarth->AutoOrientNormalsOn();
normalsCalcEarth->Update();
//cout<<"Prepare the ray using point and normal"<<endl;
// Extract the normal-vector data at the sun's cells
vtkDataArray* normalsSun=normalsCalcSun->GetOutput()->GetCellData()->GetNormals();
// Extract the normal-vector data at the earth's cells
vtkDataArray* normalsEarth=normalsCalcEarth->GetOutput()->GetCellData()->GetNormals();
// Create a dummy 'vtkPoints' to act as a container for the point coordinates where intersections are found
vtkNew<vtkPoints> intersectPoints;
vtkNew<vtkIdList> intersectCells;
// Create a dummy 'vtkList' to act as a container for the normal vectors where intersections are found.
vtkNew<vtkDoubleArray> normalsVector;
normalsVector->SetNumberOfComponents(3);
// Create a dummy 'vtkPolyData' to store points and normals
vtkNew<vtkPoints> dummy_points;
vtkNew<vtkPolyData> dummy_polydata;
for(vtkIdType i=0;i<pointsInCellCenterCalcSun->GetNumberOfPoints();++i){
// Get the coordinates of sun's cell center
double pointSun[3];
cellCenterCalcSun->GetOutput()->GetPoint(i,pointSun);
// Get normal vector at that cell
double* normals=normalsSun->GetTuple(i);
//cout<<normals[0]<<","<<normals[1]<<","<<normals[2]<<endl;
// Calculate the 'target' of the ray based on 'RayCastLength'
double pointEarth[3];
pointEarth[0]=pointSun[0]+normals[0]*RayCastLength;
pointEarth[1]=pointSun[1]+normals[1]*RayCastLength;
pointEarth[2]=pointSun[2]+normals[2]*RayCastLength;
int code=obbEarth->IntersectWithLine(pointSun,pointEarth,intersectPoints,intersectCells);
if(code<0){
continue;
}
double intersection[3];
vtkIdType cellId;
double normalIntersect[3];
vtkNew<vtkLineSource> line;
vtkNew<vtkLineSource> lineRay;
for(int i=0;i<intersectPoints->GetNumberOfPoints();++i){
if(i==1)
continue;
intersectPoints->GetPoint(i,intersection);
//cout<<"Intersection "<<i<<" : "<<intersection[0]<<","<<intersection[1]<<","<<intersection[2]<<endl;
//
cellId=intersectCells->GetId(i);
//cout<<"\tCellId "<<i<<" : "<<cellId<<endl;
normalIntersect[0]=normalsEarth->GetTuple(cellId)[0];
normalIntersect[1]=normalsEarth->GetTuple(cellId)[1];
normalIntersect[2]=normalsEarth->GetTuple(cellId)[2];
normalsVector->InsertNextTuple3(normalIntersect[0],normalIntersect[1],normalIntersect[2]);
dummy_points->InsertNextPoint(intersection[0],intersection[1],intersection[2]);
// line
line->SetPoint1(pointSun);
line->SetPoint2(intersection);
// Calculate the incident ray vector
double vecInc[3];
vecInc[0]=intersection[0]-pointSun[0];
vecInc[1]=intersection[1]-pointSun[1];
vecInc[2]=intersection[2]-pointSun[2];
// Calculate reflected ray vector
double vecRef[3];
double dot=vecInc[0]*normalIntersect[0]+vecInc[1]*normalIntersect[1]+vecInc[2]*normalIntersect[2];
vecRef[0]=vecInc[0]-2.0*dot*normalIntersect[0];
vecRef[1]=vecInc[1]-2.0*dot*normalIntersect[1];
vecRef[2]=vecInc[2]-2.0*dot*normalIntersect[2];
// Calculate the target of reflected ray
double pointRayReflectedTarget[3];
pointRayReflectedTarget[0]=intersection[0]+vecRef[0]*RayCastLength;
pointRayReflectedTarget[1]=intersection[1]+vecRef[1]*RayCastLength;
pointRayReflectedTarget[2]=intersection[2]+vecRef[2]*RayCastLength;
lineRay->SetPoint1(intersection);
lineRay->SetPoint2(pointRayReflectedTarget);
}
vtkNew<vtkPolyDataMapper> lineMapper;
vtkNew<vtkPolyDataMapper> lineRayMapper;
lineMapper->SetInputConnection(line->GetOutputPort());
lineRayMapper->SetInputConnection(lineRay->GetOutputPort());
vtkNew<vtkActor> lineActor;
vtkNew<vtkActor> lineRayActor;
lineActor->SetMapper(lineMapper);
lineActor->GetProperty()->SetColor(ColorRayHit);
lineRayActor->SetMapper(lineRayMapper);
lineRayActor->GetProperty()->SetColor(ColorRayReflected);
ren->AddActor(lineActor);
ren->AddActor(lineRayActor);
}
dummy_polydata->SetPoints(dummy_points);
dummy_polydata->GetPointData()->SetNormals(normalsVector);
vtkNew<vtkGlyph3D> glyphEarth;
glyphEarth->SetInputData(dummy_polydata);
glyphEarth->SetSourceConnection(arrow->GetOutputPort());
glyphEarth->SetVectorModeToUseNormal();
glyphEarth->SetScaleFactor(5);
vtkNew<vtkPolyDataMapper> glyphMapperEarth;
glyphMapperEarth->SetInputConnection(glyphEarth->GetOutputPort());
vtkNew<vtkActor> glyphActorEarth;
glyphActorEarth->SetMapper(glyphMapperEarth);
glyphActorEarth->GetProperty()->SetColor(ColorEarthGlyphs);
ren->AddActor(glyphActorEarth);
vtkNew<vtkCamera> camera;
camera->SetPosition(RadiusEarth,DistanceSun,RadiusEarth);
camera->SetFocalPoint(0.0,0.0,0.0);
camera->SetViewAngle(30.0);
ren->SetActiveCamera(camera);
//
vtkNew<vtkRenderWindow> renWin;
renWin->SetSize(1024,1024);
renWin->AddRenderer(ren);
vtkNew<vtkRenderWindowInteractor> iren;
iren->SetRenderWindow(renWin);
vtkNew<vtkInteractorStyleTrackballCamera> style;
iren->SetInteractorStyle(style);
iren->Initialize();
iren->Start();
return 0;
}
在代码中有详细的说明,等有空我会讲解一下。