# Physics Based Differentiable Rendering: Edge Sampling

## Introduction

Edge sampling is a method to calculate the derivative of the ray tracing result w.r.t. some scene parameters (including camera pose, scene geometry, material and light parameters).

The key idea of edge sampling is dividing the gradient integral into smooth (interior) and discontinuous (boundary) regions. For the smooth part, we use automatic differentiation. For the discontinuous, we use edge sampling to capture the changes at boundaries.

Condition:

- Focus on triangle meshes
- Assume the meshes have been preprocessed such that there is no interpenetration

- Assume no point light sources and no perfectly specular surfaces
- Approximate with area light sources and BRDFs with very low roughness

- Focus on static scenes

## Primary Visibility

Consider the 2D pixel filter integral for each pixel that integrates over the pixel filter $k$ and the radiance $L$, where the radiance itself can be another integral that integrates over light sources or the hemisphere.

The pixel color $I$ can be written as:

$$

I=\iint k(x,y)L(x,y)\mathrm dx\mathrm dy

$$

For notational convenience we will combine the pixel filter andradiance and call them scene function $f(x,y)=k(x,y)L(x,y)$. We are interested in the gradients of the integral with respect to some parameters $\Phi$ in the scene function $f(x,y;\Phi)$, such as the position of a mesh vertex:

$$

\nabla I=\nabla\iint f(x,y;\Phi)\mathrm dx\mathrm dy

$$

The integral usually does not have a closed-form solution, especially when more complex effects such as non-Lambertian BRDFs are involved. Therefore we rely on Monte Carlo integration to estimate the pixel value $I$. However, we cannot take the naive approachof applying the same Monte Carlo sampler to estimate the gradient $\nabla I$, since the scene function $f$ is not necessarily differentiable with respect to the scene parameters.

A key observation is that all the discontinuities happen at triangle edges. This allows us to explicitly integrate over the discontinuities. A 2D triangle edge splits the space into two half-spaces ($f_u$ and $f_l$ in above figure). We can model it as a Heaviside step function $\theta$:

$$

\theta(\alpha(x,y))f_u(x,y)+\theta(-\alpha(x,y))f_l(x,y)

$$

- $f_u$ represents the upper half-space
- $f_l$ represents the lower half-space
- $\alpha$ defines the edge equation formed by the triangles

For each edge with two endpoints $(a_x,a_y)$, $(b_x,b_y)$, we can construct the edge equation by forming the line $\alpha(x,y)=Ax+By+C$. If $\alpha(x,y)>0$, then the point is at upper half-space, and vice versa. For the two endpoints of the edge $\alpha(x,y)=0$. Thus by plugging inthe two endpoints we obtain:

$$

\alpha(x,y)=(a_y-b_y)x+(b_x-a_x)y+(a_xb_y-b_xa_y)

$$

We can rewrite the scene functionfas a summation of Heaviside step functions $\theta$ with edge equation $\alpha_i$ multiplied by an arbitrary function $f_i$:

$$

\iint f(x,y)\mathrm dx\mathrm dy=\sum_i\iint\theta(\alpha_i(x,y))f_i(x,y)\mathrm dx\mathrm dy

$$

- $f_i$ itself can contain Heaviside step functions, for example a triangle defines a multiplication of three Heaviside step functions
- $f_i$ caneven be an integral over light sources or the hemisphere

We want to analytically differentiate the Heaviside step function $\theta$ and explicitly integrate over its derivative – the Dirac delta function $\delta$. To do this we first swap the gradient operator inside the integral, then we use product rule to separate the integral into two:

$$

\begin{aligned}

&\nabla \iint\theta(\alpha_i(x,y))f_i(x,y)\mathrm dx\mathrm dy\\\\

&=\iint\delta(\alpha_i(x,y))\nabla\alpha_i(x,y)f_i(x,y)\mathrm dx\mathrm dy\\\\

&+\iint \nabla f_i(x,y)\theta(\alpha_i(x,y))\mathrm dx\mathrm dy

\end{aligned}

$$

Above equation shows that we can estimate the gradient using two Monte Carlo estimators. The first one estimates the integral over the edges of triangles containing the Dirac delta functions, and the second estimates the original pixel integral except the smooth function $f_i$ is replaced by its gradient, which can be computed through automatic differentiation.

To estimate the integral containing Dirac delta functions, we eliminate the Dirac function by performing variable substitution to rewrite the first term containing the Dirac delta function to an integral that integrates over the edge, that is, over the regions where $\alpha_i(x,y)=0$:

$$

\begin{aligned}

&\iint\delta(\alpha_i(x,y))\nabla\alpha_i(x,y)f_i(x,y)\mathrm dx\mathrm dy\\\\

=&\int_{\alpha_i(x,y)=0}\frac{\nabla\alpha_i(x,y)}{\|\nabla_{x,y}\alpha_i(x,y)\|}f_i(x,y)\mathrm d\sigma(x,y)

\end{aligned}

$$

- $\|\nabla_{x,y}\alpha_i(x,y)\|$ is the $L^2$ length of the gradient of the edge equation $\alpha_i$ with respect to $x$, $y$, which takes the Jacobian of the variable substitution into account.
- $\sigma(x,y)$ is the measure of the length on the edge.

The gradients of the edge equations $\alpha_i$ are:

$$

\begin{aligned}

&\|\nabla_{x,y}\alpha_i\|=\sqrt{(a_x-b_x)^2+(a_y-b_y)^2}\\\\

&\dfrac{\partial \alpha_i}{\partial a_x}=b_y-y\quad \dfrac{\partial \alpha_i}{\partial a_y}=x-b_x\\\\

&\dfrac{\partial \alpha_i}{\partial b_x}=y-a_y\quad \dfrac{\partial \alpha_i}{\partial b_y}=a_x-x\\\\

&\dfrac{\partial \alpha_i}{\partial x}=a_y-b_y\quad \dfrac{\partial \alpha_i}{\partial y}=b_x-a_x

\end{aligned}

$$

We can obtain the gradient with respect to other parameters, such as camera parameters, 3D vertex positions, or vertex normals bypropagating the derivatives from the projected triangle vertices using the chain rule:

$$

\frac{\partial \alpha}{\partial p}=\sum_{k\in{x,y}}\frac{\partial\alpha}{\partial a_k}\frac{\partial a_k}{\partial p}+\frac{\partial\alpha}{\partial b_k}\frac{\partial b_k}{\partial p}

$$

where $p$ is the desired parameter.

We use Monte Carlo sampling to estimate the Dirac integral. Recall that a triangle edge defines two half-spaces, therefore we need to compute the two values $f_l(x,y)$ and $f_u(x,y)$ on the edge. By combining above results, our Monte Carlo estimation of the Dirac integral for a single edge $E$ on a triangle can be written as:

$$

\frac{1}{N}\sum_{j=1}^N\frac{\|E\|\nabla \alpha_i(f_u(x_j,y_j)-f_l(x_j,y_j))}{P(E)\|\nabla_{x_j,y_j}\alpha_i(x_j,y_j)\|}

$$

- $\|E\|$ is the length of the edge
- $P(E)$ is the probability of selecting edge $E$

In practice, if we employ smooth shading, most of the triangleedges are in the continuous regions and the Dirac integral is zero. Only the silhouette have non-zero contribution to the gradients. We select the edges by projecting all triangle meshes to screen space and clip them against the camera frustrum. We select one silhouette edge with probability proportional to the screen space lengths. We then uniformly pick a pointon the selected edge.