MathJax.Hub.Config({tex2jax: {inlineMath: [['$','$'], ['\\(','\\)']]}}); function MyAutoRun() {    var topp=$(window).height()/2; if($(window).height()>450){ jQuery(".outline_switch_td").css({ position : "fixed", top:topp+"px" }); }  }    window.onload=MyAutoRun; $(window).resize(function(){ var bodyw=$win.width(); var _leftPaneInner_width = jQuery(".rich_html_content #leftPaneInner").width(); var _main_article_body = jQuery(".rich_html_content #main_article_body").width(); var rightw=bodyw-_leftPaneInner_width-_main_article_body-25;   var topp=$(window).height()/2; if(rightw<0||$(window).height()<455){ $("#nav-article-page").hide(); $(".outline_switch_td").hide(); }else{ $("#nav-article-page").show(); $(".outline_switch_td").show(); var topp=$(window).height()/2; jQuery(".outline_switch_td").css({ position : "fixed", top:topp+"px" }); } }); 基于密度-距离图的交互式体数据分类方法
  软件学报  2016, Vol. 27 Issue (5): 1061-1073   PDF    
基于密度-距离图的交互式体数据分类方法
周芳芳1, 高飞2, 刘勇刚2, 梁兴1, 赵颖1     
1. 中南大学信息科学与工程学院, 湖南长沙 410083;
2. 中南大学软件学院, 湖南长沙 410075
摘要: 体数据分类是体绘制中传递函数设计的核心问题.标量值-梯度模直方图作为表征体数据的一种经典二维特征空间,已被广泛应用于分类体数据.然而,大部分已有方法存在过于依赖分类算法的参数设置、运算效率低、交互复杂度高等问题.以标量值-梯度模直方图的密度分布为基础,并依据物质中心密度大且物质中心间距离远这一特性,首先快速计算每个数据点的密度及每个数据点到比其密度大的点的最小距离;然后,将所有数据点投影到密度-距离图,并以密度-距离图作为人机接口,使用户能够交互地选择多个密度中心来分类体数据并设置传递函数.通过多组实验验证,所提出的方法无需预设物质类别的数量,分割标量值-梯度模直方图的准确度较高且速度较快,所设计的密度-距离图是一个有效的人机交互接口,可以有效地引导用户完成由粗糙到精细的递进式体数据分类和可视化过程.
关键词: 体数据    传递函数    体数据分类    基于密度的聚类    
Interactive Volume Data Classification Based on Density-Distance Graph
ZHOU Fang-Fang1, GAO Fei2, LIU Yong-Gang2, LIANG Xing1, ZHAO Ying1     
1. School of Information Science and Engineering, Central South University, Changsha 410083, China;
2. School of Software, Central South University, Changsha 410075, China
Abstract: Volume data classification is a core issue of transfer function in volume rendering. Scalar-gradient magnitude histogram of volume is a classic feature space, and has been applied in volume classification for its nice result in visual extraction of boundaries between different materials. However, the design of transfer function based on scalar-gradient histogram has proven as a time-consuming and complex task which is hard for users to conduct interactions. In this paper, scalar-gradient histogram is treated as a density distribution of all voxels. This approach assumes that the density of a material center is higher than their neighbors and the distance between two material centers is far enough. By computing the minimum distance between each points and all other points with higher density in scalar-gradient histogram, a density-distance graph is constructed based on densities and minimum distances of all points. The density peaks are easily observed in the graph and can guide the users to select centers of each material as a progressive volume classification process through a set of specified interactions. Experimental results demonstrate that the presented approach does not require the prior knowledge of categories, and the volume classification is accurate with high performance.
Key words: volume data    transfer function    volume data classification    density-based clustering    

随着数据采集技术的不断发展,人们可以获取到越来越多的体数据.如何快速准确地探索体数据中有价值的信息,是一个广受关注的热点问题.体绘制技术可以形象地展示体数据,有效地揭示体数据中隐含的特征信息,是体数据分析和处理的重要手段,被普遍应用于医学诊断、工业制造、地质勘探等与人类社会发展密切相关的领域[1].

体数据分类是体绘制的核心问题,它通过分析体数据的属性,将具有相似特征的体素归为一类,并在可视化中对各类体素赋予不同的颜色和不透明度,从而区分物质和凸显用户感兴趣的区域.体数据分类的过程也被称为传递函数设计.目前已有许多传递函数设计方法,大致可以分为图像中心法和数据中心法两大类[2].在数据中心法中,研究者常常利用体数据的多维属性特征,将体数据投影到一维或二维的特征空间,然后在特征空间中实现体数据的分类和绘制.标量值-梯度模直方图[3]是一种非常有代表性的体数据二维特征空间,每个体素在直方图中的位置由该体素的标量值和梯度模决定,体数据在该特征空间中呈现漂亮的多拱形分布,拱形的桥墩区域对应不同物质的内部,拱形的弧线区域对应物质之间的边界.近年来,有许多传递函数设计方法围绕标量值-梯度模直方图展开.交互式方法提倡让用户在标量值-梯度模直方图上自主地挑选感兴趣的多边形区域[4],并依次调节颜色和不透明度.虽然交互可以更多地引入专家知识,但在缺乏合理约束时,调节理想多边形区域是一个十分繁复的过程.智能化方法主张先自动或者半自动地分割标量值-梯度模直方图,再根据分割结果设置传递函数.研究者先后探索了非参数密度估计[5]、高斯混合模型[6]、Normalized-Cut[7]等多种算法和模型,但总体而言,分割标量值-梯度模直方图的准确性过于依赖目标分类数等算法参数设置,而且算法效率普遍不高.

