**Nonlinear-density Halfspace Fog**
by James W. Walker
Original 25 May 2016
Revised 20 July 2023
# Review of Fixed-function OpenGL Fog
In the original OpenGL, fog computed a new color for a fragment by a weighted average
$$newColor = f * oldColor + (1 - f) * fogColor$$
where $f$ is a function of the distance $z$ from the camera to the fragment and depends on the
fog mode. Notice that the closer $f$ is to zero, the more the fog color replaces the
original color, so we should think of this as the transparency (not opacity) of the fog.
The **linear** fog mode defines
$$f(z) = \begin{cases}
1 & \text{if $z < start$}\\
\dfrac{end - z}{end - start} & \text{if $start \le z \le end$}\\
0 & \text{if $z > end$}
\end{cases}
$$
where $start$ and $end$ are parameters, specific distances where the fog starts and ends.
The **exponential** fog mode defines
$$f(z) = \Large e^{ - density \,* \,z }$$
where $density$ is a positive constant parameter.
Finally, **exponential squared** fog defines
$$f(z) = \Large e^{ - (density \,* \,z)^2 }\,.$$
# Hand-waving Motivation for Exponential Fog
Without doing any physics, let's see if we can motivate at least one of our fog formulae.
Notice that the general fog blending formula has the same form as the formula for alpha
blending. Hence, if we wanted to have a fixed amount of fogginess that does not depend
on distance from the camera, we could imagine placing a translucent filter of transparency
$f$ in front of the camera. If we wanted to have fogginess that starts to depend on
distance, we could use two such filters, at different distances. Fragments beyond both
filters would end up with the color
$$\begin{aligned}
newColor &= f * (f * oldColor + (1-f) * fogColor) + (1-f) * fogColor\\
&= f^2 * oldColor + (1 - f^2) * fogColor
\end{aligned}.
$$
More generally, if you used a number $n$ of identical filters at various distances, the
maximum combined fogginess would be:
$$newColor = f^n * oldColor + (1 - f^n) * fogColor$$
Now let's suppose we use lots of identical translucent filters, equally spaced at some
small distance $1/n$ where $n$ is large. The number of filters between the camera and a
fragment at distance z would be about $n * z$, so the fog transparency at that fragment
would be
$$f^{n*z} = (f^n)^z\,.$$
If you were going to use more and more filters to get progressively more smoothly
graduated fog, you wouldn't want to keep using filters of the same transparency, you
would want to increase the transparency so as to keep the cumulative transparency
approximately fixed at each distance. That is, vary $f$ so that $f^n$ is a fixed value $k$.
Therefore, we end up with a fog transparency function
$$f(z) = k^z$$
for some number $k$ between 0 and 1. We can rewrite this as
$$f(z) = \Large e^{\ln{k^z}} = e^{z \,* \,\ln{k}} \,.$$
Since we know that $k$ is between 0 and 1, $\ln{k}$ is between $-\infty$ and zero. Therefore,
if we let $d = -\ln{k}$, we will have
$$f(z) = \Large e^{-d * z}$$
for a positive parameter $d$, and this is the formula for exponential fog.
# Variable Density Exponential Fog
Suppose that as light travels from a fragment to the camera, it passes through volumes
with distinct densities $d_i$. The fog transparencies of these pieces would multiply,
so if we use exponential fog for each volume, the cumulative fog transparency would be
$$\Large \prod_i e^{-d_i * \Delta z} = e^{- \sum_i d_i * \Delta z}\,.$$
If we assume that density is a function of position, then passing to the limit will give
a fog transparency of
$$f(z) = \Large e^{-\int_0^z d(t)\, dt}\,.$$
Note that if the density function is actually a constant, then this boils back down to
the old exponential fog formula.
# Halfspace Fog
**Halfspace fog** is fog that exists only on one side of a given plane, and whose density
is a function of the distance from the fragment to the plane. Most often, the fog would
be above or below a certain height.
The discussion of halfspace
fog by Eric Lengyel[^1] which inspired this discussion uses a constant or linear density
function, but here we will extend that to a nonlinear function.
Let $\mathbf{F}$ be a 4-dimensional
vector describing a halfspace, such that a homogeneous point \(\mathbf{X}\) is on the
foggy side of the plane when \(\mathbf{F}•\mathbf{X} > 0\). We will further assume that
the \(xyz\) part of \(\mathbf{F}\) is normalized, so that \(\mathbf{F}•\mathbf{X}\) is
the signed distance between \(\mathbf{X}\) and the boundary plane. Let \(\mathbf{C}\)
be the camera position, let \(\mathbf{P}\) be the location of a fragment being rendered,
and let \(\mathbf{V} = \mathbf{C} - \mathbf{P}\), the vector from \(\mathbf{P}\) to
\(\mathbf{C}\). We can write a parametric equation of
the line segment from \(\mathbf{P}\) to
\(\mathbf{C}\) as \(\mathbf{Q}(t) = \mathbf{P} + t\,\mathbf{V}\) where \(t\) varies from
0 to 1. If density is a function $\delta(x)$ of distance from the plane, let's see what
we can say
about the density integrated along such a line segment. Let's represent that integrated
density at $\mathbf{P}$ by $g(\mathbf{P})$. Bear in mind that we will ultimately use the
fog transparency
$$f(\mathbf{P}) = \Large e^{- g(\mathbf{P})}\,.$$
## On the Foggy Side
Consider the case where \(\mathbf{P}\) and \(\mathbf{C}\) are both on the foggy side of
the plane, i.e., \(\mathbf{F}•\mathbf{P} > 0\) and \(\mathbf{F}•\mathbf{C} > 0\). Of
course that implies that the entire line segment from \(\mathbf{P}\) to \(\mathbf{C}\)
is on the foggy side.
$$g(\mathbf{P}) = \int_0^1 \delta(\mathbf{F}•\mathbf{Q}(t))\,
\left\|\mathbf{V}\right\|\, dt =
\left\|\mathbf{V}\right\| \int_0^1 \delta(\mathbf{F}•\mathbf{Q}(t))\, dt\,.
$$
Note that \(\mathbf{F}•\mathbf{Q}(t) = \mathbf{F}•\mathbf{P} + t\,\mathbf{F}•\mathbf{V}\).
For convenience, let \(p = \mathbf{F}•\mathbf{P}\) and \(v = \mathbf{F}•\mathbf{V}\), so
\(\mathbf{F}•\mathbf{Q}(t) = p + t\,v\). Thus, we need to evaluate
$$g(\mathbf{P}) = \left\|\mathbf{V}\right\|\int_{0}^{1}\delta(p + t\,v)\,dt\;.$$
Let's make the substitution \(u = p + t\,v\), so this becomes
$$g(\mathbf{P}) = \left\|\mathbf{V}\right\|\int_{p}^{p+v}\delta(u)\,\frac{1}{v}\,du =
\frac{\left\|\mathbf{V}\right\|}{v}\int_{p}^{p+v}\delta(u)\,du\;.$$
Define $c = p+v = \mathbf{F}•\mathbf{C}$, so this becomes
$$g(\mathbf{P}) =
\frac{\left\|\mathbf{V}\right\|}{c-p}\int_{p}^c\delta(u)\,du\;.$$
Now, let's suppose that our density function $\delta$ has an antiderivative $\Delta$. Our
integrated density can then be expressed as
$$g(\mathbf{P}) = \left\|\mathbf{V}\right\|
\frac{\Delta(c) - \Delta(p)}{c-p}\;.$$
## Avoiding Division by Zero
If it happens that $p = c$, then this formula
becomes zero divided by zero. That happens when $v=0$, i.e., $\mathbf{F}•\mathbf{V}=0$,
meaning that $\mathbf{P}$ and $\mathbf{C}$ are equidistant from the plane. However, our
formula does have a finite limit in that case. In fact, we can say
$$\lim_{p\to c}g(\mathbf{P}) =
\lim_{p\to c}\left\|\mathbf{V}\right\|\frac{\Delta(c) - \Delta(p)}{c-p} =
\left\|\mathbf{V}\right\| \lim_{p\to c}\frac{\Delta(c) - \Delta(p)}{c-p} =
\left\|\mathbf{V}\right\| \Delta^\prime(c) = \left\|\mathbf{V}\right\| \delta(c)\,.$$
That is of course exactly what we should expect: Integrating a constant density along a
line segment produces the density times the length of the segment.
## Crossing the Fog Plane
If \(\mathbf{P}\) and \(\mathbf{C}\) are both on the non-foggy side of
the plane, then of course the line segment encounters no fog and \(g(\mathbf{P}) = 0\).
It remains to consider the cases where \(\mathbf{P}\) and \(\mathbf{C}\) are on different
sides of the plane. In either case, the line segment encounters the plane where
\(p + t\,v = 0\), so \(t = -p/v\). One can verify that \(-p/v\) is between 0 and 1.
This parameter will take the place of one of the limits of integration.
If \(\mathbf{F}•\mathbf{P} > 0\) and \(\mathbf{F}•\mathbf{C} < 0\), the line starts in
the fog, so \(t\) ranges from 0 to \(-p/v\). As a result, the upper limit on \(u\)
becomes 0. The integral reduces to
$$g(\mathbf{P}) =
\frac{\left\|\mathbf{V}\right\|}{c-p}\int_{p}^{0}\delta(w)\,dw
= \left\|\mathbf{V}\right\| \frac{\Delta(0) - \Delta(p)}{c-p}\;.$$
Similarly, if \(\mathbf{F}•\mathbf{P} < 0\) and \(\mathbf{F}•\mathbf{C} > 0\), the lower
limit on \(u\) becomes 0, and the integral becomes
$$g(\mathbf{P}) =
\frac{\left\|\mathbf{V}\right\|}{c-p}\int_{0}^{c}\delta(w)\,dw
= \left\|\mathbf{V}\right\| \frac{\Delta(c) - \Delta(0)}{c-p}\;.$$
In these cases, there can be no division by zero, since \(c\) and \(p\) are of opposite
sign.
## Unified Formula
Let's define \(\bar{c} = \max(c, 0)\) and \(\bar{p} = \max(p, 0)\). Then we can handle
all cases as follows:
$$g(\mathbf{P}) = \begin{cases}
\left\|\mathbf{V}\right\|
\dfrac{\Delta(\bar{c}) - \Delta(\bar{p})}{c-p} & \text{if $c \ne p$}\\
\left\|\mathbf{V}\right\|\,\delta(\bar{c}) & \text{if $c = p$}
\end{cases}
\;,$$
assuming that $\delta(0) = 0$.
## Fog at an Infinite Point
If our rendering method specifies an infinite far distance, and \(\mathbf{P}\) is a point
at infinite distance, then our formula becomes undefined due to infinities. For practical
purposes, one should probably just replace the infinite point with a large but finite
point in the same direction.
## A Specific Halfspace Density Function
I want to choose a density function that starts with value zero and derivative $S$
(for sharpness) at the boundary plane, and smoothly approaches a fixed value $D$ as you
get farther from the plane. I also want the function to have an antiderivative that can
be easily computed in a shader program. There are many such functions. Here is one,
$$\delta(x) = D\,(1-{\large e^{-\frac{S}{D}x}})$$
with the antiderivative
$$\Delta(x) = D\left(x + \frac{D}{S}\,{\large e^{-\frac{S}{D}\,x}}\right)\,.$$
[^1] Lengyel, Eric, "Unified Distance Formulas for Halfspace Fog", *Journal of Graphics
Tools,* vol. 12, number 2, pp 23-32, 2007.