医学图像三维重建体绘制算法研究

贺楠楠，李向阳

(1.汝州职业技术学院，河南平顶山 467500；2.河南女子职业学院，河南郑州 450000)

1 医学图像三维重建体绘制原理

医学图像三维重建分为面绘制与体绘制。面绘制的方法有立方块法、移动立方体法、移动四面体法、剖分立方体法和表面跟踪法，最常用的是移动立方体法[1]。体绘制方法的原理是模拟光线穿过半透明的三维数据场，并计算三维数据场中的数据点对屏幕像素的投影贡献，即把三维离散的数据场按照一定的规则转换为图像设备屏幕上的二维图像[2]。

2 医学图像三维重建体绘制算法

医学图像三维重建体绘制分为从图像空间开始处理的体绘制、从物体空间开始处理的体绘制和变换域体绘制。从图像空间开始处理的体绘制是从二维图像像素开始处理的体绘制，最常用的算法是光线投射法。从物体空间开始处理的体绘制是从三维数据场中的数据点开始进行处理的体绘制。其原理为按顺序计算每个数据点对二维投影图像像素的贡献，并加以合成，形成最终的二维投影图像。常用的算法有抛雪球法、错切形变法和纹理映射法等。变换域体绘制有频域体绘制法和基于小波的体绘制法[3]。

2.1 光线投射法

光线投射法的基本原理：三维数据场为规则数据场，即三维数据点均匀分布在规则网格的网格点上。先把三维数据点按灰度值分类，再对数据点赋予相应的颜色值和不透明度值，然后重采样并合成图像。

1)数据分类

三维数据点要分为不同的类别来表示不同的物质。CT图像常使用阈值分类法把数据分为不同的类别。阈值分类法设定一个区间内的数据点为人体组织中的同一种物质，多个区间表示不同的器官组织。若三维数据场中数据点灰度值为I(x,y,z)，则满足如下条件的数据点被归为同一类。

2)颜色和不透明度赋值

二维断层图像为灰度图像，人眼对彩色图像的敏感度更高。不同的器官组织用不同的颜色标识会增加辨识度。根据数据分类的结果给不同灰度的数据点赋予不同的颜色值，灰度图像转化为伪彩色图像。

颜色赋值是把灰度值按照一定的函数关系变换为彩色值。常用的一种映射函数如图1 所示，可调整min、a1、a2、a3、max 的值，使物质的分类区间对应相应的色带，即每个色带区对应某种物质[4]。

图1 颜色映射函数

体绘制要显示三维数据场的内部结构，还要根据分类结果为每个数据点的灰度值赋予不透明度值α，α 值一般设定为跟灰度值呈线性关系，并且α值在0-1范围内。

3)重采样

按照一定的视线方向，从二维图像的每个像素点发出一条射线，穿过三维数据场。沿这条射线取n个等距的采样点，计算采样点的颜色值和不透明度值[5]。图像所处平面是二维坐标空间，而采样点所处的是三维坐标空间，因而要进行坐标系的转换。

根据视线方向从图像的像素点(a，b)(物体空间中为点t)发出的射线S为：

式中，d为采样间隔，k为采样点序号。

采样点邻近的8 个数据点为q0、q1、q2、q3、q4、q5、q6、q7，其对应的颜色值(或不透明度值)为m0、m1、m2、m3、m4、m5、m6、m7，点r11、r12、r21、r22为采样点所在平面(平行于坐标轴平面方向)截8个数据点所在立方体的分点，由线性插值法可得：

令r1、r2处的颜色值(或不透明度值)为t1、t2,由线性插值法可得：

q点的颜色值(或不透明度值)为m，由线性插值得：

图2 三线性插值

4)图像合成

图3 光线投射法体绘制

图像合成的迭代算子分为由后向前的合成和由前向后的合成。由后向前的合成是从离视点最远的采样点开始顺序合成到离视点最近的采样点，而由前向后的合成则相反。如图4所示，未参与合成的第i个采样点的颜色值和不透明度值为Cnow和αnow，将要与i合成的前边所有采样点合成后的颜色值和不透明度值为Cin和αin，与第i个采样点合成后的颜色值和不透明度值为Cout和αout。

图4 颜色和不透明度的合成

由后向前的迭代公式如下：

由前向后的迭代公式如下：

