Direct3D11でボリュームレンダリング

ボリュームレンダリングの手法はいくつかありますが、今回はもっとも単純な、下図のように長方形のスライスを何枚も重ねて実現したものを実装したいと思います。
f:id:any-programming:20171225171609p:plain

ボリュームデータ

下記ページにある、bunnyのデータを利用いたしました。
The Stanford volume data archive

もとは1ボクセルあたり16bitだったため、8bitに変換して1ファイルにまとめました。

アプリケーションコード

下図のような、マウスの位置によって回転角度の変わるアプリケーションを実装します。
f:id:any-programming:20171225195538p:plain

main関数とWndProcは下記の様に実装しました。

VolumeRenderer renderer;
float rotX = 0;
float rotY = 0;
LRESULT CALLBACK WndProc(HWND hwnd, UINT msg, WPARAM wParam, LPARAM lParam)
{
    switch (msg)
    {
    case WM_PAINT:
        renderer.Render(rotX, rotY);
        return 0;
    case WM_MOUSEMOVE:
    {
        RECT rect;
        ::GetClientRect(hwnd, &rect);
        int x = LOWORD(lParam);
        int y = HIWORD(lParam);
        rotY = 2 * 3.14f * (float)(x - rect.left) / (rect.right - rect.left);
        rotX = 2 * 3.14f * (float)(y - rect.top) / (rect.bottom - rect.top);
        InvalidateRect(hwnd, nullptr, FALSE);
        return 0;
    }
    case WM_DESTROY:
        PostQuitMessage(0);
        return 0;
    }
    return DefWindowProc(hwnd, msg, wParam, lParam);
}
int main()
{
    HWND hwnd = CreateHWND(WndProc);
    renderer = CreateVolumeRenderer(hwnd);
    Run(hwnd);

    return 0;
}


VolumeRenderer構造体とCreateVolumeRenderer関数

今回はRender関数のみをメンバとして持ったものとしました。

struct VolumeRenderer
{
    std::function<void(float rotX, float rotY)> Render;
};


HWNDを受け取りVolumeRendererを作成する関数を定義します。

複数のスライスを重ねてボリュームを描画する際、横から見てしまうとボリュームにならないため、X、Y、Z三方向のスライスセットを用意しておき随時適切なものを利用して描画するようにしています。

陰影はTextureのX、Y、Z方向の勾配を計算して、それを法線として利用しています。