标量值-梯度模直方图是一个携带密度分布信息的体数据特征空间,大量的物质内部体素会被投影到高密度的拱形桥墩区域,少量的物质边界体素则被投影到低密度的拱形弧线区域.受密度峰值快速搜索聚类算法[8]的启发,本文提出了基于密度-距离图的交互式体数据分类与可视化方法.我们首先计算标量值-梯度模直方图中每个数据点的密度值和每个数据点到比其密度更大点的最小距离,并将这些点投影到密度-最小距离投影图(简称为密度-距离图)中,然后由用户交互式地选择多个密度中心点来分割标量值-梯度模直方图.为了进一步提高密度峰值快速搜索聚类算法应用在体数据分类上的运算效率,本文重点从密度的计算、密度半径的计算和最小距离的计算这3个方面讨论了算法的实现细节.为了降低用户与标量值-梯度模直方图之间的交互复杂度,本文根据密度-距离图能够直观地表现高密度中心这一特征设计了一系列交互功能,帮助用户交互地分割标量值-梯度模直方图,选择不同物质并设置它们的颜色与不透明度.本文通过多组实验说明了密度-距离图作为人机交互接口的有效性,它可以引导用户完成由粗糙到精细的递进式体数据分类和可视化过程.

1 相关工作

体绘制的核心技术是传递函数的设计.传递函数的设计实际上是一个体数据分类的过程[9].它根据数据的特征将体数据映射为两种视觉属性:颜色和不透明度,具有相似特征的数据具有相近的颜色和不透明度.通过调整传递函数,用户可以凸显体数据中感兴趣的物质,产生更具信息量的可视化结果.传递函数的设计方法大致可以分为图像中心法和数据中心法.图像中心法通过对体绘制结果进行评估,并调节传递函数的参数直到得到满意的效果[10, 11];数据中心法则关注体数据本身的属性,在特征空间中将具有不同特征的体素映射为不同的视觉特征.

在数据中心法中,最简单的体数据特征空间是一维标量空间.体素标量值的大小本身就具有一定的含义,可以用来区分不同的物质[12].但对于复杂的体数据,不同物质的体素可能具有相似的标量值,因此,一维传递函数的分类效果往往局限性很大.更常见的研究思路是:以体素标量值为基础,同时引入体数据的其他属性或特征来提高传递函数的分类能力.比如:Correa和Ma等人[13]提出了一种size-based的传递函数设计方法 ,通过对体数据作高斯模糊来估算物质的大小,用大小信息区分同一标量值的不同物质;Fujishiro等人[14]和Zhou等人[15]关注到标量域的拓扑结构,他们考虑等值面的连通性和拓扑结构发生显著变化的关键等值面,并针对关键结构设置局部传递函数;Sereda等人[16, 17]通过沿梯度方向搜索局部最小(L)和局部最大(H)标量值,并计算LH直方图对边界面分层聚类,解决了在有噪声和模糊情况下边界面准确提取的问题;Tzeng等人[18]综合标量、梯度、三维空间、邻居等多种属性设计高维的传递函数,通过在体数据切片上标记感兴趣的区域来训练神经网络,提高数据分类的准确性;Correa等人[19]和Ruiz等人[20]在体绘制时考虑了每个体素的能见度,因为某些特征虽然被赋予了高不透明度,但在交互绘制时能见度可能很低,从而引导用户设计更好的传递函数;Guo等人[21]提出了结合维度投影和平行坐标的多变量传递函数设置方法,引领了信息可视化与体绘制相结合的研究思路.

标量值-梯度模直方图是体数据的一种经典的二维特征空间,此特征空间的X轴是体素标量值,Y轴是体素梯度模.因为物质内部体素的梯度模通常较小,物质边界体素的梯度模通常较大,并且不同物质的标量值区间通常会不一样,所以体数据中所有体素投影到该特征空间后会呈现出漂亮的多拱形现象,物质的分界面常常表现为拱形的弧形区域,物质内部则表现为拱形的桥墩区域.基于标量值-梯度模直方图的传递函数设计方法大致可以分为交互式、自动化和半自动的方法.Kniss等人[4]实现了交互工具,帮助用户在弧形区域上设置梯形或矩形来选中物质边界,但是对于多边形的位置、形状和大小的调节缺乏引导,需要反复实验,费时费力.为了形成有指导的交互,Higuera等人[22]针对拱形特征,用弧线形或抛物线形的特征选择器替代自由多边形.进一步地,Wang等人[6]先用高斯混合模型自动产生一组椭圆型传递函数,然后,用户在标量值-梯度模直方图上进一步交互地调节椭圆的大小、位置和方向.但该方法需要用户设置感兴趣特征的个数,这难以探索没有先验知识的体数据.标量值-梯度模直方图中并没有体数据的空间信息,Roettger等人[23]将体数据的空间信息引入标量值-梯度模直方图,计算图中每个点的对应体素在三维空间中的中心和半径,设计基于空间信息的传递函数.但此方法会出现过度分割的情况,着色后的标量值-梯度模直方图显得十分混乱.Huang等人[24]通过在标量值-梯度模直方图上交互地选择种子点,利用三维的区域增长算法来分割直方图,辅助用户提取感兴趣的区域特征.Maciejewski等人[5]采用非参数估计算法对标量值-梯度模直方图进行分割,并用概率密度估计来控制传递函数的设置,但是分类的准确性受到概率密度估计粒度设置的影响较大.Ip等人[7]首先使用Normalized-Cut算法自动地分割标量值-梯度模直方图,并引入信息熵指导用户在直方图上交互操作,但Normalized-Cut分割标量值-梯度模直方图的效率不高.总的来说,基于标量值-梯度模直方图的体数据分类与可视化方法一直在努力减少参数设置的影响、提升分类准确性和提高交互效率,研究者也一直在探索更新、更有效的体数据分析和处理技术.