由前向后的合成中不透明度值会慢慢增大，接近1 时，未参与合成的采样点的颜色值和不透明度值将不会对像素点有太大贡献，可忽略计算并加快绘制速度。光线投射法绘制速度比较慢、有半透明的效果并可显示内部结构，因此有较高的图像品质，但占用内存空间较大。

2.2 抛雪球法

抛雪球法把三维数据场中每个数据点看作一个能量源，当每个数据点向二维图像平面投影时，数据点的能量会扩散到二维图像的像素点上。抛雪球法先对三维数据场的数据点进行分类并赋予颜色值和不透明度值，再计算三维数据场数据点对二维图像像素的贡献。通过对图像像素有影响的采样点重构核的空间卷积域积分来计算像素值，积分的次数是重构核空间卷积域的平面投影区域内的像素数乘以采样点数。计算量比较大，所以一般使用足迹函数确定在任意观察方向三维采样点重构核空间卷积域的平面投影范围，然后查询足迹表获得采样点对屏幕像素的贡献值，来提高绘制速度。抛雪球法的绘制速度比光线投射法速度快。

2.3 错切变形法

错切变形法根据所设定的观察方向分别计算错切矩阵和变形矩阵，对构成三维数据场的断层图像序列按错切矩阵进行错切变换后投影合成，获得中间图像，然后对中间图像使用变形矩阵变换得到投影图像。错切变形法把重采样从三维降低到二维，减少了计算量，有接近实时的绘制速度，是最快的体绘制方法，并且不显著降低图像质量。

2.4 纹理映射法

纹理映射体绘制分为二维纹理映射和三维纹理映射。二维纹理映射体绘制生成三个方向(三维数据平行于三个坐标平面方向)上的二维纹理切片并合成图像。三维纹理映射体绘制先使用转换函数把体数据转换为具有一定的颜色值、光强度和不透明度值来构造三维纹理，再按照视线方向沿三个坐标轴Z轴、Y轴和X轴的顺序依次装入纹理切片，然后确定三维空间的坐标与纹理空间坐标的对应关系，最后把三维纹理切片从后向前堆叠合成图像。纹理映射中的计算由硬件实现，提高了运算速度[6]。

2.5 频域体绘制法

频域体绘制分为傅立叶体绘制、哈特里体绘制和傅立叶小波体绘制。傅立叶体绘制把三维数据通过傅立叶变换得到三维频域数据场，空间域沿任意方向θ的投影与频域三维数据场中过原点沿同一方向的切片形成傅立叶变换对，即傅立叶投影——切片定理。体绘制图像为三维数据场在二维平面上的投影。在三维空间域数据场转换后的三维频域场中截取与空间域中相对应的切片，并使用傅立叶逆变换来得到投影图像。傅立叶变换把三维数据变换到频域后，空间域的每个数据点都要用复数表示，比较耗内存。哈特里体绘制中三维数据变换到频域后是实数数据，所以计算速度更快，内存占用空间小。傅立叶小波体绘制是使用傅立叶体绘制方法得到三维频域数据和切平面提取，在频域切平面上实现小波的多分辨率分解，并在不同分辨级上重建小波系数，再通过傅立叶逆变换得到投影图像。频域体绘制中三维数据场的重采样和信号重建都是在频率域中完成的。频域算法的计算速度优于空间域算法，但图像没有半透明效果。

2.6 基于小波的体绘制法

基于小波的体绘制有小波域光线投射法和小波域足迹法两种算法。小波域光线投射法是把空间域光线投射法的光照模型在频域上求解，颜色值和不透明度值都可以在小波域求出。小波域光线投射法绘制速度慢，可得到高品质、半透明、有明暗度的图像，但内存开销比较大。小波域足迹法是把体数据分解为不同分辨率上的离散平滑逼近和离散细节信号，用小波作为重构滤波器，用平滑逼近信号重建出低分辨率级上的三维图像，然后增加细节信号逼近绘制出三维图像。小波足迹法绘制速度快，图像品质不如光线投射法和小波域光线投射法[7]。

体绘制比面绘制具有更慢的绘制速度，可以显示物体内部的结构。基于CUDA的并行计算[8]和基于GPU的并行编程平台可以加快体绘制的速度[9]。VTK 是开源的、跨平台的可视化函数库，可实现医学图像面绘制和体绘制的多种算法[10]。OpenGL是三维工具软件包，也可以实现面绘制与体绘制。