VolumeRenderer CreateVolumeRenderer(HWND hwnd) {
    // device, context, swapChain作成
    ComPtr<ID3D11Device> device;
    ComPtr<ID3D11DeviceContext> context;
    ComPtr<IDXGISwapChain> swapChain;
    std::tie(device, context, swapChain) = CreateDeviceAndSwapChain(hwnd);

    // RenderTarget、DepthBuffer作成
    ComPtr<ID3D11RenderTargetView> renderTarget = CreateRenderTarget(device.Get(), swapChain.Get());
    SIZE bufferSize = GetBufferSize(swapChain.Get());
    
    ComPtr<ID3D11DepthStencilView> depthBuffer = CreateDepthBuffer(device.Get(), bufferSize.cx, bufferSize.cy);

    // Shader作成
    std::string vertexCode =
        "cbuffer c0{                                                                                \n"
        "    float4x4 Matrix;                                                                       \n"
        "}                                                                                          \n"
        "struct VS_OUTPUT {                                                                         \n"
        "    float4 Position : SV_POSITION;                                                         \n"
        "    float3 Texcoord : TEXCOORD;                                                            \n"
        "};                                                                                         \n"
        "VS_OUTPUT main(float4 position : POSITION, float3 texcoord : TEXCOORD) {                   \n"
        "    VS_OUTPUT output = (VS_OUTPUT)0;                                                       \n"
        "    output.Position = mul(position, Matrix);                                               \n"
        "    output.Texcoord = texcoord;                                                            \n"
        "    return output;                                                                         \n"
        "}                                                                                          \n";
    std::vector<D3D11_INPUT_ELEMENT_DESC> inputDesc =
    {
        { "POSITION", 0, DXGI_FORMAT_R32G32B32_FLOAT, 0, 0,                            D3D11_INPUT_PER_VERTEX_DATA, 0 },
        { "TEXCOORD", 0, DXGI_FORMAT_R32G32B32_FLOAT, 0, D3D11_APPEND_ALIGNED_ELEMENT, D3D11_INPUT_PER_VERTEX_DATA, 0 }
    };
    ComPtr<ID3D11VertexShader> vertexShader;
    ComPtr<ID3D11InputLayout> inputLayout;
    std::tie(vertexShader, inputLayout) = CreateVertexShader(device.Get(), vertexCode, inputDesc);
    ComPtr<ID3D11Buffer> vertexConstBuffer = CreateConstantBuffer<VertexConst>(device.Get());

    std::string pixelCode =
        "Texture3D texture0 : register(t0);                                                                                                  \n"
        "SamplerState sampler0 {                                                                                                             \n"
        "    Filter = MIN_MAG_MIP_POINT;                                                                                                     \n"
        "    AddressU = Clamp;                                                                                                               \n"
        "    AddressV = Clamp;                                                                                                               \n"
        "    AddressW = Clamp;                                                                                                               \n"
        "};                                                                                                                                  \n"
        "float4 main(float4 position : SV_POSITION, float3 texcoord : TEXCOORD) : SV_TARGET {                                                \n"
        "    if (texture0.Sample(sampler0, texcoord).r < 0.5)                                                                                \n"
        "       discard;                                                                                                                     \n"
        "    float dx = texture0.Sample(sampler0, texcoord, int3(1, 0, 0)).r - texture0.Sample(sampler0, texcoord, int3(-1,  0,  0)).r;      \n"
        "    float dy = texture0.Sample(sampler0, texcoord, int3(0, 1, 0)).r - texture0.Sample(sampler0, texcoord, int3( 0, -1,  0)).r;      \n"
        "    float dz = texture0.Sample(sampler0, texcoord, int3(0, 0, 1)).r - texture0.Sample(sampler0, texcoord, int3( 0,  0, -1)).r;      \n"
        "    float3 normal = normalize(float3(dx, dy, dz));                                                                                  \n"
        "    float3 light = normalize(float3(1, 1, 1));                                                                                      \n"
        "    float gray = abs(dot(normal, light));                                                                                           \n"
        "    return float4(gray, gray, gray, 1.0);                                                                                           \n"
        "}                                                                                                                                   \n";
    ComPtr<ID3D11PixelShader> pixelShader = CreatePixelShader(device.Get(), pixelCode);

    // Texture3D作成
    std::ifstream file("bunny", std::ios::binary | std::ios::ate);
    std::vector<byte> data((size_t)file.tellg());
    file.seekg(0, std::ios::beg);
    file.read((char*)data.data(), data.size());

    ComPtr<ID3D11ShaderResourceView> volumeTexture = CreateTexture3D(device.Get(), 512, 512, data);

    // VertexBuffer, IndexBuffer作成(X、Y、Z三方向のBuffer)
    Slices slices = CreateSlices(device.Get(), 512, 512, 360);

    // culling off用 RasterrizerState作成
    ComPtr<ID3D11RasterizerState> rasterizerState = CreateRasterizerState(device.Get(), [](D3D11_RASTERIZER_DESC& desc) {
        desc.CullMode = D3D11_CULL_NONE;
    });

    // Render関数作成
    auto render = [=](float rotX, float rotY) {
        // Matrix作成
        Matrix4 mat = 
            Matrix4::Rotate(rotY, Vector3{ 0, 1, 0 }) *
            Matrix4::Rotate(rotX, Vector3{ 1, 0, 0 }) *
            Matrix4::Scale(1.5 / 512.0) *
            Matrix4::Translate(Vector3{ -512 / 2, -512 / 2, -360 / 2 });

        // RenderTarget、DepthBuffer設定
        context->OMSetRenderTargets(1, renderTarget.GetAddressOf(), depthBuffer.Get());

        // RasterizerState設定
        context->RSSetState(rasterizerState.Get());

        // Shader設定
        context->VSSetShader(vertexShader.Get(), nullptr, 0);
        context->IASetInputLayout(inputLayout.Get());
        context->PSSetShader(pixelShader.Get(), nullptr, 0);
        
        // Buffer設定
        ID3D11Buffer *vertexBuffer, *indexBuffer;
        std::tie(vertexBuffer, indexBuffer) = GetBuffers(slices, mat);        
        SetVertexBuffer<Vertex>(context.Get(), vertexBuffer);
        context->IASetIndexBuffer(indexBuffer, DXGI_FORMAT_R16_UINT, 0);

        VertexConst vertexConst = { Matrix4::OrthoD3D(2, 2, -1, 1) * mat };
        context->UpdateSubresource(vertexConstBuffer.Get(), 0, nullptr, &vertexConst, 0, 0);
        context->VSSetConstantBuffers(0, 1, vertexConstBuffer.GetAddressOf());

        // Texture3D設定
        context->PSSetShaderResources(0, 1, volumeTexture.GetAddressOf());

        // Viewport設定
        CD3D11_VIEWPORT viewport(0.0f, 0.0f, (float)bufferSize.cx, (float)bufferSize.cy);
        context->RSSetViewports(1, &viewport);

        // 描画
        float clearColor[4] = { 0.5f, 0.6f, 1.0f, 1.0f };
        context->ClearRenderTargetView(renderTarget.Get(), clearColor);
        context->ClearDepthStencilView(depthBuffer.Get(), D3D11_CLEAR_DEPTH, 1.0f, 0);

        context->IASetPrimitiveTopology(D3D11_PRIMITIVE_TOPOLOGY_TRIANGLELIST);
        context->DrawIndexed(GetBufferSize<ushort>(indexBuffer), 0, 0);

        AssertHR(swapChain->Present(0, 0));
    };

    return{ render };
}


