(VTK) DICOM

複雑なフォーマットのDICOM画像をVTKで使う。

昔フルスクラッチでDICOMを読むプログラムを書いたが非常に面倒。 マルチスライスの画像になると更に複雑になる。 vtkDICOMImageReaderを使うとあっさり画像が取得できる。

サンプル画像

公開されているDICOM画像

前のリンクは404になっていたのでMedical Image Format FAQ - Part 8のDICOM image samplesを参照。

旧情報

ネットで検索すればあちこちに転がっている。

• DICOM sample image sets

から拝借。

比較的サイズが小さいPNEUMATIX.zip)をソースのあるディレクトリで展開する。 使用する画像のあるディレクトリは

PNEUMATIX/Cardiovascular Heart-Cardiac Function/fl3d-cor

基本処理

vtkのMLで出た使い方

    #include <vtkDICOMImageReader.h>
    #include <vtkImageViewer.h>
    
    int main(){
      vtkDICOMImageReader *reader = vtkDICOMImageReader::New();
      vtkImageViewer     *viewer = vtkImageViewer::New();
      
      viewer->SetInput(reader->GetOutput());
      reader->SetDirectoryName("PNEUMATIX/Cardiovascular Heart-Cardiac Function/fl3d-cor");
      reader->Update();
      
      int VolData_Images = viewer->GetWholeZMax ();
      
      
      for (int i=0;i<=VolData_Images;i++)
        {
          viewer->SetZSlice(i);
          viewer->SetColorWindow(400);
          viewer->SetColorLevel(128);
          viewer->Render();
        }
      reader->Delete();
      viewer->Delete();
      return 0;
    }

これだけ仮にvtkDICOM.cxxとでもして保存する。 MakefileはMakefileのひな型を参考に

CC=gcc
CFLAGS = -Wno-deprecated -I /usr/X11R6/include -I /usr/local/include/vtk-5.0
LDLIBS = -L/usr/X11R6/lib -L /usr/local/lib -lvtkRendering -lvtkIO\
         -lvtkGraphics -lvtkImaging -lvtkftgl -lvtkfreetype -lGL -lXt -lSM\
         -lICE -lSM -lICE -lSM -lICE -lX11 -lXext -lX11 -lXext -lX11 -lXext\
         -lvtkFiltering -lvtkCommon -lvtksys -lpthread -lm -lvtkDICOMParser\
         -lvtkpng -lvtktiff -lvtkzlib -lvtkjpeg -lvtkexpat -lvtkMPEG2Encode

vtkDICOM:

としておく。 あとは

make

とするだけ。

DICOM画像の全スライスが次々と表示されたら成功。

VolumeRendering

上記のreader->GetOutput()さえ得られれば、VolumeRenderingは出来たも同然っぽいので、VTKのサンプルであるMedical1.cxxのコードに被せてみる。

    #include "vtkRenderer.h"
    #include "vtkRenderWindow.h"
    #include "vtkRenderWindowInteractor.h"
    #include <vtkDICOMImageReader.h>
    #include "vtkPolyDataMapper.h"
    #include "vtkActor.h"
    #include "vtkOutlineFilter.h"
    #include "vtkCamera.h"
    #include "vtkProperty.h"
    #include "vtkPolyDataNormals.h"
    #include "vtkContourFilter.h"
    //for Capture
    #include "vtkWindowToImageFilter.h"
    #include "vtkPNGWriter.h"
    //for Capture
    
    int main(){
    
      vtkRenderer *aRenderer = vtkRenderer::New();
      vtkRenderWindow *renWin = vtkRenderWindow::New();
      renWin->AddRenderer(aRenderer);
      
      vtkRenderWindowInteractor *iren = vtkRenderWindowInteractor::New();
      iren->SetRenderWindow(renWin);
      
      vtkDICOMImageReader *reader = vtkDICOMImageReader::New();
      reader->SetDirectoryName("PNEUMATIX/Cardiovascular Heart-Cardiac Function/fl3d-cor");
      reader->Update();
      
      vtkContourFilter *skinExtractor = vtkContourFilter::New();
      skinExtractor->SetInputConnection(reader->GetOutputPort());
      skinExtractor->SetValue(0, 500);
      
      vtkPolyDataNormals *skinNormals = vtkPolyDataNormals::New();
      skinNormals->SetInputConnection(skinExtractor->GetOutputPort());
      skinNormals->SetFeatureAngle(60.0);
      
      vtkPolyDataMapper *skinMapper = vtkPolyDataMapper::New();
      skinMapper->SetInputConnection(skinNormals->GetOutputPort());
      skinMapper->ScalarVisibilityOff();
      
      vtkActor *skin = vtkActor::New();
      skin->SetMapper(skinMapper);
      
      vtkOutlineFilter *outlineData = vtkOutlineFilter::New();
      outlineData->SetInputConnection(reader->GetOutputPort());
      
      vtkPolyDataMapper *mapOutline = vtkPolyDataMapper::New();
      mapOutline->SetInputConnection(outlineData->GetOutputPort());
      
      vtkActor *outline = vtkActor::New();
      outline->SetMapper(mapOutline);
      outline->GetProperty()->SetColor(0,0,0);
      
      vtkCamera *aCamera = vtkCamera::New();
      aCamera->SetViewUp (0, 0, -1);
      aCamera->SetPosition (0, 1, 0);
      aCamera->SetFocalPoint (0, 0, 0);
      aCamera->ComputeViewPlaneNormal();
      
      aRenderer->AddActor(outline);
      aRenderer->AddActor(skin);
      aRenderer->SetActiveCamera(aCamera);
      aRenderer->ResetCamera ();
      aCamera->Dolly(1.5);
    
      aRenderer->SetBackground(1,1,1);
      renWin->SetSize(640, 480);
    
      aRenderer->ResetCameraClippingRange ();
    
      iren->Initialize();
      iren->Start(); 
    
      //for Capture
      vtkWindowToImageFilter *w2if=vtkWindowToImageFilter::New();
      w2if->SetInput(renWin);
      vtkPNGWriter *wr=vtkPNGWriter::New();
      wr->SetInputConnection(w2if->GetOutputPort());
      wr->SetFileName("Capture.png)");
      wr->Write();
      w2if->Delete();
      wr->Delete();
      //for Capture
      skinExtractor->Delete();
      skinNormals->Delete();
      skinMapper->Delete();
      skin->Delete();
      outlineData->Delete();
      mapOutline->Delete();
      outline->Delete();
      aCamera->Delete();
      iren->Delete();
      renWin->Delete();
      aRenderer->Delete();
      reader->Delete();
      return 0;
    }

冗談のような話だがあっさり成功。

比較すれば分かると思うが、殆どMedical1.cxxのままである。 換えたのは読み込みの部分とv16のインスタンスを使っている部分だけの10行もない。 "//for Capture"でくくった部分は画像キャプチャ用に仕込んだコードで、"q"を押して終了した状態の画像をpng形式で保存する。