2 总体流程

本文的基本思想是,基于标量值-梯度模直方图的密度分布特性对体数据进行分类.体数据中,物质内部的体素数量要远远大于物质边界的体素数量,且物质内部体素的标量值变化不大、梯度较小,而物质边界体素的标量值变化大、梯度大.在标量值-梯度模直方图中,具有相同标量值和梯度模值的体素会投影到同一位置,因此,代表物质内部的拱形桥墩区域密度较高,而代表物质边界的拱形弧线区域密度较低.根据这一先验知识,我们可以基于密度分布来分割标量值-梯度模直方图,进而实现体数据的分类和可视化.为了准确而快速地分割标量值-梯度模直方图,本文核心方法是基于密度和最小距离的聚类算法.在标量值-梯度模直方图中,物质中心的密度通常明显高于周围点的密度,并且物质中心与比它密度高的点的最小距离会相对较大.根据这一假设,我们无需预设物质类别的数量,就可以在标量值-梯度模直方图中快速地搜索出各类物质的密度中心.同时,我们将密度-距离图作为人机交互接口,当把标量值-梯度模直方图的所有数据点投影到密度-距离图后,代表不同物质的高密度中心将凸显在投影图的中心区域,便于用户选择高密度中心点来分割标量值-梯度模直方图.

本文方法的总体流程如图1所示.图1(a)是一个牙齿体数据.图1(b)是根据文献[3]中方法得到的该数据的标量值-梯度模直方图,直方图上的4个桥墩区域分别代表牙冠、牙本质、牙髓以及包围盒4种物质,5个拱形弧区域表示这4种物质之间可能存在的边界.为了观察标量值-梯度模直方图的密度分布特征,我们用热力图对图1(b)的密度分布进行可视化,如图1(c)所示,其中红色等暖色对应密度大的区域,蓝色等冷色对应密度小的区域,可以很明显地看到,表示物质内部的拱形桥墩区域密度较高.

图1 总体流程图 Fig.1 Overview of workflow

本文借鉴了Rodriguez和Laio[8]提出的密度峰值快速搜索聚类算法寻找标量值-梯度模直方图的密度中心.该算法假设簇中心的密度具有两个特点:

(1) 簇中心的密度比周围点的密度明显较高.

(2) 簇中心点之间的距离较大.该算法计算每个数据点的密度和每个数据点到比该点密度大的点的最小 距离.

由于密度中心点的周围都是密度比它小的点,要找到比它密度更大的点,就要搜索到另外一个密度中心附近,因此密度中心点的最小距离较大;而非密度中心点附近必定有密度较大的点,因此它的最小距离较小.根据该算法的思想,在标量值-梯度模直方图中,密度大且最小距离大的点为密度中心点,对应物质的中心;密度大且最小距离小的点为非密度中心点;密度小且最小距离大的点为稀疏的孤立点.

本文采用的密度-距离图是一个X轴为密度值、Y轴为最小距离值的二维空间,用于将标量值-梯度模直方图中的点根据其密度和最小距离投影过来,帮助用户观察和挑选密度中心并用来对直方图进行分割.图1(d)图1(c)的密度-距离投影结果,投影的右侧区域有密度和最小距离都非常大的4个数据点,这4个数据点在标量值-梯度模直方图中是理想的分类中心,其他点可以根据其在直方图上的位置以及密度归入这4个分类中心对应的4类物质.图1(e)是基于密度-距离投影分割图1(c)的结果,它以4个拱形桥墩为中心形成了4个分类.我们还提供了一些交互手段,帮助用户在密度-距离图和分割后的标量值-梯度模直方图上选择感兴趣的物质,并设置相应的颜色和不透明度,从而降低设置传递函数的难度.图1(f)是只选择了两类物质并设置不同的颜色和不透明度后,牙齿体数据的体绘制效果.

3 标量值-梯度模直方图的密度中心快速搜索方法

密度峰值快速搜索聚类算法[8]的优点是无需预设目标分类数、无需迭代计算、速度快、健壮性好,而且可以找到复杂结构的聚类.但体数据的标量值分布不均且数量太大,本文将该方法直接应用于分割标量值-梯度模直方图时,遇到了密度值差异过大和计算效率不够高的问题,因此,我们根据标量值-梯度模直方图的特点,采用网格划分策略来优化密度计算、密度半径计算和最小距离计算方法.

3.1 密度的计算

我们首先定义标量值-梯度模直方图的密度模型.本文将直方图中的点pj的密度R(pj)定义为:以点pj为中心、dc为半径的圆形区域内的体素总量,公式为