VertexConstとVertexは下記の様に定義してあります。

struct VertexConst {
    Matrix4 Matrix;
};
struct Vertex {
    Vector3 Position;
    Vector3 Texcoord;
};


Slices作成

一方向のスライスセットでは横から見たときボリュームにならないため、三方向のスライスセットを用意しておきます。

Slices構造体は下記の様に定義します。

struct Slices {
    ComPtr<ID3D11Buffer> VertexBufferX;
    ComPtr<ID3D11Buffer> IndexBufferX;
    ComPtr<ID3D11Buffer> VertexBufferY;
    ComPtr<ID3D11Buffer> IndexBufferY;
    ComPtr<ID3D11Buffer> VertexBufferZ;
    ComPtr<ID3D11Buffer> IndexBufferZ;
};


Slicesを作成する関数を下記の様に実装します。

(0, 0, 0)を原点として、(width, height, depth)のボックスのSlicesを作成しています。

Slices CreateSlices(ID3D11Device* device, int width, int height, int depth) {
    float w = (float)width, h = (float)height, d = (float)depth;

    std::vector<Vertex> verticesZ;
    for (int i = 0; i < depth; ++i) {
        float t = (float)i / (depth - 1);
        verticesZ.push_back({ { 0, 0, t * d }, { 0, 0, t } });
        verticesZ.push_back({ { w, 0, t * d }, { 1, 0, t } });
        verticesZ.push_back({ { 0, h, t * d }, { 0, 1, t } });
        verticesZ.push_back({ { w, h, t * d }, { 1, 1, t } });
    }

    std::vector<Vertex> verticesY;
    for (int i = 0; i < height; ++i) {
        float t = (float)i / (height - 1);
        verticesY.push_back({ { 0, t * h, 0 }, { 0, t, 0 } });
        verticesY.push_back({ { w, t * h, 0 }, { 1, t, 0 } });
        verticesY.push_back({ { 0, t * h, d }, { 0, t, 1 } });
        verticesY.push_back({ { w, t * h, d }, { 1, t, 1 } });
    }

    std::vector<Vertex> verticesX;
    for (int i = 0; i < width; ++i) {
        float t = (float)i / (width - 1);
        verticesX.push_back({ { t * w, 0, 0 }, { t, 0, 0 } });
        verticesX.push_back({ { t * w, h, 0 }, { t, 1, 0 } });
        verticesX.push_back({ { t * w, 0, d }, { t, 0, 1 } });
        verticesX.push_back({ { t * w, h, d }, { t, 1, 1 } });
    }

    auto createIndices = [](int count) {
        std::vector<ushort> indices;
        for (int i = 0; i < count; ++i)
            for (ushort j : { 0, 2, 1, 1, 2, 3 })
                indices.push_back(4 * i + j);
        return indices;
    };

    return{
        CreateVertexBuffer(device, verticesX),
        CreateIndexBuffer(device, createIndices(width)),
        CreateVertexBuffer(device, verticesY),
        CreateIndexBuffer(device, createIndices(height)),
        CreateVertexBuffer(device, verticesZ),
        CreateIndexBuffer(device, createIndices(depth))
    };
}


