VTK 圖像基本操作_三維圖像切片提取


1.三維圖像切片提取

切片是指三維圖像中的一個切面對應的圖像。切面可以是過圖像內部一點且平行於XY、YZ、XZ平面的平面,也可以是任意的過三維圖像內部一點任意方向的平面。通過提取切片可以方便的瀏覽和分析圖像內部組織結構,是醫學圖像瀏覽軟件中的一個重要的功能。在VTK中vtkImageReslice類實現圖像切片提取功能。

下面是切片提取的代碼:

 1 #include <vtkAutoInit.h>
 2 VTK_MODULE_INIT(vtkRenderingOpenGL);  3  
 4 #include <vtkSmartPointer.h>
 5 #include <vtkImageData.h>
 6 #include <vtkMetaImageReader.h>
 7 #include <vtkMatrix4x4.h> //  8 #include <vtkImageReslice.h>
 9 #include <vtkLookupTable.h>
10 #include <vtkImageMapToColors.h>
11 #include <vtkImageActor.h>
12 #include <vtkRenderer.h>
13 #include <vtkRenderWindow.h>
14 #include <vtkRenderWindowInteractor.h>
15 #include <vtkInteractorStyleImage.h>
16  
17 int main(int argc, char* argv[]) 18 { 19     vtkSmartPointer<vtkMetaImageReader> reader =
20         vtkSmartPointer<vtkMetaImageReader>::New(); 21     reader->SetFileName("brain.mhd"); 22     reader->Update(); 23  
24     int extent[6]; 25     double spacing[3]; 26     double origin[3]; 27  
28     reader->GetOutput()->GetExtent(extent); 29     reader->GetOutput()->GetSpacing(spacing); 30     reader->GetOutput()->GetOrigin(origin); 31  
32     double center[3]; 33     center[0] = origin[0] + spacing[0] * 0.5 * (extent[0] + extent[1]); 34     center[1] = origin[1] + spacing[1] * 0.5 * (extent[2] + extent[3]); 35     center[2] = origin[2] + spacing[2] * 0.5 * (extent[4] + extent[5]); 36     //*****************************************************************// 37     static double axialElements[16] = { 38         1, 0, 0, 0, 39         0, 1, 0, 0, 40         0, 0, 1, 0, 41         0, 0, 0, 1
42  }; 43  
44     vtkSmartPointer<vtkMatrix4x4> resliceAxes =
45         vtkSmartPointer<vtkMatrix4x4>::New(); 46     resliceAxes->DeepCopy(axialElements); 47     resliceAxes->SetElement(0, 3, center[0]); 48     resliceAxes->SetElement(1, 3, center[1]); 49     resliceAxes->SetElement(2, 3, center[2]); 50  
51     vtkSmartPointer<vtkImageReslice> reslice =
52         vtkSmartPointer<vtkImageReslice>::New(); 53     reslice->SetInputConnection(reader->GetOutputPort()); 54     reslice->SetOutputDimensionality(2); 55     reslice->SetResliceAxes(resliceAxes); 56     reslice->SetInterpolationModeToLinear(); 57     //*****************************************************************// 58     vtkSmartPointer<vtkLookupTable> colorTable =
59         vtkSmartPointer<vtkLookupTable>::New(); 60     colorTable->SetRange(0, 1000); 61     colorTable->SetValueRange(0.0, 1.0); 62     colorTable->SetSaturationRange(0.0, 0.0); 63     colorTable->SetRampToLinear(); 64     colorTable->Build(); 65     vtkSmartPointer<vtkImageMapToColors> colorMap =
66         vtkSmartPointer<vtkImageMapToColors>::New(); 67     colorMap->SetLookupTable(colorTable); 68     colorMap->SetInputConnection(reslice->GetOutputPort()); 69     //*****************************************************************// 70     vtkSmartPointer<vtkImageActor> imgActor =
71         vtkSmartPointer<vtkImageActor>::New(); 72     imgActor->SetInputData(colorMap->GetOutput()); 73  
74     vtkSmartPointer<vtkRenderer> renderer =
75         vtkSmartPointer<vtkRenderer>::New(); 76     renderer->AddActor(imgActor); 77     renderer->SetBackground(1.0, 1.0, 1.0); 78  
79     vtkSmartPointer<vtkRenderWindow> renderWindow =
80         vtkSmartPointer<vtkRenderWindow>::New(); 81     renderWindow->AddRenderer(renderer); 82     renderWindow->Render(); 83     renderWindow->SetSize(640, 480); 84     renderWindow->SetWindowName("Extract3Dslice"); 85  
86     vtkSmartPointer<vtkRenderWindowInteractor> rwi =
87         vtkSmartPointer<vtkRenderWindowInteractor>::New(); 88     vtkSmartPointer<vtkInteractorStyleImage> imagestyle =
89         vtkSmartPointer<vtkInteractorStyleImage>::New(); 90     rwi->SetInteractorStyle(imagestyle); 91     rwi->SetRenderWindow(renderWindow); 92     rwi->Initialize(); 93     rwi->Start(); 94  
95     return 0; 96 }

首先通過vtkMetaImageReader讀取一張醫學三維圖像,並獲取得到圖像范圍(extent),原點和像素間隔;由這三個參數可以計算圖像的中心位置center;接下來定義了切面的變換矩陣axialElements,該矩陣的前三列分別表示x、y和z方向向量,第四列為中心點坐標;

代碼中的axialElements表示切面變換矩陣與當前坐標系一致,且切面為過中心點center,並平行於XY平面的平面???當前,定義該切面時,也可以是其他平面,甚至是任意平面,但是必須要過圖像內部點。

下面給出了一個常用的變換矩陣。

提取平行於XZ平面的切片:

1 static double coronalElements[16] = { 2  1, 0, 0, 0, 3  0, 0, 1, 0, 4  0,-1, 0, 0, 5  0, 0, 0, 1 };

提取平行於YZ平面的切片:

1 static double sagittalElements[16] = { 2 0, 0,-1, 0, 3 1, 0, 0, 0, 4 0,-1, 0, 0, 5 0, 0, 0, 1 }; 

提取斜切切片:

1 static double obliqueElements[16] = { 2 1, 0, 0, 0, 3 0, 0.866025, -0.5, 0, 4 0, 0.5, 0.866025, 0, 5 0, 0, 0, 1 };

注意使用這些變換矩陣的時候,需要將第四列替換為切片經過圖像的一個點坐標,上例中將圖像的中心添加到axialElements矩陣,並通過函數SetResliceAxes設置變換矩陣,SetOutputDimensionality(2)指定輸出的圖像為一個二維圖像; 而函數SetInterpolationModeToLinear()則指定了切面提取中的差值方式為線性差值,另外該類中還提供了其他的插值方式:
SetInterpolationModeToNearestNeighbor():最近鄰方式
SetInterpolationModeToCubic():三次線性差值
設置完畢后,執行Update()即可完成切面計算。

東靈提供的預想結果應該是:

然而,在實際運行中,卻碰到了vtkMetaImageReader不能讀取文件的問題,在之前32bit系統上,這都是正常的,暫時還不能肯定問題出現在哪里,需要進一步研究!!問題如下:


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM