李丽玲,李彬△,丁焕文,田联房
(1.华南理工大学自动化科学与工程学院,广东广州510640;2.广州军区广州总医院骨科,广东 广州510010)
脊柱矫形手术涉及到复杂的骨结构和神经组织,手术难度大,风险高,要求医生在术前制定完整的手术规划。椎弓根螺钉内固定是治疗脊柱各种病变的重要手段,该技术能提高脊柱稳定性,有利于脊柱的骨性融合[1]。但由于椎弓根局部解剖变化多端,神经血管分布复杂,置入的椎弓根螺钉如果发生偏差则有可能导致严重的并发症[2]。
计算机辅助外科手术是国内外众多学者的研究热点[3]。借助可视化技术,可以方便地实现 CT、MRI等二维图像序列的三维重建。三维模型展示了整个脊柱结构,结合人机交互技术,可以对脊柱模型进行三维姿态调整、三维切割及三维几何配准等操作。但是,以目前市面上比较成熟的3D Slicer和Mimics为例,这些系统价格昂贵且功能复杂,不适用于椎弓根螺钉植入手术的辅助设计。
针对脊柱侧凸三维评估分析及个性化手术方案设计的需求,本研究对基于移动立方体(marching cubes,MC)算法的三维重建技术、切割技术和快速成型技术展开了较深入的研究,开发了基于可视化交互的计算机辅助椎弓根螺钉植入系统[4]。本研究主要介绍该可视化系统的设计思路、系统功能模块及其关键技术,并以实例验证系统的效果。
面向脊柱矫形应用需求,本研究设计的系统结构见图1。其中医学图像数据传输模块的功能是进行DICOM格式医学图像的读取、解析与保存;二维处理与交互模块主要涉及二维图像的显示、阈值分割、区域增长及像素擦除;三维面绘制与交互模块主要实现面绘制渲染、体切割、三维轮廓线提取及椎弓根螺钉虚拟置入操作;快速成型文件接口主要是STL文件的导入和导出[5]。系统功能由各个模块相互协同实现。
图1 系统总体结构框图Fig 1 Overall system block diagram
系统读入DICOM格式医学图像后,先对图像进行阈值分割、区域增长等预处理,为提高重建模型的质量做铺垫。然后基于MC面绘制算法对脊柱进行三维重建,得到完整的脊柱模型,用于辅助医生进行脊柱病变三维评估分析。根据手术方案设计的需要,对脊柱进行三维切割,或者综合使用像素擦除功能与多平面切割功能对单独脊椎骨进行分割。基于分离得到的脊椎骨块,选择合适的剖切平面对其进行剖切,得到截面点云,并映射到二维平面,生成二维轮廓线。以二维轮廓线作为参照,完成椎弓根螺钉植入位置的设计。系统设有快速成型文件接口功能,可以导入STL文件,进行相关三维操作。也可以在生成三维模型到完成手术方案设计的过程中,导出三维模型为STL文件,并通过3D打印生成实物模型,进一步辅助医生进行手术模拟[6]。
为模拟矫正过程中各个相关脊椎骨产生的位移和旋转,需要将单个脊椎骨从脊柱中分离出来[7]。根据脊椎骨间距小、形状不规则且灰度区分度小的特点,本研究提出一种结合二维CT图像像素擦除与三维模型的多平面切割的算法对脊椎骨进行分割,分别对应于图1中二维图像像素擦除和三维模型多平面切割两个功能模块。
算法流程见图2。
基于三维连通域标记的二维图像像素擦除算法,其核心在于改变三维连通域的标记状态,从而实现体素的分割[8]。
具体步骤描述如下:
(1)对二维CT图像序列进行阈值分割和区域增长,得到最初的三维连通域∑,见图3;
(2)基于连通域∑,在目标脊椎骨与非目标脊椎骨紧密相连的像素区域上选择两个点A、B,并以此为对角点生成相应的覆盖矩阵;
(3)矩阵覆盖下的像素区域记为Γ,以此改变连通域的标记状态,见图3中蓝色区域。
在二维CT序列上进行多次擦除标记后,最初的连通域被划分为两个或者多个不同连通域。目标脊椎骨所在的连通域记为Ψ,Ψ=∑-Γ。对φ进行区域增长和三维重建,此时骨块间存在的冗余连接已被剔除,见图4箭头所指。
图3 二维图像像素擦除示意图Fig 3 Schematic diagram of Pixel Erasing based on the two-dimensional CT images
图4 擦除关键相连像素后的三维重建模型Fig 4 3D model reconstructed with erasered CT images
在脊柱矫正模拟手术中,基于三维模型的切割是一个重要的环节。目前主要有片元级模型切割和体素级切割两种方式[9-10]。由于片元级切割运算量大,效率低,所以,本系统采用体素级的切割方式。
如图5所示,体素级切割遍历三维数据场中的体素顶点,根据体素顶点与剖切平面Φ的空间几何关系对其进行标记,规定与剖切面Φ正法向量指向方向一致的点标记为1,否则标记为0;标记完毕以后,采用MC算法对标注为1的数据进行三维重建,得到切割后的三维模型[11]。
算法中的剖切平面Φ是通过在屏幕上选两点绘制出切割线,并以垂直于屏幕的方向,过切割线所在的直线生成的。根据操作需要,选择两个或者多个点,绘制出单平面或者多平面作为切割面。若是多平面切割,则根据选点的先后顺序P1→P2→P3→P4→…绘制出连续的平面,如图6所示,与所有平面正法向量一致的点标记为1,否则为0。
图5 单平面切割示意图Fig 5 Schematic diagram of single-plane cutting
图6 多平面切割示意图Fig 6 Schematic diagram of multi-plane cutting
针对脊柱侧凸疾病的个体特异性,本研究提出的是一种计算机辅助设计椎弓根螺钉植入位置的方法,对应于图一中三维轮廓线提取和椎弓根螺钉虚拟植入两个功能模块。
实现步骤描述如下:
(1)交互在脊椎骨上选取三个特征点,生成三维剖切面;
(2)遍历脊椎骨模型的三角面片,求剖切面与三角面片三条边的交点,得到截面点云,依次拟合所有交点构成脊椎骨的三维剖切轮廓模型;
(3)对三维轮廓模型进行空间几何变换,映射到xoy平面上,绘制成二维平面的轮廓;
(4)基于脊椎骨的二维轮廓线设计椎弓根螺钉的起止点。
椎弓根螺钉通常植入于脊椎骨的特征剖切面上,并以该脊椎骨的三维轮廓线作为参照坐标系[12]。基于2.2节得到的脊椎骨,在可视化窗口中对模型进行交互选点,根据图形管线逆变原理将3个屏幕坐标点变换至三维物体空间坐标,从而计算出剖切面 α的法向量为=(A,B,C)。如图 7所示,遍历构成模型的三角面片,判断三角面片的两两顶点与剖切面的位置关系,若两顶点位于剖切面的不同侧,则该两点构成的三角边与剖切面有交点。根据空间几何原理,求出交点坐标,拟合得到表面模型位于该剖切平面的三维轮廓线。
图7 三维轮廓线提取示意图Fig 7 Schematic diagram of three-dimensional contour extraction
由于在三维的空间中不易于选点,难以保证两个椎弓根螺钉位于同一平面,因此,需要把三维轮廓映射到二维平面[13],即通过空间几何变换使剖切面Φ的法向量与 z轴重合,见图8。
图8 空间几何变换示意图Fig 8 Schematic diagram of geometric transformation
具体的变换步骤如下:
(1)取模型上交互选取的第一个点 O(x0,y0,z0)为剖切面法向量=(A,B,C)的起始点。
对构成三维轮廓的点集进行上述坐标变换后,即可实现从三维空间到二维平面的映射。因为屏幕的大小有限,在绘制二维轮廓之前先对变换后的x、y左边取整,然后求出点集中x、y坐标的最大值和最小值之差|xmax-xmin|、|ymax-ymin|,按比例对二维坐标进行缩放调整使两个差值都不超过屏幕的最大显示长度。
基于得到的二维轮廓线,通过屏幕上点击选取二维轮廓线的两点为螺钉的首末两端,通过空间几何变换,把二维空间坐标映射到最初的三维空间,变换过程是三维轮廓到二维轮廓的逆过程。然后选择螺钉的半径,即可对虚拟螺钉的三维效果进行渲染绘制。
为测试本系统临床应用效果及验证系统各模块的有效性,采用广州军区广州总医院提供的某脊柱侧弯畸形患者的CT切片图像序列作为体数据源,使用本系统进行脊柱侧凸三维评估分析及椎弓根螺钉内固定手术方案设计。每组实验数据的分辨率为512×512×214。实验环境:HP Z800 Workstation;Windows 7 64位操作系统;CPU Intel Xeon(R)X5690,主频 3.47GHz;内存 192GB DDR3 1333 MHz;显卡Nvidia Quadro 6000 4G显存。
采用MC算法进行三维重建之前先对CT序列进行预处理操作。首先通过阈值分割获取骨骼部分,以人体骨骼CT值近似范围[200,4000]Hu为标准选择阈值上下限,分割得到所有骨骼数据。接着进行区域增长,去掉多余的碎骨。最后进行三维重建,得到完整的脊柱模型,见图9。
图9 MC算法三维重建效果Fig 9 3D reconstruction renderings by MC algorithm
重建得到的脊柱表面模型较为完整、光滑,没有碎骨等多余组织,同时模型渲染速度快,满足系统实时交互操作的要求。
为了进一步观察脊柱侧凸部位的病变情况,只选择该段的CT图片进行重建,结果见图10,此时根据三维模型初步判断病情为L1半椎体畸形(箭头所指),局部发生侧弯畸形,凸向右侧。
以5.1节得到的半段脊柱作为实验对象,先用基于三维模型的多平面切割初步分离L1半椎体,切割效果见图11。图11(a)为单平面切割效果图,(b)为多平面切割效果图,(c)为角形切割效果图。三种切割方式综合使用,可以满足任意角度,任意形状的快速切割。通过基于三维模型的切割,基本上可以初步分割出目标脊椎骨。
图10 侧凸部位MC算法三维重建效果Fig 10 Scoliosis site 3D reconstruction renderings by MC algorithm
图11 三维模型切割效果图(a)单平面切割;(b)多平面切割;(c)角形切割Fig 11 Three-dimensional model cutting renderings
再应用基于二维CT图像的像素擦除可以把脊椎骨之间紧密相扣的多余部分擦除,结合上述两个功能,可以快速完成脊椎骨的精确分割。
为了更好地验证本系统脊椎骨分割功能的准确性,使用本系统对L2半椎体进行分割实验,结果记为L2-1。同时医生使用Mimics系统分割同一脊椎骨,结果记为L2-2,并以此作为参考标准。然后使用Mimics系统的3D测量功能,分别测量两个脊椎骨模型的两个指标(横径IPD,矢径SD)以及两对特征点(左横突点与右横突点,棘突点与椎体上边缘点)的 3D距离,测量过程见图 12,结果见表1[14-15]。分析实验数据发现:分割得到的模型完整光滑,没有碎骨,没有非目标脊椎骨的多余部分。多次测量,平均偏差率为5.84‰,满足精确性的要求。
图12 测量过程示意图Fig 12 Schematic diagram of measure process
表1 L2半椎体测量结果Table 1 L2 hemi vertebra measurement results
以L2半椎体为实验对象,先通过选点确定剖切面,从而获得截面点云,拟合为三维轮廓线,见图13。
图13 L2半椎体三维轮廓线Fig 13 L2 hemi vertebra 3D contour
把三维轮廓线映射到二维平面,基于二维轮廓线选择椎弓根的起止点,设计植入位置和半径,见图14。图14a表示在二维轮廓线上选点确定椎弓根螺钉的植入位置,b为植入螺钉后整个模型的三维透视效果,c是同一观察角度的非透视效果。二维空间和三维空间的一一映射关系进一步确保了模拟手术的准确性。
图14 基于二维轮廓线椎弓根螺钉设计效果a.二维效果;b.三维透视效果;c.植入螺钉的三维模型Fig 14 Pedicle screw design effect based on the two-dimensional contour
本研究介绍了基于可视化交互的计算机辅助椎弓根螺钉植入系统的设计结构、核心技术和功能模块实现,并以实际病例验证了该系统的性能。不同于其他的计算机辅助系统,本系统侧重于脊椎骨的快速分割和椎弓根螺钉的模拟设计。实验结果表明,该系统满足脊柱矫形手术的直观性和交互性需求,可辅助医生进行脊柱侧凸三维评估分析及手术方案设计。本系统在开发中采用功能模块化设计思想,可针对其他手术要求进行功能方面的扩展,因此,本系统也具有较好的可扩展性与工程参考意义。