変換行列から適切なスライスセットを返す関数も定義します。

最も視線方向に向いているスライスセットを選択しています。

std::tuple<ID3D11Buffer*, ID3D11Buffer*> GetBuffers(const Slices& slices, const Matrix4& m) {
    float dotX = std::abs(m.Transform3x3(Vector3{ 1, 0, 0 }).Z);
    float dotY = std::abs(m.Transform3x3(Vector3{ 0, 1, 0 }).Z);
    float dotZ = std::abs(m.Transform3x3(Vector3{ 0, 0, 1 }).Z);
    if (dotX < dotZ && dotY < dotZ)
        return{ slices.VertexBufferZ.Get(), slices.IndexBufferZ.Get() };
    else if (dotX < dotY)
        return{ slices.VertexBufferY.Get(), slices.IndexBufferY.Get() };
    else
        return{ slices.VertexBufferX.Get(), slices.IndexBufferX.Get() };
}


Direct3D11用関数群

コード内で利用したDirect3D11関連の関数は下記の様に実装してあります。

std::tuple<ComPtr<ID3D11Device>, ComPtr<ID3D11DeviceContext>, ComPtr<IDXGISwapChain>> CreateDeviceAndSwapChain(HWND hwnd) {
#ifdef _DEBUG
    DWORD createDeviceFlag = D3D11_CREATE_DEVICE_DEBUG;
#else
    DWORD createDeviceFlag = 0;
#endif
   
    DXGI_SWAP_CHAIN_DESC desc = { 0 };
    desc.BufferDesc.Width = 0; // 現在のWindowサイズが自動で設定される
    desc.BufferDesc.Height = 0;
    desc.BufferDesc.Format = DXGI_FORMAT_B8G8R8A8_UNORM;
    desc.BufferDesc.RefreshRate.Numerator = 60;
    desc.BufferDesc.RefreshRate.Denominator = 1;
    desc.BufferDesc.Scaling = DXGI_MODE_SCALING_UNSPECIFIED;
    desc.BufferDesc.ScanlineOrdering = DXGI_MODE_SCANLINE_ORDER_UNSPECIFIED;
    desc.SampleDesc.Count = 1;
    desc.SampleDesc.Quality = 0;
    desc.BufferUsage = DXGI_USAGE_RENDER_TARGET_OUTPUT;
    desc.BufferCount = 1;
    desc.OutputWindow = hwnd;
    desc.Windowed = TRUE;
    desc.SwapEffect = DXGI_SWAP_EFFECT_DISCARD;
    desc.Flags = 0;

    ComPtr<ID3D11Device> device;
    ComPtr<ID3D11DeviceContext> context;
    ComPtr<IDXGISwapChain> swapChain;
    AssertHR(D3D11CreateDeviceAndSwapChain(
        nullptr,
        D3D_DRIVER_TYPE_HARDWARE,
        nullptr,
        createDeviceFlag,
        nullptr,
        0,
        D3D11_SDK_VERSION,
        &desc,
        &swapChain,
        &device,
        nullptr,
        &context
    ));

    return { device, context, swapChain };
}
ComPtr<ID3D11RenderTargetView> CreateRenderTarget(ID3D11Device* device, IDXGISwapChain* swapChain) {
    ComPtr<ID3D11Texture2D> backBuffer;
    AssertHR(swapChain->GetBuffer(0, IID_PPV_ARGS(&backBuffer)));
    
    ComPtr<ID3D11RenderTargetView> view;
    AssertHR(device->CreateRenderTargetView(backBuffer.Get(), nullptr, &view));
    
    return view;
}
SIZE GetBufferSize(IDXGISwapChain* swapChain) {
    DXGI_SWAP_CHAIN_DESC desc;
    AssertHR(swapChain->GetDesc(&desc));    
    return{ (LONG)desc.BufferDesc.Width, (LONG)desc.BufferDesc.Height };
}
ComPtr<ID3D11DepthStencilView> CreateDepthBuffer(ID3D11Device* device, int width, int height) {
    D3D11_TEXTURE2D_DESC textureDesc = { 0 };
    textureDesc.Width = width;
    textureDesc.Height = height;
    textureDesc.Format = DXGI_FORMAT_D16_UNORM;
    textureDesc.MipLevels = 1;
    textureDesc.ArraySize = 1;
    textureDesc.SampleDesc.Count = 1;
    textureDesc.SampleDesc.Quality = 0;
    textureDesc.Usage = D3D11_USAGE_DEFAULT;
    textureDesc.BindFlags = D3D11_BIND_DEPTH_STENCIL;
    textureDesc.CPUAccessFlags = 0;

    ComPtr<ID3D11Texture2D> texture;
    AssertHR(device->CreateTexture2D(&textureDesc, nullptr, &texture));

    CD3D11_DEPTH_STENCIL_VIEW_DESC viewDesc(D3D11_DSV_DIMENSION_TEXTURE2D);
    ComPtr<ID3D11DepthStencilView> view;
    AssertHR(device->CreateDepthStencilView(texture.Get(), &viewDesc, &view));

    return view;
}
ComPtr<ID3DBlob> CompileShader(std::string code, std::string target) {
#ifdef _DEBUG
    UINT flags1 = D3DCOMPILE_DEBUG | D3DCOMPILE_ENABLE_STRICTNESS | D3DCOMPILE_WARNINGS_ARE_ERRORS;
#else
    UINT flags1 = D3DCOMPILE_OPTIMIZATION_LEVEL3 | D3DCOMPILE_ENABLE_STRICTNESS | D3DCOMPILE_WARNINGS_ARE_ERRORS;
#endif
    ComPtr<ID3DBlob> compiled, errorMessage;
    D3DCompile(
        code.c_str(),
        code.size(),
        nullptr,
        nullptr,
        D3D_COMPILE_STANDARD_FILE_INCLUDE,
        "main",
        target.c_str(),
        flags1,
        0,
        &compiled,
        &errorMessage
    );
    if (errorMessage) 
        throw std::exception((char*)errorMessage->GetBufferPointer());

    return compiled;
}
std::tuple<ComPtr<ID3D11VertexShader>, ComPtr<ID3D11InputLayout>> CreateVertexShader(
    ID3D11Device* device, const std::vector<byte>& cso, const std::vector<D3D11_INPUT_ELEMENT_DESC>& inputDesc) {
    ComPtr<ID3D11VertexShader> shader;
    AssertHR(device->CreateVertexShader(cso.data(), cso.size(), NULL, &shader));

    ComPtr<ID3D11InputLayout> inputLayout;
    AssertHR(device->CreateInputLayout(inputDesc.data(), inputDesc.size(), cso.data(), cso.size(), &inputLayout));

    return { shader, inputLayout };
}
ComPtr<ID3D11PixelShader> CreatePixelShader(ID3D11Device* device, const std::vector<byte>& cso) {
    ComPtr<ID3D11PixelShader> shader;
    AssertHR(device->CreatePixelShader(cso.data(), cso.size(), NULL, &shader));
    return shader;
}
std::tuple<ComPtr<ID3D11VertexShader>, ComPtr<ID3D11InputLayout>> CreateVertexShader(
    ID3D11Device* device, std::string code, const std::vector<D3D11_INPUT_ELEMENT_DESC>& inputDesc) {
    auto blob = CompileShader(code, "vs_4_0");
    auto p = (byte*)blob->GetBufferPointer();
    std::vector<byte> cso(p, p + blob->GetBufferSize());
    return CreateVertexShader(device, cso, inputDesc);
}
ComPtr<ID3D11PixelShader> CreatePixelShader(ID3D11Device* device, std::string code) {
    auto blob = CompileShader(code, "ps_4_0");
    auto p = (byte*)blob->GetBufferPointer();
    std::vector<byte> cso(p, p + blob->GetBufferSize());
    return CreatePixelShader(device, cso);
}
template <class T>
ComPtr<ID3D11Buffer> CreateVertexBuffer(ID3D11Device* device, const std::vector<T>& vertices) {
    D3D11_BUFFER_DESC desc = { 0 };
    desc.Usage = D3D11_USAGE_DEFAULT;
    desc.ByteWidth = vertices.size() * sizeof(T);
    desc.BindFlags = D3D11_BIND_VERTEX_BUFFER;
    desc.CPUAccessFlags = 0;

    D3D11_SUBRESOURCE_DATA data = { 0 };
    data.pSysMem = vertices.data();

    ComPtr<ID3D11Buffer> buffer;
    AssertHR(device->CreateBuffer(&desc, &data, &buffer));

    return buffer;
}
template <class T>
void SetVertexBuffer(ID3D11DeviceContext* context, ID3D11Buffer* buffer) {
    UINT stride = sizeof(T);
    UINT offset = 0;
    context->IASetVertexBuffers(0, 1, &buffer, &stride, &offset);
}
template <class T>
ComPtr<ID3D11Buffer> CreateIndexBuffer(ID3D11Device* device, const std::vector<T>& indices) {
    D3D11_BUFFER_DESC desc = { 0 };
    desc.Usage = D3D11_USAGE_DEFAULT;
    desc.ByteWidth = indices.size() * sizeof(T);
    desc.BindFlags = D3D11_BIND_INDEX_BUFFER;
    desc.CPUAccessFlags = 0;

    D3D11_SUBRESOURCE_DATA data = { 0 };
    data.pSysMem = indices.data();

    ComPtr<ID3D11Buffer> buffer;
    AssertHR(device->CreateBuffer(&desc, &data, &buffer));

    return buffer;
}
template<class T>
UINT GetBufferSize(ID3D11Buffer* buffer) {
    D3D11_BUFFER_DESC desc;
    buffer->GetDesc(&desc);
    return desc.ByteWidth / sizeof(T);
}
template <class T>
ComPtr<ID3D11Buffer> CreateConstantBuffer(ID3D11Device* device) {
    D3D11_BUFFER_DESC desc = { 0 };
    desc.Usage = D3D11_USAGE_DEFAULT;
    desc.ByteWidth = sizeof(T);
    desc.BindFlags = D3D11_BIND_CONSTANT_BUFFER;
    desc.CPUAccessFlags = 0;

    ComPtr<ID3D11Buffer> buffer;
    AssertHR(device->CreateBuffer(&desc, nullptr, &buffer));

    return buffer;
}
ComPtr<ID3D11ShaderResourceView> CreateTexture3D(ID3D11Device* device, int width, int height, const std::vector<byte>& voxels) {
    D3D11_TEXTURE3D_DESC desc = { 0 };
    desc.Width = width;
    desc.Height = height;
    desc.Depth = voxels.size() / width / height;
    desc.Format = DXGI_FORMAT_R8_UNORM;
    desc.MipLevels = 1;
    desc.Usage = D3D11_USAGE_DEFAULT;
    desc.BindFlags = D3D11_BIND_SHADER_RESOURCE;
    desc.CPUAccessFlags = 0;
    desc.MiscFlags = 0;

    D3D11_SUBRESOURCE_DATA data = { 0 };
    data.pSysMem = voxels.data();
    data.SysMemPitch = width;
    data.SysMemSlicePitch = width * height;

    ComPtr<ID3D11Texture3D> texture;
    AssertHR(device->CreateTexture3D(&desc, &data, &texture));

    ComPtr<ID3D11ShaderResourceView> view;
    AssertHR(device->CreateShaderResourceView(texture.Get(), nullptr, &view));
    
    return view;
}
ComPtr<ID3D11RasterizerState> CreateRasterizerState(ID3D11Device* device, std::function<void(D3D11_RASTERIZER_DESC&)> f) {
    CD3D11_RASTERIZER_DESC desc = CD3D11_RASTERIZER_DESC(CD3D11_DEFAULT());
    f(desc);
    ComPtr<ID3D11RasterizerState> state;
    AssertHR(device->CreateRasterizerState(&desc, &state));
    return state;
}


