1
0
Fork 0

HDR pipeline: even better atmospheric scattering

Implementation of "A Scalable and Production Ready Sky and Atmosphere Rendering Technique" by Sébastien Hillaire (2020).
This commit is contained in:
Fernando García Liñán 2021-07-28 09:40:04 +02:00
parent 016585e176
commit 9b9ae5cf38
16 changed files with 647 additions and 238 deletions

View file

@ -65,6 +65,28 @@
</buffer>
<!-- Atmosphere LUTs -->
<buffer>
<name>transmittance</name>
<type>2d</type>
<width>256</width>
<height>64</height>
<format>rgb16f</format>
<min-filter>linear</min-filter>
<mag-filter>linear</mag-filter>
<wrap-s>clamp-to-edge</wrap-s>
<wrap-t>clamp-to-edge</wrap-t>
</buffer>
<buffer>
<name>multiple-scattering</name>
<type>2d</type>
<width>32</width>
<height>32</height>
<format>rgb16f</format>
<min-filter>linear</min-filter>
<mag-filter>linear</mag-filter>
<wrap-s>clamp-to-edge</wrap-s>
<wrap-t>clamp-to-edge</wrap-t>
</buffer>
<buffer>
<name>sky-view</name>
<type>2d</type>
@ -77,18 +99,7 @@
<wrap-t>clamp-to-edge</wrap-t>
</buffer>
<buffer>
<name>aerial-inscatter</name>
<type>2d</type>
<width>512</width>
<height>32</height>
<format>rgba16f</format>
<min-filter>linear</min-filter>
<mag-filter>linear</mag-filter>
<wrap-s>clamp-to-edge</wrap-s>
<wrap-t>clamp-to-edge</wrap-t>
</buffer>
<buffer>
<name>aerial-transmittance</name>
<name>aerial-perspective</name>
<type>2d</type>
<width>512</width>
<height>32</height>
@ -213,28 +224,70 @@
<!--=======================================================================-->
<!-- Sky-View -->
<!--
Generate the atmospheric scattering related LUTs
This is an implementation of Sébastien Hillaire's "A Scalable and
Production Ready Sky and Atmosphere Rendering Technique".
As a small summary, the sky-view and aerial perspective LUTs are used by
later passes to apply atmospheric scattering to the skydome and
opaque/transparent objects respectively. These two LUTs also rely on
another pair of LUTs (transmittance and multiple scattering) to speed up
the computations.
-->
<pass>
<name>atmos-transmittance</name>
<type>quad</type>
<effect>Effects/HDR/atmos-transmittance</effect>
<attachment>
<component>color</component>
<buffer>transmittance</buffer>
</attachment>
</pass>
<pass>
<name>atmos-multiple-scattering</name>
<type>quad</type>
<effect>Effects/HDR/atmos-multiple-scattering</effect>
<binding>
<unit>0</unit>
<buffer>transmittance</buffer>
</binding>
<attachment>
<component>color</component>
<buffer>multiple-scattering</buffer>
</attachment>
</pass>
<pass>
<name>atmos-sky-view</name>
<type>quad</type>
<effect>Effects/HDR/atmos-sky-view</effect>
<binding>
<unit>0</unit>
<buffer>transmittance</buffer>
</binding>
<binding>
<unit>1</unit>
<buffer>multiple-scattering</buffer>
</binding>
<attachment>
<component>color</component>
<buffer>sky-view</buffer>
</attachment>
</pass>
<!-- Aerial Perspective -->
<pass>
<name>atmos-aerial-perspective</name>
<type>quad</type>
<effect>Effects/HDR/atmos-aerial-perspective</effect>
<binding>
<unit>0</unit>
<buffer>transmittance</buffer>
</binding>
<binding>
<unit>1</unit>
<buffer>multiple-scattering</buffer>
</binding>
<attachment>
<component>color0</component>
<buffer>aerial-inscatter</buffer>
</attachment>
<attachment>
<component>color1</component>
<buffer>aerial-transmittance</buffer>
<component>color</component>
<buffer>aerial-perspective</buffer>
</attachment>
</pass>
@ -617,11 +670,7 @@
</binding>
<binding>
<unit>11</unit>
<buffer>aerial-inscatter</buffer>
</binding>
<binding>
<unit>12</unit>
<buffer>aerial-transmittance</buffer>
<buffer>aerial-perspective</buffer>
</binding>
<attachment>
<component>color0</component>
@ -648,6 +697,10 @@
<unit>11</unit>
<buffer>sky-view</buffer>
</binding>
<binding>
<unit>13</unit>
<buffer>transmittance</buffer>
</binding>
<attachment>
<component>color</component>
<buffer>hdr-result</buffer>

View file

@ -6,9 +6,19 @@
<program>
<vertex-shader>Shaders/HDR/trivial.vert</vertex-shader>
<fragment-shader>Shaders/HDR/atmos-aerial-perspective.frag</fragment-shader>
<fragment-shader>Shaders/HDR/atmosphere-include.frag</fragment-shader>
<fragment-shader>Shaders/HDR/atmos-include.frag</fragment-shader>
<fragment-shader>Shaders/HDR/gbuffer-include.frag</fragment-shader>
</program>
<uniform>
<name>transmittance_lut</name>
<type>sampler-2d</type>
<value type="int">0</value>
</uniform>
<uniform>
<name>multiscattering_lut</name>
<type>sampler-2d</type>
<value type="int">1</value>
</uniform>
</pass>
</technique>
</PropertyList>

