WebGL2.0でボリュームレンダリング(TypeScript)

下記記事にて、Direct3D11を用いてボリュームレンダリングを実装してみました。
Direct3D11でボリュームレンダリング - 何でもプログラミング

今回はWebGL2.0を利用して同様のものを実装してみたいと思います。

ボリュームレンダリングの手法やデータは上記記事を参照してみてください。

※現時点ではWebGL2.0はChromeFirefoxでしかサポートされていません。

HTML

表示のため、下記のhtmlを用意しました。

canvasWebGLレンタリング用で、file inputはデータを読み込むためのものです。

最終的に下図のような表示になります。
f:id:any-programming:20171226133944p:plain

<!DOCTYPE html>
<html lang="ja">
  <head>
    <meta charset="utf-8">
    <title>webGL</title>
  </head>
  <body>
    <canvas id="gl" width="200" height="200"></canvas>
    <input id="file" type="file">
    <script src="script.js"></script>
  </body>
</html>


Canvasのサイズ

htmlのcanvas要素のwidthとheightはバックバッファのサイズ、cssで設定するwidthとheightは実際の表示サイズとして利用されます。

WebGL2.d.ts

WebGL2の定義ファイルはデフォルトで含まれていないため、下記のものをダウンロードして利用しました。
GitHub - MaxGraey/WebGL2-TypeScript: WebGL2 bindings for TypeScript

いくつかWebGLRenderingContextの定義とかち合う関数があるので、適宜コメントアウトして利用してください。

main関数

Direct3Dの時と同様に、三方向のスライスセットを準備し適宜切り替えています。

また陰影もボリューム値の勾配を法線として算出しています。

function main() : void {
    // context取得
    const canvas = <HTMLCanvasElement>document.getElementById("gl");
    const gl = canvas.getContext("webgl2")!;

    // Shader作成
    const vertexCode = `#version 300 es
        layout(location = 0) in vec4 position;
        layout(location = 1) in vec3 texcoord;
        uniform mat4 matrix;
        out vec3 v_texcoord;
        void main() {
            gl_Position = matrix * position;
            v_texcoord = texcoord;
        }    
    `;

    const fragmentCode = `#version 300 es
        precision mediump float;
        precision mediump sampler3D;
        in vec3 v_texcoord;
        uniform sampler3D tex;
        out vec4 fragColor;
        void main() {
            if (texture(tex, v_texcoord).r < 0.5)
                discard;
            
            float dx = textureOffset(tex, v_texcoord, ivec3(1, 0, 0)).r - textureOffset(tex, v_texcoord, ivec3(-1,  0,  0)).r;
            float dy = textureOffset(tex, v_texcoord, ivec3(0, 1, 0)).r - textureOffset(tex, v_texcoord, ivec3( 0, -1,  0)).r;
            float dz = textureOffset(tex, v_texcoord, ivec3(0, 0, 1)).r - textureOffset(tex, v_texcoord, ivec3( 0,  0, -1)).r;
            vec3 normal = normalize(vec3(dx, dy, dz));

            vec3 light = normalize(vec3(1, 1, 1));
            float gray = abs(dot(normal, light));
            fragColor = vec4(gray, gray, gray, 1);
        }
    `;

    const program = createProgram(gl, vertexCode, fragmentCode);

    // Texture3D、VertexBuffer、IndexBuffer作成(Fileが読み込まれた際)
    let texture : WebGLTexture | undefined;
    let slices : VolumeSlices | undefined;
    const fileInput = <HTMLInputElement>document.getElementById("file");
    fileInput.onchange = ev =>
    {
        const file = (<HTMLInputElement>ev.target).files![0];
        const reader = new FileReader();
        reader.onload = _ev =>
        {
            const buf = new Uint8Array(<ArrayBuffer>reader.result);
            texture = createTexture3D(gl, 512, 512, 360, buf);
            slices = new VolumeSlices(gl, 512, 512, 360);
        };
        reader.readAsArrayBuffer(file);
    };

    // 固定の設定
    gl.useProgram(program);
    gl.enableVertexAttribArray(0);
    gl.enableVertexAttribArray(1);
    const matrixUniform = gl.getUniformLocation(program, "matrix")!;
    gl.enable(gl.DEPTH_TEST);

    // 描画関数
    function draw(rotX : number, rotY : number) : void {
        gl.clearColor(0.5, 0.6, 1.0, 1);
        gl.clear(gl.COLOR_BUFFER_BIT | gl.DEPTH_BUFFER_BIT);

        if (slices == undefined || texture == undefined)
            return;

        // 行列設定
        const mat = Matrix4.rotate(rotY, new Vector3(0, 1, 0))
            .mul(Matrix4.rotate(rotX, new Vector3(1, 0, 0)))
            .mul(Matrix4.scale(1.5 / 512))
            .mul(Matrix4.translate(new Vector3(-512.0 / 2, -512.0 / 2, -360.0 / 2)));
        gl.uniformMatrix4fv(matrixUniform, false, mat.transpose().data);

        // VertexBuffer、IndexBuffer設定
        const { vertexBuffer, indexBuffer } = slices.getBuffer(mat);
        gl.bindBuffer(gl.ARRAY_BUFFER, vertexBuffer);
        gl.vertexAttribPointer(0, 3, gl.FLOAT, false, 24, 0);
        gl.vertexAttribPointer(1, 3, gl.FLOAT, false, 24, 12);
        gl.bindBuffer(gl.ELEMENT_ARRAY_BUFFER, indexBuffer);

        // 描画
        const indexCount = gl.getBufferParameter(gl.ELEMENT_ARRAY_BUFFER, gl.BUFFER_SIZE) / 2;
        gl.drawElements(gl.TRIANGLES, indexCount, gl.UNSIGNED_SHORT, 0);
    }

    // MouseMoveで描画
    canvas.onmousemove = ev => {
        const rotX = 2 * Math.PI * ev.y / canvas.height;
        const rotY = 2 * Math.PI * ev.x / canvas.width;
        draw(rotX, rotY);
    };
}


