diff --git a/DepthSample.ipynb b/DepthSample.ipynb new file mode 100644 index 0000000..aeaff14 --- /dev/null +++ b/DepthSample.ipynb @@ -0,0 +1,410 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "5ce67ea1-4610-4ef0-96f5-cb1b9e399b3b", + "metadata": {}, + "source": [ + "# Using depth samples\n", + "\n", + "I want to take a sample from a depth buffer and find its original view position. It's not as simple as multiplying by the inverse of the projection matrix due to the pesky perspective divide." + ] + }, + { + "cell_type": "markdown", + "id": "82c5c3c7-5777-4db4-96db-3d10355f017b", + "metadata": {}, + "source": [ + "## Definitions\n", + "\n", + "**$P$ is my projection matrix**. I'm making some simplyfing assumptions by throwing some zeros in. If the projection matrix I use doesn't fit the assumptions, I'll need to decompose it first. Anyway:\n", + "\n", + "$$\n", + "P = \\begin{bmatrix}\n", + "a & 0 & c & d\\\\\n", + "0 & f & g & h\\\\\n", + "0 & 0 & k & l\\\\\n", + "0 & 0 & o & p\n", + "\\end{bmatrix}\n", + "$$\n", + "\n", + "**$v$ is a view space vertex** that I'll be transforming.\n", + "\n", + "$$\n", + "v = \\begin{bmatrix}\n", + "v_x\\\\\n", + "v_y\\\\\n", + "v_z\\\\\n", + "1\n", + "\\end{bmatrix}\n", + "$$\n", + "\n", + "**$v'$ is the transformed vertex**, $P\\cdot v$.\n", + "\n", + "$$\n", + "v' = P v = \\begin{bmatrix}\n", + "a & 0 & c & d\\\\\n", + "0 & f & g & h\\\\\n", + "0 & 0 & k & l\\\\\n", + "0 & 0 & o & p\n", + "\\end{bmatrix}\n", + "\\begin{bmatrix}\n", + "v_x\\\\\n", + "v_y\\\\\n", + "v_z\\\\\n", + "1\n", + "\\end{bmatrix}\n", + "=\n", + "\\begin{bmatrix}\n", + "av_x + cv_z + d\\\\\n", + "fv_y + gv_z + h\\\\\n", + "kv_z + l\\\\\n", + "ov_z + p\n", + "\\end{bmatrix}\n", + "$$\n", + "\n", + "**$v''$ is the vertex after perspective divide**.\n", + "\n", + "$$\n", + "v'' = \\begin{bmatrix}\n", + "\\frac{v'_x}{v'_w}\\\\\n", + "\\frac{v'_y}{v'_w}\\\\\n", + "\\frac{v'_z}{v'_w}\\\\\n", + "v'_w\n", + "\\end{bmatrix}\n", + "=\n", + "\\begin{bmatrix}\n", + "\\frac{av_x + cv_z + d}{ov_z + p}\\\\\n", + "\\frac{fv_y + gv_z + h}{ov_z + p}\\\\\n", + "\\frac{kv_z + l}{ov_z + p}\\\\\n", + "ov_z + p\n", + "\\end{bmatrix}\n", + "$$\n", + "\n", + "**$v'''$ is the vertex in NDC space**.\n", + "$$\n", + "v''' = \\begin{bmatrix}\n", + "\\frac{1}{2} & 0 & 0 & \\frac{1}{2}\\\\\n", + "0 & \\frac{1}{2} & 0 & \\frac{1}{2}\\\\\n", + "0 & 0 & \\frac{1}{2} & \\frac{1}{2}\\\\\n", + "0 & 0 & 0 & 1\n", + "\\end{bmatrix}\n", + "\\begin{bmatrix}\n", + "v''_x\\\\\n", + "v''_y\\\\\n", + "v''_z\\\\\n", + "1\n", + "\\end{bmatrix}\n", + "=\n", + "\\begin{bmatrix}\n", + "\\frac{v''_x + 1}{2}\\\\\n", + "\\frac{v''_y + 1}{2}\\\\\n", + "\\frac{v''_z + 1}{2}\\\\\n", + "1\n", + "\\end{bmatrix}\n", + "=\n", + "\\begin{bmatrix}\n", + "\\frac{a v_{x} + c v_{z} + d + o v_{z} + p}{2 \\left(o v_{z} + p\\right)}\\\\\n", + "\\frac{f v_{y} + g v_{z} + h + o v_{z} + p}{2 \\left(o v_{z} + p\\right)}\\\\\n", + "\\frac{k v_{z} + l + o v_{z} + p}{2 \\left(o v_{z} + p\\right)}\\\\\n", + "1\n", + "\\end{bmatrix}\n", + "$$\n", + "\n", + "**$s_u, s_v$ is my UV sample position**, and **$s_z$ is the sampled depth value** at that position. Conveniently,\n", + "\n", + "$$\n", + "v''' = \\begin{bmatrix}\n", + "s_u\\\\\n", + "s_v\\\\\n", + "s_z\\\\\n", + "1\n", + "\\end{bmatrix}\n", + "$$" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "id": "b102ea74-6c0d-4d3f-a6cc-accd43dcff97", + "metadata": {}, + "outputs": [], + "source": [ + "# I'm going to set up the definitions in code:\n", + "\n", + "from sympy import *\n", + "from IPython.display import Markdown, Latex\n", + "init_printing(use_unicode=True)\n", + "a, c, d, f, g, h, k, l, o, p = symbols(\"a, c, d, f, g, h, k, l, o, p\")\n", + "P = Matrix([[a, 0, c, d], [0, f, g, h], [0, 0, k, l], [0, 0, o, p]])\n", + "\n", + "v_x, v_y, v_z = symbols(\"v_x, v_y, v_z\")\n", + "v = Matrix([[v_x], [v_y], [v_z], [1]])\n", + "\n", + "v1 = P * v\n", + "v1_x, v1_y, v1_z, v1_w = v1[0,0], v1[1,0], v1[2,0], v1[3,0]\n", + "\n", + "v2 = Matrix([[v1_x/v1_w], [v1_y/v1_w], [v1_z/v1_w], [v1_w]])\n", + "v2_x, v2_y, v2_z, v2_w = v2[0,0], v2[1,0], v2[2,0], v2[3,0]\n", + "\n", + "half = Rational(0.5)\n", + "NDC = Matrix([[half, 0, 0, half], [0, half, 0, half], [0, 0, half, half], [0, 0, 0, 1]])\n", + "v3 = simplify(NDC * Matrix([[v2_x], [v2_y], [v2_z], [1]]))\n", + "v3_x, v3_y, v3_z, v3_w = v3[0,0], v3[1,0], v3[2,0], v3[3,0]" + ] + }, + { + "cell_type": "markdown", + "id": "0f2bfae7-a3fc-4f7c-9767-da87ad4d923a", + "metadata": {}, + "source": [ + "## Can I recover $v'_w$?\n", + "\n", + "I want to undo the perspective divide. If I could determine $v'_w$ from $v'''$, then I could just multiply it back on.\n", + "\n", + "$$\n", + "v'''_z = \\frac{k v_{z} + l + o v_{z} + p}{2 \\left(o v_{z} + p\\right)}\\\\\n", + "v'_w = ov_z + p\n", + "$$" + ] + }, + { + "cell_type": "markdown", + "id": "2be99309-80e4-467d-9907-703b29cf974d", + "metadata": {}, + "source": [ + "This is the definition of $v'_w$:" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "id": "a3075494-d8af-4d1a-8e87-403c6c037526", + "metadata": {}, + "outputs": [ + { + "data": { + "text/markdown": [ + "$$v'_w = o v_{z} + p$$" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "display(Markdown(f\"$$v'_w = {latex(v1_w)}$$\"))" + ] + }, + { + "cell_type": "markdown", + "id": "15a6eac4-85a1-4592-aa4e-6b6251305511", + "metadata": {}, + "source": [ + "Solving for $v_z$ I get $v_z$ in terms of $v'_w$:" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "id": "a1bf3471-988a-4a08-816f-e2439c6711d0", + "metadata": {}, + "outputs": [ + { + "data": { + "text/markdown": [ + "$$v_z = \\frac{- p + v'_{w}}{o}$$" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "v1_w_tmp = Symbol(\"v'_w\")\n", + "v_z_tmp = solve(Eq(v1_w_tmp, v1_w), v_z)[0]\n", + "display(Markdown(f\"$$v_z = {latex(v_z_tmp)}$$\"))" + ] + }, + { + "cell_type": "markdown", + "id": "005d46a0-2e1f-41f4-b99c-4beb44ff2c24", + "metadata": {}, + "source": [ + "This is the definition of $v'''_z$:" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "id": "ad1b695a-5baf-48a1-a9ed-ca83a314b196", + "metadata": {}, + "outputs": [ + { + "data": { + "text/markdown": [ + "$$v'''_z = \\frac{k v_{z} + l + o v_{z} + p}{2 \\left(o v_{z} + p\\right)}$$" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "display(Markdown(f\"$$v'''_z = {latex(v3_z)}$$\"))" + ] + }, + { + "cell_type": "markdown", + "id": "63de3df7-4dbf-495a-a6d5-a4056ac6d736", + "metadata": {}, + "source": [ + "I will substitute my new $v_z$ into $v'''_z$:" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "id": "cfa19efc-f0d6-4cd2-9ac5-774e84d84e1e", + "metadata": {}, + "outputs": [ + { + "data": { + "text/markdown": [ + "$$v'''_z = \\frac{- k \\left(p - v'_{w}\\right) + o \\left(l + v'_{w}\\right)}{2 o v'_{w}}$$" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "v3_z_tmp = simplify(v3_z.subs(v_z, v_z_tmp))\n", + "display(Markdown(f\"$$v'''_z = {latex(v3_z_tmp)}$$\"))" + ] + }, + { + "cell_type": "markdown", + "id": "3d8f8ad3-c738-4af7-844c-508b5f24b223", + "metadata": {}, + "source": [ + "Now I can just solve for $v'_w$ in terms of $v'''_z$ (which I'll have access to after sampling):" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "id": "1192cb25-1030-41f6-96c1-d22e5312a9ce", + "metadata": {}, + "outputs": [ + { + "data": { + "text/markdown": [ + "$$v'_w = \\frac{k p - l o}{k - 2 o v'''_{z} + o}$$" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "v3_z_sym = Symbol(\"v'''_z\")\n", + "v1_w_tmp = solve(Eq(v3_z_sym, v3_z_tmp), v1_w_tmp)[0]\n", + "display(Markdown(f\"$$v'_w = {latex(v1_w_tmp)}$$\"))" + ] + }, + { + "cell_type": "markdown", + "id": "cbb75cd7-062e-4932-9891-4db98f7b2d74", + "metadata": {}, + "source": [ + "As proof that this is right, I should be able to substitute the original definition of $v'''_z$ back into this equation and find the original version of $v'_w$" + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "id": "ab2c3c2e-66f0-4191-b866-9bba2663cc97", + "metadata": {}, + "outputs": [ + { + "data": { + "text/markdown": [ + "$$v'_w = o v_{z} + p$$" + ], + "text/plain": [ + "" + ] + }, + "metadata": {}, + "output_type": "display_data" + } + ], + "source": [ + "display(Markdown(f\"$$v'_w = {latex(simplify(v1_w_tmp.subs(v3_z_sym, v3_z)))}$$\"))" + ] + }, + { + "cell_type": "markdown", + "id": "7a4782b3-87b8-4b49-a216-fadc4283d1a9", + "metadata": {}, + "source": [ + "That looks perfect. Now to write a glsl function to compute $v'_w$ given a mat4:\n", + "\n", + "```glsl\n", + "float get_w(float depth_sample, mat4 proj) {\n", + " float k = proj[2][2];\n", + " float l = proj[3][2];\n", + " float o = proj[2][3];\n", + " float p = proj[3][3];\n", + " return (k * p - l * o) / (k - 2.0 * o * depth_sample + o);\n", + "}\n", + "```" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "54011004-bb74-415a-a7e5-c5870c6d5ac7", + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.8.10" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/Perspective.ipynb b/Perspective.ipynb index 41e3cfa..024c961 100644 --- a/Perspective.ipynb +++ b/Perspective.ipynb @@ -595,6 +595,138 @@ "\n", "Two things to note in comparison: their $k$ is the negated version of mine, and their $o$ is -1. TODO: figure out exactly why that is (probably coordinate system handedness)" ] + }, + { + "cell_type": "markdown", + "id": "4714bb39-341c-4187-b01b-fc54925a5983", + "metadata": {}, + "source": [ + "## Going backwards\n", + "\n", + "Something else I need to know is how to go backwards from $v''_z$ to $v_z$. It's not as straight-forward as multiplying by the inverse because of the pesky w divide. In cases where I want to do this, I won't have access to $v'_w$ to undivide it.\n", + "\n", + "An important note: if I'm sampling depth from a depth buffer in OpenGL, I'll call it $v'''_z$, then $v'''_z$ = $\\frac{v''_z+1}{2}$. Conversely that means $v''_z = 2v'''_z - 1$.\n", + "\n", + "I can say:\n", + "\n", + "$$\n", + "2v'''_z - 1 = \\frac{iv_x + jv_y + kv_z + l}{mv_x + nv_y + ov_z + p}\n", + "$$\n", + "\n", + "I wonder if I can just solve for $v_z$:" + ] + }, + { + "cell_type": "code", + "execution_count": 17, + "id": "acfba905-b876-464d-930e-c5d24c35d506", + "metadata": {}, + "outputs": [ + { + "data": { + "text/latex": [ + "$\\displaystyle v_{z} = \\frac{- i v_{x} - j v_{y} - l + 2 m v'''_{z} v_{x} - m v_{x} + 2 n v'''_{z} v_{y} - n v_{y} + 2 p v'''_{z} - p}{k - 2 o v'''_{z} + o}$" + ], + "text/plain": [ + " -i⋅vₓ - j⋅v_y - l + 2⋅m⋅v'''_z⋅vₓ - m⋅vₓ + 2⋅n⋅v'''_z⋅v_y - n⋅v_y + 2⋅p⋅\n", + "v_z = ────────────────────────────────────────────────────────────────────────\n", + " k - 2⋅o⋅v'''_z + o \n", + "\n", + "v'''_z - p\n", + "──────────\n", + " " + ] + }, + "execution_count": 17, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "i, j, k, l, m, n, o, p, v_x, v_y, v3_z = symbols('i, j, k, l, m, n, o, p, v_x, v_y, v\\'\\'\\'_z')\n", + "Eq(v_z, solve(Eq(\n", + " 2*v3_z - 1,\n", + " (i*v_x + j*v_y + k*v_z + l) / (m*v_x + n*v_y + o*v_z + p)\n", + "), v_z)[0])" + ] + }, + { + "cell_type": "markdown", + "id": "3adef650-c33a-4c51-ac6e-8c9a428b08e7", + "metadata": {}, + "source": [ + "That's pretty gnarly.\n", + "\n", + "I'm going to avoid making assumptions about the values of $o$ and $p$ because an orthographic projection matrix won't have the same values here as the perspective matrix. I *am* going to make the simplifying assumption that $i$, $j$, $m$ & $n$ are all zero and try again." + ] + }, + { + "cell_type": "code", + "execution_count": 27, + "id": "0895f989-5472-4a2b-8219-d4a3de365ff4", + "metadata": {}, + "outputs": [ + { + "data": { + "text/latex": [ + "$\\displaystyle v_{z} = \\frac{- l + 2 p v'''_{z} - p}{k - 2 o v'''_{z} + o}$" + ], + "text/plain": [ + " -l + 2⋅p⋅v'''_z - p\n", + "v_z = ───────────────────\n", + " k - 2⋅o⋅v'''_z + o" + ] + }, + "execution_count": 27, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "i = j = m = n = 0\n", + "expr = solve(Eq(\n", + " 2*v3_z - 1,\n", + " (i*v_x + j*v_y + k*v_z + l) / (m*v_x + n*v_y + o*v_z + p)\n", + "), v_z)[0]\n", + "Eq(v_z, expr)" + ] + }, + { + "cell_type": "code", + "execution_count": 31, + "id": "2963d3bc-dc5f-4786-a530-b05b8d3a7158", + "metadata": {}, + "outputs": [ + { + "data": { + "text/latex": [ + "$\\displaystyle \\left( \\frac{v'''_{z} \\left(- k t - l + 2 o t v'''_{z} - o t + 2 p v'''_{z} - p\\right)}{k - 2 o v'''_{z} + o}, \\ t\\right)$" + ], + "text/plain": [ + "⎛v'''_z⋅(-k⋅t - l + 2⋅o⋅t⋅v'''_z - o⋅t + 2⋅p⋅v'''_z - p) ⎞\n", + "⎜───────────────────────────────────────────────────────, t⎟\n", + "⎝ k - 2⋅o⋅v'''_z + o ⎠" + ] + }, + "execution_count": 31, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "solve(Eq(\n", + " s*v3_z**(-1) + t,\n", + " expr\n", + "), s, t)[0]" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "fa126eba-d46f-42de-9ecb-6f9893b92819", + "metadata": {}, + "outputs": [], + "source": [] } ], "metadata": {