View file

@ -0,0 +1,18 @@
<?xml version="1.0" encoding="utf-8"?>
<PropertyList>
<name>Effects/HDR/atmos-multiple-scattering</name>
<technique n="1">
<pass>
<program>
<vertex-shader>Shaders/HDR/trivial.vert</vertex-shader>
<fragment-shader>Shaders/HDR/atmos-multiple-scattering.frag</fragment-shader>
<fragment-shader>Shaders/HDR/atmos-include.frag</fragment-shader>
</program>
<uniform>
<name>transmittance_lut</name>
<type>sampler-2d</type>
<value type="int">0</value>
</uniform>
</pass>
</technique>
</PropertyList>

View file

@ -6,8 +6,18 @@
<program>
<vertex-shader>Shaders/HDR/trivial.vert</vertex-shader>
<fragment-shader>Shaders/HDR/atmos-sky-view.frag</fragment-shader>
<fragment-shader>Shaders/HDR/atmosphere-include.frag</fragment-shader>
<fragment-shader>Shaders/HDR/atmos-include.frag</fragment-shader>
</program>
<uniform>
<name>transmittance_lut</name>
<type>sampler-2d</type>
<value type="int">0</value>
</uniform>
<uniform>
<name>multiscattering_lut</name>
<type>sampler-2d</type>
<value type="int">1</value>
</uniform>
</pass>
</technique>
</PropertyList>

View file

@ -0,0 +1,13 @@
<?xml version="1.0" encoding="utf-8"?>
<PropertyList>
<name>Effects/HDR/atmos-transmittance</name>
<technique n="1">
<pass>
<program>
<vertex-shader>Shaders/HDR/trivial.vert</vertex-shader>
<fragment-shader>Shaders/HDR/atmos-transmittance.frag</fragment-shader>
<fragment-shader>Shaders/HDR/atmos-include.frag</fragment-shader>
</program>
</pass>
</technique>
</PropertyList>

View file

@ -68,15 +68,10 @@
<value type="int">10</value>
</uniform>
<uniform>
<name>aerial_inscatter_lut</name>
<name>aerial_perspective_lut</name>
<type>sampler-2d</type>
<value type="int">11</value>
</uniform>
<uniform>
<name>aerial_transmittance_lut</name>
<type>sampler-2d</type>
<value type="int">12</value>
</uniform>
</pass>
</technique>
</PropertyList>

View file

@ -319,6 +319,11 @@
<type>sampler-2d</type>
<value type="int">11</value>
</uniform>
<uniform>
<name>transmittance_lut</name>
<type>sampler-2d</type>
<value type="int">13</value>
</uniform>
</pass>
</technique>
<technique n="139">
@ -340,6 +345,11 @@
<type>sampler-2d</type>
<value type="int">11</value>
</uniform>
<uniform>
<name>transmittance_lut</name>
<type>sampler-2d</type>
<value type="int">13</value>
</uniform>
</pass>
</technique>

View file