Slices

三方向のスライスセットを作成、適切なものを取得できるクラスになります。

class VolumeSlices {
    private vertexBufferX : WebGLBuffer;
    private indexBufferX  : WebGLBuffer;
    private vertexBufferY : WebGLBuffer;
    private indexBufferY  : WebGLBuffer;
    private vertexBufferZ : WebGLBuffer;
    private indexBufferZ  : WebGLBuffer;
    constructor(gl : WebGLRenderingContext, width : number, height : number, depth : number) {
        const verticesZ = normalizedRange(depth, 0, 1)
            .concatMap(t => [
                0,     0,      t * depth, 0, 0, t,
                width, 0,      t * depth, 1, 0, t,
                0,     height, t * depth, 0, 1, t,
                width, height, t * depth, 1, 1, t,
            ]);

        const verticesY = normalizedRange(height, 0, 1)
            .concatMap(t => [
                0,     t * height, 0,     0, t, 0,
                width, t * height, 0,     1, t, 0,
                0,     t * height, depth, 0, t, 1,
                width, t * height, depth, 1, t, 1,
            ]);

        const verticesX = normalizedRange(width, 0, 1)
            .concatMap(t => [
                t * width, 0,      0,     t, 0, 0,
                t * width, height, 0,     t, 1, 0,
                t * width, 0,      depth, t, 0, 1,
                t * width, height, depth, t, 1, 1,
            ]);

        const createIndices = (count : number) =>
            range(count).concatMap<number>(i => [0, 2, 1, 1, 2, 3].map<number>(x => x + 4 * i));

        this.vertexBufferX = createVertexBuffer(gl, new Float32Array(verticesX));
        this.vertexBufferY = createVertexBuffer(gl, new Float32Array(verticesY));
        this.vertexBufferZ = createVertexBuffer(gl, new Float32Array(verticesZ));
        this.indexBufferX = createIndexBuffer(gl, new Uint16Array(createIndices(width)));
        this.indexBufferY = createIndexBuffer(gl, new Uint16Array(createIndices(height)));
        this.indexBufferZ = createIndexBuffer(gl, new Uint16Array(createIndices(depth)));
    }
    getBuffer(m : Matrix4) : { vertexBuffer : WebGLBuffer, indexBuffer : WebGLBuffer } {
        const dotX = Math.abs(m.transform3x3(Vector3.ex).z);
        const dotY = Math.abs(m.transform3x3(Vector3.ey).z);
        const dotZ = Math.abs(m.transform3x3(Vector3.ez).z);
        if (dotX < dotZ && dotY < dotZ)
            return { vertexBuffer : this.vertexBufferZ, indexBuffer : this.indexBufferZ };
        else if (dotX < dotY)
            return { vertexBuffer : this.vertexBufferY, indexBuffer : this.indexBufferY };
        else
            return { vertexBuffer : this.vertexBufferX, indexBuffer : this.indexBufferX };
    }
}


