利用VTK对数据进行可视化的实例

使用Visualization Toolkit 对数据进行可视化十分方便，程序的结构也正如我们所设计的可视化流程一样，十分简洁明了。下面我们将利用Visualization Toolkit，就可视化应用研究中一个十分热门的领域——计算机断层扫描CT， 简单剖析一下利用VisualizationToolkit 进行数据可视化的实现。
        我们将要进行处理的CT 数据是一个人体头部的切片数据，共有93 个切片。切片的间距是1.5mm ，每个切片由有12 个灰度级的间距为0.8 毫米的256*256 像素构成。我们打算由这些切片数据恢复出皮肤和骨骼的表面。为此我们首先应选取合适的算法，考虑到我们所处理的数据量是非常庞大的，超过了12 兆比特，而且我们只打算对表面进行重建，所以我们选择基于表面重建的经典算法Marching Cubes 算法。而对于Marching Cubes算法，Visualization Toolkit 中已经有封装好的vtkMarchingCubes 给予支持，这样就进一步简化了我们的工作。

       Marching Cube 算法是Lorensen 等人于1987 提出的，是三维数据场等值面生成的经典算法，是体素单元内等值面抽取技术的代表，并一直沿用至今。
       Marching Cube 算法的基本思想是在数据体中将位于两个相临切片上的8 个相临的体素构成一个立方体（cube），逐个处理数据场中的立方体，分类出与等值面相交的立方体，采用插值计算出等值面与立方体边的交点。根据立方体每一顶点与等值面的相对位置，将等值面与立方体边的交点按一定方式连接生成等值面，作为等值面在该立方体内的一个逼近表示。
4.1 读取数据
       首先，我们要做的事情是读取切片数据，并将其转换为我们的开发工具VisualizationToolkit 所支持的一种数据表达形式；然后根据其物理结构建立起相应的模型，我们给CT数据建立的是比较抽象的等值面模型；最后将物理组件与抽象的模型结合在一起来建立对CT 数据的可视化，以帮助用户正确理解数据。我们所要进行处理的是有结构点阵数据，其拓扑和几何都是隐含知道的，所以我们只需要知道数据的维数、数据源和数据空间。利用Visualization Toolkit 中的vtkVolume16Reader 我们可以很方便的读取切片数据，只需要告诉读取数据对象我们的CT 数据的一些参数，如切片之间的间距、切片上像素之间的间距以及所读取切片的起始段（如从第1 个切片到45 个切片），读取数据的代码如下所示：

 

vtkVolume16Reader  * Reader  =  vtkVolume16Reader::New();              // 建立一个读取对象
       Reader -> SetDataDimensions( 256 , 256 );                                  // 设置数据的维数
       Reader -> SetDataByteOrderToLittleEndian ();
       Reader -> SetFilePrefix ( " .. /.. /.. /vtkdata/headsq/quarter " );       // 设置所读取切片数据文件的路径
      Reader -> SetImageRange( 1 ,  93 );                                            // 设置读取切片的起始段
      Reader -> SetDataSpacing ( 0.8 ,  0.8 ,  1.5 );                                // 设置切片之间的间距和像素之间的间距 

 

 

 

 

4.2 提取等值面
       下面我们就可以用Marching Cubes 算法对所读取的数据进行处理了。首先利用vtkMarchingCubes 类来提取出某一CT 值的等值面，但这时的等值面其实仍只是一些三角面片，还必须由vtkStripper 类将其拼接起来形成连续的等值面。这样就把读取的原始数据经过处理转换为应用数据，也即由原始的点阵数据转换为多边形数据然后由 vtkPolyDataMapper 将其映射为几何数据，并将其属性赋给窗口中代表它的演员，将结果显示出来。

 

