CSharpGL(25)一个用raycast实现体渲染VolumeRender的例子
本文涉及的VolumeRendering相关的C#代码是从（https://github.com/toolchainX/Volume_Rendering_Using_GLSL）的C++代码转换来的。
效果图
下载
CSharpGL已在GitHub开源，欢迎对OpenGL有兴趣的同学加入（https://github.com/bitzhuwei/CSharpGL）
实现思路
raycast
用一个3D纹理存储整个模型数据。
如下图所示，从画布上的一个像素点出发，垂直于画布的方向上，一条射线穿过3D纹理，每隔一小段距离采样一次颜色，累加起来就是此像素点应有的颜色。
起始点/终点
从射线接触3D纹理的第一个点开始（起始点），到射线离开3D纹理的位置（终点），这段距离就是要采样的范围。
终点
那么如何获取终点的位置？
办法是：渲染一个恰好包围3D纹理的立方体，且只渲染此立方体的背面（用glCullFace(GL_FRONT);）。因为背面就是终点啊。另外，要把这个渲染结果弄到一个2D纹理上。这就需要一个与画布大小相同的2D纹理来记录终点的位置，即需要一个新的FrameBuffer。（详情可参考http://www.opengl-tutorial.org/intermediate-tutorials/tutorial-14-render-to-texture，以后有时间我会把这个翻译一下）
渲染背面的shader非常简单。
1 // for raycasting 2 #version 400 3 4 layout(location = 0) in vec3 position; 5 layout(location = 1) in vec3 color; 6 7 out vec3 passColor; 8 9 uniform mat4 MVP; 10 11 12 void main() 13 { 14 passColor = color; 15 //passColor = vec4(1, 1, 1, 1); 16 gl_Position = MVP * vec4(position, 1.0); 17 }
1 // for raycasting 2 #version 400 3 4 in vec3 passColor; 5 layout (location = 0) out vec4 FragColor; 6 7 8 void main() 9 { 10 FragColor = vec4(passColor, 1.0); 11 }
渲染时则需要注意启用新的framebuffer。这样就会渲染到指定的2D纹理上。
1 // render to texture 2 glBindFramebufferEXT(OpenGL.GL_FRAMEBUFFER_EXT, frameBuffer[0]); 3 OpenGL.Clear(OpenGL.GL_COLOR_BUFFER_BIT | OpenGL.GL_DEPTH_BUFFER_BIT); 4 this.backfaceRenderer.Render(arg); 5 glBindFramebufferEXT(OpenGL.GL_FRAMEBUFFER_EXT, 0);
注意，上图中左侧是立方体的前面，右侧是立方体的背面。
起始点/raycast
终点有了，起点也就知道怎么找了：渲染一个恰好包围3D纹理的立方体，且只渲染此立方体的前面（用glCullFace(GL_BACK);）。
起始点和终点都有了，就可以通过累加颜色来计算一个像素点的颜色值了。
1 #version 400 2 3 in vec3 EntryPoint; 4 in vec4 ExitPointCoord; 5 6 uniform sampler2D ExitPoints; 7 uniform sampler3D VolumeTex; 8 uniform sampler1D TransferFunc; 9 uniform float StepSize = 0.001f; 10 uniform vec2 ScreenSize; 11 uniform vec4 backgroundColor = vec4(0, 0, 0, 0);// value in glClearColor(value); 12 layout (location = 0) out vec4 FragColor; 13 14 void main() 15 { 16 // ExitPointCoord is normalized device coordinate 17 vec3 exitPoint = texture(ExitPoints, gl_FragCoord.st / ScreenSize).xyz; 18 // that will actually give you clip-space coordinates rather than 19 // normalised device coordinates, since you're not performing the perspective 20 // division which happens during the rasterisation process (between the vertex 21 // shader and fragment shader 22 // vec2 exitFragCoord = (ExitPointCoord.xy / ExitPointCoord.w + 1.0)/2.0; 23 // vec3 exitPoint = texture(ExitPoints, exitFragCoord).xyz; 24 25 //background need no raycasting 26 if (EntryPoint == exitPoint) { discard; } 27 28 vec3 direction = exitPoint - EntryPoint; 29 float directionLength = length(direction); // the length from front to back is calculated and used to terminate the ray 30 vec3 deltaDirection = direction * (StepSize / directionLength); 31 32 vec3 voxelCoord = EntryPoint; 33 vec3 colorAccumulator = vec3(0.0); // The dest color 34 float alphaAccumulator = 0.0f; 35 float lengthAccumulator = 0.0; 36 float intensity; 37 vec4 colorSample; // The src color 38 39 for(int i = 0; i < 1600; i++) 40 { 41 // get scaler value in the volume data 42 intensity = texture(VolumeTex, voxelCoord).x; 43 // get mapped color from 1-D texture 44 colorSample = texture(TransferFunc, intensity); 45 // modulate the value of colorSample.a 46 // front-to-back integration 47 if (colorSample.a > 0.0) { 48 // accomodate for variable sampling rates (base interval defined by mod_compositing.frag) 49 colorSample.a = 1.0 - pow(1.0 - colorSample.a, StepSize * 200.0f); 50 colorAccumulator += (1.0 - alphaAccumulator) * colorSample.rgb * colorSample.a; 51 alphaAccumulator += (1.0 - alphaAccumulator) * colorSample.a; 52 } 53 voxelCoord += deltaDirection; 54 lengthAccumulator += StepSize; 55 if (lengthAccumulator >= directionLength) 56 { 57 colorAccumulator = colorAccumulator * alphaAccumulator 58 + (1 - alphaAccumulator) * backgroundColor.rgb; 59 break; // terminate if opacity > 1 or the ray is outside the volume 60 } 61 else if (alphaAccumulator > 1.0) 62 { 63 alphaAccumulator = 1.0; 64 break; 65 } 66 } 67 FragColor = vec4(colorAccumulator, alphaAccumulator); 68 }
Raycast所需的vertex shader和backface.vert几乎一样。
1 #version 400 2 3 4 layout (location = 0) in vec3 position; 5 // have to use this variable!!!, or it will be very hard to debug for AMD video card 6 layout (location = 1) in vec3 color; 7 8 9 out vec3 EntryPoint; 10 out vec4 ExitPointCoord; 11 12 uniform mat4 MVP; 13 14 void main() 15 { 16 EntryPoint = color; 17 gl_Position = MVP * vec4(position,1.0); 18 ExitPointCoord = gl_Position; 19 }
而（在CSharpGL中）所需的渲染指令也只需一句话。
this.raycastRenderer.Render(arg);
Miscellaneous
在实现上述过程之前，必须初始化很多东西：3D纹理，附带了2D纹理的frameBuffer，用于渲染背面的shader和立方体模型，用于渲染正面/raycast的shader和立方体模型，从float类型的intensity值到vec3类型的颜色值的转换功能（1D纹理），设置uniform变量。
1 protected override void DoInitialize() 2 { 3 InitBackfaceRenderer(); 4 5 InitRaycastRenderer(); 6 7 initTFF1DTex(@"10RaycastVolumeRender ff.dat"); 8 int[] viewport = OpenGL.GetViewport(); 9 initFace2DTex(viewport[2], viewport[3]); 10 initVol3DTex(@"10RaycastVolumeRenderhead256.raw", 256, 256, 225); 11 initFrameBuffer(viewport[2], viewport[3]); 12 13 //this.depthTest = new DepthTestSwitch(); 14 15 RaycastingSetupUniforms(); 16 } 17 18 private void RaycastingSetupUniforms() 19 { 20 // setting uniforms such as 21 // ScreenSize 22 // StepSize 23 // TransferFunc 24 // ExitPoints i.e. the backface, the backface hold the ExitPoints of ray casting 25 // VolumeTex the texture that hold the volume data i.e. head256.raw 26 int[] viewport = OpenGL.GetViewport(); 27 this.raycastRenderer.SetUniform("ScreenSize", new vec2(viewport[2], viewport[3])); 28 this.raycastRenderer.SetUniform("StepSize", g_stepSize); 29 this.raycastRenderer.SetUniform("TransferFunc", new samplerValue(BindTextureTarget.Texture1D, transferFunc1DTexObj[0], OpenGL.GL_TEXTURE0)); 30 this.raycastRenderer.SetUniform("ExitPoints", new samplerValue(BindTextureTarget.Texture2D, backface2DTexObj[0], OpenGL.GL_TEXTURE1)); 31 this.raycastRenderer.SetUniform("VolumeTex", new samplerValue(BindTextureTarget.Texture3D, vol3DTexObj[0], OpenGL.GL_TEXTURE2)); 32 var clearColor = new float[4]; 33 OpenGL.GetFloat(GetTarget.ColorClearValue, clearColor); 34 this.raycastRenderer.SetUniform("backgroundColor", clearColor.ToVec4()); 35 } 36 37 private void initFrameBuffer(int texWidth, int texHeight) 38 { 39 // create a depth buffer for our framebuffer 40 var depthBuffer = new uint[1]; 41 OpenGL.GetDelegateFor<OpenGL.glGenRenderbuffersEXT>()(1, depthBuffer); 42 OpenGL.GetDelegateFor<OpenGL.glBindRenderbufferEXT>()(OpenGL.GL_RENDERBUFFER, depthBuffer[0]); 43 OpenGL.GetDelegateFor<OpenGL.glRenderbufferStorageEXT>()(OpenGL.GL_RENDERBUFFER, OpenGL.GL_DEPTH_COMPONENT, texWidth, texHeight); 44 45 // attach the texture and the depth buffer to the framebuffer 46 OpenGL.GetDelegateFor<OpenGL.glGenFramebuffersEXT>()(1, frameBuffer); 47 OpenGL.GetDelegateFor<OpenGL.glBindFramebufferEXT>()(OpenGL.GL_FRAMEBUFFER_EXT, frameBuffer[0]); 48 OpenGL.GetDelegateFor<OpenGL.glFramebufferTexture2DEXT>()(OpenGL.GL_FRAMEBUFFER_EXT, OpenGL.GL_COLOR_ATTACHMENT0_EXT, OpenGL.GL_TEXTURE_2D, backface2DTexObj[0], 0); 49 OpenGL.GetDelegateFor<OpenGL.glFramebufferRenderbufferEXT>()(OpenGL.GL_FRAMEBUFFER_EXT, OpenGL.GL_DEPTH_ATTACHMENT_EXT, OpenGL.GL_RENDERBUFFER, depthBuffer[0]); 50 checkFramebufferStatus(); 51 //OpenGL.Enable(GL_DEPTH_TEST); 52 } 53 54 private void checkFramebufferStatus() 55 { 56 uint complete = OpenGL.GetDelegateFor<OpenGL.glCheckFramebufferStatusEXT>()(OpenGL.GL_FRAMEBUFFER_EXT); 57 if (complete != OpenGL.GL_FRAMEBUFFER_COMPLETE_EXT) 58 { 59 throw new Exception("framebuffer is not complete"); 60 } 61 } 62 63 private void initVol3DTex(string filename, int width, int height, int depth) 64 { 65 var data = new UnmanagedArray<byte>(width * height * depth); 66 unsafe 67 { 68 int index = 0; 69 int readCount = 0; 70 byte* array = (byte*)data.Header.ToPointer(); 71 using (var fs = new FileStream(filename, FileMode.Open, FileAccess.Read)) 72 using (var br = new BinaryReader(fs)) 73 { 74 int unReadCount = (int)fs.Length; 75 const int cacheSize = 1024 * 1024; 76 do 77 { 78 int min = Math.Min(cacheSize, unReadCount); 79 var cache = new byte[min]; 80 readCount = br.Read(cache, 0, min); 81 if (readCount != min) 82 { throw new Exception(); } 83 84 for (int i = 0; i < readCount; i++) 85 { 86 array[index++] = cache[i]; 87 } 88 unReadCount -= readCount; 89 } while (readCount > 0); 90 } 91 } 92 93 OpenGL.GenTextures(1, vol3DTexObj); 94 // bind 3D texture target 95 OpenGL.BindTexture(OpenGL.GL_TEXTURE_3D, vol3DTexObj[0]); 96 OpenGL.TexParameteri(OpenGL.GL_TEXTURE_3D, OpenGL.GL_TEXTURE_MAG_FILTER, (int)OpenGL.GL_LINEAR); 97 OpenGL.TexParameteri(OpenGL.GL_TEXTURE_3D, OpenGL.GL_TEXTURE_MIN_FILTER, (int)OpenGL.GL_LINEAR); 98 OpenGL.TexParameteri(OpenGL.GL_TEXTURE_3D, OpenGL.GL_TEXTURE_WRAP_S, (int)OpenGL.GL_REPEAT); 99 OpenGL.TexParameteri(OpenGL.GL_TEXTURE_3D, OpenGL.GL_TEXTURE_WRAP_T, (int)OpenGL.GL_REPEAT); 100 OpenGL.TexParameteri(OpenGL.GL_TEXTURE_3D, OpenGL.GL_TEXTURE_WRAP_R, (int)OpenGL.GL_REPEAT); 101 // pixel transfer happens here from client to OpenGL server 102 OpenGL.PixelStorei(OpenGL.GL_UNPACK_ALIGNMENT, 1); 103 OpenGL.TexImage3D(OpenGL.GL_TEXTURE_3D, 0, (int)OpenGL.GL_INTENSITY, 104 width, height, depth, 0, 105 OpenGL.GL_LUMINANCE, OpenGL.GL_UNSIGNED_BYTE, data.Header); 106 data.Dispose(); 107 } 108 109 private void initFace2DTex(int width, int height) 110 { 111 OpenGL.GenTextures(1, backface2DTexObj); 112 OpenGL.BindTexture(OpenGL.GL_TEXTURE_2D, backface2DTexObj[0]); 113 OpenGL.TexParameteri(OpenGL.GL_TEXTURE_2D, OpenGL.GL_TEXTURE_WRAP_S, (int)OpenGL.GL_REPEAT); 114 OpenGL.TexParameteri(OpenGL.GL_TEXTURE_2D, OpenGL.GL_TEXTURE_WRAP_T, (int)OpenGL.GL_REPEAT); 115 OpenGL.TexParameteri(OpenGL.GL_TEXTURE_2D, OpenGL.GL_TEXTURE_MIN_FILTER, (int)OpenGL.GL_NEAREST); 116 OpenGL.TexParameteri(OpenGL.GL_TEXTURE_2D, OpenGL.GL_TEXTURE_MAG_FILTER, (int)OpenGL.GL_NEAREST); 117 OpenGL.TexImage2D(OpenGL.GL_TEXTURE_2D, 0, OpenGL.GL_RGBA16F, width, height, 0, OpenGL.GL_RGBA, OpenGL.GL_FLOAT, IntPtr.Zero); 118 } 119 120 private void initTFF1DTex(string filename) 121 { 122 // read in the user defined data of transfer function 123 byte[] tff; 124 using (var fs = new FileStream(filename, FileMode.Open, FileAccess.Read)) 125 using (var br = new BinaryReader(fs)) 126 { 127 tff = br.ReadBytes((int)fs.Length); 128 } 129 OpenGL.GenTextures(1, transferFunc1DTexObj); 130 OpenGL.BindTexture(OpenGL.GL_TEXTURE_1D, transferFunc1DTexObj[0]); 131 OpenGL.TexParameteri(OpenGL.GL_TEXTURE_1D, OpenGL.GL_TEXTURE_WRAP_S, (int)OpenGL.GL_REPEAT); 132 OpenGL.TexParameteri(OpenGL.GL_TEXTURE_1D, OpenGL.GL_TEXTURE_MIN_FILTER, (int)OpenGL.GL_NEAREST); 133 OpenGL.TexParameteri(OpenGL.GL_TEXTURE_1D, OpenGL.GL_TEXTURE_MAG_FILTER, (int)OpenGL.GL_NEAREST); 134 OpenGL.PixelStorei(OpenGL.GL_UNPACK_ALIGNMENT, 1); 135 OpenGL.TexImage1D(OpenGL.GL_TEXTURE_1D, 0, OpenGL.GL_RGBA8, 256, 0, OpenGL.GL_RGBA, OpenGL.GL_UNSIGNED_BYTE, tff); 136 } 137 138 private void InitRaycastRenderer() 139 { 140 var shaderCodes = new ShaderCode[2]; 141 shaderCodes[0] = new ShaderCode(File.ReadAllText(@"10RaycastVolumeRender aycasting.vert"), ShaderType.VertexShader); 142 shaderCodes[1] = new ShaderCode(File.ReadAllText(@"10RaycastVolumeRender aycasting.frag"), ShaderType.FragmentShader); 143 var map = new PropertyNameMap(); 144 map.Add("position", "position"); 145 map.Add("color", "color"); 146 this.raycastRenderer = new Renderer(model, shaderCodes, map); 147 this.raycastRenderer.Initialize(); 148 this.raycastRenderer.SwitchList.Add(new CullFaceSwitch(CullFaceMode.Back, true)); 149 } 150 151 private void InitBackfaceRenderer() 152 { 153 var shaderCodes = new ShaderCode[2]; 154 shaderCodes[0] = new ShaderCode(File.ReadAllText(@"10RaycastVolumeRenderackface.vert"), ShaderType.VertexShader); 155 shaderCodes[1] = new ShaderCode(File.ReadAllText(@"10RaycastVolumeRenderackface.frag"), ShaderType.FragmentShader); 156 var map = new PropertyNameMap(); 157 map.Add("position", "position"); 158 map.Add("color", "color"); 159 this.backfaceRenderer = new Renderer(model, shaderCodes, map); 160 this.backfaceRenderer.Initialize(); 161 this.backfaceRenderer.SwitchList.Add(new CullFaceSwitch(CullFaceMode.Front, true)); 162 }
2018-01-16
最近终于解决了在某些电脑上不显示VR的情况。原来是我没有把3个Texture分别绑定到不同的Unit上。
总结
当然，也可以先渲染出起始点，然后再找到终点的时候计算各个像素点的颜色值。
raycast做volume rendering的这个例子中，最耗空间的是3D纹理。但是这是无法避免的。其他空间和时间耗费都是极少的。
欢迎对OpenGL有兴趣的同学关注（https://github.com/bitzhuwei/CSharpGL）