in maths, shader, tutorial, Unity3D

Atmospheric Scattering Shader

Share Button

This tutorial finally concludes our journey to simulate Rayleigh Scattering for a planet’s atmosphere.

The next (and final) part will show how to change the shader to also include an additional type of scattering, known as Mie Scattering.

You can find all the post in this series here:

You can refer to the Atmospheric Scattering Cheatsheet for a complete reference of all the equations used.

You can download the Unity package for this tutorial at the bottom of the page.

Sampling the View Ray

Let’s recall the equation for the atmospheric scattering that we have recently derived:

    \[I= I_S  \sum_{P \in \overline{AB}} {S\left(\lambda, \theta, h\right)  T\left(\overline{CP}\right)  T\left(\overline{PA}\right) ds }\]

The amount of light that we receive is equal to the amount of light emitted from the sun, I_S, multiplied by the sum of the individual contributions from each point P in the segment \overline{AB}.

We could go one and implement that function directly in our shader. However, there are few optimisations that can be done. It was hinted, in a previous tutorial, that the expression could be simplified even further. The first step we can take is to decompose the scattering function into its two basic components:

    \[S \left(\lambda, \theta, h\right ) = \beta \left(\lambda, h \right ) \gamma\left(\theta\right) = \beta \left(\lambda\right )\rho\left(h\right) \gamma\left(\theta\right)\]

The phase function \gamma\left(\theta\right) and the scattering coefficient at sea level \beta \left(\lambda\right ) are constant with respects to the summation, since the angle \theta and the wavelength \lambda do not depend on the sampled point. Hence, they can be factored out:

    \[I = I_S  \, \beta \left(\lambda\right )  \gamma\left(\theta\right) \sum_{P \in \overline{AB}} {   T\left(\overline{CP}\right)  T\left(\overline{PA}\right) \rho\left(h\right) ds }\]

This new expression is mathematically equivalent to the previous one, but is more efficient to calculate since some of the most heavy parts have been taken out of the summation.

We are not ready to start implementing it. There are infinitely many points P we should take into account. A reasonable approximation to I is to split \overline{AB} into several smaller segments of length ds, and to accumulate the contribution of each individual segment. In doing so, we assume that each segment is small enough to have a constant density. Generally speaking, that is not the case, but if the ds is small enough, we can still reach a reasonably good approximation.

The number of segments in  \overline{AB} is referred to as the view samples, since all segments lie on the view ray. In the shader, this will be the _ViewSamples property. By having it as a property, it is accessible from the material inspector. This allows us to reduce the precision of the shader, in favour of its performance.

The following piece of code allows looping through all the segments in the atmosphere.

The variable time is used to keep track of how far we are from the origin point O, and is incremented by ds after each iteration.

Optical Depth PA

Each point along the view ray \overline{AB} brings its own contribution to the final colour of the pixel that we are drawing. That contribution is, mathematically speaking, the quantity inside the summation:

    \[I = I_S  \, \beta \left(\lambda\right )  \gamma\left(\theta\right) \sum_{P \in \overline{AB}}   \underset{\text{light contribution of}\,L\left(P\right)}{\underbrace{T\left(\overline{CP}\right)  T\left(\overline{PA}\right) \rho\left(h\right) ds}}\]

Like we did in the previous paragraph, let’s try to simplify it. We can expand the expression above even further, by replacing T with its actual definition:

    \[T\left(\overline{XY}\right) =\exp\left\{ - \beta\left(\lambda\right) D\left(\overline{XY}\right) \right\}\]

The product of the transmittance over \overline{CP} and \overline{PA} becomes:

    \[T\left(\overline{CP}\right) T\left(\overline{PA}\right)=\]

    \[=\underset{T\left(\overline{CP}\right) }{\underbrace{ \exp\left\{- \beta\left(\lambda\right) D\left(\overline{CP}\right)\right \} }} \, \underset{T\left(\overline{PA}\right) }{\underbrace{ \exp\left\{- \beta\left(\lambda\right) D\left(\overline{PA}\right) \right \} }}=\]

    \[= \exp\left\{- \beta\left(\lambda\right) \left( D\left(\overline{CP}\right) + D\left(\overline{PA}\right) \right) \right \}\]

The combined transmittance is modelled as an exponential decay with which coefficient is the sum of the optical depths over the path travelled by the light (\overline{CP} and \overline{PA}), multiplied by the scattering coefficient at sea level (\beta with h=0).

The first quantity that we start calculating is the optical depth for the segment \overline{PA}, which goes from the point of entry through the atmosphere, up to the point we are currently sampling in the for loop. Let’s recall the definition of optical depth:

    \[D\left( \overline{PA}\right)=\sum_{Q \in \overline{PA}} { \exp\left\{-\frac{h_Q}{H}\right\} } \, ds\]

If one had to implement this naively, we would create a function called opticalDepth that samples points between P and A in a loop. This is possible, but very inefficient. As a matter of fact, D\left( \overline{PA}\right) is the optical depth of the segment that we are already analysing the outermost for loop. We can save a lot of calculations if we calculate the optical depth the current segmenet centred at P ( opticalDepthSegment), and keep accumulating it in the for loop ( opticalDepthPA).

Light Sampling

If we look back at the expression for the light contribution of P, we see that the only quantity needed is the optical depth of the segment \overline{CP},

    \[L\left(P\right) = \underset{\text{combined transmittance}}{\underbrace{\exp\left\{- \beta\left(\lambda\right) \left( D\left(\overline{CP}\right) + D\left(\overline{PA}\right) \right) \right \}}} \, \underset{\text{optical depth of}\,ds}{\underbrace{\rho\left(h\right) ds}}\]

We will move the code that calculates the optical depth of the segment \overline{CP} into a function called lightSampling. The name comes from the light ray, which is the segment that starts at P and is directed towards the sun. We have called the point where it exits the atmosphere is C.

The  lightSampling function, however, will not just calculate the optical depth of \overline{CP}. So far we have only considered the contribution of the atmosphere, ignoring the role of the actual planet. Our equations do not account for the possibility that the light ray that goes from P towards the sun might hit the planet. If this happens, all calculations done so far must be discarded, as no light would actually reach the camera.

In the diagram above, it’s easy to see that the light contribution of P_0 should be ignored, since no sun light is reaching P_0. While looping through the points in between P and C, the  lightSampling function will also check if the planet has been hit. This can be done by checking if the altitude of a point is negative.

The function above first calculates the point C using rayInstersect. It then divides the segment \overline{PA} into _LightSamples segments of length ds. The calculation for the optical depth is the same used in the outermost loop.

The function returns false if the planet has been hit. We can use this to update the missing code the outermost loop, replacing the  ....

Now that we have taken into account all the elements, our shader is complete.

Coming Next…

This post (finally!) completes the volumetric shader that simulates atmospheric scattering. So far, we have only taken into account the contribution of the Rayleigh scattering. There are many optical phenomena that cannot be explained with Rayleigh scattering alone. The next post will introduce a second type of scattering, which is known as Mie scattering.

You can find all the post in this series here:

You can refer to the Atmospheric Scattering Cheatsheet for a complete reference of all the equations used.


You can download all the assets necessary to reproduce the volumetric atmospheric scattering presented in this tutorial.

Feature Standard Premium
Volumetric Shader
Clouds Support
Day & Night Cycle
Earth 8K
Mars 8K
Venus 8K
Neptune 2K
Download Standard Premium
📧 Stay updated

A new tutorial is released every week.

💖 Support this blog

This websites exists thanks to the contribution of patrons on Patreon. If you think these posts have either helped or inspired you, please consider supporting this blog.


Write a Comment