\[R({{p}_{j}})=\sum\nolimits_{i=1}^{n}{L({{p}_{i}})h(||{{p}_{i}}-{{p}_{j}}||-{{d}_{c}})}\] (1)
其中,n表示圆形区域内点的数量;h(x)是一个阶跃函数,满足:如果x<0,h(x)=1;否则,h(x)=0;||pi-pj||表示pipj之间的距离;dc是密度半径;L(pi)与投影到点pi的体素个数有关,是点pi的权值函数,表示为
L(pi)=logS(pi) (2)
其中,S(pi)表示投影到点pi处的体素个数.在标量值-梯度模直方图中,具有相同标量值和梯度模的体素会投影在同一点pi上.通过计算发现,每个点pi包含的体素个数S(pi)差异很大,如物质内部的点可能包含上万个体素,而边界点则可能仅包含几个体素.这将导致不同点的密度差异很大,容易忽略密度值小的点.本文对体素的个数取对数,不但保持包含体素个数的大小关系,而且能够减小密度计算中的体素权重之间的差距.

如果标量值-梯度模直方图中每个点的标量值和梯度模为连续的浮点数,那么计算一个点的密度需要遍历整个样本空间,计算该点与所有样本点之间的距离,才能确定距离小于密度半径dc的点.计算密度的时间复杂度是O(n2),n是样本点的个数.考虑到体数据的标量值是离散的,为了提高密度计算的效率,本文将标量值-梯度模直方图划分为一个二维的均匀网格,如图2所示,样本点的位置都位于网格交点上.在计算密度时,只需遍历该点周围的正方形内的点,从而大大减少了密度计算的时间复杂度.所以,根据公式(1)计算点pj的密度搜索以点pj为中心、以dc为半径的圆内的样本点,则转换为在一个均匀的网格上统计以pj为中心、以H为边长的正方形内所有的点,H=2h+1且h为不小于dc的最小正整数,此时,密度计算的时间复杂度减少为O(n×H2).

图2 网格化的标量值-梯度模直方图 Fig.2 Grided scalars gradient magnititudes histogram
3.2 密度半径的计算

密度半径dc是密度峰值快速搜索聚类算法的核心参数,Rodriguez等人[8]总结了半径dc的设置经验,即,对聚类空间中所有点,其半径范围内点的个数的平均值占全样本空间点个数的1%~2%.但在实验中我们发现,对于具有较大样本点数量的标量值-梯度模直方图,dc的设置仍需要进一步讨论:

· 首先,如果样本空间内点的数量较少或分散稀疏时,直接使用经验值很难获得准确的分类结果;

· 其次,如果使用迭代的方法计算dc,即,给定一个dc初始值并迭代调节dc的大小,直到满足平均邻居点数在1%~2%范围内这个条件,那么计算半径dc成为一项耗时的工作.

本文提出一种改进的dc快速计算方法,它根据样本点分布的疏密程度自动地调节dc的值,流程如图3所示.

图3 计算密度半径dc的流程图 Fig.3 Flow diagram of calculating radius dc

Step 1. 在标量值-梯度模直方图中,随机地选择一些点作为计算密度半径dc的样本,采样比例为ratio1.

Step 2. 计算所有样本点包含的体素权重之和T,并估计全部点的分布情况Ttotal,其中,m是样本的个数:

\[T=\sum\nolimits_{i=1}^{m}{L({{p}_{i}})}\] (3)
\[{{T}_{total}}=\frac{T}{rati{{o}_{1}}}\] (4)

Step 3. 根据Ttotal设定阈值Ttotalxratio2,计算最小半径ri,使得每个样本邻域内点的个数不小于这个阈值.

Step 4. 计算m个样本点的半径ri的平均值,设为dc的值:

\[{{d}_{c}}=\frac{\sum\nolimits_{i=1}^{m}{{{r}_{i}}}}{m}\] (5)

本方法通过随机选取标量值-梯度模直方图中的点,并计算样本的平均半径来获取dc,避免了迭代地调整dc的大小,提高了计算效率.另外,公式(3)中使用体素数量的对数不仅可以表现点数量的大小,还能反映标量值-梯度模直方图中点分布的稀疏情况.通过多组实验我们发现:当采样比例ratio1设为2%,而邻域点阈值比例ratio2设为0.1%时,可以获得较好的标量值-梯度模直方图分割效率和分割结果.

3.3 最小距离的计算

在密度峰值快速搜索聚类算法中,最耗时的阶段是最小距离的计算.最小距离d是当前点pc到密度更大点的最短距离,即

\[d({{p}_{c}})={{\min }_{R(p)>R({{p}_{c}})}}(||p-{{p}_{c}}||)\] (6)

我们统计了实验中5个体数据的标量值-梯度模直方图中样本点的最小距离d,我们发现,绝大部分点的最小距离d通常都小于3.因此,针对标量值-梯度模直方图,我们提出了由近及远的扩散式最小距离计算方法.在计算某个点的最小距离时,以该点为中心,由近及远地扩散式地搜索,当搜索到比点密度更大的点时,立即停止搜索,从而避免了遍历整个样本空间,提高了算法运算效率.

具体步骤如图4所示.

图4 计算点pc的最短距离的流程图 Fig.4 Flow diagram of computing the minimum distance of the point pc

Step 1. 首先令H=3,H为以当前样本点为中心的正方形的边长(即,边包含的点数).

Step 2. 以当前样本点为中心、H为边长的正方形的最外层(当H=3时,就是除了当前点以外的正方形中的其他样本点)搜索密度大于当前点密度的点.