WebGL関数群

main関数内で利用されているWebGL関連の関数は下記のように実装しています。

function compileShader(gl : WebGLRenderingContext, shader : WebGLShader, code : string) : void {
    gl.shaderSource(shader, code);
    gl.compileShader(shader);

    console.log(gl.getShaderInfoLog(shader));
    if (!gl.getShaderParameter(shader, gl.COMPILE_STATUS))
        throw new Error("compile error");
}

function createProgram(gl : WebGLRenderingContext, vertexCode : string, fragmentCode : string) : WebGLProgram {
    const vertexShader = gl.createShader(gl.VERTEX_SHADER)!;
    compileShader(gl, vertexShader, vertexCode);

    const fragmentShader = gl.createShader(gl.FRAGMENT_SHADER)!;
    compileShader(gl, fragmentShader, fragmentCode);

    const program = gl.createProgram()!;
    gl.attachShader(program, vertexShader);
    gl.attachShader(program, fragmentShader);
    gl.linkProgram(program);

    console.log(gl.getProgramInfoLog(program));
    if (!gl.getProgramParameter(program, gl.LINK_STATUS))
        throw new Error("program error");

    return program;
}
function createVertexBuffer(gl : WebGLRenderingContext, vertices : Float32Array) : WebGLBuffer {
    const buffer = gl.createBuffer()!;

    gl.bindBuffer(gl.ARRAY_BUFFER, buffer);
    gl.bufferData(gl.ARRAY_BUFFER, vertices, gl.STATIC_DRAW);

    return buffer;
}
function createIndexBuffer(gl : WebGLRenderingContext, indices : Uint16Array) : WebGLBuffer {
    const buffer = gl.createBuffer()!;

    gl.bindBuffer(gl.ELEMENT_ARRAY_BUFFER, buffer);
    gl.bufferData(gl.ELEMENT_ARRAY_BUFFER, indices, gl.STATIC_DRAW);

    return buffer;
}
function createTexture3D(gl : WebGL2RenderingContext, width : number, height : number, depth : number, voxels : Uint8Array) : WebGLTexture {
    const texture = gl.createTexture()!;

    gl.bindTexture(gl.TEXTURE_3D, texture);
    gl.pixelStorei(gl.UNPACK_ALIGNMENT, 1);
    gl.texImage3D(gl.TEXTURE_3D, 0, gl.R8, width, height, depth, 0, gl.RED, gl.UNSIGNED_BYTE, voxels);
    gl.texParameteri(gl.TEXTURE_3D, gl.TEXTURE_MAG_FILTER, gl.NEAREST);
    gl.texParameteri(gl.TEXTURE_3D, gl.TEXTURE_MIN_FILTER, gl.NEAREST);
    gl.texParameteri(gl.TEXTURE_3D, gl.TEXTURE_WRAP_S, gl.CLAMP_TO_EDGE);
    gl.texParameteri(gl.TEXTURE_3D, gl.TEXTURE_WRAP_T, gl.CLAMP_TO_EDGE);
    gl.texParameteri(gl.TEXTURE_3D, gl.TEXTURE_WRAP_R, gl.CLAMP_TO_EDGE);

    return texture;
}