@ -1,7 +1,18 @@
// An implementation of Sébastien Hillaire's "A Scalable and Production Ready
// Sky and Atmosphere Rendering Technique".
//
// This shader generates the aerial perspective LUT. This LUT is used by opaque
// and transparent objects to apply atmospheric scattering. In-scattering is
// stored in the RGB channels, while transmittance is stored in the alpha
// channel.
// Unlike the paper, we are using a tiled 2D texture instead of a true 3D
// texture. For some reason the overhead of rendering to a texture a lot of
// times (the depth of the 3D texture) seems to be too high, probably because
// OSG is not sharing state between those passes.
#version 330 core
layout(location = 0) out vec3 out_inscatter;
layout(location = 1) out vec3 out_transmittance;
out vec4 fragColor;
in vec2 texCoord;
@ -9,39 +20,111 @@ uniform mat4 fg_ViewMatrixInverse;
uniform vec3 fg_CameraPositionCart;
uniform vec3 fg_CameraPositionGeod;
uniform vec3 fg_SunDirectionWorld;
uniform vec2 fg_NearFar;
uniform sampler2D transmittance_lut;
uniform sampler2D multiscattering_lut;
const float PI = 3.141592653;
const float ATMOSPHERE_RADIUS = 6471e3;
const float TOTAL_SLICES = 16.0;
const float DEPTH_RANGE = 32000.0;
const int AERIAL_PERSPECTIVE_SAMPLES = 20;
const vec3 ONE_OVER_THREE = vec3(1.0 / 3.0);
vec3 positionFromDepth(vec2 pos, float depth);
void calculateScattering(in vec3 rayOrigin,
in vec3 rayDir,
in float tmax,
in vec3 lightDir,
in float earthRadius,
out vec3 inscatter,
out vec3 transmittance);
float raySphereIntersection(vec3 ro, vec3 rd, float radius);
vec3 sampleMedium(in float height,
out float mieScattering, out float mieAbsorption,
out vec3 rayleighScattering, out vec3 ozoneAbsorption);
float miePhaseFunction(float cosTheta);
float rayleighPhaseFunction(float cosTheta);
vec3 getValueFromLUT(sampler2D lut, float sunCosTheta, float normalizedHeight);
void main()
{
vec3 up = normalize(fg_CameraPositionCart);
float sunCosTheta = dot(fg_SunDirectionWorld, up);
// Account for the depth layer we are currently in
// FIXME: We should probably be writing the pixel center
float x = texCoord.x * TOTAL_SLICES;
vec2 coord = vec2(fract(x), texCoord.y);
float depth = mix(fg_NearFar.x, DEPTH_RANGE, ceil(x) / TOTAL_SLICES);
// Depth goes from the 0 to DEPTH_RANGE in a squared distribution.
// The first slice is not at 0 since that would waste a slice.
float w = ceil(x) / TOTAL_SLICES;
w *= w;
float depth = w * DEPTH_RANGE;
vec3 fragPos = positionFromDepth(coord * 2.0 - 1.0, 0.0);
vec3 rayDir = vec4(fg_ViewMatrixInverse * vec4(normalize(fragPos), 0.0)).xyz;
float earthRadius = length(fg_CameraPositionCart) - fg_CameraPositionGeod.z;
float cameraHeight = length(fg_CameraPositionCart);
float earthRadius = cameraHeight - fg_CameraPositionGeod.z;
vec3 inscatter, transmittance;
calculateScattering(fg_CameraPositionCart,
rayDir,
depth,
fg_SunDirectionWorld,
earthRadius,
inscatter, transmittance);
out_inscatter = inscatter;
out_transmittance = transmittance;
vec3 rayOrigin = fg_CameraPositionCart;
float atmosDist = raySphereIntersection(rayOrigin, rayDir, ATMOSPHERE_RADIUS);
float groundDist = raySphereIntersection(rayOrigin, rayDir, earthRadius);
float tmax;
if (cameraHeight < ATMOSPHERE_RADIUS) {
// We are inside the atmosphere
if (groundDist < 0.0) {
// No ground collision, use the distance to the outer atmosphere
tmax = atmosDist;
} else {
// Use the distance to the ground
tmax = groundDist;
}
} else {
// We are in outer space, skip
fragColor = vec4(0.0, 0.0, 0.0, 1.0);
return;
}
// Clip the max distance to the depth of this slice
tmax = min(tmax, depth);
float cosTheta = dot(rayDir, fg_SunDirectionWorld);
float miePhase = miePhaseFunction(cosTheta);
float rayleighPhase = rayleighPhaseFunction(cosTheta);
vec3 L = vec3(0.0);
vec3 throughput = vec3(1.0);
float t = 0.0;
for (int i = 0; i < AERIAL_PERSPECTIVE_SAMPLES; ++i) {
float newT = ((float(i) + 0.3) / AERIAL_PERSPECTIVE_SAMPLES) * tmax;
float dt = newT - t;
t = newT;
vec3 samplePos = rayOrigin + rayDir * t;
float height = length(samplePos) - earthRadius;
float normalizedHeight = height / (ATMOSPHERE_RADIUS - earthRadius);
float mieScattering, mieAbsorption;
vec3 rayleighScattering, ozoneAbsorption;
vec3 extinction = sampleMedium(height, mieScattering, mieAbsorption,
rayleighScattering, ozoneAbsorption);
vec3 sampleTransmittance = exp(-dt*extinction);
vec3 sunTransmittance = getValueFromLUT(
transmittance_lut, sunCosTheta, normalizedHeight);
vec3 multiscattering = getValueFromLUT(
multiscattering_lut, sunCosTheta, normalizedHeight);
vec3 S =
rayleighScattering * (rayleighPhase * sunTransmittance + multiscattering) +
mieScattering * (miePhase * sunTransmittance + multiscattering);
vec3 Sint = (S - S * sampleTransmittance) / extinction;
L += throughput * Sint;
throughput *= sampleTransmittance;
}
// Instead of storing an entire vec3, store the mean of its components
float transmittance = dot(throughput, ONE_OVER_THREE);
fragColor = vec4(L, transmittance);
}

View file