Step 3. 判断Step 2中是否搜索到密度更大的点:若搜索到,则执行Step 4;反之,则令H=H+2,执行Step 2.

Step 4. 在最外层中搜索到密度更大的样本点中距离最小的点,获得d的值.

4 密度-距离图的布局与交互设计方法

密度峰值快速搜索聚类算法需要计算每个数据点的密度值和每个数据点到比其密度大的点的最小距离,密度中心点通常具有较大密度值和较大最小距离值,而且密度中心点的数量远小于非密度中心体素的数量.因此,当投影所有点到一个以密度为X轴、最小距离为Y轴的二维空间时,则可得到如图1(d)所示的密度-距离图.在该投影图中,所有密度中心点将位于中心偏右的区域,并且附近没有其他噪声点干扰,用户可以快速获得密度中心点.密度-距离图的这一特征,使它成为一个很理想的人机交互接口,用户可以灵活地选择感兴趣的密度中心对标量值-梯度模直方图实现由粗到细的递进式分割,用户也可以快速地设置感兴趣物质的颜色和不透明度来获得不同的体绘制效果.本节我们将具体讨论密度-距离图的布局和交互功能设计.

4.1 密度-距离图的布局设计

由于标量值-梯度模直方图中大量的非密度中心点会投影到密度-距离图靠近两个坐标轴的区域中,造成数据点分布严重不均衡,并且部分区域数据点高度重叠等问题.为了提高交互体验和视觉效果,在投影前我们进行了取对数和归一化处理,公式如下:

\[{{x}_{i}}=Width\frac{\log {{r}_{i}}-\min (\log r)}{\max (\log r)-\min (\log r)}\] (7)
\[{{y}_{i}}=Hight\frac{\log {{d}_{i}}-\min (\log d)}{\max (\log d)-\min (\log d)}\] (8)
其中,ridi分别是标量值-梯度模直方图中当前数据点i的密度与最小距离,xiyi是数据点i在投影图中的位置,WidthHeight分别是投影图的宽度和高度.取对数和归一化运算使得相互距离远的点之间差距减小,而距离近的点之间的差距变大,从而有助于整体上减少投影点互相覆盖重叠,提高屏幕利用率.

在密度-距离图中,用户交互时更关心密度和最小距离都较大的密度中心点.大量重叠出现在坐标轴附近的点不但会干扰用户的交互操作,而且会延长投影图的布局时间.因此,我们提供了阈值过滤功能,用户可以分别为密度和最小距离设置阈值,让标量值-梯度模直方图中小于阈值的数据点不参与投影.

4.2 密度-距离图的交互设计

在密度-距离图中,每个数据点用一个圆表示,圆的初始颜色是随机分配的,圆半径采用min(r,d),即,取当前数据点的密度r与最小距离d中较小的值.如果圆点半径大,说明密度和最小距离都大,那么它是标量值-梯度模直方图的密度中心的可信度就越大.当用户不断将感兴趣的密度中心点加入分类中心集合时,标量值-梯度模直方图中每个数据点将归类于离自己最近的分类中心,从而实现对标量值-梯度模直方图的交互式分割.

我们在密度-距离图中提供的交互功能包括:

(1) 添加分类中心,用户可将密度-距离图中任意一个数据点当作密度中心,并加入分类中心集合,用于分割标量值-梯度模直方图.

(2) 自动改变形状,在密度-距离图中作为分类中心的数据点会自动由圆形变为正方形.

(3) 删除分类中心,用户可将任意一个分类中心移出分类中心集合,并重新分割标量值-梯度模直方图.

(4) 更改分类的颜色,用户可以在密度-距离图中更改每个分类中心的颜色,在标量值-梯度模直方图中,该分类对应的所有数据点的颜色也随之改变.

(5) 微调位置,用户可以移动密度-距离图中任意数据点,改变其初始投影位置,便于与其他数据点区分 开来.

基于密度-距离图交互地分割标量值-梯度模直方图是一个从粗糙到精细的划分过程,下面我们以牙齿体数据为例,说明如何通过交互操作来递进地分割标量值-梯度模直方图.如图5所示,在密度-距离图的中心偏右区域中,有4个明显的密度高且最小距离大的密度中心点,它们可以作为候选分类中心.在Step 1中,我们选择最大的紫色密度中心点作为分类中心,这时,圆点变为正方形,表示该点将被用来分割标量值-梯度模直方图.由于只有一个分类中心,根据自顶向下的分类策略,标量值-梯度模直方图上所有点都变成了紫色.在Step 2中,我们将密度-距离图中第二大的密度中心点也作为分类中心,并将其颜色设置为深绿色.由于在计算密度-距离投影过程中,我们保留了直方图上每个点的最近距离的大密度点,所以对直方图上的点进行分类时,收敛到同一分类中心的所有点归为一类.随着第2个分类中心的加入,根据上述分类方法,标量值-梯度模直方图上所有点被分为紫色和深绿色两类.在Step 3中,我们将另外两个密度中心点也加入进来,4个分类中心将标量值-梯度模直方图分割为4个连通区域,我们可以注意到,4个拱形桥墩刚好被划分为4类,分别代表牙齿体数据的4种不同物质.

图5 基于密度-距离图交互式分割标量值-梯度模直方图的示意图 Fig.5 Illustration of interactive segmentation of scalars-gradient histogram based on density-distance graph

