用到的算法/数据结构多如牛毛:
各种树、图为主,其他如栈、队列、散列表、并查集。。。
贪心、回溯、动态规划、遗传算法、矩阵变换。。
在一个问题下很难回答好。。 先简单介绍一下和图相关的。
1. 和什么图打交道
CFG(Control Flow Graph)
控制流图是对程序中分支跳转关系的抽象,描述程序所有可能执行路径
节点是语句集合(basic block);
每个basic block有唯一入口和出口;
如果A到B有边,表示A执行完后可能执行B
PDG(Program Dependence Graph)
PDG在编译器中用得不多,常见于软件工程/安全相关的应用(程序切片、安全信息流等)
SSA(Single Static Assignment)
SSA简化了很多数据流分析问题。
其他图
DJ Graph, Loop Nesting Forest, Program Structure Tree等等。
可参考:IR for Program Analysis。下面主要介绍CFG
2. CFG初步处理
CFG构造
dominator树生成
在CFG中,如果A是B的dominator,则从程序入口执行到B的任意路径一定经过A
控制依赖分析
根据dominator和post-dominator分析依赖关系。数据依赖、控制依赖信息在自动并行化中尤其重要(如果循环的每次迭代都没有依赖,那么可以并行处理)
控制流图化简
在复杂度相同的情况下,CFG的规模影响算法的效果。如果一个CFG仅通过如下变换能化简为一个节点,则它是可化简的:
如果节点n有唯一的前驱,那么将其和其前驱合并为一个节点
如果节点存在到自身的边,那么将该边删除
构造SSA
SSA可以由CFG构造。
3. CFG与数据流分析
下面才进入主题。。
一般的文献介绍DFA(Data flow analysis),都会用几个基础的分析为例:Constant Propagation,Range propagation,Avaliable expressions,Reaching Definition。而Reaching Definition的一个应用,就是大家喜闻乐见的“跳转到定义处”(真要做到“智能”跳转并不简单)
这部分涉及东西较多,一些算法也和”图“并不直接相关,不再展开。
PS,很多DFA问题可以用graph reachability统一建模,强烈推荐此文:
Program analysis via graph reachability
软件调试技术包括:
1、分析和推理;
设计人员和开发人员根据软件缺陷问题的信息,分析和推理调试软件。
根据软件程序架构自顶向下缩小定位范围,确定可能发生问题的软件组件。
根据软件功能,软件运行时序定位软件问题。
根据算法原理,分析和确定缺陷问题发生的根源。
2、归纳类比法;
归纳法是一种从特殊推断一般的系统化思考方法,归纳法调试的基本思想是:从一些线索(错误征兆)着手,通过分析它们之间的关系来找出错误。该方法主要是根据积累的工作经验和案例处理调试工作。
根据工作经验和比对程序设计中类似问题的处理方式进行调试工作。
咨询相关部门和有经验的相关人员。
查找相关文档和案例,为处理问题提供思路和方法。在软件开发过程中,通常对每个缺陷问题进行跟踪管理,将解决问题的方案和过程详细记录。
收集出错的信息,列出数据,包括输入,输出,归纳整理,发现规律,从线索除法,寻找线索之间的联系。也就意味着:从特殊到一般。
3、跟踪回朔;
在小程序中常用的一种有效的调试方法,一旦发现了错误,人们先分析错误的征兆,确定最先发现“症状“的位置然后,人工沿程序的控制流程,向回追踪源程序代码,直到找到错误根源或确定错误产生的范围。
例如,程序中发现错误处是某个打印语句,通过输出值可推断程序在这一点上变量的值,再从这一点出发,回溯程序的执行过程,反复思考:“如果程序在这一点上的状态(变量的值)是这样,那么程序在上一点的状态一定是这样···“直到找到错误所在。
在软件开发通常采用基线与版本管理。基线为程序代码开发提供统一的开发基点,基线的建立有助于分清楚各个阶段存在的问题,便于对缺陷问题定位。软件版本在软件产品的开发过程中生成了一个版本树。软件产品实际上是某个软件版本,新产品的开发通常是在某个软件版本的基础上进行开发。
开发过程中发现有问题,可以回退至版本树上的稳定版本,查找问题根源。
通过基线版本序列可以追踪产品的各种问题,可以重新建立基于某个版本的配置,可以重现软件开发过程中的软件缺陷和各种问题,进行定位并查找问题根源。
4、增量调试;
软件开发大多采用软件配置管理和持续集成技术。开发人员每天将评代码提交到版本库。持续集成人员完成集成构建工作。
可以通过控制持续集成的粒度(构建时间间隔),控制开发人员提交到版本库的程序代码量,从而便于对缺陷问题定位。
通常每天晚上进行持续集成工作,发现问题时,开发人员实际上只需要调试处理当天编写的代码。
5、写出能重现问题的最短代码;
采用程序切片和插桩技术写出能重现问题的最短代码调试软件模块。
程序切片程序切片是通过在特定位置消除那些不影响表达式计算的所有语句,把程序减少到最小化形式,并仍能产生给定的行为。
使用切片技术,可以把一个规模较大并且较复杂的软件模块转换成多个切片程序。这些切片程序相对原来的程序,简单并且易于调试和测试。
程序插桩程序插桩方法是在被测程序中插入某些语句或者程序段来获取各种信息。通过这些信息进一步了解执行过程中程序的一些动态特性。一个软件组件的独立调试和测试需要采用插桩技术,该组件调用或运行需要桩模块。在软件模块的调试过程中程序切片和程序插桩可以结合起来使用。
6、日志追踪技术;
日志是一种记录机制,软件模块持续集成构建过程中,日志文件记录了有用信息。若构建失败,通过查看日志文件,将信息反馈给相关人员进行软件调试。
7、调试和测试融合的技术;
测试驱动开发。
测试驱动开发是一种不同于传统软件开发流程的开发方法。在编写某个功能的代码之前先编写测试代码,然后编写测试通过的功能代码,这有助于编写简洁可用和高质量的代码。
开发与测试融合。
程序开发人员除了进行程序代码的编写,白盒测试,也要完成基本的功能测试设计和执行。这样有助于程序开发人员更好地开展调试工作。
程序开发人员可以通过交叉测试来解决测试心理学的问题(不能自己测试自己)。采用这种模式测试人员的数量会减少,专业的测试人员去做其他复杂的测试工作。
研发中的很多低级缺陷会尽早在开发过程中被发现,从而减少缺陷后期发现的成本。
8、强行排错;
这种调试方法目前使用较多,效率较低,它不需要过多的思考,比较省脑筋。例如:
通过内存全部打印来调试,在这大量的数据中寻找出错的位置。
在程序特定位置设置打印语句,把打印语句插在出错的源程序的各个关键变量改变部位,重要分支部位,子程序调用部位,跟踪程序的执行,监视重要变量的变化
自动调用工具,利用某些程序语言的调试功能或专门的交互式调试工具,分析程序的动态过程,而不必修改程序。
应用以上任一种方法之前,都应当对错误的征兆进行全面彻底的分析,得出对出错位置及错误性质的推测,再使用一种适当的调试方法来检验推测的正确性。
9、演绎法调试;
演绎法是一种从一般原理或前提出发,经过排除和精华的过程来推导出结论的思考方法,演绎法排错是测试人员首先根据已有的测试用例,设想及枚举出所有可能出错的原因作为假设,然后再用原始测试数据或新的测试,从中逐个排除不可能正确的假设,最后,再用测试数据验证余下的假设确是出错的原因。
列举所有可能出错原因的假设,把所有可能的错误原因列成表,通过它们,可以组织,分析现有数据。
利用已有的测试数据,排除不正确的假设。
仔细分析已有的数据,寻找矛盾,力求排除前一步列出所有原因,如果所有原因都被排除了,则需要补充一些数据(测试用例),以建立新的假设。
改进余下的假设;
利用已知的线索,进一步改进余下的假设,使之更具体化,以便可以精确地确定出错位置;
证明余下的假设。
扩展资料:
软件调试技术的内容:
CPU的调试支持,包括异常、断点、单步执行、分支监视、JTAG、MCE等。
Windows *** 作系统中的调试设施,包括内核调试引擎、用户态调试予系统、验证器、Dr.Watson、WER、ETW、故障转储、WHEA等。
VisualC/C++编译器的调试支持,包括编译期检查、运行期检查,以及调试符号。
WinDBG调试器的发展历史、模块结构、工作模型、使用方法、主要调试功能的实现细节,以及遍布全书的应用实例。
内核调试、用户态调试、JIT调试、远程调试的原理、实现和用法。异常的概念、分发方法、处理方法(SEH、VEH、CppEH),未处理异常,以及编译器编译异常处理代码的方法。
调试符号的作用、产生过程、存储格式和使用方法。栈和堆的结构布局、工作原理和有关的软件问题,包括栈的自动增长和溢出,缓;中区溢出,溢出攻击,内存泄漏,堆崩溃等。
参考资料:百度百科-软件调试
Marching Cubes算法(医学图像三维绘制中的面绘制)2007-08-16 00:50建议要看的资料[1] Lorensen W E, Cline H E .Marching cubes: a high-resoulution 3D suface construction algorithm [J], Computer Graphics,1987, 21(4):163~169
[2]集成化医学影像算法平台理论与实践田捷,赵明昌,何晖光 清华大学出版社2004年10月
[3]Polygonising a scalar field Also known as: "3D Contouring", "Marching Cubes", "Surface Reconstruction"
http://local.wasp.uwa.edu.au/~pbourke/geometry/polygonise/Marching Cubes;
[4]www.3dmed.net
Marching Cubes算法工作原理
Marching Cubes算法是三维数据场等值面生成的经典算法,是体素单元内等值面抽取技术的代表。
等值面是空间中所有具有某个相同值的点的集合。它可以表示为, ,c是常数。则称F(f)为体数据f中的等值面。
在MC算法中,假定原始数据是离散的三维空间规则数据场。用于医疗诊断的断层扫描(CT)及核磁共振成像(MRI) 等产生的图像均属于这一类型。MC算法的基本思想是逐个处理数据场中的体素,分类出与等值面相交的体素,采用插值计算出等值面与体素棱边的交点。根据体素中每一顶点与等值面的相对位置,将等值面与立方体边的交点按一定方式连接生成等值面,作为等值面在该立方体内的一个逼近表示。在计算出关于体数据场内等值面的有关参数后山常用的图形软件包或硬件提供的面绘制功能绘制出等值面。
图2.1 离散的三维空间规则数据场中的一个体素
2.1.1 MC算法的主要步骤
1. 确定包含等值面的体素
离散的三维空间规则数据场中的一个体素可以用图2.1表示。8个数据点分别位于该体素的8个角点上。MC算法的基本假设是:沿着体素的边其数据场呈局部连续线性变化,根据这个假设,可认为,如果两个相邻采样点一个为正点,一个为负点,则它们连成的边上一定存在且仅有一个等值点 (设等值面值为c)。如果得到了体素各条边上的等值点,就可以以这些点为顶点,用一系列的三角形拟合出该体素中的等值面。因此确定立方体体素中等值面的分布是该算法的基础。
首先对体素的8个顶点进行分类,以判断其顶点是位于等值面之外,还是位于等值面之内。再根据8个顶点的状态,确定等值面的剖分模式。顶点分类规则为:
1. 如体素顶点的数据值大于或等于等值面的值,则定义该顶点位于等值面之外, 记为正点,即“1“
2. 如体素顶点的数据值小于等值面的值,则定义该顶点位于等值面之内,记为负点, 即“0"
由于每个体素共有8个顶点,且每个顶点有正负两种状态,所以等值面可能以 =256种方式与一个体素相交。通过列举出这256种情况,就能创建一张表格,利用它可以查出任意体素中的等值面的三角面片表示。如果考虑互补对称性,将等值面的值和8个角点的函数值的大小关系颠倒过来,即将体素的顶点标记置反(0变为1, 1变为0),这样做不会影响该体素的8个角点和等值面之间的拓扑结构,可将256种方式简化成128种。其次,再利用旋转对称性,可将这128种构型进一步简化成15种。图3.2给出了这15种基本构型[131其中黑点标记为“1”的角点。
图2.2 分布状态表
图2.3 体素角点分布不同情况
基于上面的分析,MC算法中用一个字节的空间构造了一个体素状态表,如图2.2所示,该状态表中的每一位可表示出该体元中的一个角点的0或1的状态。根据这一状态表,就可知道当前体素属于图2.3中哪一种情况,以及等值面将与哪一条边相交。
2.求等值面与体元边界的交点
在确定体素内三角剖分模式后,就要计算三角片顶点位置。当三维离散数据场的密度较高时,即当体素很小时,可以假定函数值沿体素边界呈线性变化,这就是MC算法的基本假设。因此,根据这一基本假设,可以直接用线性插值计算等值面与体素边的交点。
对于当前被处理体素的某一条边,如果其两顶点 , 的标记值不同,那么等值面一定与此边相交,且仅有一个交点。交点为 其中P代表等值点坐标, , 代表两个端点的坐标, , 代表两个端点的灰度值,v代表域值。求出等值面与体素棱边的交点以后,就可以将这些交点连接成三角形或多边形,形成等值面的一部分。
3.等值面的法向量计算
为了利用图形硬件显示等值面图象,必须给出形成等值面的各三角面片的法向分量,选择适当的局部面光照模型进行光照计算,以生成真实感图形。
对于等值面上的每一点,其沿面的切线方向的梯度分量应该是零,因此,该点的梯度矢量的方向也就代表了等值面在该点的法向量,当梯度值非零。所幸的是等值面往往是由两种具有不同密度的物质的分解面,因此其上的每点的梯度矢量均不为零,即
Mc算法采用中心差分方法求采样点p〔m ,n, k ) 处的梯度矢量,公式如下:
Gx=〔g(i+1,j,k)-g(i-1,j,k)〕/2dx
Gy=〔g(i,j+1,k)-g(i,j-1,k)〕/2dy
Gz=〔g(i,j,k+1)-g(i,j,k-1)〕/2dz
其中D(i,j ,k)是切片k在像素(i,j)的密度, , , 是立方体边的长度。对g进行归一化,得到(gx/|g|,gy/|g|,gz/|g|)作为(i,j,k)上的单位法向量。然后,对体素八个顶点上法向量进行线性插值就可得到位于体素棱边上的三角片的各个顶点上的法向量。设计算得到的某个三角片的三个顶点上的单位法向量分别为( , 和 ),这个三角片的几何重心为 ,则该三角片的法向量起始于 ,终止于 。代入Gourand光照模型公式,就可计算出小三角片表面的光强(灰度)。将其投影在某个特定的二维平面上进行显示,从而显示出物体富有光感的整个表面形态。其中我们在内存中保留四个切片来计算立方体中所有顶点梯度。
2.1.2 MC算法流程
1、将三维离散规则数据场分层读入内存
2、扫描两层数据,逐个构造体素,每个体素中的8个角点取自相邻的两层
3、将体素每个角点的函数值与给定的等值面值c做比较,根据比较结果,构造
该体素的状态表
4、根据状态表,得出将与等值面有交点的边界体素
5、通过线性插值方法计算出体素棱边与等值面的交点
6、利用中心差分方法,求出体素各角点处的法向量,再通过线性插值方法,求出三角面片各顶点处的法向
7,根据各三角面片上各顶点的坐标及法向量绘制等值面图像。
========================================================
MC代码
MarchingCubes(float lowThreshold,float highThreshold,float XSpace,float YSpace,float ZSpace)
{
//记录生成的顶点数和面数初始时应该为0
m_vNumber=0
m_fNumber=0
//当前Cube中生成的顶点和面数
int vertPos,facePos
//包围盒的尺寸用于绘制程序计算整个场景的包围盒,用于调整观察位置,以使整个场景尽可能占满整个窗口。
float min[3],max[3]
min[0]=min[1]=min[2]=max[0]=max[1]=max[2]=0//初始化
//当前扫描层的切片数据和一个临时的切片数据
short *pSliceA,*pSliceB,*pSliceC,*pSliceD,*tempSlice
pSliceA=pSliceB=pSliceC=tempSlice=NULL
int imageWidth,imageHeight,imageSize,sliceNumber
imageWidth=imageHeight=512//我们是512×512的数据
imageSize=imageWidth*imageHeight
sliceNumber=m_FileCount-1
if((highThreshold*lowThreshold)==0)
{
return 0
}
pSliceD =new short [imageSize]
//因为等值面是每相邻两层切片为单位进行提取的,所以在处理后两层切片时难免生成前两层切片已经生成的顶点,这时候就用下面的数组记录哪些边上的顶点已经生成了,如果遇到已经生成的顶点就不再重复计算而是直接使用记录的索引,否则就生成新的顶点。
long *bottomXEdge=new long[imageSize]
long *bottomYEdge=new long[imageSize]
long *topXEdge=new long[imageSize]
long *topYEdge=new long[imageSize]
long *zEdge=new long[imageSize]
tempSlice=new short [imageSize]
if(bottomXEdge==NULL||bottomYEdge==NULL||
topXEdge==NULL||topYEdge==NULL||
zEdge==NULL||tempSlice==NULL)
{
return 0//错误
}
//初始化数据
memset(bottomXEdge,-1,sizeof(long)*imageSize)
memset(bottomYEdge,-1,sizeof(long)*imageSize)
memset(topXEdge,-1,sizeof(long)*imageSize)
memset(topYEdge,-1,sizeof(long)*imageSize)
memset(zEdge,-1,sizeof(long)*imageSize)
memset(tempSlice,0,sizeof(short)*imageSize)
//计算某一层顶点和三角时所需要的一些变量
//一些循环变量
int i,j,k,w,r
//cube类型
unsigned char cubeType(0)
//计算法向量
float dx[8],dy[8],dz[8],squaroot
//记录某个Cube生成
float vertPoint[12][6]
int cellVerts[12]//what use
int triIndex[5][3]//每个cube最多生成5条边
//用于记录已生成顶点索引的临时变量
int offset
//当前cube8个顶点的灰度值
short cubegrid[8]
long *edgeGroup
//得到数据
pSliceD=m_volumeData
pSliceB=tempSlice
pSliceA=tempSlice
int tt,tt1
//扫描4层切片的顺序
/*
-----------------------D |
-----------------------B |
-----------------------C |
-----------------------A |
V
*/
//marching cubes 算法开始实行 ?第一次循环时,只读入一个切片?
for(i=0i<=(sliceNumber)i++)
{
pSliceC=pSliceA
pSliceA=pSliceB
pSliceB=pSliceD
if(i>=sliceNumber-2)
{
pSliceD=tempSlice
}
else
{
pSliceD+=imageSize
}
for(j=0j<imageHeight-1++j)
for(k=0k<imageWidth-1++k)
/* for(j=10j<imageHeight-5j++)//调试用
for(k=10k<imageWidth-5k++)*/
{
//得到八个顶点的灰度值step0
cubegrid[0]=pSliceA[j*imageWidth+k]
cubegrid[1]=pSliceA[j*imageWidth+k+1]
cubegrid[2]=pSliceA[(j+1)*imageWidth+k+1]
cubegrid[3]=pSliceA[(j+1)*imageWidth+k]
cubegrid[4]=pSliceB[j*imageWidth+k]
cubegrid[5]=pSliceB[j*imageWidth+k+1]
cubegrid[6]=pSliceB[(j+1)*imageWidth+k+1]
cubegrid[7]=pSliceB[(j+1)*imageWidth+k]
//计算cube的类型
cubeType=0
for(w=0w<8w++)
{
if((cubegrid[w]>lowThreshold)&&(cubegrid[w]<highThreshold))//需要画的点
{
cubeType|=(1<<w)
}
}//end for计算cube的类型
if(cubeType==0||cubeType==255)
{
continue
}
for(w=0w<12w++) //初始化cellVerts表到零
{
cellVerts[w]=-1
}
//计算6个方向相邻点的象素差值(用于计算法向量)
if(k==0)
{
dx[0]=pSliceA[j*imageWidth+1]
dx[3]=pSliceA[(j+1)*imageWidth+1]
dx[4]=pSliceB[j*imageWidth+1]
dx[7]=pSliceB[(j+1)*imageWidth+1]
}
else
{
dx[0]=pSliceA[j*imageWidth+k+1]
-pSliceA[j*imageWidth+k-1]
dx[3]=pSliceA[(j+1)*imageWidth+k+1]
-pSliceA[(j+1)*imageWidth+k-1]
dx[4]=pSliceB[j*imageWidth+k+1]
-pSliceB[j*imageWidth+k-1]
dx[7]=pSliceB[(j+1)*imageWidth+k+1]
-pSliceB[(j+1)*imageWidth+k-1]
}
if(k==imageWidth-2)
{
dx[1]=-pSliceA[j*imageWidth+imageWidth-2]
dx[2]=-pSliceA[(j+1)*imageWidth+imageWidth-2]
dx[5]=-pSliceB[j*imageWidth+imageWidth-2]
dx[6]=-pSliceB[(j+1)*imageWidth+imageWidth-2]
}
else
{
dx[1]=pSliceA[j*imageWidth+k+2]
-pSliceA[j*imageWidth+k]
dx[2]=pSliceA[(j+1)*imageWidth+k+2]
-pSliceA[(j+1)*imageWidth+k]
dx[5]=pSliceB[j*imageWidth+k+2]
-pSliceB[j*imageWidth+k]
dx[6]=pSliceB[(j+1)*imageWidth+k+2]
-pSliceB[(j+1)*imageWidth+k]
}
if(j==0)
{
dy[0]=pSliceA[imageWidth+k]
dy[1]=pSliceA[imageWidth+k+1]
dy[4]=pSliceB[imageWidth+k]
dy[5]=pSliceB[imageWidth+k+1]
}
else
{
dy[0]=pSliceA[(j+1)*imageWidth+k]
-pSliceA[(j-1)*imageWidth+k]
dy[1]=pSliceA[(j+1)*imageWidth+k+1]
-pSliceA[(j-1)*imageWidth+k+1]
dy[4]=pSliceB[(j+1)*imageWidth+k]
-pSliceB[(j-1)*imageWidth+k]
dy[5]=pSliceB[(j+1)*imageWidth+k+1]
-pSliceB[(j-1)*imageWidth+k+1]
}
if(j==imageHeight-2)
{
dy[2]=-pSliceA[(imageHeight-2)*imageWidth+k+1]
dy[3]=-pSliceA[(imageHeight-2)*imageWidth+k]
dy[6]=-pSliceB[(imageHeight-2)*imageWidth+k+1]
dy[7]=-pSliceB[(imageHeight-2)*imageWidth+k]
}
else
{
dy[2]=pSliceA[(j+2)*imageWidth+k+1]-pSliceA[j*imageWidth+k+1]
dy[3]=pSliceA[(j+2)*imageWidth+k]-pSliceA[j*imageWidth+k]
dy[6]=pSliceB[(j+2)*imageWidth+k+1]-pSliceB[j*imageWidth+k+1]
dy[7]=pSliceB[(j+2)*imageWidth+k]-pSliceB[j*imageWidth+k]
}
dz[0]=pSliceB[j*imageWidth+k]
-pSliceC[j*imageWidth+k]
dz[1]=pSliceB[j*imageWidth+k+1]
-pSliceC[j*imageWidth+k+1]
dz[2]=pSliceB[(j+1)*imageWidth+k+1]
-pSliceC[(j+1)*imageWidth+k+1]
dz[3]=pSliceB[(j+1)*imageWidth+k]
-pSliceC[(j+1)*imageWidth+k]
dz[4]=pSliceD[j*imageWidth+k]
-pSliceA[j*imageWidth+k]
dz[5]=pSliceD[j*imageWidth+k+1]
-pSliceA[j*imageWidth+k+1]
dz[6]=pSliceD[(j+1)*imageWidth+k+1]
-pSliceA[(j+1)*imageWidth+k+1]
dz[7]=pSliceD[(j+1)*imageWidth+k]
-pSliceA[(j+1)*imageWidth+k]
//计算三角形顶点的坐标和梯度
vertPos=0
facePos=0
for(w=0w<12w++)
{
if(g_EdgeTable[cubeType]&(1<<w)) //what …..
{
//根据g_edgeTable[256]对应值判断cube的那一条边与等值面有交点
switch(w)
{
case 0:
offset=j*imageWidth+k
edgeGroup=bottomXEdge
break
case 1:
offset=j*imageWidth+k+1
edgeGroup=bottomYEdge
break
case 2:
offset=(j+1)*imageWidth+k
edgeGroup=bottomXEdge
break
case 3:
offset=j*imageWidth+k
edgeGroup=bottomYEdge
break
case 4:
offset=j*imageWidth+k
edgeGroup=topXEdge
break
case 5:
offset=j*imageWidth+k+1
edgeGroup=topYEdge
break
case 6:
offset=(j+1)*imageWidth+k
edgeGroup=topXEdge
break
case 7:
offset=j*imageWidth+k
edgeGroup=topYEdge
break
case 8:
offset=j*imageWidth+k
edgeGroup=zEdge
break
case 9:
offset=j*imageWidth+k+1
edgeGroup=zEdge
break
case 10:
offset=(j+1)*imageWidth+k+1
edgeGroup=zEdge
break
case 11:
offset=(j+1)*imageWidth+k
edgeGroup=zEdge
break
}//对应switch的{。。。end for//根据g_EdgeTable对应值判断cube的那一条边与等值面有交点
//该边上的顶点是否已经在上一层中生成
if(edgeGroup[offset]==-1)
{
int index1,index2
short s1,s2,s
float x1,y1,z1,nx1,ny1,nz1
float x2,y2,z2,nx2,ny2,nz2
//得到该边两端点的索引进而得到两点的灰度值
index1=g_CoordTable[w][3]
index2=g_CoordTable[w][4]
s1=cubegrid[index1]
s2=cubegrid[index2]
if(s1<highThreshold&&s1>lowThreshold)
{
if(s2>=highThreshold)
{
s=highThreshold
}
else if(s2<=lowThreshold)
{
s=lowThreshold
}
}
else if(s2<highThreshold&&s2>lowThreshold)
{
if(s1>=highThreshold)
{
s=highThreshold
}
else if(s1<=lowThreshold)
{
s=lowThreshold
}
}
//计算两端点实际坐标
x1=(k+g_CoordVertex[index1][0])*XSpace
y1=(j+g_CoordVertex[index1][1])*YSpace
z1=(i+g_CoordVertex[index1][2])*ZSpace
x2=(k+g_CoordVertex[index2][0])*XSpace
y2=(j+g_CoordVertex[index2][1])*YSpace
z2=(i+g_CoordVertex[index2][2])*ZSpace
//计算两端点的法向量
nx1=dx[index1]/XSpace
ny1=dy[index1]/YSpace
nz1=dz[index1]/ZSpace
nx2=dx[index2]/XSpace
ny2=dy[index2]/YSpace
nz2=dz[index2]/ZSpace
float factor=((float)(s-s1))/((float)(s2-s1))
//插值计算交点坐标
vertPoint[vertPos][0]=factor*(x2-x1)+x1
vertPoint[vertPos][1]=factor*(y2-y1)+y1
vertPoint[vertPos][2]=factor*(z2-z1)+z1
//计算法向量
vertPoint[vertPos][3]=factor*(nx1-nx2)-nx1
vertPoint[vertPos][4]=factor*(ny1-ny2)-ny1
vertPoint[vertPos][5]=factor*(nz1-nz2)-nz1
//法向量归一化
squaroot=sqrt(vertPoint[vertPos][3]*vertPoint[vertPos][3]+vertPoint[vertPos][4]*vertPoint[vertPos][4]
+vertPoint[vertPos][5]*vertPoint[vertPos][5])
if(squaroot<=0)squaroot=1.0
vertPoint[vertPos][3]/=squaroot
vertPoint[vertPos][4]/=squaroot
vertPoint[vertPos][5]/=squaroot
//更新包围盒数据
if(min[0]>vertPoint[vertPos][0])
{
min[0]=vertPoint[vertPos][0]
}
if(min[1]>vertPoint[vertPos][1])
{
min[1]=vertPoint[vertPos][1]
}
if(min[2]>vertPoint[vertPos][2])
{
min[2]=vertPoint[vertPos][2]
}
if(max[0]<vertPoint[vertPos][0])
{
max[0]=vertPoint[vertPos][0]
}
if(max[1]<vertPoint[vertPos][1])
{
max[1]=vertPoint[vertPos][1]
}
if(max[2]<vertPoint[vertPos][2])
{
max[2]=vertPoint[vertPos][2]
}
//记录新生成的顶点索引
cellVerts[w]=m_vNumber
edgeGroup[offset]=cellVerts[w]
m_vNumber++
vertPos++
} //end if(edgeGroup[offset]==-1) ////
else
{
//若该点已经在上一层生成,则直接得到其索引
cellVerts[w]=edgeGroup[offset]
}
} // end对应if(g_EdgeTable[cubeType]&(1<<w)) //
} //对应for(w=0w<12w++)
//保存当前cubes 顶点和法向量
tt1=m_vNumber-vertPos
for(tt=0tt<vertPostt++)
{
vPointNomal[tt1+tt][0]=vertPoint[tt][0]
vPointNomal[tt1+tt][1]=vertPoint[tt][1]
vPointNomal[tt1+tt][2]=vertPoint[tt][2]
vPointNomal[tt1+tt][3]=vertPoint[tt][3]
vPointNomal[tt1+tt][4]=vertPoint[tt][4]
vPointNomal[tt1+tt][5]=vertPoint[tt][5]
}
// memcpy(vPointNomal+6*(m_vNumber-vertPos) ,vertPoint,sizeof(float)*6*vertPos)
//记录新生成的三角面片数据
w=0
while (g_TriTable[cubeType][w]!=-1)
{
for(r=0r<3r++)
{
triIndex[facePos][r]=cellVerts[g_TriTable[cubeType][w++]]
if(triIndex[facePos][r]<0)
{
AfxMessageBox("有问题",MB_OK,0)
}
}
facePos++
m_fNumber++
} //end 对应while (g_TriTable[cubeType][w]!=-1)
//保存面数据
tt1=m_fNumber-facePos
for(tt=0tt<facePostt++)
{
pFace[tt1+tt][0]=triIndex[tt][0]
pFace[tt1+tt][1]=triIndex[tt][1]
pFace[tt1+tt][2]=triIndex[tt][2]
}
// memcpy(pFace+3*(m_fNumber-facePos)*sizeof(long),triIndex,sizeof(int)*3*facePos)
} memcpy(bottomXEdge,topXEdge,sizeof(short)*imageSize)
memcpy(bottomYEdge,topYEdge,sizeof(short)*imageSize)
memset(topXEdge,-1,sizeof(short)*imageSize)
memset(topYEdge,-1,sizeof(short)*imageSize)
memset(zEdge,-1,sizeof(short)*imageSize)
}
delete []tempSlice
delete []bottomXEdge
delete []bottomYEdge
delete []topXEdge
delete []topYEdge
delete []zEdge
return 1
}
在OnDraw
glBegin(GL_TRIANGLES)
for(i=0i<pDoc->m_fNumberi++)
{
glNormal3fv(&(pDoc->vPointNomal[pDoc->pFace[ i][0]][3]))
glVertex3fv(&(pDoc->vPointNomal[pDoc->pFace[ i][0]][0]))
glNormal3fv(&(pDoc->vPointNomal[pDoc->pFace[ i ][1]][3]))
glVertex3fv(&(pDoc->vPointNomal[pDoc->pFace[ i ][1]][0]))
glNormal3fv(&(pDoc->vPointNomal[pDoc->pFace[ i ][2]][3]))
glVertex3fv(&(pDoc->vPointNomal[pDoc->pFace[ i ][2]][0]))
}
glEnd()
以上代码只用于理解,未测试
欢迎分享,转载请注明来源:内存溢出
评论列表(0条)