While re-deriving some equations relevant to Monte-Carlo path tracing (specifically, the probability distribution of sampling a specific light direction from Sampling the GGX Distribution of Visible Normals), I stumbled upon the need to compute the determinant of the Jacobian of the operator that reflects a vector $V$ from the surface orthogonal to a vector $H$, with respect to the vector $H$:
$$L(H,V) = 2(H\cdot V)H – V$$
$$\det \frac{\partial L}{\partial H} = ?$$
as a function on the unit sphere.
To make it more malleable, I rewrote it as
$$L(H,V) = 2(H\cdot V)H – (H\cdot H)V$$
which is the same function if $\|H\|=1$, but is also a homogeneous function of $H$ of order 2, which means that it is a product of radial and spherical components, and thus the determinant of the spherical part of the Jacobian is just the determinant of this function viewed as
$$L(\cdot,V) : \mathbb R^3\rightarrow \mathbb R^3$$
divided by the derivative of the radial part (which is the derivative of $R^2$ evaluated at $R=1$, which is simply $2$).
Now, to compute the 3D Jacobian, I transformed the function into matrix form
$$L(H,V) = 2HH^TV – VH^TH$$
and applied the usual matrix calculus to obtain
$$\frac{\partial L}{\partial H} = 2(H^TV\cdot I+HV^T-VH^T)$$
This matrix looks quite special: neglecting the factor of 2, its diagonal elements are just $H\cdot V$, while the off-diagonal part is antisymmetric (and seems to be a matrix form of the bivector $H\wedge V$).
However, computing its determinant by hand proved to be quite tedious. Using sympy I was able to obtain the formula
$$\det \frac{\partial L}{\partial H} = 8\cdot H^2V^2(H\cdot V) \quad\quad(1)$$
After substituting $\|H\|=\|V\|=1$ and dividing by 2 (the radial part of the Jacobian), we get $4(H\cdot V)$, which coincides with the result from the paper linked above.
Now, my question is: how can one obtain equation (1) without fully expanding the Jacobian in coordinates?
My thoughts:
-
It seems that one can treat the $HV^T-VH^T$ as two rank-one updates for the matrix $H^T V\cdot I$, and use the rank-one determinant update formula twice to compute the determinant. However, given that the problem is so geometric in nature, is doesn't feel like the right approach. Also, since the rank-one update formula requires the adjugate of the updated matrix, this hardly seems easier than computing the determinant entirely from scratch.
-
Using the Clifford algebra $\operatorname{Cl}(\mathbb R^3)$ we can express our reflection using the Clifford product simply as $L=H\cdot V\cdot H$, and the derivative (as a linear operator on the tangent space of the algebra) is $\frac{\partial L}{\partial H}(\delta) = H\cdot V\cdot \delta + \delta \cdot V \cdot H$. The Clifford algebra formulation seems to be somehow related to the $H\cdot V + H\wedge V$ expansion we've seen earlier. However, I don't know how this can help in computing the determinant.
Best Answer
Define $$ A = (H \cdot V)I + H V^T - V H^T. $$ Let $H,V$ be unit vectors. Let $Q=H \times V$ be a vector orthogonal to them. Consider $A$ in the $H,V,Q$ frame: $$ A H = 2(H \cdot V)H - V ,$$ $$ A V = H, $$ $$ A Q = (H \cdot V)Q. $$ It is now clear that $A$ only performs a scaling in the $Q$ direction and some rotation in the $H,V$ plane. Thus, $$\det A = H \cdot V.$$