int Dim[3] = { 200,200,200 };//体数据维度大小
int* Data = (int*)malloc(sizeof(int) * Dim[0] * Dim[1] * Dim[2]);
GLubyte CData[200][200][200][4];//存储颜色和不透明度
glm::vec4 smallCubeC = glm::vec4(1.0, 1.0, 0.0, 1.0);//小立方体颜色
glm::vec4 middleSphereC = glm::vec4(1.0, 0.0, 0.0, 1.0);//中间球体颜色
glm::vec4 largeCubeC = glm::vec4(1.0, 1.0, 1.0, 1.0);//大立方体颜色
float smallCubeD = 0.05;//小立方体不透明度
float middleSphereD = 0.015;//中间球体不透明度
float largeCubeD = 0.018;//大立方体不透明度
void GenCube(int x, int y, int z, int side, int density, int* Data, int* Dim)
{
	int max_x = x + side, max_y = y + side, max_z = z + side;
	int Dimxy = Dim[0] * Dim[1];
	for (int k = z; k < max_z; k++)
	{
		for (int j = y; j < max_y; j++)
		{
			for (int i = x; i < max_x; i++)
			{
				Data[k * Dimxy + j * Dim[0] + i] = density;
			}
		}
	}
}


void GenSphere(int x, int y, int z, int radius, int density, int* Data, int* Dim)
{
	int radius2 = radius * radius;
	int Dimxy = Dim[0] * Dim[1];
	for (int k = 0; k < Dim[2]; k++)
	{
		for (int j = 0; j < Dim[1]; j++)
		{
			for (int i = 0; i < Dim[0]; i++)
			{
				if ((i - x) * (i - x) + (j - y) * (j - y) + (k - z) * (k - z) <= radius2)
				{
					Data[k * Dimxy + j * Dim[0] + i] = density;
				}
			}
		}
	}
}

void Classify(GLubyte CData[200][200][200][4], int* Data, int* Dim)//按照所在位置为每个体数据点赋值，颜色和不透明度
{
	int* LinePS = Data;
	for (int k = 0; k < Dim[2]; k++)
	{
		for (int j = 0; j < Dim[1]; j++)
		{
			for (int i = 0; i < Dim[0]; i++)
			{
				if (LinePS[0] <= 100)
				{
					//白色
					CData[i][j][k][0] = 255.0 * largeCubeC[0];
					CData[i][j][k][1] = 255.0 * largeCubeC[1];
					CData[i][j][k][2] = 255.0 * largeCubeC[2];
					CData[i][j][k][3] = largeCubeD*255.0;
				}
				else if (LinePS[0] <= 200)
				{
					//红色
					CData[i][j][k][0] = 255.0 * middleSphereC[0];
					CData[i][j][k][1] = 255.0 * middleSphereC[1];
					CData[i][j][k][2] = 255.0 * middleSphereC[2];
					CData[i][j][k][3] = middleSphereD*255.0;
				}
				else
				{
					//黄色
					CData[i][j][k][0] = 255.0 * smallCubeC[0];
					CData[i][j][k][1] = 255.0 * smallCubeC[1];
					CData[i][j][k][2] = 255.0 * smallCubeC[2];
					CData[i][j][k][3] = smallCubeD*255.0;
				}
				LinePS++;
			}
		}
	}
	//return CDdata[200][200][200][4];
}

void GenerateVolume(int* Data, int* Dim)
{
	GenCube(0, 0, 0, 200, 100, Data, Dim);//大正方体
	GenSphere(100, 100, 100, 80, 200, Data, Dim);//球体
	GenCube(70, 70, 70, 60, 300, Data, Dim);//小正方体
}

                          
                        

                          
                            //volume texture ID
GLuint textureID;

bool LoadVolume() {
	GenerateVolume(Data, Dim);//生成原始体数据
	Classify(CData, Data, Dim);//对体数据分类赋予对应颜色值和不透明度

	//generate OpenGL texture
	glGenTextures(1, &textureID);
	glBindTexture(GL_TEXTURE_3D, textureID);
	// set the texture parameters
	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_S, GL_CLAMP);
	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_T, GL_CLAMP);
	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_R, GL_CLAMP);
	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MAG_FILTER, GL_LINEAR);
	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MIN_FILTER, GL_LINEAR_MIPMAP_LINEAR);
	//set the mipmap levels (base and max)
	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_BASE_LEVEL, 0);
	glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_MAX_LEVEL, 4);
	glTexImage3D(GL_TEXTURE_3D, 0, GL_RGBA, 200, 200, 200, 0, GL_RGBA, GL_UNSIGNED_BYTE, CData);//将体数据存入3D纹理
	GL_CHECK_ERRORS

	//generate mipmaps
	glGenerateMipmap(GL_TEXTURE_3D);
	return true;
}

                          
                        

                          
                            #version 330 core
layout(location = 0) in vec3 vVertex; //object space vertex position
//uniform
uniform mat4 MVP;   
smooth out vec3 vUV; 
void main()
{  
	gl_Position = MVP*vec4(vVertex.xyz,1);
	vUV = vVertex + vec3(0.5);
}

                          
                        

                          
                            #version 330 core
layout(location = 0) out vec4 vFragColor;	//fragment shader output
smooth in vec3 vUV;				//用于采样3D纹理的坐标
								
uniform sampler3D	volume;		//体数据纹理
uniform vec3		camPos;		//相机位置
uniform vec3		step_size;	//采样步长
//constants
const int MAX_SAMPLES = 300;	
const vec3 texMin = vec3(0);	//最小纹理坐标
const vec3 texMax = vec3(1);	//最大纹理坐标

void main()
{ 
	vec3 dataPos = vUV;    //光线投射起始点坐标
	vec3 geomDir;          //光线步进方向
	if(abs(camPos.x)<=0.5&&abs(camPos.y)<=0.5&&abs(camPos.z)<=0.5)//当相机也就是视点位于体渲染对象内部时，起始点应该改为相机视点的位置坐标作为起始点
	{
		dataPos=camPos+vec3(0.5);
	}
	geomDir = normalize((vUV-vec3(0.5)) - camPos); //由视点坐标和起始点坐标相减得到沿视线方向步进的方向的向量
	vec3 dirStep = geomDir * step_size; 
	bool stop = false; 
	vec4 cumc=vec4(0);
	//沿射线方向采样累积颜色和不透明度
	for (int i = 0; i < MAX_SAMPLES; i++) {
		dataPos = dataPos + dirStep;		
		stop = dot(sign(dataPos-texMin),sign(texMax-dataPos)) < 3.0;
		if (stop) 
			break;
		vec4 samplec=texture(volume, dataPos).rgba;//获取采样点颜色值和不透明度
		cumc[0]+=samplec.r*samplec[3]*(1-cumc[3]);
		cumc[1]+=samplec.g*samplec[3]*(1-cumc[3]);
		cumc[2]+=samplec.b*samplec[3]*(1-cumc[3]);
		cumc[3]+=samplec.a*(1-cumc[3]);
		if( cumc[3]>0.99)
			break;
	} 
	vFragColor=cumc.rgba;
}