Vector3、Matrix4

動作確認レベルの実装ですので、お好みの行列ライブラリで置き換えてください。

class Vector3 {
    constructor(
        readonly x : number,
        readonly y : number,
        readonly z : number,
    ) {}
    add(v : Vector3) : Vector3 {
        return new Vector3(this.x + v.x, this.y + v.y, this.z + v.z);
    }
    sub(v : Vector3) : Vector3 {
        return new Vector3(this.x - v.x, this.y - v.y, this.z - v.z);
    }
    mul(s : number) : Vector3 {
        return new Vector3(this.x * s, this.y * s, this.z * s);
    }
    div(s : number) : Vector3 {
        return new Vector3(this.x / s, this.y / s, this.z / s);
    }
    length() : number {
        return Math.sqrt(this.x * this.x + this.y * this.y + this.z * this.z);
    }
    normalize() : Vector3 {
        return this.div(this.length());
    }
    dot(v : Vector3) : number {
        return this.x * v.x + this.y * v.y + this.z * v.z;
    }
    cross(v : Vector3) : Vector3 {
        return new Vector3(this.y * v.z - this.z * v.y, this.z * v.x - this.x * v.z, this.x * v.y - this.y * v.x);
    }
    static readonly ex = new Vector3(1, 0, 0);
    static readonly ey = new Vector3(0, 1, 0);
    static readonly ez = new Vector3(0, 0, 1);
}
class Matrix4 {
    constructor(readonly data: number[])
    { }
    mul(m : Matrix4) : Matrix4 {
        return new Matrix4([
            this.data[0] * m.data[0] + this.data[1] * m.data[4] + this.data[2] * m.data[8]  + this.data[3] * m.data[12],
            this.data[0] * m.data[1] + this.data[1] * m.data[5] + this.data[2] * m.data[9]  + this.data[3] * m.data[13],
            this.data[0] * m.data[2] + this.data[1] * m.data[6] + this.data[2] * m.data[10] + this.data[3] * m.data[14],
            this.data[0] * m.data[3] + this.data[1] * m.data[7] + this.data[2] * m.data[11] + this.data[3] * m.data[15],

            this.data[4] * m.data[0] + this.data[5] * m.data[4] + this.data[6] * m.data[8]  + this.data[7] * m.data[12],
            this.data[4] * m.data[1] + this.data[5] * m.data[5] + this.data[6] * m.data[9]  + this.data[7] * m.data[13],
            this.data[4] * m.data[2] + this.data[5] * m.data[6] + this.data[6] * m.data[10] + this.data[7] * m.data[14],
            this.data[4] * m.data[3] + this.data[5] * m.data[7] + this.data[6] * m.data[11] + this.data[7] * m.data[15],

            this.data[8] * m.data[0] + this.data[9] * m.data[4] + this.data[10] * m.data[8]  + this.data[11] * m.data[12],
            this.data[8] * m.data[1] + this.data[9] * m.data[5] + this.data[10] * m.data[9]  + this.data[11] * m.data[13],
            this.data[8] * m.data[2] + this.data[9] * m.data[6] + this.data[10] * m.data[10] + this.data[11] * m.data[14],
            this.data[8] * m.data[3] + this.data[9] * m.data[7] + this.data[10] * m.data[11] + this.data[11] * m.data[15],

            this.data[12] * m.data[0] + this.data[13] * m.data[4] + this.data[14] * m.data[8]  + this.data[15] * m.data[12],
            this.data[12] * m.data[1] + this.data[13] * m.data[5] + this.data[14] * m.data[9]  + this.data[15] * m.data[13],
            this.data[12] * m.data[2] + this.data[13] * m.data[6] + this.data[14] * m.data[10] + this.data[15] * m.data[14],
            this.data[12] * m.data[3] + this.data[13] * m.data[7] + this.data[14] * m.data[11] + this.data[15] * m.data[15],
        ]);
    }
    transform(v : Vector3) : Vector3 {
        const x = this.data[0]  * v.x + this.data[1]  * v.y + this.data[2]  * v.z + this.data[3];
        const y = this.data[4]  * v.x + this.data[5]  * v.y + this.data[6]  * v.z + this.data[7];
        const z = this.data[8]  * v.x + this.data[9]  * v.y + this.data[10] * v.z + this.data[11];
        const w = this.data[12] * v.x + this.data[13] * v.y + this.data[14] * v.z + this.data[15];
        return new Vector3(x / w, y / w, z / w);
    }
    transform3x3(v : Vector3) : Vector3 {
        return new Vector3(
            this.data[0] * v.x + this.data[1] * v.y + this.data[2]  * v.z,
            this.data[4] * v.x + this.data[5] * v.y + this.data[6]  * v.z,
            this.data[8] * v.x + this.data[9] * v.y + this.data[10] * v.z,
        );
    }
    transpose() : Matrix4 {
        return new Matrix4([
            this.data[0], this.data[4], this.data[8],  this.data[12],
            this.data[1], this.data[5], this.data[9],  this.data[13],
            this.data[2], this.data[6], this.data[10], this.data[14],
            this.data[3], this.data[7], this.data[11], this.data[15],
        ]);
    }
    getMatrix3Data() : number[] {
        return [
            this.data[0], this.data[1], this.data[2],
            this.data[4], this.data[5], this.data[6],
            this.data[8], this.data[9], this.data[10],
        ];
    }
    static readonly identity =
        new Matrix4([
            1, 0, 0, 0,
            0, 1, 0, 0,
            0, 0, 1, 0,
            0, 0, 0, 1
        ]);
    static translate(v : Vector3) : Matrix4 {
        return new Matrix4([
            1, 0, 0, v.x,
            0, 1, 0, v.y,
            0, 0, 1, v.z,
            0, 0, 0, 1
        ]);
    }
    static rotate(radian : number, axis : Vector3) : Matrix4 {
        const x = axis.x, y = axis.y, z = axis.z;
        const s = Math.sin(radian);
        const c = Math.cos(radian);
        return new Matrix4([
            x * x * (1 - c) + c,     x * y * (1 - c) - z * s, z * x * (1 - c) + y * s, 0,
            x * y * (1 - c) + z * s, y * y * (1 - c) + c,     y * z * (1 - c) - x * s, 0,
            z * x * (1 - c) - y * s, y * z * (1 - c) + x * s, z * z * (1 - c) + c,     0,
            0,                       0,                       0,                       1
        ]);
    }
    static scale(s : number) : Matrix4 {
        return new Matrix4([
            s, 0, 0, 0,
            0, s, 0, 0,
            0, 0, s, 0,
            0, 0, 0, 1
        ]);
    }
    static ortho(width : number, height : number, near : number, far : number) : Matrix4 {
        return new Matrix4([
            2.0 / width, 0,            0,                   0,
            0,           2.0 / height, 0,                   0,
            0,           0,            -2.0 / (far - near), -(far + near) / (far - near),
            0,           0,            0,                   1
        ]);
    }
    static perspective(fovy : number, aspect : number, near : number, far : number) : Matrix4 {
        const f = 1.0 / Math.tan(fovy / 2.0);
        return new Matrix4([
            f / aspect, 0, 0,                           0,
            0,          f, 0,                           0,
            0,          0, (far + near) / (near - far), 2 * far * near / (near - far),
            0,          0, -1,                          0
        ]);
    }
}


Array関連の拡張

Slices作成内で利用しているArray関連の関数は下記のように実装してあります。

interface Array<T> {
    concatMap<U>(f : (value : T) => U[]) : U[];
}
Array.prototype.concatMap = function(f : any) {
    return this.reduce((dst, x) => dst.concat(f(x)), []);
};

function range(count : number) : number[] {
    const ary = Array<number>(count);
    for (let i = 0; i < count; ++i)
        ary[i] = i;
    return ary;
}
function normalizedRange(count : number, min : number, max : number) : number[] {
    return range(count).map<number>(i => ((count - 1 - i) * min + i * max) / (count - 1));
}