被分割后的标量值-梯度模直方图会形成多个连续的色块,每个色块代表一个分类,每个分类对应一种体数据中可能存在的物质.用户不但可以通过色块的大小和点的密集程度清晰地分辨每个类所对应的体素数量,还可以直接根据分类结果和色彩设置生成传递函数,并输出绘制结果.同时,我们还提供了一些交互功能,用来操作分割后的标量值-梯度模直方图.标量值-梯度模直方图的交互以色块为单位,用户既可以选择感兴趣的色块并调节对应物质在绘制时的颜色和不透明度,也可以直接删除色块或者将色块对应的物质设置为不可见,这些交互操作都会及时地反映到体绘制结果中.

5 实验与结果分析

在实验环节,我们首先用Pig和Engine两个体数据来演示本文提出并实现的交互式体数据分类方法;然后,我们用几组对比实验来说明本文方法在分类准确性和分类速度上的优势.实验环境为处理器Intel Core i5- 4200M,主频2.5GHz,内存4G和Windows7操作系统.

· Pig体数据的交互式分类与可视化

Pig体数据是一个512x512x134的猪形存钱罐,如图6所示.

图6 Pig体数据的交互式分类与可视化示意图 Fig.6 Illustration of interactive classification and visualization of Pig volume data

Step 1的左边是Pig体数据在未分类时的标量值-梯度模直方图,中间是初始体绘制效果图.未经分类的体数据呈现出混沌的绘制效果,难以判断Pig体数据包含哪些物质.Step 1右边是密度-距离图,dc为2.0,4个明显的大圆点是标量值-梯度模直方图的最佳候选密度中心,它将引导我们完成体数据分类并获得更好的绘制效果.在Step 2中,我们选择了密度-距离图最靠右的绿色和蓝色两个候选密度中心点,并将其加入分类中心集合,这时,标量值-梯度模直方图被相对应地分割成绿色和蓝色两个区域,这两个区域也同时将Pig体数据分为两大类物质.绿色物质是体数据的包围盒,当我们降低绿色包围盒的不透明度时,猪形存钱罐的轮廓从绿色包围盒中隐现出来,同时还可以隐约看到存钱罐的支撑物.在Step 3中,我们添加另两个候选密度中心到分类中心集合,这时,Step 2的标量值-梯度模直方图的绿色色块被进一步分割出两个较小的色块,在对应的体绘制结果中,存钱罐的底座已经非常明显了,而且存钱罐内也可以隐约看到硬币的存在.为了能够更加清楚地观察存钱罐及其内部物体,在Step 4中我们重新设置了4个色块的颜色,将背景环境对应的色块直接设为不可见,并在不透明度配置面板中用高斯型函数设定不透明度,这时,已经可以清晰地看到存钱罐本体、存钱罐底座和存钱罐内的硬币这3种不同物质.

· Engine体数据的交互式分类与可视化

Engine是一个256x256x256的引擎体数据.如图7所示,dc设置为2.0,引擎的密度-距离图中有两个半径明显很大的点,它们在标量值-梯度模直方图中相对应的两个密度中心是理想的分类中心,可以将引擎体数据分为两大类物质.于是,我们在Step 1中用这两个分类中心将标量值-梯度模直方图分割为深蓝色和浅蓝色两个色块,在体绘制结果中可以发现:深蓝色对应引擎的一些内部部件,浅蓝色对应引擎壳.在引擎的密度-距离图中还有一些半径较大点,如果将这些点加入分类中心集合,可以将引擎进一步细分.在Step 2中,我们分析密度-距离图中密度大而最小距离小的点对分类结果的影响.在标量值-梯度模直方图上,这种点会靠近某个大密度中心,它通常与某个大密度中心是同一类物质,如果将其加入分类中心集合,则很可能将某大物质的一些体素误分为新物质.在Step 2的密度-距离图中,一个密度较大而最小距离很小的橙色点被加入分类中心集合,随后,标量值-梯度模直方图上的浅蓝色块被划分出一小块橙色区域.在体绘制时,我们只将这类橙色体素显示出来,可以发现,橙色体素遍布了原来浅蓝色引擎壳大部分位置,这说明这次细分并没有形成新物质类型.在Step 3中,我们将密度-距离图中一个密度较小的红色点加入分类中心集合,在绘制结果中,浅蓝色引擎壳中少量体素被分类为红色,这些红色体素零星地分散在引擎壳上.这说明,如果将密度较小的点作为分类中心,它会将其附近少量密度更小的点作为一个新类.

图7 Engine体数据的交互式分类与可视化示意图 Fig.7 Illustration of interactive classification and visualization of Engine volume data

· 本文方法的计算效率统计

计算标量值-梯度模直方图中各数据点的密度和最小距离,是整个算法中最耗时的部分.为了避免计算密度和最小距离时采用全局搜索策略,加快计算速度,我们根据体数据和标量值-梯度模直方图的特点,在第3节详细讨论了如何改进计算策略.表1是对5个体数据进行对比实验的结果,密度半径dc都为2.

表1 密度和最小距离的计算效率对比 Table 1 Computing efficiency comparison of density and minimum distance

实验结果表明,改进的计算策略比全局搜索策略要快很多,5个体数据的计算都可以在1s内完成.这非常有利于提高交互效率和改善用户体验.

· 与其他方法的对比实验