@ -0,0 +1,78 @@
// An implementation of Sébastien Hillaire's "A Scalable and Production Ready
// Sky and Atmosphere Rendering Technique".
#version 330 core
const float PI = 3.141592653;
// Atmosphere parameters
// Section 2.1 of [Bruneton08], units are inverse meters
const float mie_density_height_scale = 8.33333e-4; // Hm=1.2km
const float rayleigh_density_height_scale = 1.25e-4; // Hr=8km
const float mie_scattering = 3.996e-6;
const float mie_absorption = 4.4e-6;
const vec3 rayleigh_scattering = vec3(5.802, 13.558, 33.1) * vec3(1e-6);
const vec3 ozone_absorption = vec3(0.650, 1.881, 0.085) * vec3(1e-6);
const float g = 0.8;
const float gg = g*g;
const float mie_phase_scale = 3.0/(8.0*PI);
const float rayleigh_phase_scale = 3.0/(16.0*PI);
// Returns the distance between ro and the first intersection with the sphere
// or -1.0 if there is no intersection.
// -1.0 is also returned if the ray is pointing away from the sphere.
float raySphereIntersection(vec3 ro, vec3 rd, float radius)
{
float b = dot(ro, rd);
float c = dot(ro, ro) - radius*radius;
if (c > 0.0 && b > 0.0) return -1.0;
float d = b*b - c;
if (d < 0.0) return -1.0;
if (d > b*b) return (-b+sqrt(d));
return (-b-sqrt(d));
}
// Sample the sky medium properties at a height in meters, 0 being the ground
// and ~100km being the top of the atmosphere.
// Returns the total atmospheric extinction.
vec3 sampleMedium(in float height,
out float mieScattering, out float mieAbsorption,
out vec3 rayleighScattering, out vec3 ozoneAbsorption)
{
float densityMie = exp(-mie_density_height_scale * height);
float densityRayleigh = exp(-rayleigh_density_height_scale * height);
float densityOzone = max(0.0, 1.0 - abs(height-25.0e3)/15.0e3);
mieScattering = mie_scattering * densityMie;
mieAbsorption = mie_absorption * densityMie;
rayleighScattering = rayleigh_scattering * densityRayleigh;
ozoneAbsorption = ozone_absorption * densityOzone;
return mieScattering + mieAbsorption + rayleighScattering + ozoneAbsorption;
}
// Approximation of the Mie phase function with the Cornette-Shanks phase function
float miePhaseFunction(float cosTheta)
{
float num = (1.0 - gg) * (1.0 + cosTheta*cosTheta);
float den = (2.0 + gg) * pow((1.0 + gg - 2.0 * g * cosTheta), 1.5);
return mie_phase_scale * num / den;
}
float rayleighPhaseFunction(float cosTheta)
{
return rayleigh_phase_scale * (1.0 + cosTheta*cosTheta);
}
// Sample one of the LUTs (transmittance or multiple scattering) for a given
// normalized height inside the atmosphere [0,1] and the cosine of the Sun
// zenith angle.
vec3 getValueFromLUT(sampler2D lut, float sunCosTheta, float normHeight)
{
float x = clamp(sunCosTheta * 0.5 + 0.5, 0.0, 1.0);
float y = clamp(normHeight, 0.0, 1.0);
return texture(lut, vec2(x, y)).rgb;
}

View file

@ -0,0 +1,135 @@
// An implementation of Sébastien Hillaire's "A Scalable and Production Ready
// Sky and Atmosphere Rendering Technique".
//
// This shader generates the multiple scattering LUT.
#version 330 core
out vec3 fragColor;
in vec2 texCoord;
uniform vec3 fg_CameraPositionCart;
uniform vec3 fg_CameraPositionGeod;
uniform sampler2D transmittance_lut;
const float PI = 3.141592653;
const float ATMOSPHERE_RADIUS = 6471e3;
const int SQRT_SAMPLES = 4;
const float INV_SAMPLES = 1.0 / float(SQRT_SAMPLES*SQRT_SAMPLES);
const int MULTIPLE_SCATTERING_SAMPLES = 20;
const vec3 ground_albedo = vec3(0.3);
float raySphereIntersection(vec3 ro, vec3 rd, float radius);
vec3 sampleMedium(in float height,
out float mieScattering, out float mieAbsorption,
out vec3 rayleighScattering, out vec3 ozoneAbsorption);
float miePhaseFunction(float cosTheta);
float rayleighPhaseFunction(float cosTheta);
vec3 getValueFromLUT(sampler2D lut, float sunCosTheta, float normalizedHeight);
vec3 generateRayDir(float theta, float phi)
{
float cosPhi = cos(phi);
float sinPhi = sin(phi);
float cosTheta = cos(theta);
float sinTheta = sin(theta);
return vec3(cosTheta * sinPhi, sinTheta * sinPhi, cosPhi);
}
void main()
{
float sunCosTheta = texCoord.x * 2.0 - 1.0;
vec3 sunDir = vec3(-sqrt(1.0 - sunCosTheta*sunCosTheta), 0.0, sunCosTheta);
float earthRadius = length(fg_CameraPositionCart) - fg_CameraPositionGeod.z;
float altitude = mix(earthRadius, ATMOSPHERE_RADIUS, texCoord.y);
vec3 rayOrigin = vec3(0.0, 0.0, altitude);
vec3 Ltotal = vec3(0.0);
vec3 LMStotal = vec3(0.0);
for (int i = 0; i < SQRT_SAMPLES; ++i) {
for (int j = 0; j < SQRT_SAMPLES; ++j) {
float theta = 2.0 * PI * (float(i) + 0.5) / float(SQRT_SAMPLES);
float phi = PI * (float(j) + 0.5) / float(SQRT_SAMPLES);
vec3 rayDir = generateRayDir(theta, phi);
float atmosDist = raySphereIntersection(rayOrigin, rayDir, ATMOSPHERE_RADIUS);
float groundDist = raySphereIntersection(rayOrigin, rayDir, earthRadius);
float tmax;
if (groundDist < 0.0) {
// No ground collision, use the distance to the outer atmosphere
tmax = atmosDist;
} else {
// Use the distance to the ground
tmax = groundDist;
}
float cosTheta = dot(rayDir, sunDir);
float miePhase = miePhaseFunction(cosTheta);
float rayleighPhase = rayleighPhaseFunction(-cosTheta);
vec3 L = vec3(0.0);
vec3 LMS = vec3(0.0);
vec3 throughput = vec3(1.0);
float t = 0.0;
for (int k = 0; k < MULTIPLE_SCATTERING_SAMPLES; ++k) {
float newT = ((float(k) + 0.3) / MULTIPLE_SCATTERING_SAMPLES) * tmax;
float dt = newT - t;
t = newT;
vec3 samplePos = rayOrigin + rayDir * t;
float height = length(samplePos) - earthRadius;
float normalizedHeight = height / (ATMOSPHERE_RADIUS - earthRadius);
float mieScattering, mieAbsorption;
vec3 rayleighScattering, ozoneAbsorption;
vec3 extinction = sampleMedium(height, mieScattering, mieAbsorption,
rayleighScattering, ozoneAbsorption);
vec3 sampleTransmittance = exp(-dt*extinction);
vec3 sunTransmittance = getValueFromLUT(
transmittance_lut, sunCosTheta, normalizedHeight);
vec3 S = (rayleighScattering * rayleighPhase +
mieScattering * miePhase) * sunTransmittance;
// Not using the power serie
vec3 MS = mieScattering + rayleighScattering;
vec3 MSint = (MS - MS * sampleTransmittance) / extinction;
LMS += throughput * MSint;
vec3 Sint = (S - S * sampleTransmittance) / extinction;
L += throughput * Sint;
throughput *= sampleTransmittance;
}
if (groundDist >= 0.0) {
// Account for bounced light off the Earth
vec3 p = rayOrigin + rayDir * groundDist;
float pHeight = length(p);
vec3 up = p / pHeight;
float normHeight = (pHeight - earthRadius)
/ (ATMOSPHERE_RADIUS - earthRadius);
float sunZenithCosTheta = dot(sunDir, up);
vec3 transmittanceFromGround = getValueFromLUT(
transmittance_lut, sunZenithCosTheta, normHeight);
L += transmittanceFromGround * throughput
* clamp(sunZenithCosTheta, 0.0, 1.0) * ground_albedo / PI;
}
Ltotal += L * INV_SAMPLES;
LMStotal += LMS * INV_SAMPLES;
}
}
fragColor = Ltotal / (1.0 - LMStotal);
}