Vector3、Matrix4

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

struct Vector3 {
    float X, Y, Z;
    Vector3 operator - () const { 
        return{ -X, -Y, -Z };
    }
    Vector3 operator + (const Vector3& v) const { 
        return{ X + v.X, Y + v.Y, Z + v.Z }; 
    }
    Vector3 operator - (const Vector3& v) const { 
        return{ X - v.X, Y - v.Y, Z - v.Z }; 
    }
    Vector3 operator * (float s) const { 
        return{ s * X, s * Y, s * Z }; 
    }
    Vector3 operator / (float s) const { 
        return{ X / s, Y / s, Z / s }; 
    }
    friend Vector3 operator * (float s, const Vector3& v) { 
        return v * s; 
    }
    float Length() const { 
        return std::sqrt(X * X + Y * Y + Z * Z); 
    }
    Vector3 Normalize() const {
        float length = this->Length();
        return{ X / length, Y / length, Z / length };
    }
    float Dot(const Vector3& v) const { 
        return X * v.X + Y * v.Y + Z * v.Z; 
    }
    Vector3 Cross(const Vector3& v) const { 
        return{ Y * v.Z - Z * v.Y, Z * v.X - X * v.Z, X * v.Y - Y * v.X }; 
    }
};
struct Matrix4 {
    float M00, M01, M02, M03;
    float M10, M11, M12, M13;
    float M20, M21, M22, M23;
    float M30, M31, M32, M33;
    Matrix4 operator * (const Matrix4& m) const {
        return{
            M00 * m.M00 + M01 * m.M10 + M02 * m.M20 + M03 * m.M30,
            M00 * m.M01 + M01 * m.M11 + M02 * m.M21 + M03 * m.M31,
            M00 * m.M02 + M01 * m.M12 + M02 * m.M22 + M03 * m.M32,
            M00 * m.M03 + M01 * m.M13 + M02 * m.M23 + M03 * m.M33,            

            M10 * m.M00 + M11 * m.M10 + M12 * m.M20 + M13 * m.M30,
            M10 * m.M01 + M11 * m.M11 + M12 * m.M21 + M13 * m.M31,
            M10 * m.M02 + M11 * m.M12 + M12 * m.M22 + M13 * m.M32,
            M10 * m.M03 + M11 * m.M13 + M12 * m.M23 + M13 * m.M33,            

            M20 * m.M00 + M21 * m.M10 + M22 * m.M20 + M23 * m.M30,
            M20 * m.M01 + M21 * m.M11 + M22 * m.M21 + M23 * m.M31,
            M20 * m.M02 + M21 * m.M12 + M22 * m.M22 + M23 * m.M32,
            M20 * m.M03 + M21 * m.M13 + M22 * m.M23 + M23 * m.M33,            

            M30 * m.M00 + M31 * m.M10 + M32 * m.M20 + M33 * m.M30,
            M30 * m.M01 + M31 * m.M11 + M32 * m.M21 + M33 * m.M31,
            M30 * m.M02 + M31 * m.M12 + M32 * m.M22 + M33 * m.M32,
            M30 * m.M03 + M31 * m.M13 + M32 * m.M23 + M33 * m.M33
        };
    }
    Vector3 Transform(const Vector3& v) const {
        float x = M00 * v.X + M01 * v.Y + M02 * v.Z + M03;
        float y = M10 * v.X + M11 * v.Y + M12 * v.Z + M13;
        float z = M20 * v.X + M21 * v.Y + M22 * v.Z + M23;
        float w = M30 * v.X + M31 * v.Y + M32 * v.Z + M33;
        return Vector3{x, y, z} / w;
    }
    Vector3 Transform3x3(const Vector3& v) const {
        return Vector3{
            M00 * v.X + M01 * v.Y + M02 * v.Z,
            M10 * v.X + M11 * v.Y + M12 * v.Z,
            M20 * v.X + M21 * v.Y + M22 * v.Z,
        };
    }
    static Matrix4 Identity() {
        return{ 
            1, 0, 0, 0,
            0, 1, 0, 0,
            0, 0, 1, 0,
            0, 0, 0, 1 
        };
    }
    static Matrix4 Translate(const Vector3& v) {
        return{
            1, 0, 0, v.X,
            0, 1, 0, v.Y,
            0, 0, 1, v.Z,
            0, 0, 0, 1
        };
    }
    static Matrix4 Rotate(float radian, const Vector3& axis) {
        float x = axis.X, y = axis.Y, z = axis.Z;
        float s = std::sin(radian);
        float c = std::cos(radian);
        return{
            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 Matrix4 Scale(float s) {
        return{
            s, 0, 0, 0,
            0, s, 0, 0,
            0, 0, s, 0,
            0, 0, 0, 1
        };
    }
    static Matrix4 OrthoD3D(float width, float height, float nearz, float farz) {
        return{
            2 / width, 0,          0,                  0,
            0,         2 / height, 0,                  0,
            0,         0,          1 / (farz - nearz), -nearz / (farz - nearz),
            0,         0,          0,                  1
        };
    }
    static Matrix4 PerspectiveD3D(float fovy, float aspect, float nearz, float farz) {
        float f = 1 / std::tan(fovy / 2.0f);
        return{
            f / aspect, 0, 0,                     0,
            0,          f, 0,                     0,
            0,          0, farz / (farz - nearz), -nearz * farz / (farz - nearz),
            0,          0, 1,                     0
        };
    }
};