均值漂移也是一种常见的基于密度的聚类算法[25].它通过迭代搜索密度峰值来实现数据聚类,常被用于分割图像.下面我们比较本文方法与均值漂移的实验结果.实验对象是128x128x128的Sphere体数据,密度半径dc为2.0.如图8(a)所示,首先,我们选中密度-距离图的5个密度中心点,将标量值-梯度模直方图分割成5个色块,然后用其中3个色块输出绘制结果.在结果中我们可以发现,最小的球和中等的球被分为一类,这是由于体素较少的物质在标量值-梯度模直方图中密度非常低,因为本文方法会忽略密度过小的物质.当均值漂移的带宽设置为16时,可以较好地分割标量值-梯度模直方图,如图8(b)所示,最小的球被正确地分为一类,但在右上的稀疏区域出现了许多被过度分割的小色块,这会造成同一类物质被错分为许多类.

图8 本文方法和均值漂移对Sphere体数据分割并可视化结果的对比 Fig.8 Comparison of segmented and visual results on Sphere volume data by our method and mean shift

我们继续讨论改变两种算法的核心参数设置时对标量值-梯度模直方图分割结果的影响.均值漂移算法的核心参数是搜索带宽,本文方法的核心参数是密度半径dc,我们在多个体数据上进行了多组对比实验.图9显示了3组Pig体数据的标量值-梯度模直方图的分割结果,其中,本文方法的dc设置为2.0时的结果可参见图6.

图9 不同参数设置对Pig体数据分割结果的对比 Fig.9 Comparison of segmented results on Pig volume data with different parameter settings

实验结果表明,分割标量值-梯度模直方图时均值漂移受到核心参数的影响更大.另外,本文方法采用交互式分割,用户可以参与分割过程,从而可以进一步保持在不同dc设置情况下分割结果的稳定性.

从计算效率来看,在表2所示的5个实验结果中,均值漂移的带宽均设为16,本文方法的密度计算半径dc均设为2,本文的方法分割标量值-梯度模直方图明显快于均值漂移.另外,在文献[7]中提到,Normalized-Cut分割这5个体数据的标量值-梯度模直方图平均需要约15s.

表2 本文方法和均值漂移的计算效率对比 Table 2 Computing efficiency comparison of our method and mean shift
6 总结与展望

本文实现了一种基于密度-距离图的交互式体数据分类方法.该方法以体数据的标量值-梯度模直方图的密度分布为基础,首先计算每个数据点的密度和最小距离,并构造密度-距离图;然后,用户交互地在密度-距离图中选择高密度中心点对标量值-梯度模直方图进行分割,进而设置传递函数并完成体绘制.通过多组实验表明:本文方法分割标量值-梯度模直方图的准确率受算法核心参数设置的影响较小,分割速度快,且无需预设物质类别的数量;同时,本文可交互的密度-距离图是一个有效的人机接口,能够简化传递函数设置过程,引导用户实现由粗到细的递进式体数据分类和可视化.

本文方法还存在一些有待改进之处:首先,本文算法假设分类中心点是密度很高的点,但是当体素数量不多的小物质投影到标量值-梯度模直方图时,对应的数据点的密度也比较小,这将导致小物质在分类时会被忽略;其次,根据体数据投影到标量值-梯度模直方图的原理,多种物质或多个边界可能被投影到直方图的相近甚至相同区域中,这也是影响本文方法分类结果的重要因素,一个潜在的解决思路是,引入体数据的空间信息来提高分类的准确性;最后,本文方法将标量值-梯度模直方图上所有点都归入相应物质类型,没有考虑某些点可能不属于当前添加的所有分类中心的情况,也没有专门提取物质边界和分析边界的不确定性,但物质边界的提取和分析有助于理解体数据和提高体绘制效果,因此,后续研究应该注重体数据分类与边界分析相结合的思路.

