From 04f15656c846ea38ba2c31d0bd8235cfb8c97db8 Mon Sep 17 00:00:00 2001 From: Douglas Dwyer Date: Fri, 25 Mar 2022 22:48:00 -0400 Subject: [PATCH] Add optimized isosurface rendering; add the ability to view cross-sectional slices of geometry --- mcxcloud/frontend/index.html | 280 +++++++++++++++++++++++++---------- 1 file changed, 200 insertions(+), 80 deletions(-) diff --git a/mcxcloud/frontend/index.html b/mcxcloud/frontend/index.html index fd359225..0cc0727e 100644 --- a/mcxcloud/frontend/index.html +++ b/mcxcloud/frontend/index.html @@ -351,8 +351,36 @@

Schema

Preview

- Colormap - +
+ Display Mode
+ +
+ + +
+
+ Colormap
+
+ +
+
+ Cross Sections
+
+ X + + +
+
+ Y + + +
+
+ Z + + +
+
Clip Plane @@ -1480,48 +1508,8 @@

Backend

} } - -var VolumeRenderShader1 = { - uniforms: { - 'u_size': { value: new THREE.Vector3( 1, 1, 1 ) }, - 'u_renderstyle': { value: 0 }, - 'u_renderthreshold': { value: 0.5 }, - 'u_clim': { value: new THREE.Vector2( 1, 1 ) }, - 'u_data': { value: null }, - 'u_cmdata': { value: null } - }, - vertexShader: [ - ' varying vec4 v_nearpos;', - ' varying vec4 v_farpos;', - ' varying vec3 v_position;', - - ' void main() {', - // Prepare transforms to map to "camera view". See also: - // https://threejs.org/docs/#api/renderers/webgl/WebGLProgram - ' mat4 viewtransformf = modelViewMatrix;', - ' mat4 viewtransformi = inverse(modelViewMatrix);', - - // Project local vertex coordinate to camera position. Then do a step - // backward (in cam coords) to the near clipping plane, and project back. Do - // the same for the far clipping plane. This gives us all the information we - // need to calculate the ray and truncate it to the viewing cone. - ' vec4 position4 = vec4(position, 1.0);', - ' vec4 pos_in_cam = viewtransformf * position4;', - - // Intersection of ray and near clipping plane (z = -1 in clip coords) - ' pos_in_cam.z = -pos_in_cam.w;', - ' v_nearpos = viewtransformi * pos_in_cam;', - - // Intersection of ray and far clipping plane (z = +1 in clip coords) - ' pos_in_cam.z = pos_in_cam.w;', - ' v_farpos = viewtransformi * pos_in_cam;', - - // Set varyings and output pos - ' v_position = position;', - ' gl_Position = projectionMatrix * viewMatrix * modelMatrix * position4;', - ' }', - ].join( '\n' ), - fragmentShader: [ +function createFragmentShader(mode) { + return [ ' precision highp float;', ' precision mediump sampler3D;', @@ -1532,6 +1520,8 @@

Backend

' uniform sampler3D u_data;', ' uniform sampler2D u_cmdata;', + ' uniform vec3 u_minslice;', + ' uniform vec3 u_maxslice;', ' varying vec3 v_position;', ' varying vec4 v_nearpos;', @@ -1561,14 +1551,13 @@

Backend

// Calculate unit vector pointing in the view direction through this fragment. ' vec3 view_ray = normalize(nearpos.xyz - farpos.xyz);', - // Compute the (negative) distance to the front surface or near clipping plane. // v_position is the back face of the cuboid, so the initial distance calculated in the dot // product below is the distance from near clip plane to the back of the cuboid ' float distance = dot(nearpos - v_position, view_ray);', ' vec3 cmp = (-vec3(0.5) - v_position) / view_ray;', - ' vec3 cmpu = u_size / view_ray + cmp;', + ' vec3 cmpu = (u_size * vec3(1.0, 1.0, 1.0)) / view_ray + cmp;', ' cmp = min(cmp, cmpu);', ' distance = max(distance, max(cmp.x, max(cmp.y, cmp.z)));', @@ -1582,11 +1571,32 @@

Backend

' vec3 step = ((v_position - front) / u_size) / float(nsteps);', ' vec3 start_loc = front / u_size;', + // Clip starting position and step count within the bounding box defined by u_minslice and u_maxslice + // allowing for cross-sectional views. + ` vec3 nstepminbound = mix(vec3(float(0)), (u_maxslice - start_loc) / step, greaterThan(start_loc, u_maxslice)); + float skips = max(nstepminbound.x, max(nstepminbound.y, nstepminbound.z)); + nstepminbound = mix(vec3(float(0)), (u_minslice - start_loc) / step, lessThan(start_loc, u_minslice)); + skips = max(skips, max(nstepminbound.x, max(nstepminbound.y, nstepminbound.z))); + start_loc += skips * step; + nsteps -= int(skips + 0.5);`, + + ` vec3 nstepmaxbound = mix(vec3(float(nsteps)), ceil((u_maxslice - start_loc) / step), greaterThan((start_loc + float(nsteps) * step), u_maxslice)); + nstepmaxbound = min(nstepmaxbound, mix(vec3(float(nsteps)), ceil((u_minslice - start_loc) / step), lessThan((start_loc + float(nsteps) * step), u_minslice))); + nsteps = int(min(nstepmaxbound.x, min(nstepmaxbound.y, nstepmaxbound.z)) + 0.5);`, + + // For testing: show the number of steps. This helps to establish // whether the rays are correctly oriented //'gl_FragColor = vec4(0.0, float(nsteps) / 1.0 / u_size.x, 1.0, 1.0);', //'return;', - ' cast_mip(start_loc, step, nsteps, view_ray);', + ` + if(` + mode + `) { + cast_mip(start_loc, step, nsteps, view_ray); + } + else { + cast_iso(start_loc, step, nsteps, view_ray); + } + `, ' }', @@ -1597,14 +1607,14 @@

Backend

' vec4 apply_colormap(float val) {', - ' val = (val - u_clim[0]) / (u_clim[1] - u_clim[0]);', + ' val = (val - u_clim.x) / (u_clim.y - u_clim.x);', ' return texture2D(u_cmdata, vec2(val, 0.5));', ' }', ' void cast_mip(vec3 start_loc, vec3 step, int nsteps, vec3 view_ray) {', - ' float max_val = -1e6;', + ' float max_val = 0.0;', ' vec3 loc = start_loc;', // Enter the raycasting loop. In WebGL 1 the loop index cannot be compared with @@ -1619,58 +1629,54 @@

Backend

' loc += step;', ' }', - // Refine location, gives crispier images - ' vec3 iloc = loc + step * (-0.5);', - ' vec3 istep = step / float(REFINEMENT_STEPS);', - ' for (int i=0; i= nsteps)', - ' break;', - + isWebGL2Available() ? 'for (int iter=0; iter<=nsteps; iter++) {' : 'for (int iter=0; iter<=MAX_STEPS; iter++) {\nif(iter >= nsteps) break;', // Sample from the 3D texture - ' float val = sample1(loc);', + ' val = sample1(loc);', - ' if (val > low_threshold) {', + ' if (val > u_renderthreshold) {', // Take the last interval in smaller steps - ' vec3 iloc = loc - 0.5 * step;', - ' vec3 istep = step / float(REFINEMENT_STEPS);', - ' for (int i=0; i u_renderthreshold) {', - ' gl_FragColor = add_lighting(val, iloc, dstep, view_ray);', - ' return;', - ' }', - ' iloc += istep;', - ' }', + ' break;', ' }', // Advance location deeper into the volume ' loc += step;', ' }', + + ' if (val > u_renderthreshold) {', + // Take the last interval in smaller steps + ' vec3 iloc = loc - 0.5 * step;', + ' vec3 istep = step / float(REFINEMENT_STEPS);', + ' for (int i=0; i u_renderthreshold) {', + ' gl_FragColor = add_lighting(val, iloc, dstep, view_ray);', + ' return;', + ' }', + ' iloc += istep;', + ' }', + ' gl_FragColor = add_lighting(val, iloc, dstep, view_ray);', + ' }', + ' else {', + ' gl_FragColor = vec4(0.0);', + ' }', ' }', - ' vec4 add_lighting(float val, vec3 loc, vec3 step, vec3 view_ray)', ' {', // Calculate color by incorporating lighting @@ -1735,7 +1741,58 @@

Backend

' final_color.a = color.a;', ' return final_color;', ' }', - ].join( '\n' ) + ].join( '\n' ); +} + +var MipRenderShader = { + uniforms: { + 'u_size': { value: new THREE.Vector3( 1, 1, 1 ) }, + 'u_renderstyle': { value: 0 }, + 'u_renderthreshold': { value: 0.5 }, + 'u_clim': { value: new THREE.Vector2( 1, 1 ) }, + 'u_data': { value: null }, + 'u_cmdata': { value: null }, + 'u_minslice': { value: new THREE.Vector3( 0, 0, 0 ) }, + 'u_maxslice': { value: new THREE.Vector3( 1, 1, 1 ) } + }, + vertexShader: [ + ' varying vec4 v_nearpos;', + ' varying vec4 v_farpos;', + ' varying vec3 v_position;', + + ' void main() {', + // Prepare transforms to map to "camera view". See also: + // https://threejs.org/docs/#api/renderers/webgl/WebGLProgram + ' mat4 viewtransformf = modelViewMatrix;', + ' mat4 viewtransformi = inverse(modelViewMatrix);', + + // Project local vertex coordinate to camera position. Then do a step + // backward (in cam coords) to the near clipping plane, and project back. Do + // the same for the far clipping plane. This gives us all the information we + // need to calculate the ray and truncate it to the viewing cone. + ' vec4 position4 = vec4(position, 1.0);', + ' vec4 pos_in_cam = viewtransformf * position4;', + + // Intersection of ray and near clipping plane (z = -1 in clip coords) + ' pos_in_cam.z = -pos_in_cam.w;', + ' v_nearpos = viewtransformi * pos_in_cam;', + + // Intersection of ray and far clipping plane (z = +1 in clip coords) + ' pos_in_cam.z = pos_in_cam.w;', + ' v_farpos = viewtransformi * pos_in_cam;', + + // Set varyings and output pos + ' v_position = position;', + ' gl_Position = projectionMatrix * viewMatrix * modelMatrix * position4;', + ' }', + ].join( '\n' ), + fragmentShader: createFragmentShader(true) +}; + +var IsoRenderShader = { + uniforms: MipRenderShader.uniforms, + vertexShader: MipRenderShader.vertexShader, + fragmentShader: createFragmentShader(false) }; // parse url -> merge options -> refreshUI() -> initJsoneditor() -> direct link @@ -2612,16 +2669,18 @@

Backend

}; // Material - const shader = VolumeRenderShader1; + const shader = document.getElementById('mip-radio-button').checked ? MipRenderShader : IsoRenderShader; const uniforms = THREE.UniformsUtils.clone( shader.uniforms ); uniforms[ "u_data" ].value = texture; uniforms[ "u_size" ].value.set( dim[0], dim[1], dim[2] ); uniforms[ "u_clim" ].value.set( volume.min(), volume.max() ); - uniforms[ "u_renderstyle" ].value = 0; + uniforms[ "u_renderstyle" ].value = document.getElementById('mip-radio-button').checked ? 0 : 1; uniforms[ "u_renderthreshold" ].value = 0.2; uniforms[ "u_cmdata" ].value = cmtextures[ "viridis" ]; + uniforms[ "u_minslice" ].value.set( parseFloat($("#cross-x-low").val()), parseFloat($("#cross-y-low").val()), parseFloat($("#cross-z-low").val()) ); + uniforms[ "u_maxslice" ].value.set( parseFloat($("#cross-x-hi").val()), parseFloat($("#cross-y-hi").val()), parseFloat($("#cross-z-hi").val()) ); lastclim=uniforms[ "u_clim" ].value; $("#clim-low").prop( "disabled", false ); @@ -2843,6 +2902,14 @@

Backend

storage.userinfo = JSON.stringify({ fullname: $("#fullname").val(), inst: $("#inst").val(), email: $("#email").val(), netname: $("#netname").val()}); } +function setcrosssectionsizes() { + if(lastvolume !== null){ + lastvolume.material.uniforms[ "u_minslice" ].value.set( parseFloat($("#cross-x-low").val()), parseFloat($("#cross-y-low").val()), parseFloat($("#cross-z-low").val()) ); + lastvolume.material.uniforms[ "u_maxslice" ].value.set( parseFloat($("#cross-x-hi").val()), parseFloat($("#cross-y-hi").val()), parseFloat($("#cross-z-hi").val()) ); + renderer.updateComplete = false; + } +} + var canvas = $("#canvas"); var bbxsize=[0,0,0]; @@ -2976,6 +3043,59 @@