View file

@ -1,3 +1,12 @@
// An implementation of Sébastien Hillaire's "A Scalable and Production Ready
// Sky and Atmosphere Rendering Technique".
//
// This shader generates the sky-view texture. Since the sky generally has low
// frequency detail, it's possible to pre-compute it on a small texture and
// sample it later when rendering the skydome. This effectively bypasses the
// need for raymarching on screen-sized textures, which is specially costly on
// larger resolutions like 4K.
#version 330 core
out vec3 fragColor;
@ -8,22 +17,29 @@ uniform vec3 fg_CameraPositionCart;
uniform vec3 fg_CameraPositionGeod;
uniform vec3 fg_SunDirectionWorld;
uniform sampler2D transmittance_lut;
uniform sampler2D multiscattering_lut;
const float PI = 3.141592653;
void calculateScattering(in vec3 rayOrigin,
in vec3 rayDir,
in float tmax,
in vec3 lightDir,
in float earthRadius,
out vec3 inscatter,
out vec3 transmittance);
const float ATMOSPHERE_RADIUS = 6471e3;
const int SCATTERING_SAMPLES = 32;
float raySphereIntersection(vec3 ro, vec3 rd, float radius);
vec3 sampleMedium(in float height,
out float mieScattering, out float mieAbsorption,
out vec3 rayleighScattering, out vec3 ozoneAbsorption);
float miePhaseFunction(float cosTheta);
float rayleighPhaseFunction(float cosTheta);
vec3 getValueFromLUT(sampler2D lut, float sunCosTheta, float normalizedHeight);
void main()
{
// Always leave the sun right in the middle of the texture as the skydome
// model is already being rotated.
float sunCosTheta = dot(fg_SunDirectionWorld, normalize(fg_CameraPositionCart));
vec3 lightDir = vec3(-sqrt(1.0 - sunCosTheta*sunCosTheta), 0.0, sunCosTheta);
vec3 up = normalize(fg_CameraPositionCart);
float sunCosTheta = dot(fg_SunDirectionWorld, up);
vec3 sunDir = vec3(-sqrt(1.0 - sunCosTheta*sunCosTheta), 0.0, sunCosTheta);
float azimuth = 2.0 * PI * texCoord.x; // [0, 2pi]
// Apply a non-linear transformation to the elevation to dedicate more
@ -32,20 +48,67 @@ void main()
float elev = l*l * sign(l) * PI * 0.5; // [-pi/2, pi/2]
vec3 rayDir = vec3(cos(elev) * cos(azimuth), cos(elev) * sin(azimuth), sin(elev));
float cameraPosLength = length(fg_CameraPositionCart);
// Since FG internally uses WG84 coordinates, we use the current Earth
// radius under the camera, which varies with the latitude. For practical
// purposes we model the Earth as a perfect sphere with this radius.
float earthRadius = cameraPosLength - fg_CameraPositionGeod.z;
float cameraHeight = length(fg_CameraPositionCart);
float earthRadius = cameraHeight - fg_CameraPositionGeod.z;
vec3 rayOrigin = vec3(0.0, 0.0, cameraPosLength);
vec3 rayOrigin = vec3(0.0, 0.0, cameraHeight);
vec3 inscatter, transmittance;
calculateScattering(rayOrigin,
rayDir,
9.0e8,
lightDir,
earthRadius,
inscatter, transmittance);
fragColor = inscatter;
float atmosDist = raySphereIntersection(rayOrigin, rayDir, ATMOSPHERE_RADIUS);
float groundDist = raySphereIntersection(rayOrigin, rayDir, earthRadius);
float tmax;
if (cameraHeight < ATMOSPHERE_RADIUS) {
// We are inside the atmosphere
if (groundDist < 0.0) {
// No ground collision, use the distance to the outer atmosphere
tmax = atmosDist;
} else {
// Use the distance to the ground
tmax = groundDist;
}
} else {
// We are in outer space, skip
fragColor = vec3(0.0);
return;
}
float cosTheta = dot(rayDir, sunDir);
float miePhase = miePhaseFunction(cosTheta);
float rayleighPhase = rayleighPhaseFunction(-cosTheta);
vec3 L = vec3(0.0);
vec3 throughput = vec3(1.0);
float t = 0.0;
for (int i = 0; i < SCATTERING_SAMPLES; ++i) {
float newT = ((float(i) + 0.3) / SCATTERING_SAMPLES) * tmax;
float dt = newT - t;
t = newT;
vec3 samplePos = rayOrigin + rayDir * t;
float height = length(samplePos) - earthRadius;
float normalizedHeight = height / (ATMOSPHERE_RADIUS - earthRadius);
float mieScattering, mieAbsorption;
vec3 rayleighScattering, ozoneAbsorption;
vec3 extinction = sampleMedium(height, mieScattering, mieAbsorption,
rayleighScattering, ozoneAbsorption);
vec3 sampleTransmittance = exp(-dt*extinction);
vec3 sunTransmittance = getValueFromLUT(
transmittance_lut, sunCosTheta, normalizedHeight);
vec3 multiscattering = getValueFromLUT(
multiscattering_lut, sunCosTheta, normalizedHeight);
vec3 S =
rayleighScattering * (rayleighPhase * sunTransmittance + multiscattering) +
mieScattering * (miePhase * sunTransmittance + multiscattering);
vec3 Sint = (S - S * sampleTransmittance) / extinction;
L += throughput * Sint;
throughput *= sampleTransmittance;
}
fragColor = L;
}