参考文献
[1] Pfister H, Lorensen B, Bajaj C, Kindlmann G, Schroeder W, Avila LS, Martin K, Machiraju R, Lee J. The transfer function bakeoff. IEEE Computer Graphics and Applications, 2001,21(3):16-22.[doi:10.1109/38.920623]
[2] Guo HQ, Yuan XR. Survey on transfer functions in volume visualization. Journal of Computer-Aided Design & Computer Graphics, 2012,24(10):1249-1258(in Chinese with English abstract).
[3] Kindlmann G, Durkin JW. Semi-Automatic generation of transfer functions for direct volume rendering. In:Proc. of the 1998 Symp. of Volume Visualization. New York:ACM Press, 1998. 79-86.[doi:10.1145/288126.288167]
[4] Kniss J, Kindlmann G, Hansen C. Multidimensional transfer functions for interactive volume rendering. IEEE Trans. on Visualization and Computer Graphics, 2002,8(3):270-285.[doi:10.1109/TVCG.2002.1021579]
[5] Maciejewski R, Woo I, Chen W, Ebert DS. Structuring feature space:A non-parametric method for volumetric transfer function generation. IEEE Trans. on Visualization & Computer Graphics, 2009,15(6):1473-1480.[doi:10.1109/TVCG.2009.185]
[6] Wang YH, Chen W, Zhang J, Dong TX, Shan GH, Chi XB. Efficient volume exploration using the Gaussian mixture model. IEEE Trans. on Visualization and Computer Graphics, 2011,17(11):1560-1573.[doi:10.1109/TVCG.2011.97]
[7] Ip CY, Varshney A, JaJa J. Hierarchical exploration of volumes using multilevel segmentation of the intensity-gradient histograms. IEEE Trans. on Visualization and Computer Graphics, 2012,18(12):2355-2363.[doi:10.1109/TVCG.2012.231]
[8] Rodriguez A, Laio A. Clustering by fast search and find of density peaks. Science, 2014,344(6191):1492-1496.[doi:10.1126/science.1242072]
[9] Zhou FF, Zhao Y, Ma KL. Parallel mean shift for interactive volume segmentation. In:Proc. of the 1st Int'l Conf. on Machine Learning in Medical Imaging. LNCS 6357, Berlin, Heidelberg:Springer-Verlag, 2010. 67-75.[doi:10.1007/978-3-642-15948-0_9]
[10] Marks J, Andalman B, Beardsley PA, Freeman W, Gibson S, Hodgins J, Kang T, Mirtich B, Pfister H, Ruml W, Ryall K, Seims J, Shieber S. Design galleries:A general approach to setting parameters for computer graphics and animation. In:Proc. of the 24th Annual Conf. on Computer Graphics and Interactive Techniques. Los Angeles:Association for Computing Machinery, 1997. 389-400.[doi:10.1145/258734.258887]
[11] Tory M, Potts S, Moller T. A parallel coordinates style interface for exploratory volume visualization. IEEE Trans. on Visualization and Computer Graphics, 2005,11(1):71-80.[doi:10.1109/TVCG.2005.2]
[12] Zhou FF, Fan XP, Yang B. Prospects and current studies on designing transfer function in volume rendering. Journal of Image and Graphics, 2008,13(6):1034-1047(in Chinese with English abstract).
[13] Carlos D, Ma KL. Size-Based transfer functions:A new volume exploration technique. IEEE Trans. on Visualization and Computer Graphics, 2008,14(6):1380-1387.[doi:10.1109/TVCG.2008.162]
[14] Fujishiro I, Azuma T, Takeshima Y. Automating transfer function design for comprehensible volume rendering based on 3D field topology analysis. In:Proc. of the 10th IEEE Visualization 1999 Conf. (VIS'99). Washington:IEEE Computer Society, 1999. 467-563.[doi:10.1109/VISUAL.1999.809932]
[15] Zhou JL, Takatsuka M. Automatic transfer function generation using contour tree controlled residue flow model and color harmonics. IEEE Trans. on Visualization and Computer Graphics, 2009,15(6):1481-1488.[doi:10.1109/TVCG.2009.120]
[16] Sereda P, Bartroli AV, Serlie IWO, Gerritsen FA. Visualization of boundaries in volumetric datasets using LH histograms. IEEE Trans. on Visualization and Computer Graphics, 2006,12(2):208-218.[doi:10.1109/TVCG.2006.39]
[17] Sereda P, Vilanova A, Gerritsen FA. Automating transfer function design for volume rendering using hierarchical clustering of material boundaries. In:Proc. of the 8th Joint Eurographics/IEEE VGTC Conf. on Visualization. EuroGraphics Association, 2006. 243-250.[doi:10.2312/VisSym/EuroVis06/243-250]
[18] Tzeng FY, Lum EB, Ma KL. An intelligent system approach to higher-dimensional classification of volume data. IEEE Trans. on Visualization and Computer Graphics, 2005,11(3):273-284.[doi:10.1109/TVCG.2005.38]
[19] Correa CD, Ma KL. Visibility histograms and visibility-driven transfer functions. IEEE Trans. on Visualization and Computer Graphics, 2011,17(2):192-204.[doi:10.1109/TVCG.2010.35]
[20] Ruiz M, Bardera A, Boada I, Viola I, Feixas M, Sbert M. Automatic transfer functions based on informational divergence. IEEE Trans. on Visualization and Computer Graphics, 2011,17(12):1932-1941.[doi:10.1109/TVCG.2011.173]
[21] Guo HQ, Xiao H, Yuan XR. Scalable multivariate volume visualization and analysis based on dimension projection and parallel coordinates. IEEE Trans. on Visualization and Computer Graphics, 2012,18(9):1397-1410.[doi:10.1109/TVCG.2012.80]
[22] Higuera FV, Sauber N, Tomandl B, Nimsky C, Greiner G, Hastreiter P. Automatic adjustment of bidimensional transfer functions for direct volume visualization of intracranial aneurysms. In:Proc. of the SPIE Medical Imaging. SPIE, 2004. 275-284.[doi:10. 1117/12.535534]
[23] Roettger S, Bauer M, Stamminger M. Spatialized transfer functions. In:Proc. of the 7th Joint Eurographics/IEEE VGTC Conf. on Visualization. Eurographics Association, 2005. 271-278.[doi:10.2312/VisSym/EuroVis05/271-278]
[24] Huang RZ, Ma KL. RGVis:Region growing based techniques for volume visualization. In:Proc. of the Pacific Conf. on Computer Graphics and Applications. Washington:IEEE Computer Society, 2003. 355-363.[doi:10.1109/PCCGA.2003.1238277]
[25] Zhou FF, Huang W, Li JC. Extending dimensions in radviz based on mean shift. In:Proc. of the 2015 IEEE Pacific Visualization Symp. Hangzhou:IEEE, 2015. 111-115.[doi:10.1109/PACIFICVIS.2015.7156365]
[2] 郭翰琦,袁晓如.体数据可视化传递函数研究.计算机辅助设计与图形学学报,2012,24(10):1249-1258.
[12] 周芳芳,樊晓平,杨斌.体绘制中传递函数设计的研究现状与展望.中国图像图形学报,2008,13(6):1034-1047.