maxmustermann/armory
1
0
Fork 0
Code Issues Pull requests Projects Releases Wiki Activity

Use 2D LUT for Nishita skies

Browse source
This commit is contained in:
Moritz Brückner 2021-03-25 23:25:41 +01:00
parent 05c14238f2
commit 8d812548c4
4 changed files with 159 additions and 66 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

63
Shaders/std/sky.glsl
View file

@ -1,12 +1,13 @@
/* Various sky functions
* =====================
*
* Nishita model is based on https://github.com/wwwtyro/glsl-atmosphere(Unlicense License)
* Nishita model is based on https://github.com/wwwtyro/glsl-atmosphere (Unlicense License)
*
* Changes to the original implementation:
* - r and pSun parameters of nishita_atmosphere() are already normalized
* - Some original parameters of nishita_atmosphere() are replaced with pre-defined values
* - Implemented air, dust and ozone density node parameters (see Blender source)
* - Replaced the inner integral calculation with a LUT lookup
*
* Reference for the sun's limb darkening and ozone calculations:
* [Hill] Sebastien Hillaire. Physically Based Sky, Atmosphere and Cloud Rendering in Frostbite
@ -19,15 +20,18 @@
#ifndef _SKY_GLSL_
#define _SKY_GLSL_
// OpenGl ES doesn't support 1D textures so we use a 1 px height sampler2D here...
uniform sampler2D nishitaLUT;
#define PI 3.141592
#ifndef PI
#define PI 3.141592
#endif
#ifndef HALF_PI
#define HALF_PI 1.570796
#endif
#define nishita_iSteps 16
#define nishita_jSteps 8
// The values here are taken from Cycles code if they
// These values are taken from Cycles code if they
// exist there, otherwise they are taken from the example
// in the glsl-atmosphere repo
#define nishita_sun_intensity 22.0
@ -49,7 +53,17 @@ uniform sampler2D nishitaLUT;
// Values from [Hill: 60]
#define sun_limb_darkening_col vec3(0.397, 0.503, 0.652)
#define heightToLUT(h) (textureLod(nishitaLUT, vec2(clamp(h * (1 / 60000.0), 0.0, 1.0), 0.0), 0.0).xyz * 10.0)
float random(vec2 coords) {
return fract(sin(dot(coords.xy, vec2(12.9898,78.233))) * 43758.5453);
}
vec3 nishita_lookupLUT(const float height, const float sunTheta) {
vec2 coords = vec2(
sqrt(height * (1 / nishita_atmo_radius)),
0.5 + 0.5 * sign(sunTheta - HALF_PI) * sqrt(abs(sunTheta * (1 / HALF_PI) - 1))
);
return textureLod(nishitaLUT, coords, 0.0).rgb;
}
/* Approximates the density of ozone for a given sample height. Values taken from Cycles code. */
float nishita_density_ozone(const float height) {
@ -112,40 +126,21 @@ vec3 nishita_atmosphere(const vec3 r, const vec3 r0, const vec3 pSun, const floa
// Calculate the height of the sample.
float iHeight = length(iPos) - rPlanet;
// Calculate the optical depth of the Rayleigh and Mie scattering for this step.
vec3 iLookup = heightToLUT(iHeight);
float odStepRlh = iLookup.x * iStepSize;
float odStepMie = iLookup.y * iStepSize;
// Calculate the optical depth of the Rayleigh and Mie scattering for this step
float odStepRlh = exp(-iHeight / nishita_rayleigh_scale) * density.x * iStepSize;
float odStepMie = exp(-iHeight / nishita_mie_scale) * density.y * iStepSize;
// Accumulate optical depth.
iOdRlh += odStepRlh;
iOdMie += odStepMie;
// Calculate the step size of the secondary ray.
float jStepSize = nishita_rsi(iPos, pSun, nishita_atmo_radius).y / float(nishita_jSteps);
// Idea behind this: "Rotate" everything by iPos (-> iPos is the new zenith) and then all calculations for the
// inner integral only depend on the sample height (iHeight) and sunTheta (angle between sun and new zenith).
float sunTheta = acos(dot(normalize(iPos), normalize(pSun)));
vec3 jODepth = nishita_lookupLUT(iHeight, sunTheta);// * vec3(14000000 / 255, 14000000 / 255, 2000000 / 255);
// Initialize the secondary ray time.
float jTime = 0.0;
// Initialize optical depth accumulators for the secondary ray.
vec3 jODepth = vec3(0.0); // (Rayleigh, Mie, ozone)
// Sample the secondary ray.
for (int j = 0; j < nishita_jSteps; j++) {
// Calculate the secondary ray sample position.
vec3 jPos = iPos + pSun * (jTime + jStepSize * 0.5);
// Calculate the height of the sample.
float jHeight = length(jPos) - rPlanet;
// Accumulate the optical depth.
vec3 jLookup = heightToLUT(jHeight);
jODepth += jLookup * jStepSize;
// Increment the secondary ray time.
jTime += jStepSize;
}
// Apply dithering to reduce visible banding
jODepth += mix(-1000, 1000, random(r.xy));
// Calculate attenuation.
vec3 attn = exp(-(

2
Sources/armory/object/Uniforms.hx
View file

@ -21,7 +21,7 @@ class Uniforms {
public static function textureLink(object: Object, mat: MaterialData, link: String): kha.Image {
if (link == "_nishitaLUT") {
if (armory.renderpath.Nishita.data == null) armory.renderpath.Nishita.recompute(Scene.active.world);
return armory.renderpath.Nishita.data.optDepthLUT;
return armory.renderpath.Nishita.data.lut;
}
#if arm_ltc
else if (link == "_ltcMat") {

158
Sources/armory/renderpath/Nishita.hx
View file

@ -5,36 +5,77 @@ import kha.graphics4.TextureFormat;
import kha.graphics4.Usage;
import iron.data.WorldData;
import iron.math.Vec2;
import iron.math.Vec3;
import armory.math.Helper;
/**
Utility class to control the Nishita sky model.
**/
class Nishita {
public static var data: NishitaData = null;
/**
Recompute the nishita lookup table. Call this function after updating
the sky density settings.
**/
public static function recompute(world: WorldData) {
if (world == null || world.raw.sun_direction == null) return;
if (data == null) data = new NishitaData();
// TODO
data.recompute(1.0, 1.0, 1.0);
data.computeLUT(new Vec3(1.0, 1.0, 1.0));
}
}
/**
This class holds the precalculated result of the inner scattering integral
of the Nishita sky model. The outer integral is calculated in
[`armory/Shaders/std/sky.glsl`](https://github.com/armory3d/armory/blob/master/Shaders/std/sky.glsl).
@see `armory.renderpath.Nishita`
**/
class NishitaData {
static inline var LUT_WIDTH = 16;
/** Maximum ray height as defined by Cycles **/
static inline var MAX_HEIGHT = 60000;
public var lut: kha.Image;
static inline var RAYLEIGH_SCALE = 8e3;
static inline var MIE_SCALE = 1.2e3;
/**
The amount of individual sample heights stored in the LUT (and the width
of the LUT image).
**/
public static var lutHeightSteps = 128;
/**
The amount of individual sun angle steps stored in the LUT (and the
height of the LUT image).
**/
public static var lutAngleSteps = 128;
public var optDepthLUT: kha.Image;
/**
Amount of steps for calculating the inner scattering integral. Heigher
values are more precise but take longer to compute.
**/
public static var jSteps = 8;
/** Radius of the atmosphere in meters. **/
public static var radiusAtmo = 6420000;
/**
Radius of the planet in meters. The default value is the earth radius as
defined in Cycles.
**/
public static var radiusPlanet = 6360000;
/** Rayleigh scattering scale parameter. **/
public static var rayleighScale = 8e3;
/** Mie scattering scale parameter. **/
public static var mieScale = 1.2e3;
public function new() {}
/** Approximates the density of ozone for a given sample height. **/
function getOzoneDensity(height: FastFloat): FastFloat {
// Values are taken from Cycles code
if (height < 10000.0 || height >= 40000.0) {
return 0.0;
}
@ -45,36 +86,91 @@ class NishitaData {
}
/**
The RGBA texture layout looks as follows:
R = Rayleigh optical depth at height \in [0, 60000]
G = Mie optical depth at height \in [0, 60000]
B = Ozone optical depth at height \in [0, 60000]
A = Unused
Ray-sphere intersection test that assumes the sphere is centered at the
origin. There is no intersection when result.x > result.y. Otherwise
this function returns the distances to the two intersection points,
which might be equal.
**/
public function recompute(densityFacAir: FastFloat, densityFacDust: FastFloat, densityFacOzone: FastFloat) {
optDepthLUT = kha.Image.create(LUT_WIDTH, 1, TextureFormat.RGBA32, Usage.StaticUsage);
function raySphereIntersection(rayOrigin: Vec3, rayDirection: Vec3, sphereRadius: Int): Vec2 {
// Algorithm is described here: https://en.wikipedia.org/wiki/Line%E2%80%93sphere_intersection
var a = rayDirection.dot(rayDirection);
var b = 2.0 * rayDirection.dot(rayOrigin);
var c = rayOrigin.dot(rayOrigin) - (sphereRadius * sphereRadius);
var d = (b * b) - 4.0 * a * c;
var textureData = optDepthLUT.lock();
for (i in 0...LUT_WIDTH) {
// Get the height for each LUT pixel i (in [-1, 1] range)
var height = (i / LUT_WIDTH) * 2 - 1;
// Ray does not intersect the sphere
if (d < 0.0) return new Vec2(1e5, -1e5);
return new Vec2(
(-b - Math.sqrt(d)) / (2.0 * a),
(-b + Math.sqrt(d)) / (2.0 * a)
);
}
/**
Computes the LUT texture for the given density values.
@param density 3D vector of air density, dust density, ozone density
**/
public function computeLUT(density: Vec3) {
var imageData = new haxe.io.Float32Array(lutHeightSteps * lutAngleSteps * 4);
for (x in 0...lutHeightSteps) {
var height = (x / (lutHeightSteps - 1));
// Use quadratic height for better horizon precision
// See https://sebh.github.io/publications/egsr2020.pdf (5.3)
height = 0.5 + 0.5 * Helper.sign(height) * Math.sqrt(Math.abs(height));
height *= MAX_HEIGHT; // Denormalize
height *= height;
height *= radiusAtmo; // Denormalize
// Make sure we use 32 bit floats
var optDepthRayleigh: FastFloat = Math.exp(-height / RAYLEIGH_SCALE) * densityFacAir;
var optDepthMie: FastFloat = Math.exp(-height / MIE_SCALE) * densityFacDust;
var optDepthOzone: FastFloat = getOzoneDensity(height) * densityFacOzone;
for (y in 0...lutAngleSteps) {
var sunTheta = y / (lutAngleSteps - 1) * 2 - 1;
// 10 is the maximum density, so we divide by it to be able to use normalized values
textureData.set(i * 4 + 0, Std.int(optDepthRayleigh * 255 / 10));
textureData.set(i * 4 + 1, Std.int(optDepthMie * 255 / 10));
textureData.set(i * 4 + 2, Std.int(optDepthOzone * 255 / 10));
textureData.set(i * 4 + 3, 255); // Unused
// Improve horizon precision
// See https://sebh.github.io/publications/egsr2020.pdf (5.3)
sunTheta = Helper.sign(sunTheta) * sunTheta * sunTheta;
sunTheta = sunTheta * Math.PI / 2 + Math.PI / 2; // Denormalize
var jODepth = sampleSecondaryRay(height, sunTheta, density);
var pixelIndex = (x + y * lutHeightSteps) * 4;
imageData[pixelIndex + 0] = jODepth.x;
imageData[pixelIndex + 1] = jODepth.y;
imageData[pixelIndex + 2] = jODepth.z;
imageData[pixelIndex + 3] = 1.0; // Unused
}
}
optDepthLUT.unlock();
lut = kha.Image.fromBytes(imageData.view.buffer, lutHeightSteps, lutAngleSteps, TextureFormat.RGBA128, Usage.StaticUsage);
}
/**
Calculates the integral for the secondary ray.
**/
public function sampleSecondaryRay(height: FastFloat, sunTheta: FastFloat, density: Vec3): Vec3 {
// Reconstruct values from the shader
var iPos = new Vec3(0, 0, height + radiusPlanet);
var pSun = new Vec3(0.0, Math.sin(sunTheta), Math.cos(sunTheta)).normalize();
var jTime: FastFloat = 0.0;
var jStepSize: FastFloat = raySphereIntersection(iPos, pSun, radiusAtmo).y / jSteps;
// Optical depth accumulators for the secondary ray (Rayleigh, Mie, ozone)
var jODepth = new Vec3();
for (i in 0...jSteps) {
// Calculate the secondary ray sample position and height
var jPos = iPos.clone().add(pSun.clone().mult(jTime + jStepSize * 0.5));
var jHeight = jPos.length() - radiusPlanet;
// Accumulate optical depth
var optDepthRayleigh = Math.exp(-jHeight / rayleighScale) * density.x;
var optDepthMie = Math.exp(-jHeight / mieScale) * density.y;
var optDepthOzone = getOzoneDensity(jHeight) * density.z;
jODepth.addf(optDepthRayleigh, optDepthMie, optDepthOzone);
jTime += jStepSize;
}
return jODepth.mult(jStepSize);
}
}

2
blender/arm/material/cycles_nodes/nodes_texture.py
View file

@ -373,6 +373,8 @@ def parse_sky_nishita(node: bpy.types.ShaderNodeTexSky, state: ParserState) -> v
curshader = state.curshader
curshader.add_include('std/sky.glsl')
curshader.add_uniform('vec3 sunDir', link='_sunDirection')
curshader.add_uniform('sampler2D nishitaLUT', link='_nishitaLUT', included=True,
tex_addr_u='clamp', tex_addr_v='clamp')
planet_radius = 6360e3 # Earth radius used in Blender
ray_origin_z = planet_radius + node.altitude