View file

@ -0,0 +1,55 @@
// An implementation of Sébastien Hillaire's "A Scalable and Production Ready
// Sky and Atmosphere Rendering Technique".
//
// This shader generates the transmittance LUT. It stores the transmittance to
// the Sun through the atmosphere for a given Sun zenith angle and a height
// inside the atmosphere (0 being the ground).
#version 330 core
out vec3 fragColor;
in vec2 texCoord;
uniform vec3 fg_CameraPositionCart;
uniform vec3 fg_CameraPositionGeod;
const float ATMOSPHERE_RADIUS = 6471e3;
const int TRANSMITTANCE_STEPS = 40;
float raySphereIntersection(vec3 ro, vec3 rd, float radius);
vec3 sampleMedium(in float height,
out float mieScattering, out float mieAbsorption,
out vec3 rayleighScattering, out vec3 ozoneAbsorption);
void main()
{
float sunCosTheta = texCoord.x * 2.0 - 1.0;
vec3 sunDir = vec3(-sqrt(1.0 - sunCosTheta*sunCosTheta), 0.0, sunCosTheta);
float earthRadius = length(fg_CameraPositionCart) - fg_CameraPositionGeod.z;
float altitude = mix(earthRadius, ATMOSPHERE_RADIUS, texCoord.y);
vec3 rayOrigin = vec3(0.0, 0.0, altitude);
float dist = raySphereIntersection(rayOrigin, sunDir, ATMOSPHERE_RADIUS);
float t = 0.0;
vec3 transmittance = vec3(1.0);
for (int i = 0; i < TRANSMITTANCE_STEPS; ++i) {
float newT = ((float(i) + 0.3) / TRANSMITTANCE_STEPS) * dist;
float dt = newT - t;
t = newT;
vec3 samplePos = rayOrigin + sunDir * t;
float height = length(samplePos) - earthRadius;
float mieScattering, mieAbsorption;
vec3 rayleighScattering, ozoneAbsorption;
vec3 extinction = sampleMedium(height, mieScattering, mieAbsorption,
rayleighScattering, ozoneAbsorption);
transmittance *= exp(-dt * extinction);
}
fragColor = transmittance;
}

View file