//  从切片数据中提取出皮肤
     vtkMarchingCubes  * skinExtractor  =  vtkMarchingCubes::New();        // 建立一个Marching Cubes 算法的对象
     skinExtractor -> SetInput(Reader -> GetOutput());                           // 获得所读取的CT 数据
     skinExtractor -> SetValue( 0 ,  500 );                                               // 提取出CT 值为500 的皮肤
     vtkStripper  * skinStripper  =  vtkStripper::New();                            // 建立三角带对象
     skinStripper -> SetInput(skinExtractor -> GetOutput());                    // 将生成的三角片连接成三角带
     vtkPolyDataMapper  * skinMapper  =  vtkPolyDataMapper::New();      // 建立一个数据映射对象
     skinMapper -> SetInput(skinStripper -> GetOutput());                      // 将三角带映射为几何数据
     vtkActor  * skin  =  vtkActor::New();                                            // 建立一个代表皮肤的演员
     skin -> SetMapper(skinMapper);                                                 // 获得皮肤几何数据的属性
     skin -> GetProperty() -> SetDiffuseColor( 1 , . 49 , . 25 );                      // 设置皮肤颜色的属性
     skin -> GetProperty() -> SetSpecular(. 3 );                                     // 设置反射率
     skin -> GetProperty() -> SetSpecularPower( 20 );                            // 设置反射光强度 

 

 

 

 

 

利用同样的方法，我们也可以提取出骨骼的等值面。只是骨骼的CT 值是1150 左右而已。所以只要在SetValue()方法中将参数设置为1150 就可以了。而且Visualization Toolkit支持多表面重建，所以在实际应用中我们可以设置多个参数值，提取出多个等值面并同时显示出来。在这个应用实例中我们只对皮肤和骨骼地的等值面进行了重建。
4.3 显示结果
       通过前面这些工作，我们基本上已经完成了对数据的读取处理映射等步骤，下面我们就要对数据进行显示了。

     vtkRenderer  * ren  =  vtkRenderer::New();                                                // 建立绘制者
     vtkRenderWindow  * renWindow  =  vtkRenderWindow::New();                      // 建立绘制窗口
     renWindow -> AddRenderer(ren);                                                           // 将绘制者加入绘制窗口
     vtkRenderWindowInteractor  * iren  =  vtkRenderWindowInteractor::New();    // 对绘制结果进行交互操作
     iren -> SetRenderWindow (renWindow);                                                  //  告诉绘制者将要在绘制窗口中进行显示的演员
     ren -> AddActor(skin);                                                                        // 皮肤
     ren -> AddActor(bone);                                                                      // 骨骼 

 

 

 

 

 

 

 

//在我的VTK 5.0.2中，相应的测试程序代码如下：

#include "vtkVolume16Reader.h"
       #include "vtkRenderWindowInteractor.h"
       #include "vtkRenderer.h"
       #include "vtkRenderWindow.h"
      #include "vtkMarchingCubes.h"
      #include "vtkStripper.h"
      #include "vtkActor.h"
      #include "vtkPolyDataMapper.h"
      #include "vtkProperty.h"
//读取RAW文件，提取等值面。

int main()
{
      vtkVolume16Reader *reader=vtkVolume16Reader ::New();
      reader->SetDataDimensions(64,64);
      reader->SetDataByteOrderToLittleEndian();
      reader->SetFilePrefix("D://headsq//quarter");
      reader->SetImageRange(1,93);
      reader->SetDataSpacing(3.2,3.2,1.5);

vtkMarchingCubes *boneExtractor=vtkMarchingCubes::New();
      boneExtractor->SetInput((vtkDataObject *)reader->GetOutput());
      boneExtractor->SetValue(0,500);
      vtkStripper *boneStripper=vtkStripper::New();
      boneStripper->SetInput(boneExtractor->GetOutput());
 
       vtkPolyDataMapper *boneMapper=vtkPolyDataMapper::New();
       boneMapper->SetInput(boneStripper->GetOutput());

vtkActor *bone=vtkActor::New();
       bone->SetMapper(boneMapper);

bone->GetProperty()->SetDiffuseColor(.1,.94,.52);
       bone->GetProperty()->SetSpecular(.3);
       bone->GetProperty()->SetSpecularPower(20);
       vtkRenderer *ren=vtkRenderer::New();
       vtkRenderWindow *renWindow=vtkRenderWindow::New();
       renWindow->AddRenderer(ren);
       vtkRenderWindowInteractor *iren=vtkRenderWindowInteractor::New();
       ren->AddActor(bone);
       iren->Initialize();
       iren->Start();
       reader->Delete();
       iren->Delete();

       return 0;
}