体绘制
体绘制中光线投射的主要步骤如下：
（1）计算光线上的点与体数据交点，并且将所有的交点进行累加
vec4 DirectVolumeRendering(vec3 front, vec3 back, float upperThreshold, float threshold, float offset)
{
//1.0计算视点的位置和投影的方向
vec3 dir = back - front;
float len = length(dir);
vec3 norm_dir = normalize(dir);
vec4 result = vec4(0,0,0,0);//RGBA
float raylength = 0;
//2.0计算光线穿过体数据累加投影得到的RGBA值
while(result.a < 1.0 && raylength < len )
{
//2.1得到光线在raylength长度处与体数据相交的位置pos
vec3 pos = front+raylength*norm_dir;
//2.2得到体数据中pos位置的灰度值
float valuef = texture3D(blockPool, pos).a;
//2.3根据灰度值，对照颜色表，得到对应的RGBA值
vec4 texValue = texture1D(tranf, valuef);
//2.4根据光线传输的距离，对透明值做加权处理
texValue.a = texValue.a *stepsize *sizeX;
//2.5根据pos位置处RGBA值，以及之前计算得到的result值，继续累加得到result
result.rgb += (1.0-result.a)*texValue.a*texValue.rgb;
result.a += (1.0-result.a)*texValue.a;
//2.6光线长度继续向外延伸
raylength += stepsize;
}
return result;
}
（2）添加阴影，如果光线上有一点的值大于阈值Outsideisovalue,则认为该光线需要添加阴影
float Shadow(vec3 pos, vec3 raydir)
{
//SHADOWS
vec4 MV_pos = gl_ModelViewMatrix*vec4(pos, 1.0);
//vec3 shadowdir = normalize((MV_pos.xyz/MV_pos.w) - (MV_light.xyz/MV_light.w));
vec3 shadowdir = -normalize(lightposition - pos);
vec3 delta_shadow = shadowdir*0.002;
//pos = pos-normalize(raydir)*stepsize;
float shadowsample;
for(int i = 0; i < 1400; i++)
{
shadowsample = texture3D(volumetexture,pos).x;
if(shadowsample > Outsideisovalue+0.02)
{
return 1.0;
}
pos += delta_shadow;
if (!(pos.x>0 && pos.x<1 && pos.y>0 && pos.y<1 && pos.z>0 && pos.z<1))
break;
}
return 0.0;
}
（3）添加光照模型
vec3 PhongShading(vec3 pos, vec3 diffuseColor)
{
vec3 shadedValue=vec3(0,0,0);
vec3 G;
G.x = (texture3D(volumetexture, pos+vec3(stepsize,0,0) ).x - texture3D(volumetexture, pos+vec3(-stepsize,0,0) ).x);
G.y = (texture3D(volumetexture, pos+vec3(0,stepsize,0) ).y - texture3D(volumetexture, pos+vec3(0,-stepsize,0) ).y);
G.z = (texture3D(volumetexture, pos+vec3(0,0,stepsize) ).z - texture3D(volumetexture, pos+vec3(0,0,-stepsize) ).z);
G = gl_NormalMatrix * normalize(G);
G = normalize(G);
vec4 MV_light = (gl_ModelViewMatrix*vec4(lightposition.xyz,1.0));
vec3 L = normalize(MV_light.xyz/MV_light.w -eposition);
vec3 V = -normalize(eposition.xyz);
//specularcolor
vec3 H = normalize(V+L);
float GdotH = pow(max(dot(G, H), 0.0), floor(shininess));
float GdotL = max(dot(G, L), 0.0);
vec3 specular = ks * vec3(1.0,1.0,1.0) * GdotH;
vec3 diffuse = kd * diffuseColor.rgb * GdotL;
vec3 ambient = ka * diffuseColor.rgb;
shadedValue = specular + diffuse + ambient;
return shadedValue;
}
实现效果如下：
实现过程中的几点说明：
在OpenGL中，将模型显示在窗口中需要经过这样几个步骤：
（1）建立模型，得到模型顶点的坐标和法向量；对于复杂模型，还需要给出顶点顺序。
（2）建立世界坐标系，将模型放入到世界坐标系中，即将模型的顶点坐标和法向量等信息转为世界坐标系中的坐标。
（3）建立视图坐标系，即设定观察坐标系，通过确定视点，观察点和向上向量，三个向量确定视图坐标系。（lookat方法）
（4）设定透视方式，在视图坐标系的基础上选择裁剪方式，可以想象为相机调节焦距来控制显示的范围。
（5）设定画布的大小，可以设定最后输出的窗口大小。
这里主要是说明一下LookAt的使用，该方法有9个参数，分别是视点坐标、目标位置坐标和向上坐标方向。实际上就是确定了一个视图坐标系，也可以通过设定当前矩阵来替代LookAt方法。
这里前三列向量表示该坐标系的三个坐标轴，第四列实际上是原点的位置。下面给出我实现的方法：
public double[] getCoordinate(double []vector){
double []angle = new double[2];
double len = Math.sqrt(vector[0]*vector[0] +vector[1]*vector[1] +vector[2] * vector[2]);
if(Math.abs(len)<0.001) {
len = 1;
}
for(int i = 0; i < 3; i ++){
vector[i] = vector[i]/len;
}
angle[0] = Math.asin(vector[2]);
angle[1] = Math.acos(vector[0]/Math.cos(angle[0]));
return angle;
}
public void NormalVector(double[] vector) {
double len = Math.sqrt(vector[0] * vector[0] + vector[1] * vector[1] + vector[2] * vector[2]);
if (Math.abs(len) < 0.001) {
len = 1;
}
for (int i = 0; i < 3; i++) {
vector[i] = vector[i] / len;
}
}
public double[] VectorCMulVector(double []vector1,double []vector2){
double[] vector = {0,0,0};
vector[0] = vector1[1] * vector2[2] - vector1[2] * vector2[1];
vector[1] = vector1[2] * vector2[0] - vector1[0] * vector2[2];
vector[2] = vector1[0] * vector2[1] - vector1[1] * vector2[0];
return vector;
}
public void gluLookAt(GL2 gl,double eye0 , double eye1, double eye2, double ref0 , double ref1 ,double ref2 , double up_dir0, double up_dir1, double up_dir2) {
double [] directMat = new double[16];
for(int i = 0; i < 16; i ++){
directMat[i] = 0;
}
directMat[15]= 1;
double []fvDirect = {ref0-eye0,ref1-eye1,ref2-eye2};
NormalVector(fvDirect);
double []fvUpD = {up_dir0,up_dir1,up_dir2,};
NormalVector(fvUpD);
double []fvC = VectorCMulVector(fvDirect, fvUpD);
NormalVector(fvC);
double []fvUp = VectorCMulVector(fvC, fvDirect);
NormalVector(fvUp);
fvDirect[0] = -fvDirect[0];
fvDirect[1] = -fvDirect[1];
fvDirect[2] = -fvDirect[2];
directMat[0] = fvC[0];
directMat[4] = fvC[1];
directMat[8] = fvC[2];
directMat[1] = fvUp[0];
directMat[5] = fvUp[1];
directMat[9] = fvUp[2];
directMat[2] = fvDirect[0];
directMat[6] = fvDirect[1];
directMat[10] = fvDirect[2];
directMat[3] = (float)eye0;
directMat[7] = (float)eye1;
directMat[11] = (float)eye2;
for(int i = 0; i < 16; i ++){
modelView[i] = (float)directMat[i];
}
}
这里的modelView就是当前坐标矩阵了。
在漫游时LookAt方法比较有效，能够很方便的进行调节。在体绘制中，需要得到视图矩阵，根据视图矩阵来进行光线投射。通过上面的方法，可以直接获取到视图矩阵，完全不依赖于OpenGL，自己来实现。
通过这个方法，想明白OpenGL坐标系的关系，对OpenGL做深入的研究，可以说矩阵变换是基础，不可不深入了解。
具体参考了http://blog.csdn.net/liuzhidejll/article/details/7209640博客
在体绘制的过程中主要参考了下面一些资料：
[1]康玉之,GPU 编程与 CG语言之阳春白雪下里巴人
[2]Advanced Illumination Techniques for GPU-Based Volume Raycasting
[3]CUDA Volume Rendering 简介
这些资料在我的资源里都可免费下载到，希望对体绘制的实现有一些帮助。