@ -1,133 +0,0 @@
#version 330 core
uniform vec3 beta_rayleigh = vec3(5.5e-6, 13.0e-6, 22.4e-6);
uniform float beta_mie = 21e-6;
uniform vec3 beta_absortion = vec3(2.04e-5, 4.97e-5, 1.95e-6);
uniform float beta_ambient = 0.0;
uniform float rayleigh_scale_height = 8e3;
uniform float mie_scale_height = 1.2e3;
uniform float absortion_scale_height = 30e3;
uniform float absortion_falloff = 3e3;
uniform int num_samples = 32;
uniform int num_light_samples = 4;
const float PI = 3.141592653;
const float ATMOSPHERE_RADIUS = 6471e3;
const vec3 SUN_INTENSITY = vec3(20.0);
vec2 raySphereIntersection(vec3 ro, vec3 rd, float radius)
{
vec3 tc = -ro;
float b = dot(tc, rd);
float d = b*b - dot(tc, tc) + radius*radius;
if (d < 0.0) return vec2(-1.0);
float s = sqrt(d);
return vec2(b-s, b+s);
}
void calculateScattering(in vec3 rayOrigin,
in vec3 rayDir,
in float tmax,
in vec3 lightDir,
in float earthRadius,
out vec3 inscatter,
out vec3 transmittance)
{
vec2 hit = raySphereIntersection(rayOrigin, rayDir, ATMOSPHERE_RADIUS);
vec2 hitEarth = raySphereIntersection(rayOrigin, rayDir, earthRadius - 1.0);
if (hitEarth.y > 0.0)
tmax = max(0.0, hitEarth.x);
float tmin = max(hit.x, 0.0);
tmax = min(hit.y, tmax);
if (tmax < 0.0)
discard;
float stepSize = (tmax - tmin) / float(num_samples);
const float g = 0.758; // Mie scattering direction
const float gg = g*g;
float mu = dot(rayDir, lightDir);
float mumu = mu*mu;
float phaseRayleigh = 3.0 / (50.2654824574 /* 16*PI */) * (1.0 + mumu);
float phaseMie = 3.0 / (25.1327412287 /* 8*PI */) * ((1.0 - gg) * (mumu + 1.0)) /
(pow(1.0 + gg - 2.0 * mu * g, 1.5) * (2.0 + gg));
float opticalDepthRayleigh = 0.0;
float opticalDepthMie = 0.0;
float opticalDepthAbsortion = 0.0;
float primaryTime = tmin;
vec3 extinctionFactor = vec3(0.0);
vec3 totalRayleigh = vec3(0.0);
vec3 totalMie = vec3(0.0);
for (int i = 0; i < num_samples; ++i) {
vec3 samplePoint = rayOrigin + rayDir * (primaryTime + stepSize * 0.5);
float altitude = length(samplePoint) - earthRadius;
float densityRayleigh = exp(-altitude / rayleigh_scale_height);
float densityMie = exp(-altitude / mie_scale_height);
float densityAbsortion = clamp((1.0 / cosh((absortion_scale_height - altitude) / absortion_falloff)) * densityRayleigh, 0.0, 1.0);
float stepOpticalDepthRayleigh = densityRayleigh * stepSize;
float stepOpticalDepthMie = densityMie * stepSize;
float stepOpticalDepthAbsortion = densityAbsortion * stepSize;
opticalDepthRayleigh += stepOpticalDepthRayleigh;
opticalDepthMie += stepOpticalDepthMie;
opticalDepthAbsortion += stepOpticalDepthAbsortion;
vec2 pl = raySphereIntersection(samplePoint, lightDir, ATMOSPHERE_RADIUS);
float stepSizeLight = pl.y / float(num_light_samples);
float opticalDepthLightRayleigh = 0.0;
float opticalDepthLightMie = 0.0;
float opticalDepthLightAbsortion = 0.0;
float secondaryTime = 0.0;
int j;
for (j = 0; j < num_light_samples; ++j) {
vec3 samplePointLight = samplePoint + lightDir *
(secondaryTime + stepSizeLight * 0.5);
float altitudeLight = length(samplePointLight) - earthRadius;
if (altitudeLight < 0.0)
break;
float densityLightRayleigh = exp(-altitudeLight / rayleigh_scale_height);
float densityLightMie = exp(-altitudeLight / mie_scale_height);
float densityLightAbsortion = clamp((1.0 / cosh((absortion_scale_height - altitudeLight) / absortion_falloff)) * densityLightRayleigh, 0.0, 1.0);
opticalDepthLightRayleigh += densityLightRayleigh * stepSizeLight;
opticalDepthLightMie += densityLightMie * stepSizeLight;
opticalDepthLightAbsortion += densityLightAbsortion * stepSizeLight;
secondaryTime += stepSizeLight;
}
if (j == num_light_samples) {
vec3 tau =
beta_rayleigh * (opticalDepthRayleigh + opticalDepthLightRayleigh) +
beta_mie * (opticalDepthMie + opticalDepthLightMie) +
beta_absortion * (opticalDepthAbsortion + opticalDepthLightAbsortion);
vec3 attenuation = exp(-tau);
extinctionFactor += attenuation;
totalRayleigh += stepOpticalDepthRayleigh * attenuation;
totalMie += stepOpticalDepthMie * attenuation;
}
primaryTime += stepSize;
}
transmittance = exp(-(beta_rayleigh * opticalDepthRayleigh +
beta_mie * opticalDepthMie +
beta_absortion * opticalDepthAbsortion));
inscatter = SUN_INTENSITY *
(totalRayleigh * beta_rayleigh * phaseRayleigh +
totalMie * beta_mie * phaseMie +
opticalDepthRayleigh * beta_ambient);
}

View file

