feat: Volume rendering with OIT

packages/core/examples/volume_rendering/main.ts

@@ -30,17 +30,16 @@ const volumeLayer = new VolumeLayer({
     },
     {
       visible: true,
-      color: [1, 1, 1],
+      color: [0, 0, 1],
       contrastLimits: [108, 353],
     },
     {
       visible: true,
-      color: [1, 1, 1],
+      color: [0, 1, 0],
       contrastLimits: [144, 3825],
     },
   ],
 });
-volumeLayer.opacityMultiplier = 0.001;
 const idetik = new Idetik({
   canvas: document.querySelector<HTMLCanvasElement>("#canvas")!,
   viewports: [
@@ -81,21 +80,21 @@ const opacityControls = {
     return (Math.log10(volumeLayer.opacityMultiplier) + 3) / 4;
   },
   set opacity(sliderValue: number) {
-    volumeLayer.opacityMultiplier = Math.pow(10, sliderValue * 4 - 3);
+    volumeLayer.opacityMultiplier = 10 ** (sliderValue * 4 - 3);
   },
 };
 volumeFolder
   .add(opacityControls, "opacity", 0, 1, 0.01)
   .name("Opacity")
   .decimals(2);
 volumeFolder
-  .add(volumeLayer, "earlyTerminationAlpha", 0.8, 1.0, 0.01)
+  .add(volumeLayer, "earlyTerminationAlpha", 0.01, 1.0, 0.01)
   .name("Early termination threshold");
 
 const overlaysFolder = gui.addFolder("Debug");
 overlaysFolder
   .add(volumeLayer, "debugShowWireframes")
   .name("Show tile wireframes");
 overlaysFolder
-  .add(volumeLayer, "debugShowDegenerateRays")
-  .name("Show degenerate rays (length 0)");
+  .add(volumeLayer, "debugSetOITWeightOne")
+  .name("Set OIT weight to 1 (no depth)");

packages/core/src/core/renderable_object.ts

@@ -33,6 +33,10 @@ export abstract class RenderableObject extends Node {
     return stale;
   }
 
+  public get visible() {
+    return true;
+  }
+
   public get geometry() {
     return this.geometry_;
   }

packages/core/src/layers/volume_layer.ts

@@ -1,14 +1,20 @@
-import { Chunk, ChunkSource, SliceCoordinates } from "../data/chunk";
-import { Layer, RenderContext } from "../core/layer";
+import type { Chunk, ChunkSource, SliceCoordinates } from "../data/chunk";
+import { Layer, type RenderContext } from "../core/layer";
 import { VolumeRenderable } from "../objects/renderable/volume_renderable";
-import { IdetikContext } from "../idetik";
-import { ChunkStoreView, INTERNAL_POLICY_KEY } from "../core/chunk_store_view";
-import { ImageSourcePolicy } from "../core/image_source_policy";
+import type { IdetikContext } from "../idetik";
+import {
+  type ChunkStoreView,
+  INTERNAL_POLICY_KEY,
+} from "../core/chunk_store_view";
+import type { ImageSourcePolicy } from "../core/image_source_policy";
 import { Texture3D } from "../objects/textures/texture_3d";
 import { RenderablePool } from "../utilities/renderable_pool";
 import { vec3 } from "gl-matrix";
 import { sortFrontToBack } from "../math/sort_by_distance";
-import { ChannelProps, ChannelsEnabled } from "../objects/textures/channel";
+import type {
+  ChannelProps,
+  ChannelsEnabled,
+} from "../objects/textures/channel";
 
 export type VolumeLayerProps = {
   source: ChunkSource;
@@ -38,9 +44,9 @@ export class VolumeLayer extends Layer implements ChannelsEnabled {
 
   // TODO: Make a debug config object to manage debug options
   private debugShowWireframes_ = false;
-  public debugShowDegenerateRays = false;
+  public debugSetOITWeightOne = false;
   public relativeStepSize = 1.0;
-  public opacityMultiplier = 0.1;
+  public opacityMultiplier = 0.5;
   public earlyTerminationAlpha = 0.99;
 
   public get debugShowWireframes() {
@@ -102,6 +108,7 @@ export class VolumeLayer extends Layer implements ChannelsEnabled {
   private createVolume(chunk: Chunk) {
     const volume = new VolumeRenderable(
       Texture3D.createWithChunk(chunk),
+      chunk.chunkIndex.c,
       this.channelProps_
     );
     this.updateVolumeChunk(volume, chunk);
@@ -238,7 +245,7 @@ export class VolumeLayer extends Layer implements ChannelsEnabled {
 
   public getUniforms(): Record<string, unknown> {
     return {
-      DebugShowDegenerateRays: Number(this.debugShowDegenerateRays),
+      DebugSetOITWeightOne: Number(this.debugSetOITWeightOne),
       RelativeStepSize: this.relativeStepSize * this.interactiveStepSizeScale_,
       OpacityMultiplier: this.opacityMultiplier,
       EarlyTerminationAlpha: this.earlyTerminationAlpha,

packages/core/src/objects/renderable/volume_renderable.ts

@@ -5,40 +5,51 @@ import type { TextureDataType } from "../textures/texture";
 import type { Texture3D } from "../textures/texture_3d";
 import { vec3 } from "gl-matrix";
 import {
-  Channel,
-  ChannelProps,
+  type Channel,
+  type ChannelProps,
   validateChannel,
   validateChannels,
 } from "../textures/channel";
 
 export class VolumeRenderable extends RenderableObject {
   public voxelScale: vec3 = vec3.fromValues(1, 1, 1);
+  public channelIndex: number = 0;
 
   private channels_: Required<Channel>[];
 
-  constructor(texture: Texture3D, channels: ChannelProps[] = []) {
+  constructor(
+    texture: Texture3D,
+    channelIndex: number = 0,
+    channels: ChannelProps[] = []
+  ) {
     super();
     this.geometry = new BoxGeometry(1, 1, 1, 1, 1, 1);
     this.setTexture(0, texture);
     this.programName = dataTypeToVolumeShader(texture.dataType);
     this.cullFaceMode = "front";
     this.depthTest = false;
-    // TODO handle visibility property of channels
+    this.channelIndex = channelIndex;
     this.channels_ = validateChannels(texture, channels);
   }
 
+  public get visible() {
+    const channel = this.getChannelOrDefault(this.channelIndex);
+    return channel.visible;
+  }
+
   public get type() {
     return "VolumeRenderable";
   }
 
   public override getUniforms(): Record<string, unknown> {
-    const channel = this.channels_[0] ?? validateChannel(this.textures[0], {});
+    const channel = this.getChannelOrDefault(this.channelIndex);
     const { color, contrastLimits } = channel;
     return {
       VoxelScale: this.voxelScale,
       Color: color.rgb,
       ValueOffset: -contrastLimits[0],
       ValueScale: 1 / (contrastLimits[1] - contrastLimits[0]),
+      DebugShowChunkBoundaries: Number(this.wireframeEnabled),
     };
   }
 
@@ -58,6 +69,12 @@ export class VolumeRenderable extends RenderableObject {
 
     this.channels_[channelIndex] = newChannel;
   }
+
+  private getChannelOrDefault(channelIndex: number): Required<Channel> {
+    return (
+      this.channels_[channelIndex] ?? validateChannel(this.textures[0], {})
+    );
+  }
 }
 
 function dataTypeToVolumeShader(dataType: TextureDataType): Shader {

packages/core/src/renderers/shaders/index.ts

@@ -10,6 +10,8 @@ import wireframeFragmentShader from "./wireframe_frag.glsl";
 import volumeVertexShader from "./volume_vert.glsl";
 import volumeFragmentShader from "./volume_frag.glsl";
 import labelImage from "./label_image_frag.glsl";
+import transparentCompositeVertexShader from "./transparent_composite_vert.glsl";
+import transparentCompositeFragmentShader from "./transparent_composite_frag.glsl";
 
 export type Shader =
   | "floatScalarImage"
@@ -24,6 +26,7 @@ export type Shader =
   | "uintScalarImage"
   | "uintScalarImageArray"
   | "uintVolume"
+  | "transparentComposite"
   | "wireframe";
 
 type ShaderCode = {
@@ -92,4 +95,8 @@ export const shaderCode: Record<Shader, ShaderCode> = {
     fragment: volumeFragmentShader,
     fragmentDefines: new Map([["TEXTURE_DATA_TYPE_UINT", "1"]]),
   },
+  transparentComposite: {
+    vertex: transparentCompositeVertexShader,
+    fragment: transparentCompositeFragmentShader,
+  },
 };

packages/core/src/renderers/shaders/transparent_composite_frag.glsl

@@ -0,0 +1,16 @@
+#version 300 es
+
+precision highp float;
+uniform sampler2D Sampler0;
+uniform sampler2D Sampler1;
+
+in vec2 TexCoords;
+layout (location = 0) out vec4 fragColor;
+
+void main() {
+    vec4 tex0 = texture(Sampler0, TexCoords);
+    vec4 tex1 = texture(Sampler1, TexCoords);
+    vec4 accum = vec4(tex0.rgb, tex1.r);
+    float revealage = tex0.a;
+    fragColor = vec4(accum.rgb / clamp(accum.a, 1e-4, 5e4), revealage);
+}

packages/core/src/renderers/shaders/transparent_composite_vert.glsl

@@ -0,0 +1,16 @@
+#version 300 es
+
+out vec2 TexCoords;
+
+void main() {
+    vec2 pos;
+    if (gl_VertexID == 0) {
+        pos = vec2(-1.0, -1.0);
+    } else if (gl_VertexID == 1) {
+        pos = vec2(3.0, -1.0);
+    } else {
+        pos = vec2(-1.0, 3.0);
+    }
+    TexCoords = pos * 0.5 + 0.5;
+    gl_Position = vec4(pos, 0.0, 1.0);
+}

packages/core/src/renderers/shaders/volume_frag.glsl

@@ -3,7 +3,8 @@
 
 precision highp float;
 
-layout (location = 0) out vec4 fragColor;
+layout(location = 0) out vec4 FragData0;
+layout(location = 1) out float FragData1;
 
 #if defined TEXTURE_DATA_TYPE_INT
 uniform mediump isampler3D ImageSampler;
@@ -14,21 +15,30 @@ uniform mediump sampler3D ImageSampler;
 #endif
 
 uniform highp vec3 CameraPositionModel;
+uniform mat4 Projection;
+uniform mat4 ModelView;
+uniform mat4 ModelViewProjection;
 in highp vec3 PositionModel;
 
 // The bounding box in model space is normalized to -0.5 to 0.5
 vec3 boundingboxMin = vec3(-0.50);
 vec3 boundingboxMax = vec3(0.50);
 
 // Volume rendering parameters
-uniform bool DebugShowDegenerateRays;
 uniform float OpacityMultiplier;
 uniform float EarlyTerminationAlpha;
 uniform float RelativeStepSize;
 uniform float ValueScale;
 uniform float ValueOffset;
 uniform vec3 Color;
 uniform vec3 VoxelScale;
+uniform bool DebugShowChunkBoundaries;
+uniform bool DebugSetOITWeightOne;
+
+float computeOITWeight(float alpha, float depth) {
+    float d = (1.0 - depth);
+    return alpha * max(1e-2, 3e3 * d * d * d);
+}
 
 vec2 findBoxIntersectionsAlongRay(vec3 rayOrigin, vec3 rayDir, vec3 boxMin, vec3 boxMax) {
     vec3 reciprocalRayDir = 1.0 / rayDir;
@@ -47,6 +57,25 @@ vec2 findBoxIntersectionsAlongRay(vec3 rayOrigin, vec3 rayDir, vec3 boxMin, vec3
 }
 
 void main() {
+    // Initialize outputs to defaults for early exit cases
+    FragData0 = vec4(0.0, 0.0, 0.0, 0.0);
+    FragData1 = 0.0;
+
+    if (DebugShowChunkBoundaries) {
+        vec3 distToMin = abs(PositionModel - boundingboxMin);
+        vec3 distToMax = abs(PositionModel - boundingboxMax);
+        bvec3 nearMin = lessThan(distToMin, vec3(0.01));
+        bvec3 nearMax = lessThan(distToMax, vec3(0.01));
+        bvec3 nearBoundary = bvec3(nearMin.x || nearMax.x, nearMin.y || nearMax.y, nearMin.z || nearMax.z);
+        int numDimsOnBoundary = int(nearBoundary.x) + int(nearBoundary.y) + int(nearBoundary.z);
+        bool onEdge = numDimsOnBoundary >= 2;
+        if (onEdge) {
+            FragData0 = vec4(1.0, 1.0, 1.0, 0.4);
+            FragData1 = 1.0;
+        }
+        return;
+    }
+
     // Step 1 - calculate where the ray enters and exits the volume
 
     // The ray in model space goes from the camera to the point on the back face
@@ -56,11 +85,6 @@ void main() {
     float tEnter = rayIntersections.x;
     float tExit = rayIntersections.y;
 
-    if (DebugShowDegenerateRays && (tExit == tEnter)) {
-        fragColor = vec4(1.0, 0.0, 0.0, 1.0);
-        return;
-    }
-
     vec3 entryPoint = CameraPositionModel + RayDirModel * tEnter;
     entryPoint = clamp(entryPoint + 0.5, 0.0, 1.0);
     vec3 exitPoint = CameraPositionModel + RayDirModel * tExit;
@@ -88,25 +112,38 @@ void main() {
 
     // Step 3 - perform the ray marching and compositing in front to back order
     vec3 position = entryPoint;
+    vec4 clipPosition;
     vec4 accumulatedColor = vec4(0.0);
-    float sampledData, sampleAlpha, blendedSampleAlpha;
+    vec3 sampleColor;
+    float revealage = 1.0;
+    float sampledData, sampleAlpha, blendedSampleAlpha, rayDepth, weight, scaledData;
 
     // Scale intensity to opacity.
-    // OpacityMultiplier and MaxIntensity control the transfer function.
+    // ValueScale is the
     // worldSpaceStepSize corrects for anisotropic voxels.
-    float intensityScale = OpacityMultiplier * worldSpaceStepSize * ValueScale;
+    float intensityScale = worldSpaceStepSize * ValueScale;
 
     // March until we reach the number of samples or accumulate enough opacity
-    for (int i = 0; i < numSamples && accumulatedColor.a < EarlyTerminationAlpha; i++) {
+    for (int i = 0; i < numSamples && revealage > (1.0 - EarlyTerminationAlpha); i++) {
+        // Sample the volume data and convert to color and opacity
         sampledData = vec4(texture(ImageSampler, position)).r;
-        sampleAlpha = (sampledData + ValueOffset) * intensityScale;
-        blendedSampleAlpha = (1.0 - accumulatedColor.a) * sampleAlpha;
+        scaledData = (sampledData + ValueOffset) * intensityScale;
+        sampleAlpha = clamp(scaledData * OpacityMultiplier, 0.0, 1.0);
+        sampleColor = Color * sampleAlpha;
+
+        // Weighted blended OIT
+        clipPosition = ModelViewProjection * vec4(position, 1.0);
+        rayDepth = (clipPosition.z / clipPosition.w) * 0.5 + 0.5;
+        weight = computeOITWeight(sampleAlpha, rayDepth);
+        if (DebugSetOITWeightOne) {
+            weight = 1.0;
+        }
+        accumulatedColor += vec4(sampleColor, sampleAlpha) * weight;
+        revealage *= clamp(1.0 - sampleAlpha, 0.0, 1.0);
 
-        // Front-to-back compositing
-        accumulatedColor.a += blendedSampleAlpha;
-        accumulatedColor.rgb += Color * blendedSampleAlpha;
         position += stepIncrement;
     }
 
-    fragColor = accumulatedColor;
+    FragData0 = vec4(accumulatedColor.rgb, 1.0 - revealage);
+    FragData1 = accumulatedColor.a;
 }