Backend

$("#credits").show(); }); +$("#mip-radio-button").on('change', function() { + if(lastvolume !== null){ + const unfs = lastvolume.material.uniforms; + lastvolume.material = new THREE.ShaderMaterial( { + uniforms: THREE.UniformsUtils.clone( MipRenderShader.uniforms ), + vertexShader: MipRenderShader.vertexShader, + fragmentShader: MipRenderShader.fragmentShader, + side: THREE.BackSide + } ); + lastvolume.material.uniforms = unfs; + + renderer.updateComplete = false; + } +}); + +$("#iso-radio-button").on('change', function() { + if(lastvolume !== null){ + const unfs = lastvolume.material.uniforms; + lastvolume.material = new THREE.ShaderMaterial( { + uniforms: THREE.UniformsUtils.clone( IsoRenderShader.uniforms ), + vertexShader: IsoRenderShader.vertexShader, + fragmentShader: IsoRenderShader.fragmentShader, + side: THREE.BackSide + } ); + lastvolume.material.uniforms = unfs; + renderer.updateComplete = false; + } +}); + +$("#cross-x-low").on('input', function() { + setcrosssectionsizes(); +}); + +$("#cross-y-low").on('input', function() { + setcrosssectionsizes(); +}); + +$("#cross-z-low").on('input', function() { + setcrosssectionsizes(); +}); + +$("#cross-x-hi").on('input', function() { + setcrosssectionsizes(); +}); + +$("#cross-y-hi").on('input', function() { + setcrosssectionsizes(); +}); + +$("#cross-z-hi").on('input', function() { + setcrosssectionsizes(); +}); + /* -------------------------------------------------------- main function */ parseUrl()