@ -12,14 +12,12 @@ uniform sampler2D ao_tex;
uniform samplerCube prefiltered_envmap;
uniform sampler2DShadow shadow_tex;
uniform sampler2D dfg_lut;
uniform sampler2D aerial_inscatter_lut;
uniform sampler2D aerial_transmittance_lut;
uniform sampler2D aerial_perspective_lut;
uniform mat4 fg_ViewMatrix;
uniform mat4 fg_ViewMatrixInverse;
uniform vec3 fg_SunDirection;
uniform vec3 fg_CameraPositionCart;
uniform vec2 fg_NearFar;
uniform mat4 fg_LightMatrix_csm0;
uniform mat4 fg_LightMatrix_csm1;
@ -329,30 +327,34 @@ vec3 BRDF(in vec3 albedo, in float metalness, in float roughness,
}
//------------------------------------------------------------------------------
// Atmospheric scattering
float map(float value, float min1, float max1, float min2, float max2) {
return min2 + (value - min1) * (max2 - min2) / (max1 - min1);
}
vec3 sampleAerialPerspectiveSlice(sampler2D tex, int slice)
vec4 sampleAerialPerspectiveSlice(int slice)
{
// Sample at the pixel center
float offset = slice * AERIAL_LUT_TILE_SIZE + AERIAL_LUT_TEXEL_SIZE * 0.5;
float x = texCoord.x * (AERIAL_LUT_TILE_SIZE - AERIAL_LUT_TEXEL_SIZE) + offset;
return texture(tex, vec2(x, texCoord.y)).rgb;
return texture(aerial_perspective_lut, vec2(x, texCoord.y));
}
vec3 sampleAerialPerspective(sampler2D tex, vec3 zero, float depth)
vec4 sampleAerialPerspective(float depth)
{
vec3 color;
depth = min(abs(depth), AERIAL_MAX_DEPTH);
float d = map(depth, fg_NearFar.x, AERIAL_MAX_DEPTH, 0.0, AERIAL_SLICES);
if (d <= 1.0) {
color = mix(zero, sampleAerialPerspectiveSlice(tex, 0), d);
vec4 color;
// Map to [0,1]
float w = depth / AERIAL_MAX_DEPTH;
// Squared distribution
w = sqrt(clamp(w, 0.0, 1.0));
// Remap to [0,16] to sample the right tile
w *= AERIAL_SLICES;
if (w <= 1.0) {
// Handle special case of fragments behind the first slice
color = mix(vec4(0.0, 0.0, 0.0, 1.0), sampleAerialPerspectiveSlice(0), w);
} else {
d -= 1.0;
color = mix(sampleAerialPerspectiveSlice(tex, int(floor(d))),
sampleAerialPerspectiveSlice(tex, int(ceil(d))),
fract(d));
w -= 1.0; // [0,15]
// Manually linearly interpolate between slices
color = mix(sampleAerialPerspectiveSlice(int(floor(w))),
sampleAerialPerspectiveSlice(int(ceil(w))),
fract(w));
}
return color;
}
@ -403,10 +405,8 @@ void main()
float shadowFactor = getShadowing(pos, n, NdotL);
vec3 color = ambient + brdf * sunIlluminance * shadowFactor;
vec3 inscatter = sampleAerialPerspective(
aerial_inscatter_lut, vec3(0.0), length(pos));
vec3 transmittance = sampleAerialPerspective(
aerial_transmittance_lut, vec3(1.0), length(pos));
vec4 aerialPerspective = sampleAerialPerspective(length(pos));
color = color * aerialPerspective.a + aerialPerspective.rgb * SUN_INTENSITY;
fragHdrColor = color * transmittance + inscatter;
fragHdrColor = color;
}

View file

@ -2,19 +2,35 @@
out vec4 fragColor;
in vec3 rayDir;
in vec3 vRayDir;
in vec3 vRayDirView;
uniform sampler2D sky_view_lut;
uniform sampler2D transmittance_lut;
uniform vec3 fg_SunDirection;
const float PI = 3.141592653;
const vec3 SUN_INTENSITY = vec3(20.0);
const float sun_solid_angle = 0.53*PI/180.0; // ~half a degree
const float sun_cos_solid_angle = cos(sun_solid_angle);
void main()
{
vec3 rayDir = normalize(vRayDir);
float azimuth = atan(rayDir.y, rayDir.x) / PI * 0.5 + 0.5;
// Undo the non-linear transformation from the sky-view LUT
float l = asin(rayDir.z);
float elev = sqrt(abs(l) / (PI * 0.5)) * sign(l) * 0.5 + 0.5;
vec3 color = texture(sky_view_lut, vec2(azimuth, elev)).rgb;
color *= SUN_INTENSITY;
// Render the Sun disk
// XXX: Apply transmittance
float cosTheta = dot(normalize(vRayDirView), fg_SunDirection);
if (cosTheta >= sun_cos_solid_angle)
color += SUN_INTENSITY;
fragColor = vec4(color, 1.0);
}

View file

@ -2,12 +2,15 @@
layout(location = 0) in vec4 pos;
out vec3 rayDir;
out vec3 vRayDir;
out vec3 vRayDirView;
uniform mat4 osg_ModelViewMatrix;
uniform mat4 osg_ModelViewProjectionMatrix;
void main()
{
gl_Position = osg_ModelViewProjectionMatrix * pos;
rayDir = normalize(pos.xyz);
vRayDir = normalize(pos.xyz);
vRayDirView = (osg_ModelViewMatrix * vec4(vRayDir, 0.0)).xyz;
}