matlab可视化作图

第五章图像重建§5.5三维图像重建的可视化研究[6]医学影像的三维可视化技术，是利用一系列二维切片进行边界识别那样的不同技术处理，重新还原出被检物体的三维图像，然后进行定性和定量分析，可以从两维图像中获取三维结构信息，为医生提供逼真的显示手段和定量分析工具，在辅助医生诊断、手术仿真、引导治疗等方面发挥重要作用。两种方法：面绘制技术和体绘制技术。5.5.1立体步进法（MarchingCubes）－面绘制面绘制：由三维空间均匀数据场构造中间几何图元，如小三角形、小曲面等，然后再用传统的计算机图形学技术实现面绘制，加上光照模型、阴影处理，使得重建的三维图像产生真实感。立体步进法（MarchingCubes，MC）的要旨是：将三维数据网格分成许多体元，根据物体表面特征的信息，给出物体等值面的相关参数值，再逐个测试体元的8个角点是否在等值面上，通过线性插值得出体元中的哪些点在等值面上，用连接这些点得到的三角形或多边形来代替立方体，由这些全部的三角形或多边形得到三维数据场的三维表面信息。可用下图所示的两维数据网格简单理解。顺次扫描各个正方形，看每个正方形是否都在实心圆内，从而得出这些顶点的状态。若一个正方形两个相邻顶点中只有一个顶点在黄心圆外、一个顶点在黄心圆内时，就能确定边界穿过了由这两顶点形成的边，这样可以画出一个交点。如此扫描全部小正方形，就可以得到边界与网格线相交的全部1第五章图像重建交点，相连后就可以得到该圆的轮廓信息。当然应使轮廓光滑，其要点有：1）减少网格尺寸大小；2）用非线性插值得到更多中间点，再进行相连。对三维数据场，类似要确定物体表面轮廓信息。主要步骤是：5.5.2确定包含表面轮廓面的体元设表面轮廓面的某个特征参数为，实际上相同参数所构成的面就是等值面。在三维数据场中，等值面与体元相交的是三角形面或多边形面。如下图所示体元中，等值面与AB、AD和AE边相交。0CA(1)B(0)F(0)E(0)C(0)G(0)D(0)H(0)由于每个体元有8个角点，每个角点有0，1两个状态，每个体元2第五章图像重建的8个角点有个状态，要从中求出每个体元的等值面是困难的。但可简化为15种基本状态，然后利用不同的对称性如：8225=61）绕3条坐标轴中的任一轴旋转；2）绕3条坐标轴中的任一轴镜像反转；可以组合成8个角点的全部256种状态，从8个角的状态可以判断出所求等值面将与哪些体元相交，穿过体元的等值面的形状等信息。1．等值面与体元相交的边界交点如何求相交点的位置？当三维离散数据场的密度很高、体元很小时，用各角点的函数值在体元边界上线性内插，可以求得等值面与体元边界的交点。当三维离散数据场密度稀疏、体元较大时，可以通过三次线性插值求得等值面与边界的交点，等值面不是简单的两个角点的直线连接，而是一条曲线。如下图所示，PP11P12P21P22P1P2•82256=点函数值可以通过三次线性内插为：()01234567,,fxyzaaxayazaxyayzazxaxyz=+++++++3第五章图像重建其中系数(0,1,2,...7)iai=由体元8个角点处的函数值决定，如果用户等值面的值为，则等值面方程为0C()0,,fxyzC=由此等值面方程可以方便地求出某等值面与体元边界面的交线方程。考虑某边界面所在的平面方程0zz=，代入前式，可得如下的一条双曲线：()()()()030160250470aazaazxaaz0yaazxyC+++++++=。2．求等值面的法向矢量各三角面片的法向矢量：因为重建图像表面的法向矢量代表了局部的弯曲性，也决定了镜面反射的方向，是采用光照模型和消除各三角面片之间的明暗度不连续性变化所必须知道的量。我们一般只知道多边形的近似表示，只能针对每一个三角形面片或多边形面片，根据其所在的平面方程的系数决定该小片的法线，且取外法线方向。如果每一个多边形的平面方程已经知道，则多边形的顶点的法线，取包围该顶点的各多边形的法线的平均值。5.5.3体绘制这里不需要构造中间几何图元，而是由三维数据场直接产生屏幕上的二维图像，由体元直接绘制出来。体绘制可以在空间域或在频域进行。4第五章图像重建体直接绘制的实质：将已经采集到的三维数据场重新采样，按照一定的规则转换为图像显示缓存中二维离散信号。1．图像空间扫描的体绘制方法这是从图像空间（显示图像的屏幕）到物体空间（三维离散数据场）的过程：由显示图像屏幕上的每个像素点位置向物体空间发出射线，该射线与物体空间相交许多点－这些交点即为物体空间上新的采样点。选择适当的重构元素，对离散的三维数据场进行卷积运算，重构原始的图像信号，以求得新采样点的数据。根据重采样的Nyquist频率极限对三维连续数据场进行重新采样，得出重采样点的灰度值。在重采样的基础上，进行图像合成，即计算全体采样点对屏幕像素的贡献，得到屏幕上每一个像素点的光强度值。例、光线投射体绘制法主要步骤：1）对物体空间的三维离散数据场进行预处理，包括各层图像的噪声和冗余数据。根据数据值不同进行组织分类，对不同属性的数据赋予不同的颜色和不同的透明度。2）从显示屏幕上的每一个像素点按照观察方向发出一条投影光线，这条射线穿过三维数据场，沿着这条射线选择K个等距离的采样点，当采样点落在原来物体空间的网格点上，就用原来的灰度值来代替；当采样点不落在网格点上时，就用该点临近8个网格点的值进5第五章图像重建行插值得出新采样点不透明度值和颜色值。透射光线物体空间二维图像•3）采用表面明暗计算增强三维逼真效果，突出显示不同组织的边界面。明暗计算通常是基于物体表面的法向信息。可用各采样点的梯度值来代替法向值。设三维数据场中某采样点的函数值以(),,fijk表示，则采用中心差分方法求出该采样点的梯度值，即()()()()()()1,,1,,2,1,,1,,,1,,12xyzGradfijkfijkGradfijkfijkGradfijkfijk=+−−⎡⎤⎣⎦=+−−⎡⎤⎣⎦=+−−⎡⎤⎣⎦2得到各采样点的梯度值之后，可以用光照模型计算出各采样点处的漫反射分量，更加突出地显示出体数据中的边界面。4）计算每条射线对屏幕像素点的贡献。就是沿每条射线方向从前往后根据各采样点的颜色值和不透明值，按照一定的规则加以合成，6第五章图像重建得出屏幕像素点的昀终灰度值或颜色值。将屏幕上各像素点的颜色值都计算出来之后，形成了一幅图像。2．从物体空间到图像空间的方法逐层、逐行、逐个地计算物体空间的三维空间网格点上的数据对屏幕像素的贡献，并加以合成，形成昀后的图像。关键：如何计算物体空间三维网格的数据对二维屏幕像素的贡献，需要确定物体空间的每一个采样点对屏幕上多个点的贡献。把全部采样点遍历一次之后，把屏幕上同一像素点的多次贡献叠加起来可得昀终贡献。著名的算法有L.Westover的足迹表法和Jwilhelms的体元投射法。5.5.4MATLAB中的三维图形举例[4]MATLAB中三维图像数据的绘制和图形辅助分析工具相当丰富。准备了多种图形修饰和色彩设置函数：设置图形视角光影效果消隐和透视镂空和裁切可绘制三维线性图和各种三维曲面；做一些特殊三维图形：三维直方图三维柄图7第五章图像重建三维饼图三维箭图„三维线性图形绘制t=(0:0.02:2)*pi;x=sin(t);y=cos(t);z=cos(2*t);plot3(x,y,z,'b-',x,y,z,'bd');%Plotlinesandpointsin3-Dspace.这里实际上已经使用了多个命令，%将多个曲线重叠。view([-70,58]);%VIEW(AZ,EL)andVIEW([AZ,EL])settheangleoftheview%fromwhichanobserverseesthecurrent3-Dplot.AZistheazimuth%orhorizontalrotationandEListheverticalelevation(bothindegrees).%Azimuthrevolvesaboutthez-axis,withpositivevaluesindicatingcounter-%clockwiserotationoftheviewpoint.Positivevaluesofelevation%correspondtomovingabovetheobject;negativevaluesmovebelow.%VIEW([XYZ])setstheviewangleinCartesiancoordinates.The%magnitudeofvectorX,Y,Zisignored.boxon;legend('chain','gem');8第五章图像重建„三维曲面－椭圆抛物面x=-6:0.3:6;y=x;[X,Y]=meshgrid(x,y);％产生二维的网格数据%[X,Y]=MESHGRID(x,y)transformsthedomainspecifiedbyvectors%xandyintoarraysXandYthatcanbeusedfortheevaluationoffunctionsof%twovariablesand3-Dsurfaceplots.TherowsoftheoutputarrayXarecopies%ofthevectorxandthecolumnsoftheoutputarrayYarecopiesofthevectory.Z=X.^2/9+Y.^2/9;%tocalculatethefunctionvaluesofthegridpointsubplot(1,2,1);mesh(X,Y,Z);%MESH(X,Y,Z,C)plotsthecoloredparametricmeshsurfacedefined%byfourmatrixarguments.TheviewpointisspecifiedbyVIEW.Theaxis9第五章图像重建%labelsaredeterminedbytherangeofX,YandZ,orbythecurrentsettingof%AXIS.ThecolorscalingisdeterminedbytherangeofC,orbythecurrent%settingofCAXIS.Thescaledcolorvaluesareusedasindicesintothecurrent%COLORMAP.title('椭圆抛物面网线图');subplot(1,2,2);surf(X,Y,Z);%SURF(X,Y,Z,C)plotsthecoloredparametricsurfacedefinedby%fourmatrixarguments.TheviewpointisspecifiedbyVIEW.Theaxislabels%aredeterminedbytherangeofX,YandZ,orbythecurrentsettingof%AXIS.ThecolorscalingisdeterminedbytherangeofC,orbythecurrent%settingofCAXIS.Thescaledcolorvaluesareusedasindicesintothecurrent%COLORMAP.TheshadingmodelissetbySHADING.%SURF(X,Y,Z)usesC=Z,socolorisproportionaltosurfaceheight.title('椭圆抛物面网面图')10第五章图像重建„三维曲面－阔边帽面x=-6.5:0.5:6.5;y=x;%产生x和y两个阵列[X,Y]=meshgrid(x,y);%toproducethegriddataintwodimensionsR=sqrt(X.^2+Y.^2)+eps;%eps可以避免R在分母趋于0时无法定义Z=sin(R)./R;subplot(1,2,1);mesh(X,Y,Z);title('阔边帽面网线图')subplot(1,2,2);sur