16
0
mirror of https://github.com/MaSzyna-EU07/maszyna.git synced 2026-07-22 17:29:18 +02:00

Some renaming

This commit is contained in:
2025-04-15 01:32:56 +02:00
parent 12d6a1578d
commit a39d972126
864 changed files with 29 additions and 34 deletions

View File

@@ -0,0 +1,6 @@
#ifndef AERIAL_PERSPECTIVE_HLSLI
#define AERIAL_PERSPECTIVE_HLSLI
SamplerState g_AerialPerspectiveLutSampler : register(s13);
#endif

View File

@@ -0,0 +1,8 @@
#ifndef ALPHA_MASK_HLSLI
#define ALPHA_MASK_HLSLI
void AlphaMask(in float alpha) {
if(g_DrawConstants.m_AlphaThreshold >= 0. ? (alpha < g_DrawConstants.m_AlphaThreshold) : (alpha >= -g_DrawConstants.m_AlphaThreshold)) discard;
}
#endif

File diff suppressed because it is too large Load Diff

File diff suppressed because it is too large Load Diff

View File

@@ -0,0 +1,7 @@
RWTexture2D<float4> g_ColorImage : register(u0);
Texture2D<float2> g_AutoExposureTexture : register(t0);
[numthreads(8, 8, 1)]
void CS_ApplyAutoExposure(uint3 PixCoord : SV_DispatchThreadID) {
g_ColorImage[PixCoord.xy] = max(1.e-6, g_ColorImage[PixCoord.xy] * pow(2., g_AutoExposureTexture[uint2(0, 0)].x));
}

View File

@@ -0,0 +1,90 @@
Texture2D<float3> g_Source : register(t0);
Texture2D<float3> g_SourceAdd : register(t1);
RWTexture2D<float4> g_Output : register(u0);
SamplerState g_SamplerLinearClamp : register(s0);
cbuffer BloomConstants : register(b0){
float4 g_PrefilterVector;
float g_MixFactor;
}
float3 Prefilter(float3 c) {
float brightness = max(c.r, max(c.g, c.b));
float soft = brightness - g_PrefilterVector.y;
soft = clamp(soft, 0, g_PrefilterVector.z);
soft = soft * soft * g_PrefilterVector.w;
float contribution = max(soft, brightness - g_PrefilterVector.x);
contribution /= max(brightness, 0.00001);
return c * contribution;
}
[numthreads(8, 8, 1)]
void CS_BloomPrefilter(uint3 PixCoord : SV_DispatchThreadID) {
g_Output[PixCoord.xy].xyz = Prefilter(g_Source[PixCoord.xy]);
}
[numthreads(8, 8, 1)]
void CS_BloomDownsample(uint3 PixCoord : SV_DispatchThreadID) {
uint2 pixel = PixCoord.xy;
uint2 source_size;
g_Source.GetDimensions(source_size.x, source_size.y);
uint2 output_size;
g_Output.GetDimensions(output_size.x, output_size.y);
float2 uv = (pixel + .5) / output_size;
float2 offset = 1. / source_size;
float3 a = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(-2, 2), 0.);
float3 b = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(0, 2), 0.);
float3 c = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(2, 2), 0.);
float3 d = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(-2, 0), 0.);
float3 e = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(0, 0), 0.);
float3 f = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(2, 0), 0.);
float3 g = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(-2, -2), 0.);
float3 h = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(0, -2), 0.);
float3 i = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(2, -2), 0.);
float3 j = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(-1, 1), 0.);
float3 k = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(1, 1), 0.);
float3 l = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(-1, -1), 0.);
float3 m = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(1, -1), 0.);
g_Output[pixel].xyz = e*0.125+(a+c+g+i)*0.03125+(b+d+f+h)*0.0625+(j+k+l+m)*0.125;
}
float3 SampleTent(Texture2D<float3> source, SamplerState source_sampler, float2 uv) {
uint2 source_size;
source.GetDimensions(source_size.x, source_size.y);
float2 offset = 1. / source_size;
float3 a = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(-1, 1), 0.);
float3 b = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(0, 1), 0.);
float3 c = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(1, 1), 0.);
float3 d = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(-1, 0), 0.);
float3 e = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(0, 0), 0.);
float3 f = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(1, 0), 0.);
float3 g = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(-1, -1), 0.);
float3 h = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(0, -1), 0.);
float3 i = g_Source.SampleLevel(g_SamplerLinearClamp, uv + offset * float2(1, -1), 0.);
return e*0.25+(a+c+g+i)*0.0625+(b+d+f+h)*0.125;
}
[numthreads(8, 8, 1)]
void CS_BloomUpsample(uint3 PixCoord : SV_DispatchThreadID) {
uint2 output_size;
g_Output.GetDimensions(output_size.x, output_size.y);
float2 uv = (PixCoord.xy + .5) / output_size;
g_Output[PixCoord.xy].xyz += SampleTent(g_Source, g_SamplerLinearClamp, uv);
}
[numthreads(8, 8, 1)]
void CS_BloomApply(uint3 PixCoord : SV_DispatchThreadID) {
uint2 output_size;
g_Output.GetDimensions(output_size.x, output_size.y);
float2 uv = (PixCoord.xy + .5) / output_size;
g_Output[PixCoord.xy].xyz = g_SourceAdd[PixCoord.xy] + SampleTent(g_Source, g_SamplerLinearClamp, uv) * g_MixFactor;
}

View File

@@ -0,0 +1,61 @@
#ifndef COLOR_TRANSFORM_HLSLI
#define COLOR_TRANSFORM_HLSLI
float3 ACEScg_to_XYZ(in float3 c);
float3 REC709_to_XYZ(in float3 c);
float3 AdobeRGB_to_XYZ(in float3 c);
float3 XYZ_to_REC2020(in float3 c);
float3 XYZ_to_REC709(in float3 c);
float3 ACEScg_to_REC709(in float3 c);
float3 ACEScg_to_XYZ(in float3 c){
return
c.r * float3(0.66245418, 0.27222872, -0.00557465) +
c.g * float3(0.13400421, 0.67408177, 0.00406073) +
c.b * float3(0.15618769, 0.05368952, 1.0103391);
}
float3 REC709_to_XYZ(in float3 c) {
return
c.r * float3(0.4123908, 0.21263901, 0.01933082) +
c.g * float3(0.35758434, 0.71516868, 0.11919478) +
c.b * float3(0.18048079, 0.07219232, 0.95053215);
}
float3 AdobeRGB_to_XYZ(in float3 c) {
return
c.r * float3(0.57667, 0.29734, 0.02703) +
c.g * float3(0.18556, 0.62736, 0.07069) +
c.b * float3(0.18823, 0.07529, 0.99134);
}
float3 XYZ_to_REC2020(in float3 c) {
return
c.r * float3(1.71665119, -0.66668435, 0.01763986) +
c.g * float3(-0.35567078, 1.61648124, -0.04277061) +
c.b * float3(-0.25336628, 0.01576855, 0.94210312);
}
float3 XYZ_to_REC709(in float3 c) {
return
c.r * float3(3.24096994, -0.96924364, 0.05563008) +
c.g * float3(-1.53738318, 1.8759675, -0.20397696) +
c.b * float3(-0.49861076, 0.04155506, 1.05697151);
}
float3 ACEScg_to_REC709(in float3 c) {
return
c.r * float3(1.70507964, -0.12970053, -0.02416634) +
c.g * float3(-0.62423346, 1.13846855, -0.12461416) +
c.b * float3(-0.08084618, -0.00876802, 1.1487805);
}
float REC709_to_Luma(in float3 c) {
return dot(c, float3(.2126, .7152, .0722));
}
float XYZ_to_Luma(in float3 c) {
return REC709_to_Luma(XYZ_to_REC709(c));
}
#endif

View File

@@ -0,0 +1,150 @@
#include "color_transform.hlsli"
/* ---------------------------------------------------------------------------------------------- */
/* Resource bindings */
/* ---------------------------------------------------------------------------------------------- */
Texture2D<float4> g_ColorImage : register(t0);
globallycoherent RWTexture2D<float2> g_AutoExposureTexture : register(u0);
globallycoherent RWStructuredBuffer<uint> g_SpdAtomicCounter : register(u2);
globallycoherent RWTexture2D<float> g_ImgMip6 : register(u3);
/* --------------------------------------- Push constants --------------------------------------- */
struct PushConstantsCompLuminance {
uint m_num_mips;
uint m_num_workgroups;
float m_delta_time;
float m_min_luminance_ev;
float m_max_luminance_ev;
};
#ifdef SPIRV
[[vk::push_constant]] ConstantBuffer<PushConstantsCompLuminance> g_PushConstants;
#else
cbuffer g_Const : register(b1) { PushConstantsCompLuminance g_PushConstants; }
#endif
/* ---------------------------------------------------------------------------------------------- */
/* Function interface */
/* ---------------------------------------------------------------------------------------------- */
uint GetNumMips() {
return g_PushConstants.m_num_mips;
}
uint GetNumWorkGroups() {
return g_PushConstants.m_num_workgroups;
}
float GetDeltaTime() {
return g_PushConstants.m_delta_time;
}
float ComputeAutoExposureFromLavg(float l_avg) {
l_avg = exp(l_avg);
const float s = 100.; //ISO arithmetic speed
const float k = 12.5;
float exposure_iso100 = log2((l_avg * s) / k);
const float q = .65;
float l_max = (78. / (q * s)) * pow(2., exposure_iso100);
return 1. / l_max;
}
void OutputLuminance(float4 out_value) {
float prev = g_AutoExposureTexture[uint2(0, 0)].y;
float result = out_value.r;
result = clamp(
result,
log(.125 * exp2(g_PushConstants.m_min_luminance_ev)),
log(.125 * exp2(g_PushConstants.m_max_luminance_ev)));
if (prev < 1.e8) // Compare Lavg, so small or negative values
{
result = prev + (result - prev) * (1. - exp(-GetDeltaTime()));
}
g_AutoExposureTexture[uint2(0, 0)] = float2(ComputeAutoExposureFromLavg(result), result);
}
/* ---------------------------------------------------------------------------------------------- */
/* SPD implementation */
/* ---------------------------------------------------------------------------------------------- */
#define A_GPU
#define A_HLSL
#include "amd/ffx_a.h"
/* ---------------------------------------- Shared memory --------------------------------------- */
groupshared AU1 g_SpdCounter;
groupshared AF4 g_SpdIntermediate[16][16];
/* ------------------------------------ Image load-store ops ------------------------------------ */
AF4 SpdLoadSourceImage(ASU2 p, AU1 slice){
return AF4(log(max(1.e-6, XYZ_to_Luma(g_ColorImage[p].rgb))), 0., 0., 0.);
}
AF4 SpdLoad(ASU2 p, AU1 slice) {
return float4(g_ImgMip6[p], 0., 0., 0.);
}
void SpdStore(ASU2 p, AF4 value, AU1 mip, AU1 slice) {
if(mip == 5) {
g_ImgMip6[p].r = value;
}
else if(mip == GetNumMips() - 1 && all(p == ASU2(0, 0))) {
OutputLuminance(value);
}
}
/* -------------------------------- Shared memory load-store ops -------------------------------- */
AF4 SpdLoadIntermediate(ASU1 x, ASU1 y) {
return g_SpdIntermediate[x][y];
}
void SpdStoreIntermediate(ASU1 x, ASU1 y, AF4 value) {
g_SpdIntermediate[x][y] = value;
}
/* ------------------------------------- Atomic counter ops ------------------------------------- */
void SpdIncreaseAtomicCounter(AU1 slice) {
InterlockedAdd(g_SpdAtomicCounter[0], 1, g_SpdCounter);
}
AU1 SpdGetAtomicCounter() {
return g_SpdCounter;
}
void SpdResetAtomicCounter(AU1 slice) {
g_SpdAtomicCounter[0] = 0;
}
/* ------------------------------------- Reduction function ------------------------------------- */
AF4 SpdReduce4(AF4 v0, AF4 v1, AF4 v2, AF4 v3) {
return (v0 + v1 + v2 + v3) * .25;
}
/* ----------------------------------------- SPD wrapper ---------------------------------------- */
#include "amd/ffx_spd.h"
[numthreads(256,1,1)]
void CS_ComputeAvgLuminance(uint3 WorkGroupId : SV_GroupID, uint LocalThreadIndex : SV_GroupIndex) {
SpdDownsample(
WorkGroupId.xy,
LocalThreadIndex,
GetNumMips(),
GetNumWorkGroups(),
0);
};

View File

@@ -0,0 +1,31 @@
#ifndef DRAW_CONSTANTS_HLSLI
#define DRAW_CONSTANTS_HLSLI
struct VertexPushConstants {
float4x3 m_ModelView;
float4x3 m_ModelViewHistory;
float3 m_Diffuse;
float m_AlphaThreshold;
float m_AlphaMult;
float m_SelfIllum;
};
#ifdef SPIRV
[[vk::push_constant]] ConstantBuffer<VertexPushConstants> g_DrawConstants;
#else
cbuffer g_Const : register(b1) { VertexPushConstants g_DrawConstants; }
#endif
#endif
float4x3 GetModelView() {
return g_DrawConstants.m_ModelView;
}
float4x3 GetModelViewHistory() {
return g_DrawConstants.m_ModelViewHistory;
}

View File

@@ -0,0 +1,19 @@
#ifndef DRAW_CONSTANTS_SHADOW_HLSLI
#define DRAW_CONSTANTS_SHADOW_HLSLI
struct VertexPushConstantsShadow {
float4x4 m_ModelViewProjection;
float m_AlphaThreshold;
};
#ifdef SPIRV
[[vk::push_constant]] ConstantBuffer<VertexPushConstantsShadow> g_DrawConstants;
#else
cbuffer g_Const : register(b1) { VertexPushConstantsShadow g_DrawConstants; }
#endif
#endif

View File

@@ -0,0 +1,436 @@
// Shamelessly stolen from https://www.3dgep.com/forward-plus/
#include "light.hlsli"
#include "primitives.hlsli"
#define FORWARD_PLUS_TILE_SIZE 8
#define FORWARD_PLUS_LIGHTS_PER_TILE 64
struct ComputeShaderInput
{
uint3 groupID : SV_GroupID;
uint3 groupThreadID : SV_GroupThreadID;
uint3 dispatchThreadID : SV_DispatchThreadID;
uint groupIndex : SV_GroupIndex;
};
/* ---------------------------------------------------------------------------------------------- */
/* Constant Buffers */
/* ---------------------------------------------------------------------------------------------- */
// Global variables
cbuffer DispatchParams : register(b4)
{
// Number of groups dispatched. (This parameter is not available as an HLSL system value!)
uint3 numThreadGroups;
// Total number of threads dispatched. (Also not available as an HLSL system value!)
// Note: This value may be less than the actual number of threads executed
// if the screen size is not evenly divisible by the block size.
uint3 numThreads;
};
// Parameters required to convert screen space coordinates to view space.
cbuffer ScreenToViewParams : register( b3 )
{
float4x4 InverseProjection;
float4x4 WorldToViewSpace;
float2 ScreenDimensions;
uint NumLights;
}
/* ---------------------------------------------------------------------------------------------- */
/* Shader Resources */
/* ---------------------------------------------------------------------------------------------- */
// The depth from the screen space texture.
Texture2D<float> DepthTextureVS : register(t3);
StructuredBuffer<PackedLight> Lights : register(t8);
// Precomputed frustums for the grid.
StructuredBuffer<Frustum> in_Frustums : register(t9);
/* ---------------------------------------------------------------------------------------------- */
/* Unordered Access views */
/* ---------------------------------------------------------------------------------------------- */
// View space frustums for the grid cells.
RWStructuredBuffer<Frustum> out_Frustums: register(u0);
// Global counter for current index into the light index list.
// "o_" prefix indicates light lists for opaque geometry while
// "t_" prefix indicates light lists for transparent geometry.
globallycoherent RWStructuredBuffer<uint> o_LightIndexCounter: register(u1);
globallycoherent RWStructuredBuffer<uint> t_LightIndexCounter: register(u2);
// Light index lists and light grids.
globallycoherent RWStructuredBuffer<uint> o_LightIndexList : register(u3);
globallycoherent RWStructuredBuffer<uint> t_LightIndexList : register(u4);
globallycoherent RWTexture2D<uint2> o_LightGrid : register(u5);
globallycoherent RWTexture2D<uint2> t_LightGrid : register(u6);
/* ---------------------------------------------------------------------------------------------- */
/* Group shared globals */
/* ---------------------------------------------------------------------------------------------- */
groupshared uint uMinDepth;
groupshared uint uMaxDepth;
groupshared Frustum GroupFrustum;
// Opaque geometry light lists.
groupshared uint o_LightCount;
groupshared uint o_LightIndexStartOffset;
groupshared uint o_LightList[FORWARD_PLUS_LIGHTS_PER_TILE];
// Transparent geometry light lists.
groupshared uint t_LightCount;
groupshared uint t_LightIndexStartOffset;
groupshared uint t_LightList[FORWARD_PLUS_LIGHTS_PER_TILE];
/* ---------------------------------------------------------------------------------------------- */
/* Utility functions */
/* ---------------------------------------------------------------------------------------------- */
// Add the light to the visible light list for opaque geometry.
void o_AppendLight(uint lightIndex);
// Add the light to the visible light list for transparent geometry.
void t_AppendLight(uint lightIndex);
// Convert clip space coordinates to view space
float4 ClipToView(float4 clip);
// Convert screen space coordinates to view space.
float4 ScreenToView(float4 screen);
/* ---------------------------------------------------------------------------------------------- */
/* Frustum computation */
/* ---------------------------------------------------------------------------------------------- */
[numthreads(FORWARD_PLUS_TILE_SIZE, FORWARD_PLUS_TILE_SIZE, 1)]
void CS_ComputeFrustums(in ComputeShaderInput input) {
// View space eye position is always at the origin.
const float3 eyePos = (float3)0.;
// Compute the 4 corner points on the far clipping plane to use as the
// frustum vertices.
float4 screenSpace[5];
// Top left point
screenSpace[0] = float4(input.dispatchThreadID.xy * FORWARD_PLUS_TILE_SIZE, -1.0f, 1.0f);
// Top right point
screenSpace[1] = float4(float2(input.dispatchThreadID.x + 1, input.dispatchThreadID.y) * FORWARD_PLUS_TILE_SIZE, -1.0f, 1.0f);
// Bottom left point
screenSpace[2] = float4(float2(input.dispatchThreadID.x, input.dispatchThreadID.y + 1) * FORWARD_PLUS_TILE_SIZE, -1.0f, 1.0f);
// Bottom right point
screenSpace[3] = float4(float2(input.dispatchThreadID.x + 1, input.dispatchThreadID.y + 1) * FORWARD_PLUS_TILE_SIZE, -1.0f, 1.0f);
// View ray
screenSpace[4] = float4(float2(input.dispatchThreadID.x + .5, input.dispatchThreadID.y + .5) * FORWARD_PLUS_TILE_SIZE, -1.0f, 1.0f);
float3 viewSpace[5];
// Now convert the screen space points to view space
for (int i = 0; i < 5; i++)
{
viewSpace[i] = -ScreenToView(screenSpace[i]).xyz;
}
// Now build the frustum planes from the view space points
Frustum frustum;
// Left plane
frustum.planes[0] = ComputePlane(eyePos, viewSpace[2], viewSpace[0]);
// Right plane
frustum.planes[1] = ComputePlane(eyePos, viewSpace[1], viewSpace[3]);
// Top plane
frustum.planes[2] = ComputePlane(eyePos, viewSpace[0], viewSpace[1]);
// Bottom plane
frustum.planes[3] = ComputePlane(eyePos, viewSpace[3], viewSpace[2]);
frustum.m_view_ray = normalize(viewSpace[4]);
float cos_frustum_angle = dot(normalize(viewSpace[0]), normalize(viewSpace[4]));
float sin_frustum_angle = sqrt(1. - saturate(cos_frustum_angle * cos_frustum_angle));
frustum.m_tan_frustum_angle = sin_frustum_angle / cos_frustum_angle;
//frustum.cone = ComputeConeFrustum((float3)0., normalize(viewSpace[4]), abs());
// Store the computed frustum in global memory (if our thread ID is in bounds of the grid).
if (input.dispatchThreadID.x < numThreads.x && input.dispatchThreadID.y < numThreads.y)
{
uint index = input.dispatchThreadID.x + (input.dispatchThreadID.y * numThreads.x);
out_Frustums[index] = frustum;
}
}
/* ---------------------------------------------------------------------------------------------- */
/* Light culling */
/* ---------------------------------------------------------------------------------------------- */
float sdSphere( float3 pos, float s )
{
return length(pos)-s;
}
float sdSolidAngle(float3 pos, float2 c, float ra)
{
float2 p = float2( length(pos.xz), pos.y );
//float2 p = float2( length(pos.xz), pos.y );
float l = length(p) - ra;
float m = length(p - c*clamp(dot(p,c),0.0,ra) );
return max(l,m*sign(c.y*p.x-c.x*p.y));
}
float sdSolidAngleOutside(float3 pos, float2 c, float ra)
{
float2 p = float2( length(pos.xz), pos.y );
//float2 p = float2( length(pos.xz), pos.y );
float l = length(p) - ra;
float m = length(p - c*clamp(dot(p,c),0.0,ra) );
return max(l,-m*sign(c.y*p.x-c.x*p.y));
}
[numthreads(FORWARD_PLUS_TILE_SIZE, FORWARD_PLUS_TILE_SIZE, 1)]
void CS_CullLights(in ComputeShaderInput input) {
// Calculate min & max depth in threadgroup / tile.
int2 texCoord = input.dispatchThreadID.xy;
float fDepth = DepthTextureVS[texCoord].r;
uint uDepth = asuint(fDepth);
if (!input.groupIndex) // Avoid contention by other threads in the group.
{
uMinDepth = 0xffffffff;
uMaxDepth = 0;
o_LightCount = 0;
t_LightCount = 0;
GroupFrustum = in_Frustums[input.groupID.x + (input.groupID.y * numThreadGroups.x)];
}
GroupMemoryBarrierWithGroupSync();
InterlockedMin(uMinDepth, uDepth);
InterlockedMax(uMaxDepth, uDepth);
GroupMemoryBarrierWithGroupSync();
float fMinDepth = asfloat(uMaxDepth);
float fMaxDepth = asfloat(uMinDepth);
// Convert depth values to view space.
float minDepthVS = ClipToView(float4(0, 0, fMinDepth, 1)).z;
float maxDepthVS = ClipToView(float4(0, 0, fMaxDepth, 1)).z;
float nearClipVS = ClipToView(float4(0, 0, 1, 1)).z;
// Clipping plane for minimum depth value
// (used for testing lights within the bounds of opaque geometry).
Plane minPlane = { float3(0, 0, -1), -minDepthVS };
Plane farPlane = { float3(0, 0, 1), maxDepthVS };
Plane nearPlane = { float3(0, 0, -1), -nearClipVS };
for (uint i = input.groupIndex; i < NumLights; i += FORWARD_PLUS_TILE_SIZE * FORWARD_PLUS_TILE_SIZE) {
Light light = UnpackLight(Lights[i]);
float3 lightPosVS = mul(WorldToViewSpace, float4(light.m_origin, 1.)).xyz;
Sphere sphere = { lightPosVS, light.m_radius };
float tan_angle = GroupFrustum.m_tan_frustum_angle;
float3 rd = GroupFrustum.m_view_ray;
float dist = nearClipVS / rd.z;
float maxDist = maxDepthVS / rd.z;
float minDistOpaque = minDepthVS / rd.z;
bool intersection = false;
bool is_opaque = false;
if(light.m_cos_outer < -1.5) {
while(dist < maxDist) {
float3 ro = dist * rd;
float distance = sdSphere(ro - lightPosVS, light.m_radius);
if(distance <= dist * tan_angle) {
intersection = true;
is_opaque = dist > minDistOpaque;
break;
}
dist += distance;
}
if(!is_opaque) {
dist = minDistOpaque - 1.e-2;
while(dist < maxDist) {
float3 ro = dist * rd;
float distance = sdSphere(ro - lightPosVS, light.m_radius);
if(distance <= dist * tan_angle) {
is_opaque = true;
break;
}
dist += distance;
}
}
} else {
float3 lightDirVS = mul((float3x3)WorldToViewSpace, light.m_direction);
float3 tang = select(lightDirVS.z < .999, float3(0., 0., 1.), float3(0., 1., 0.));
tang = normalize(tang - lightDirVS * dot(tang, lightDirVS));
float3 bitang = cross(lightDirVS, tang);
float3x3 transform = {tang, lightDirVS, bitang};
float2 angle = {sqrt(1. - light.m_cos_outer * light.m_cos_outer), light.m_cos_outer};
while(dist < maxDist) {
float3 ro = dist * rd;
float distance = sdSolidAngle(mul(transform, ro - lightPosVS), angle, light.m_radius);
if(distance <= dist * tan_angle) {
intersection = true;
is_opaque = dist > minDistOpaque;
break;
}
dist += distance;
}
if(!is_opaque) {
dist = minDistOpaque - 1.e-2;
while(dist < maxDist) {
float3 ro = dist * rd;
float distance = sdSolidAngle(mul(transform, ro - lightPosVS), angle, light.m_radius);
if(distance <= dist * tan_angle) {
is_opaque = true;
break;
}
dist += distance;
}
}
}
if(intersection) {
t_AppendLight(i);
if (is_opaque)
{
// Add light to light list for opaque geometry.
o_AppendLight(i);
}
}
//continue;
//
//
//if(!DoQueryInfiniteCone(sphere, GroupFrustum.cone)) continue;
//if (dot(light.m_color, (float3)1.)) {
// float cosOuterCone = light.m_cos_outer;
//
// if (cosOuterCone == 0. || true) { // Point light
// if(SphereInsidePlane(sphere, nearPlane) || SphereInsidePlane(sphere, farPlane))
// continue;
//
// //if (SphereInsideFrustum(sphere, GroupFrustum, nearClipVS, maxDepthVS))
// //{
// // Add light to light list for transparent geometry.
// t_AppendLight(i);
//
// if (!SphereInsidePlane(sphere, minPlane))
// {
// // Add light to light list for opaque geometry.
// o_AppendLight(i);
// }
// //}
// }
// else { // Spot light
// float sinOuterCone = sqrt(saturate(1. - cosOuterCone * cosOuterCone));
// float tanOuterCone = sinOuterCone / cosOuterCone;
// float3 lightDirVS = mul((float3x3)WorldToViewSpace, light.m_direction);
// float coneRadius = tanOuterCone * light.m_radius;
// Cone cone = { lightPosVS, light.m_radius, lightDirVS, coneRadius };
// //if (SphereInsideFrustum(sphere, GroupFrustum, nearClipVS, maxDepthVS) &&
// // ConeInsideFrustum(cone, GroupFrustum, nearClipVS, maxDepthVS))
// //{
// // Add light to light list for transparent geometry.
// t_AppendLight(i);
//
// //if (!SphereInsidePlane(sphere, minPlane)
// // && !ConeInsidePlane(cone, minPlane))
// //{
// // Add light to light list for opaque geometry.
// o_AppendLight(i);
// //}
// //}
// }
//}
}
// Wait till all threads in group have caught up.
GroupMemoryBarrierWithGroupSync();
// Update global memory with visible light buffer.
// First update the light grid (only thread 0 in group needs to do this)
if (!input.groupIndex)
{
// Update light grid for opaque geometry.
InterlockedAdd(o_LightIndexCounter[0], o_LightCount, o_LightIndexStartOffset);
o_LightGrid[input.groupID.xy] = uint2(o_LightIndexStartOffset, o_LightCount);
// Update light grid for transparent geometry.
InterlockedAdd(t_LightIndexCounter[0], t_LightCount, t_LightIndexStartOffset);
t_LightGrid[input.groupID.xy] = uint2(t_LightIndexStartOffset, t_LightCount);
}
GroupMemoryBarrierWithGroupSync();
// Now update the light index list (all threads).
// For opaque geometry.
for (uint i = input.groupIndex; i < o_LightCount; i += FORWARD_PLUS_TILE_SIZE * FORWARD_PLUS_TILE_SIZE)
{
o_LightIndexList[o_LightIndexStartOffset + i] = o_LightList[i];
}
// For transparent geometry.
for (uint i = input.groupIndex; i < t_LightCount; i += FORWARD_PLUS_TILE_SIZE * FORWARD_PLUS_TILE_SIZE)
{
t_LightIndexList[t_LightIndexStartOffset + i] = t_LightList[i];
}
}
/* ---------------------------------------------------------------------------------------------- */
/* Function implementations */
/* ---------------------------------------------------------------------------------------------- */
// Add the light to the visible light list for opaque geometry.
void o_AppendLight(uint lightIndex)
{
uint index; // Index into the visible lights array.
InterlockedAdd(o_LightCount, 1, index);
if (index < FORWARD_PLUS_LIGHTS_PER_TILE)
{
o_LightList[index] = lightIndex;
}
}
// Add the light to the visible light list for transparent geometry.
void t_AppendLight(uint lightIndex)
{
uint index; // Index into the visible lights array.
InterlockedAdd(t_LightCount, 1, index);
if (index < FORWARD_PLUS_LIGHTS_PER_TILE)
{
t_LightList[index] = lightIndex;
}
}
// Convert clip space coordinates to view space
float4 ClipToView( float4 clip )
{
// View space position.
float4 view = mul( InverseProjection, clip );
// Perspective projection.
view = view / view.w;
return view;
}
// Convert screen space coordinates to view space.
float4 ScreenToView( float4 screen )
{
// Convert to normalized texture coordinates
float2 texCoord = screen.xy / ScreenDimensions;
// Convert to clip space
float4 clip = float4( float2( texCoord.x, 1.0f - texCoord.y ) * 2.0f - 1.0f, screen.z, screen.w );
return ClipToView( clip );
}

View File

@@ -0,0 +1,34 @@
#ifndef FORWARD_PLUS_LIGHT_HLSLI
#define FORWARD_PLUS_LIGHT_HLSLI
struct Light {
float3 m_origin;
float m_radius;
float3 m_direction;
float3 m_color;
float m_cos_inner;
float m_cos_outer;
};
struct PackedLight {
float3 m_origin;
float m_radius;
uint4 m_packed_data;
};
Light UnpackLight(in PackedLight packed_light) {
Light light;
light.m_origin = packed_light.m_origin;
light.m_radius = packed_light.m_radius;
light.m_direction.x = f16tof32(packed_light.m_packed_data[0]);
light.m_direction.y = f16tof32(packed_light.m_packed_data[0] >> 16);
light.m_direction.z = f16tof32(packed_light.m_packed_data[1]);
light.m_cos_inner = f16tof32(packed_light.m_packed_data[1] >> 16);
light.m_color.x = f16tof32(packed_light.m_packed_data[2]);
light.m_color.y = f16tof32(packed_light.m_packed_data[2] >> 16);
light.m_color.z = f16tof32(packed_light.m_packed_data[3]);
light.m_cos_outer = f16tof32(packed_light.m_packed_data[3] >> 16);
return light;
}
#endif

View File

@@ -0,0 +1,259 @@
// Partially stolen from https://www.3dgep.com/forward-plus/
#ifndef FORWARD_PLUS_PRIMITIVES_HLSLI
#define FORWARD_PLUS_PRIMITIVES_HLSLI
/* ---------------------------------------------------------------------------------------------- */
/* Plane */
/* ---------------------------------------------------------------------------------------------- */
struct Plane
{
float3 N; // Plane normal.
float d; // Distance to origin.
};
// Compute a plane from 3 noncollinear points that form a triangle.
// This equation assumes a right-handed (counter-clockwise winding order)
// coordinate system to determine the direction of the plane normal.
Plane ComputePlane( float3 p0, float3 p1, float3 p2 )
{
Plane plane;
float3 v0 = p1 - p0;
float3 v2 = p2 - p0;
plane.N = normalize( cross( v0, v2 ) );
// Compute the distance to the origin using p0.
plane.d = dot( plane.N, p0 );
return plane;
}
// Check to see if a point is fully behind (inside the negative halfspace of) a plane.
bool PointInsidePlane( float3 p, Plane plane )
{
return dot( plane.N, p ) - plane.d < 0;
}
/* ---------------------------------------------------------------------------------------------- */
/* Frustum */
/* ---------------------------------------------------------------------------------------------- */
struct ConeFrustum
{
float3 m_origin;
float m_cos_angle;
float3 m_direction;
float m_sin_angle;
float m_cos_angle_sqr;
float m_inv_sin_angle;
};
ConeFrustum ComputeConeFrustum(in float3 origin, in float3 direction, in float cos_angle) {
ConeFrustum cone;
cone.m_origin = origin;
cone.m_direction = direction;
cone.m_cos_angle = cos_angle;
cone.m_cos_angle_sqr = cos_angle * cos_angle;
cone.m_sin_angle = sqrt(1. - saturate(cone.m_cos_angle_sqr));
cone.m_inv_sin_angle = 1. / cone.m_sin_angle;
return cone;
}
// Four planes of a view frustum (in view space).
// The planes are:
// * Left,
// * Right,
// * Top,
// * Bottom.
// The back and/or front planes can be computed from depth values in the
// light culling compute shader.
struct Frustum
{
Plane planes[4]; // left, right, top, bottom frustum planes.
float3 m_view_ray;
float m_tan_frustum_angle;
};
/* ---------------------------------------------------------------------------------------------- */
/* Sphere */
/* ---------------------------------------------------------------------------------------------- */
struct Sphere
{
float3 m_origin; // Center point.
float m_radius; // Radius.
};
// Check to see if a sphere is fully behind (inside the negative halfspace of) a plane.
// Source: Real-time collision detection, Christer Ericson (2005)
bool SphereInsidePlane( Sphere sphere, Plane plane )
{
return dot( plane.N, sphere.m_origin ) - plane.d < -sphere.m_radius;
}
// Check to see of a light is partially contained within the frustum.
bool SphereInsideFrustum( Sphere sphere, Frustum frustum, float zNear, float zFar )
{
bool result = true;
// First check depth
// Note: Here, the view vector points in the -Z axis so the
// far depth value will be approaching -infinity.
if ( sphere.m_origin.z - sphere.m_radius > zNear || sphere.m_origin.z + sphere.m_radius < zFar )
{
result = false;
}
// Then check frustum planes
for ( int i = 0; i < 4 && result; ++i )
{
if ( SphereInsidePlane( sphere, frustum.planes[i] ) )
{
result = false;
}
}
return result;
}
/* ---------------------------------------------------------------------------------------------- */
/* Cone */
/* ---------------------------------------------------------------------------------------------- */
struct Cone
{
float3 T; // Cone tip.
float h; // Height of the cone.
float3 d; // Direction of the cone.
float r; // bottom radius of the cone.
};
// Check to see if a cone if fully behind (inside the negative halfspace of) a plane.
// Source: Real-time collision detection, Christer Ericson (2005)
bool ConeInsidePlane( Cone cone, Plane plane )
{
// Compute the farthest point on the end of the cone to the positive space of the plane.
float3 m = cross( cross( plane.N, cone.d ), cone.d );
float3 Q = cone.T + cone.d * cone.h - m * cone.r;
// The cone is in the negative halfspace of the plane if both
// the tip of the cone and the farthest point on the end of the cone to the
// positive halfspace of the plane are both inside the negative halfspace
// of the plane.
return PointInsidePlane( cone.T, plane ) && PointInsidePlane( Q, plane );
}
bool ConeInsideFrustum( Cone cone, Frustum frustum, float zNear, float zFar )
{
bool result = true;
Plane nearPlane = { float3( 0, 0, -1 ), -zNear };
Plane farPlane = { float3( 0, 0, 1 ), zFar };
// First check the near and far clipping planes.
if ( ConeInsidePlane( cone, nearPlane ) || ConeInsidePlane( cone, farPlane ) )
{
result = false;
}
// Then check frustum planes
for ( int i = 0; i < 4 && result; ++i )
{
if ( ConeInsidePlane( cone, frustum.planes[i] ) )
{
result = false;
}
}
return result;
}
/* ---------------------------------------------------------------------------------------------- */
/* AABB */
/* ---------------------------------------------------------------------------------------------- */
struct Aabb {
float3 min;
float3 max;
};
float3 GetAabbCorner( Aabb aabb, float4x4 transform, uint index ) {
float3 t = float3(
( index & 1 ) ? 1. : 0.,
( index & 2 ) ? 1. : 0.,
( index & 4 ) ? 1. : 0. );
return mul(
transform,
float4( ( 1. - t ) * aabb.min + t * aabb.max, 1. ) ).xyz;
}
bool AabbInsidePlane( Aabb aabb, float4x4 transform, Plane plane ) {
bool result = true;
for( uint i = 0; i < 8 && result; ++i ) {
result = PointInsidePlane( GetAabbCorner( aabb, transform, i ), plane );
}
return result;
}
bool AabbInsideFrustum( Aabb aabb, float4x4 transform, Frustum frustum, float zNear, float zFar )
{
bool result = true;
Plane nearPlane = { float3( 0, 0, -1 ), -zNear };
Plane farPlane = { float3( 0, 0, 1 ), zFar };
// First check the near and far clipping planes.
if ( AabbInsidePlane( aabb, transform, nearPlane ) || AabbInsidePlane( aabb, transform, farPlane ) )
{
result = false;
}
// Then check frustum planes
for ( int i = 0; i < 4 && result; ++i )
{
if ( AabbInsidePlane( aabb, transform, frustum.planes[i] ) )
{
result = false;
}
}
return result;
}
bool DoQueryInfiniteCone(in Sphere sphere, in ConeFrustum cone)
{
float3 U = cone.m_origin - (sphere.m_radius * cone.m_inv_sin_angle) * cone.m_direction;
float3 CmU = sphere.m_origin - U;
float AdCmU = dot(cone.m_direction, CmU);
if (AdCmU > 0.)
{
float sqrLengthCmU = dot(CmU, CmU);
if (AdCmU * AdCmU >= sqrLengthCmU * cone.m_cos_angle_sqr)
{
float3 CmV = sphere.m_origin - cone.m_origin;
float AdCmV = dot(cone.m_direction, CmV);
if (AdCmV < -sphere.m_radius)
{
return false;
}
float rSinAngle = sphere.m_radius * cone.m_sin_angle;
if (AdCmV >= -rSinAngle)
{
return true;
}
float sqrLengthCmV = dot(CmV, CmV);
return sqrLengthCmV <= sphere.m_radius * sphere.m_radius;
}
}
return false;
}
#endif

View File

@@ -0,0 +1,10 @@
#ifndef GBUFFER_CONTACT_SHADOWS_HLSLI
#define GBUFFER_CONTACT_SHADOWS_HLSLI
Texture2D<float> g_ContactShadows : register(t12);
float GetContactShadows(in uint2 pixel_position) {
return g_ContactShadows[pixel_position];
}
#endif

View File

@@ -0,0 +1,106 @@
#include "math.hlsli"
#include "lighting_functions.hlsli"
#include "color_transform.hlsli"
#include "manul/gbuffer_ssao.hlsli"
struct ViewConstants {
float4x4 m_inverse_model_view;
float4x4 m_inverse_projection;
};
struct PushConstantsSpotLight {
float3 m_origin;
float m_radius_sqr;
float3 m_direction;
float m_cos_inner_cone;
float3 m_color;
float m_cos_outer_cone;
};
cbuffer g_ViewConst : register(b0) { ViewConstants g_View; }
#ifdef SPIRV
[[vk::push_constant]] ConstantBuffer<PushConstantsSpotLight> g_SpotLight;
#else
cbuffer g_SpotLightConst : register(b1) { PushConstantsSpotLight g_SpotLight; }
#endif
RWTexture2D<float4> g_OutDiffuse : register(u0);
Texture2D<float> g_GbufferDepth : register(t0);
Texture2D<float3> g_GbufferNormal : register(t1);
Texture2D<float4> g_GbufferParams : register(t2);
Texture2D<float3> g_GbufferColor : register(t3);
float2 PixelToCS(in float2 pixel, in float2 size) {
return ((pixel + .5) / size - .5) * float2(2., -2.);
}
float3 ReconstructPos(in float2 cs, in float depth) {
float4 ndc = float4(cs, depth, 1.);
ndc = mul(g_View.m_inverse_projection, ndc);
return ndc.xyz / ndc.w;
}
[numthreads(8, 8, 1)]
void CS_Clear(uint3 PixCoord : SV_DispatchThreadID, uint3 GroupID : SV_GroupID, uint GroupIndex : SV_GroupIndex) {
uint2 pixel = PixCoord.xy;
g_OutDiffuse[pixel].rgb = 0.;
}
[numthreads(8, 8, 1)]
void CS_Spotlight(uint3 PixCoord : SV_DispatchThreadID, uint3 GroupID : SV_GroupID, uint GroupIndex : SV_GroupIndex) {
uint2 gbuffer_dimensions;
uint2 pixel = PixCoord.xy;
g_OutDiffuse.GetDimensions(gbuffer_dimensions.x, gbuffer_dimensions.y);
float3 position = ReconstructPos(PixelToCS(PixCoord.xy, gbuffer_dimensions), g_GbufferDepth[pixel]);
float3 view = mul((float3x3)g_View.m_inverse_model_view, position);
float3 L_offset = g_SpotLight.m_origin - view;
float3 L = normalize(L_offset);
if(dot(L_offset, L_offset) > g_SpotLight.m_radius_sqr) return;
float3 albedo = g_GbufferColor[pixel];
float4 params = g_GbufferParams[pixel];
float3 normal = UnpackNormalXYZ(g_GbufferNormal[pixel]);
SurfaceData surface_data;
surface_data.m_Albedo = albedo;
surface_data.m_Alpha = 1.;
surface_data.m_Metalness = params.r;
surface_data.m_Roughness = params.g;
#ifdef GBUFFER_SSAO_HLSLI
surface_data.m_DiffuseOcclusion = min(params.b, GetSSAO(pixel));
#else
surface_data.m_DiffuseOcclusion = params.b;
#endif
surface_data.m_Normal = mul((float3x3)g_View.m_inverse_model_view, normal);
surface_data.m_View = -normalize(view);
surface_data.m_Reflect = reflect(-surface_data.m_View, surface_data.m_Normal);
surface_data.m_NdotV = max(dot(surface_data.m_Normal, surface_data.m_View), 0.);
surface_data.m_SpecularF0 = float3(.04, .04, .04) * params.a * 2.;
surface_data.m_SpecularF0 = lerp(surface_data.m_SpecularF0, surface_data.m_Albedo, surface_data.m_Metalness);
surface_data.m_SpecularF = FresnelSchlickRoughness(surface_data.m_NdotV, surface_data.m_SpecularF0, surface_data.m_Roughness);
surface_data.m_HorizonFading = 1.6;
surface_data.m_SpecularOcclusion = ComputeSpecOcclusion(surface_data.m_NdotV, surface_data.m_DiffuseOcclusion, surface_data.m_Roughness);
float3 lit = float3(0., 0., 0.);
DirectionalLight directionalLight;
directionalLight.m_LightVector = L;
directionalLight.m_Color = max(REC709_to_XYZ(g_SpotLight.m_color), 0.);
ApplyDirectionalLight(lit, directionalLight, surface_data, 1.);
float cone = saturate((dot(L, -g_SpotLight.m_direction) - g_SpotLight.m_cos_inner_cone) / (g_SpotLight.m_cos_outer_cone - g_SpotLight.m_cos_inner_cone));
float attenuation = 1. / dot(L_offset, L_offset) * (1. - cone);
g_OutDiffuse[pixel].rgb += lit * attenuation;
}

View File

@@ -0,0 +1,27 @@
#ifndef GBUFFER_SSAO_HLSLI
#define GBUFFER_SSAO_HLSLI
Texture2D<float> g_SSAO : register(t5);
float4 R8G8B8A8_UNORM_to_FLOAT4( uint packedInput )
{
float4 unpackedOutput;
unpackedOutput.x = (float)( packedInput & 0x000000ff ) / 255;
unpackedOutput.y = (float)( ( ( packedInput >> 8 ) & 0x000000ff ) ) / 255;
unpackedOutput.z = (float)( ( ( packedInput >> 16 ) & 0x000000ff ) ) / 255;
unpackedOutput.w = (float)( packedInput >> 24 ) / 255;
return unpackedOutput;
}
void DecodeVisibilityBentNormal( const uint packedValue, out float visibility, out float3 bentNormal )
{
float4 decoded = R8G8B8A8_UNORM_to_FLOAT4( packedValue );
bentNormal = decoded.xyz * 2.0.xxx - 1.0.xxx; // could normalize - don't want to since it's done so many times, better to do it at the final step only
visibility = decoded.w;
}
float GetSSAO(in uint2 pixel_position) {
return g_SSAO[pixel_position];
}
#endif

View File

@@ -0,0 +1,124 @@
#ifndef LIGHTING_HLSLI
#define LIGHTING_HLSLI
#include "math.hlsli"
#include "material_common.hlsli"
#include "view_data.hlsli"
#include "lighting_functions.hlsli"
#include "forward_plus/light.hlsli"
SamplerState g_SamplerLinearClamp : register(s8);
TextureCube<float3> g_DiffuseEnvmap : register(t8);
TextureCube<float3> g_SpecularEnvmap : register(t9);
Texture2D<float2> g_BrdfLUT : register(t10);
Texture2D<uint2> g_LightGrid : register(t16);
StructuredBuffer<uint> g_LightIndexBuffer : register(t17);
StructuredBuffer<PackedLight> g_LightBuffer : register(t18);
#define LIGHTING_NEEDS_PIXELPOSITION 1
#if LIGHTING_NEEDS_PIXELPOSITION
void ApplyMaterialLighting(out float4 lit, in MaterialData material, in uint2 pixel_position);
#else
void ApplyMaterialLighting(out float4 lit, in MaterialData material);
#endif
void ApplyIBL(inout float3 lit, in SurfaceData material);
void ApplyForwardPlus(inout float3 lit, in SurfaceData material, in float3 position, in uint2 pixel_position);
#if LIGHTING_NEEDS_PIXELPOSITION
void ApplyMaterialLighting(out float4 lit, in MaterialData material, in uint2 pixel_position)
#else
void ApplyMaterialLighting(out float4 lit, in MaterialData material)
#endif
{
// Camera-centric world position
float3 view = mul((float3x3)g_InverseModelView, material.m_Position);
// Convert material to surface data
SurfaceData surface_data;
surface_data.m_Albedo = clamp(material.m_MaterialAlbedoAlpha.rgb, .02, .9);
surface_data.m_Alpha = material.m_MaterialAlbedoAlpha.a;
surface_data.m_Metalness = material.m_MaterialParams.r;
surface_data.m_Roughness = material.m_MaterialParams.g;
#ifdef GBUFFER_SSAO_HLSLI
surface_data.m_DiffuseOcclusion = min(material.m_MaterialParams.b, GetSSAO(pixel_position));
#else
surface_data.m_DiffuseOcclusion = material.m_MaterialParams.b;
#endif
surface_data.m_Normal = mul((float3x3)g_InverseModelView, material.m_MaterialNormal);
surface_data.m_View = -normalize(view);
surface_data.m_Reflect = reflect(-surface_data.m_View, surface_data.m_Normal);
surface_data.m_NdotV = max(dot(surface_data.m_Normal, surface_data.m_View), 0.);
surface_data.m_SpecularF0 = float3(.04, .04, .04) * material.m_MaterialParams.a * 2.;
surface_data.m_SpecularF0 = lerp(surface_data.m_SpecularF0, surface_data.m_Albedo, surface_data.m_Metalness);
surface_data.m_SpecularF = FresnelSchlickRoughness(surface_data.m_NdotV, surface_data.m_SpecularF0, surface_data.m_Roughness);
surface_data.m_HorizonFading = 1.6;
surface_data.m_SpecularOcclusion = ComputeSpecOcclusion(surface_data.m_NdotV, surface_data.m_DiffuseOcclusion, surface_data.m_Roughness);
lit.rgb = material.m_MaterialEmission * surface_data.m_Alpha;
lit.a = 1. - (saturate((1. - surface_data.m_Alpha) * lerp(1. - dot(surface_data.m_SpecularF, 1./3.), 0., surface_data.m_Metalness)));
float shadow = 1.;
#ifdef SHADOW_HLSLI
shadow = GetShadow(view, material);
#endif
#ifdef GBUFFER_CONTACT_SHADOWS_HLSLI
shadow = min(shadow, GetContactShadows(pixel_position));
#endif
// Apply IBL cubemap
ApplyIBL(lit.rgb, surface_data);
// Apply directional light
DirectionalLight directionalLight;
directionalLight.m_LightVector = g_LightDir.xyz;
directionalLight.m_Color = g_LightColor.rgb;
ApplyDirectionalLight(lit.rgb, directionalLight, surface_data, shadow);
ApplyForwardPlus(lit.rgb, surface_data, view, pixel_position);
}
void ApplyIBL(inout float3 lit, in SurfaceData material) {
float3 kS = material.m_SpecularF;
float3 kD = 1. - kS;
kD = lerp(kD, 0., material.m_Metalness);
float2 brdf = g_BrdfLUT.SampleLevel(g_SamplerLinearClamp, float2(material.m_NdotV, material.m_Roughness), 0.);
lit += material.m_Albedo * kD * g_DiffuseEnvmap.SampleLevel(g_SamplerLinearClamp, material.m_Normal, 0.) * material.m_DiffuseOcclusion * material.m_Alpha;
lit += g_SpecularEnvmap.SampleLevel(g_SamplerLinearClamp, material.m_Reflect, sqrt(material.m_Roughness) * 5.) * (material.m_SpecularF * brdf.xxx + brdf.yyy) * material.m_SpecularOcclusion * HorizonFading(material.m_NdotV, material.m_HorizonFading);
}
void ApplyForwardPlus(inout float3 lit, in SurfaceData material, in float3 position, in uint2 pixel_position) {
uint2 light_range = g_LightGrid[pixel_position / 8];
for(int i = 0; i < light_range.y; ++i) {
Light light = UnpackLight(g_LightBuffer[g_LightIndexBuffer[light_range.x + i]]);
float3 L_offset = light.m_origin - position;
float L_dist_sqr = dot(L_offset, L_offset);
float3 L_dir = normalize(L_offset);
//lit += light.m_color;
if(L_dist_sqr > light.m_radius * light.m_radius) continue;
float cos_point_deviation = dot(-L_dir, light.m_direction);
if(light.m_cos_outer > -1.5 && cos_point_deviation < light.m_cos_outer) continue;
float cone = light.m_cos_outer > -1.5 ? saturate((cos_point_deviation - light.m_cos_inner) / (light.m_cos_outer - light.m_cos_inner)) : 0.;
float attenuation = 1. / dot(L_offset, L_offset) * (1. - cone);
DirectionalLight directionalLight;
directionalLight.m_LightVector = L_dir;
directionalLight.m_Color = light.m_color * attenuation;
ApplyDirectionalLight(lit.rgb, directionalLight, material, 1.);
}
}
#endif

View File

@@ -0,0 +1,122 @@
#ifndef LIGHTING_FUNCTIONS_HLSLI
#define LIGHTING_FUNCTIONS_HLSLI
struct SurfaceData {
float3 m_Albedo;
float m_Alpha;
float m_Metalness;
float m_Roughness;
float m_DiffuseOcclusion;
float m_SpecularOcclusion;
float3 m_View;
float3 m_Normal;
float3 m_Reflect;
float m_NdotV;
float3 m_SpecularF0;
float3 m_SpecularF;
float m_HorizonFading;
};
struct DirectionalLight {
float3 m_LightVector;
float3 m_Color;
};
void ApplyDirectionalLight(inout float3 lit, in DirectionalLight light, in SurfaceData material, in float shadow);
float3 FresnelSchlickRoughness(float cosTheta, float3 F0, float roughness);
float FresnelSchlickRoughness(float cosTheta, float F0, float roughness);
float DistributionGGX(float3 N, float3 H, float roughness);
float GeometrySchlickGGX(float NdotV, float roughness);
float GeometrySmith(float3 N, float3 V, float3 L, float roughness);
float ComputeSpecOcclusion(float NdotV, float AO, float roughness);
float HorizonFading(float NdotL, float horizonFade);
void ApplyDirectionalLight(inout float3 lit, in DirectionalLight light, in SurfaceData material, in float shadow) {
float3 H = normalize(material.m_View + light.m_LightVector);
float3 radiance = light.m_Color;
float NdotL = max(dot(material.m_Normal, light.m_LightVector), 0.);
// Cook-Torrance BRDF
float NDF = DistributionGGX(material.m_Normal, H, material.m_Roughness);
float G = GeometrySmith(material.m_Normal, material.m_View, light.m_LightVector, material.m_Roughness);
float3 numerator = NDF * G * material.m_SpecularF;
float denominator = 4. * material.m_NdotV * NdotL + 1e-4; // + 0.0001 to prevent divide by zero
float3 specular = numerator / denominator;
// kS is equal to Fresnel
float3 kS = material.m_SpecularF;
// for energy conservation, the diffuse and specular light can't
// be above 1.0 (unless the surface emits light); to preserve this
// relationship the diffuse component (kD) should equal 1.0 - kS.
float3 kD = 1. - kS;
// multiply kD by the inverse metalness such that only non-metals
// have diffuse lighting, or a linear blend if partly metal (pure metals
// have no diffuse light).
kD = lerp(kD, 0., material.m_Metalness);
// add to outgoing radiance Lo
//lit += shadow * (kD * material.m_Albedo * ONE_OVER_PI * material.m_DiffuseOcclusion * material.m_Alpha + specular * material.m_SpecularOcclusion * HorizonFading(NdotL, material.m_HorizonFading)) * radiance * NdotL; // note that we already multiplied the BRDF by the Fresnel (kS) so we won't multiply by kS again
lit += shadow * (kD * material.m_Albedo * ONE_OVER_PI * material.m_Alpha + specular * HorizonFading(NdotL, material.m_HorizonFading)) * radiance * NdotL; // note that we already multiplied the BRDF by the Fresnel (kS) so we won't multiply by kS again
}
float FresnelSchlickRoughness(float cosTheta, float F0, float roughness)
{
return F0 + (max(1. - roughness, F0) - F0) * pow(saturate(1. - cosTheta), 5.);
}
float3 FresnelSchlickRoughness(float cosTheta, float3 F0, float roughness)
{
return F0 + (max((float3)(1. - roughness), F0) - F0) * pow(saturate(1. - cosTheta), 5.);
}
float DistributionGGX(float3 N, float3 H, float roughness)
{
float a = roughness*roughness;
float a2 = a*a;
float NdotH = max(dot(N, H), 0.);
float NdotH2 = NdotH*NdotH;
float num = a2;
float denom = (NdotH2 * (a2 - 1.) + 1.);
denom = PI * denom * denom;
return num / denom;
}
float GeometrySchlickGGX(float NdotV, float roughness)
{
float r = (roughness + 1.);
float k = (r*r) / 8.;
float num = NdotV;
float denom = NdotV * (1. - k) + k;
return num / denom;
}
float GeometrySmith(float3 N, float3 V, float3 L, float roughness)
{
float NdotV = max(dot(N, V), 0.);
float NdotL = max(dot(N, L), 0.);
float ggx2 = GeometrySchlickGGX(NdotV, roughness);
float ggx1 = GeometrySchlickGGX(NdotL, roughness);
return ggx1 * ggx2;
}
float ComputeSpecOcclusion(float NdotV, float AO, float roughness)
{
return saturate(pow(NdotV + AO, exp2(-16. * roughness - 1.)) - 1. + AO);
}
float HorizonFading(float NdotL, float horizonFade)
{
float horiz = saturate(1.0 + horizonFade * NdotL);
return horiz * horiz;
}
#endif

View File

@@ -0,0 +1,113 @@
#include "color_transform.hlsli"
struct VertexInput {
float3 m_Position : Position;
float3 m_Normal : Normal;
float2 m_TexCoord : TexCoord;
float4 m_Tangent : Tangent;
};
struct VertexOutput {
float3 m_Position : Position;
};
struct PixelInput {
float3 m_Position : Position;
float3 m_Normal : Normal;
float4 m_PositionSV : SV_Position;
float4 m_PositionCS : PositionCS;
float4 m_HistoryPositionCS : HistoryPositionCS;
};
struct PixelOutput {
float4 m_Color : SV_Target0;
float4 m_Emission : SV_Target1;
float4 m_Params : SV_Target2;
float4 m_Normal : SV_Target3;
float2 m_Motion : SV_Target4;
};
cbuffer VertexConstants : register(b0) {
float4x4 g_JitteredProjection;
float4x4 g_Projection;
float4x4 g_ProjectionHistory;
};
struct LinePushConstants {
float4x3 m_ModelView;
float4x3 m_ModelViewHistory;
float3 m_Color;
float m_LineWeight;
float m_Metalness;
float m_Roughness;
};
#ifdef SPIRV
[[vk::push_constant]] ConstantBuffer<LinePushConstants> g_DrawConstants;
#else
cbuffer g_Const : register(b1) { LinePushConstants g_DrawConstants; }
#endif
VertexOutput vtx_main(in VertexInput vs_in) {
VertexOutput result;
result.m_Position = vs_in.m_Position;
return result;
}
#define NUM_LINE_VERTS 4
static const float3 g_VertexNormals[NUM_LINE_VERTS] = {
float3(-1., 0., 0.),
float3(0., 1., 0.),
float3(1., 0., 0.),
float3(0., -1., 0.) };
PixelInput ComputePoint(in float3 position, in float3 normal);
[maxvertexcount(NUM_LINE_VERTS * 6)]
void geo_main(line VertexOutput gs_in[2], inout TriangleStream<PixelInput> gs_out) {
float3 normal = normalize(gs_in[1].m_Position - gs_in[0].m_Position);
float3 up = normal.z < 0.99984769515 ? float3(0., 0., 1.) : float3(0., 1., 0.);
float3 tangent = normalize(cross(normal, up));
up = cross(tangent, normal);
float3x3 offset_matrix = float3x3(tangent, up, normal);
for(int i = 0; i < NUM_LINE_VERTS; ++i) {
float3 normal = mul(g_VertexNormals[i], offset_matrix);
float3 next_normal = mul(g_VertexNormals[(i + 1) % NUM_LINE_VERTS], offset_matrix);
float3 offset = normal * .5 * g_DrawConstants.m_LineWeight;
float3 next_offset = next_normal * .5 * g_DrawConstants.m_LineWeight;
gs_out.Append(ComputePoint(gs_in[0].m_Position + offset, normal));
gs_out.Append(ComputePoint(gs_in[0].m_Position + next_offset, next_normal));
gs_out.Append(ComputePoint(gs_in[1].m_Position + offset, normal));
gs_out.Append(ComputePoint(gs_in[1].m_Position + next_offset, next_normal));
gs_out.RestartStrip();
}
};
PixelOutput pix_main(in PixelInput ps_in) {
PixelOutput result;
result.m_Color.rgb = REC709_to_XYZ(g_DrawConstants.m_Color);
result.m_Color.a = 1.;
result.m_Emission = 0.;
result.m_Params = float4(g_DrawConstants.m_Metalness, g_DrawConstants.m_Roughness, 1., .5);
result.m_Normal.rgb = ps_in.m_Normal * .5 + .5;
result.m_Normal.a = 1.;
result.m_Motion = (ps_in.m_HistoryPositionCS.xy / ps_in.m_HistoryPositionCS.w) - (ps_in.m_PositionCS.xy / ps_in.m_PositionCS.w);
result.m_Motion = result.m_Motion * float2(.5, -.5);
return result;
}
PixelInput ComputePoint(in float3 position, in float3 normal) {
PixelInput result;
result.m_Position = mul(float4(position, 1.), g_DrawConstants.m_ModelView).xyz;
result.m_Normal = mul(float4(normal, 0.), g_DrawConstants.m_ModelView).xyz;
result.m_PositionSV = mul(g_JitteredProjection, float4(result.m_Position, 1.));
result.m_PositionCS = mul(g_Projection, float4(result.m_Position, 1.));
float3 history_position = mul(float4(position, 1.), g_DrawConstants.m_ModelViewHistory);
result.m_HistoryPositionCS = mul(g_ProjectionHistory, float4(history_position, 1.));
return result;
}

View File

@@ -0,0 +1 @@
#define DECLARE_PUSH_CONSTANTS(sig, slot)

View File

@@ -0,0 +1,102 @@
#ifndef MATERIAL_HLSLI
#define MATERIAL_HLSLI
#define MOTION_VECTORS (1 << 0)
#define FORWARD_LIGHTING (1 << 1)
#define PASS_GBUFFER (MOTION_VECTORS)
#define PASS_FORWARD (FORWARD_LIGHTING)
#define PASS_CUBEMAP (0)
struct PixelOutput {
#if PASS & FORWARD_LIGHTING
float4 m_Color : SV_Target0;
float m_CompositionMask : SV_Target1;
#if PASS & MOTION_VECTORS
float2 m_Motion : SV_Target2;
#endif
#else
float4 m_Color : SV_Target0;
float4 m_Emission : SV_Target1;
float4 m_Params : SV_Target2;
float4 m_Normal : SV_Target3;
#if PASS & MOTION_VECTORS
float2 m_Motion : SV_Target4;
#endif
#endif
};
struct PixelInput {
float3 m_Position : Position;
float3 m_Normal : Normal;
float2 m_TexCoord : TexCoord;
float4 m_Tangent : Tangent;
float4 m_PositionSV : SV_Position;
float4 m_PositionCS : PositionCS;
#if PASS & MOTION_VECTORS
float4 m_HistoryPositionCS : HistoryPositionCS;
#endif
};
#include "draw_constants.hlsli"
#include "material_common.hlsli"
#if PASS & FORWARD_LIGHTING
#include "shadow.hlsli"
#include "lighting.hlsli"
#include "sky.hlsli"
Texture2D<float> g_GbufferDepth : register(t12);
#endif
void MaterialPass(inout MaterialData material);
PixelOutput main(in PixelInput ps_in) {
float2 uv = (ps_in.m_PositionCS.xy/ps_in.m_PositionCS.w) * float2(.5, -.5) + .5;
MaterialData material;
material.m_Position = ps_in.m_Position;
material.m_Normal = ps_in.m_Normal;
material.m_Tangent = ps_in.m_Tangent.xyz;
material.m_Bitangent = ps_in.m_Tangent.w * cross(ps_in.m_Normal, ps_in.m_Tangent.xyz);
material.m_TexCoord = ps_in.m_TexCoord;
material.m_PixelCoord = ps_in.m_PositionSV.xy;
material.m_PositionNDC = ps_in.m_PositionCS / ps_in.m_PositionCS.w;
material.m_MaterialAlbedoAlpha = float4(1., 1., 1., 1.);
material.m_MaterialEmission = float3(0., 0., 0.);
material.m_MaterialParams = float4(0., .5, 1., .5); // Metalness.Roughness.Occlusion.Specular
material.m_MaterialNormal = float3(0., 0., 1.);
MaterialPass(material);
material.m_MaterialAlbedoAlpha.rgb = saturate(material.m_MaterialAlbedoAlpha.rgb);
material.m_MaterialEmission = max(material.m_MaterialEmission, 0.);
PixelOutput ps_out;
#if PASS & FORWARD_LIGHTING
#if LIGHTING_NEEDS_PIXELPOSITION
// Should not really happen?
ApplyMaterialLighting(ps_out.m_Color, material, (uint2)ps_in.m_PositionSV.xy);
#else
ApplyMaterialLighting(ps_out.m_Color, material);
#endif
float3 view = mul((float3x3)g_InverseModelView, material.m_Position);
ApplyAerialPerspective(ps_out.m_Color.rgb, ps_out.m_Color.a, normalize(view), g_LightDir, length(material.m_Position)/2500.);
//ps_out.m_Color.rgb += CalcAtmosphere(normalize(view), g_LightDir, g_Altitude, length(view)) * material.m_MaterialAlbedoAlpha.a;
ps_out.m_CompositionMask = ps_out.m_Color.a > 0.;
#else
ps_out.m_Color.rgb = material.m_MaterialAlbedoAlpha.rgb;
ps_out.m_Emission.rgb = material.m_MaterialEmission;
ps_out.m_Params = material.m_MaterialParams;
ps_out.m_Normal.rgb = material.m_MaterialNormal * .5 + .5;
#endif
#if PASS & MOTION_VECTORS
ps_out.m_Motion = (ps_in.m_HistoryPositionCS.xy / ps_in.m_HistoryPositionCS.w) - (ps_in.m_PositionCS.xy / ps_in.m_PositionCS.w);
ps_out.m_Motion = ps_out.m_Motion * float2(.5, -.5);
#endif
return ps_out;
}
#endif

View File

@@ -0,0 +1,20 @@
#ifndef MATERIAL_COMMON_HLSLI
#define MATERIAL_COMMON_HLSLI
struct MaterialData {
float3 m_Position;
float3 m_PositionDDX;
float3 m_PositionDDY;
float3 m_Tangent;
float3 m_Bitangent;
float3 m_Normal;
float2 m_TexCoord;
uint2 m_PixelCoord;
float4 m_PositionNDC;
float4 m_MaterialAlbedoAlpha;
float3 m_MaterialEmission;
float4 m_MaterialParams; // Metalness.Roughness.Occlusion.Specular
float3 m_MaterialNormal;
};
#endif

View File

@@ -0,0 +1,28 @@
#ifndef MATH_HLSLI
#define MATH_HLSLI
#define PI 3.1415926535897932384626433832795
#define TWO_PI 6.283185307179586476925286766559
#define PI_OVER_TWO 1.5707963267948966192313216916398
#define PI_OVER_THREE 1.0471975511965977461542144610932
#define PI_OVER_SIX 0.52359877559829887307710723054658
#define ONE_OVER_PI 0.31830988618379067153776752674503
#define TWO_OVER_PI 0.63661977236758134307553505349006
#define ONE_OVER_TWO_PI 0.15915494309189533576888376337251
float3 UnpackNormalXY(in float2 normal) {
float3 result;
result.xy = 2. * (normal - .5);
result.z = sqrt(saturate(1. - dot(result.xy, result.xy)));
return result;
}
float3 UnpackNormalXYZ(in float3 normal) {
return normalize(normal - .5);
}
float InvLerp(float a, float b, float v) {
return (v - a) / (b - a);
}
#endif

View File

@@ -0,0 +1,87 @@
#ifndef SHADOW_HLSLI
#define SHADOW_HLSLI
#include "view_data.hlsli"
#include "material_common.hlsli"
#define CASCADE_COUNT 8
cbuffer ShadowProjections : register(b11) {
float4x4 g_CascadeProjections[CASCADE_COUNT];
float4 g_CascadePlanes[(CASCADE_COUNT + 3) / 4];
}
SamplerComparisonState g_ShadowSampler : register(s11);
Texture2DArray<float> g_ShadowMap : register(t11);
float3 PosToShadow(in float4x4 cascade, in float3 pos) {
float4 shadow_view = mul(cascade, float4(pos, 1.));
shadow_view /= shadow_view.w;
shadow_view.xy = shadow_view.xy * float2(.5, -.5) + .5;
return shadow_view.xyz;
}
bool CalculateShadowView(in float4x4 cascade, in float plane, in MaterialData data, out float3 shadow_view, out float3 shadow_view_ddx, out float3 shadow_view_ddy) {
#ifdef DEFERRED_LIGHTING_PASS
if(dot(data.m_Position, data.m_Position) > plane) return false;
shadow_view = PosToShadow(cascade, data.m_Position);
shadow_view_ddx = PosToShadow(cascade, data.m_PositionDDX) - shadow_view;
shadow_view_ddy = PosToShadow(cascade, data.m_PositionDDY) - shadow_view;
return true;
#else
if(dot(data.m_Position, data.m_Position) > plane) return false;
shadow_view = PosToShadow(cascade, data.m_Position);
//if(any(shadow_view <= 0. || shadow_view >= 1.)) return false;
//shadow_view_ddx = PosToShadow(cascade, data.m_Position + ddx(data.m_Position)) - shadow_view;
//shadow_view_ddy = PosToShadow(cascade, data.m_Position + ddy(data.m_Position)) - shadow_view;
return true;
#endif
}
float GetShadow(in float3 view, in MaterialData data) {
float shadow = 1.;
uint3 shadowmap_size;
g_ShadowMap.GetDimensions(shadowmap_size.x, shadowmap_size.y, shadowmap_size.z);
for(int i = 0; i < shadowmap_size.z; ++i) {
float3 shadowView, vShadowTexDDX, vShadowTexDDY;
if(CalculateShadowView(mul(g_CascadeProjections[i], g_InverseModelView), g_CascadePlanes[i / 4][i % 4], data, shadowView, vShadowTexDDX, vShadowTexDDY)) {
shadow = 0.;
#ifdef DEFERRED_LIGHTING_PASS
float2x2 screenToShadow = float2x2(vShadowTexDDX.xy, vShadowTexDDY.xy);
float fDeterminant = determinant(screenToShadow);
float fInvDeterminant = 1.0f / fDeterminant;
float2x2 shadowToScreen = float2x2 (
screenToShadow._22 * fInvDeterminant,
screenToShadow._12 * -fInvDeterminant,
screenToShadow._21 * -fInvDeterminant,
screenToShadow._11 * fInvDeterminant );
float2 vRightTexelDepthRatio = mul( float2(1./shadowmap_size.x, 0.),
shadowToScreen );
float2 vUpTexelDepthRatio = mul( float2(0., 1./shadowmap_size.y),
shadowToScreen );
float fUpTexelDepthDelta =
vUpTexelDepthRatio.x * vShadowTexDDX.z
+ vUpTexelDepthRatio.y * vShadowTexDDY.z;
float fRightTexelDepthDelta =
vRightTexelDepthRatio.x * vShadowTexDDX.z
+ vRightTexelDepthRatio.y * vShadowTexDDY.z;
#endif
[unroll]
for (int y = -1; y <= 1; ++y) {
[unroll]
for (int x = -1; x <= 1; ++x) {
#ifdef DEFERRED_LIGHTING_PASS
shadow += g_ShadowMap.SampleCmpLevelZero(g_ShadowSampler, float3(shadowView.xy, float(i)), shadowView.z + x * fRightTexelDepthDelta + y * fUpTexelDepthDelta, int2(x, y)) / 9.;
#else
shadow += g_ShadowMap.SampleCmpLevelZero(g_ShadowSampler, float3(shadowView.xy, float(i)), shadowView.z, int2(x, y)) / 9.;
#endif
}
}
break;
}
}
return shadow;
}
#endif

View File

@@ -0,0 +1,429 @@
#ifndef SKY_HLSLI
#define SKY_HLSLI
#include "math.hlsli"
#define SKY_INF 1.#INF
static const float3x3 CAS2RGB = float3x3(
1.6218, -0.4493, 0.0325,
-0.0374, 1.0598, -0.0742,
-0.0283, -0.1119, 1.0491);
static const float3x3 RGB2CAS = float3x3(
6.2267e-1, 2.6392e-1, -6.2375e-4,
2.3324e-2, 9.6056e-1, 6.7215e-2,
1.9285e-2, 1.0958e-1, 9.6035e-1);
struct Atmosphere
{
float3 sunDirection;
float atmosphereRadius;
float earthRadius;
float Hr;
float Hm;
float3 betaR;
float3 betaM;
float3 irradiance;
};
void swap(inout float a, inout float b)
{
float temp = a;
a = b;
b = temp;
}
bool solveQuadratic(in float a, in float b, in float c, out float x1, out float x2)
{
if (b == 0)
{
// Handle special case where the the two vector ray.dir and V are perpendicular
// with V = ray.orig - sphere.centre
if (a == 0)
return false;
x1 = 0;
x2 = sqrt(-c / a);
return true;
}
float discr = b * b - 4 * a * c;
if (discr < 0)
return false;
float q = (b < 0.f) ? -0.5f * (b - sqrt(discr)) : -0.5f * (b + sqrt(discr));
x1 = q / a;
x2 = c / q;
return true;
}
bool solveQuadratic(in float a, in float b, in float c)
{
if (b == 0)
{
// Handle special case where the the two vector ray.dir and V are perpendicular
// with V = ray.orig - sphere.centre
return a != 0;
}
return b * b - 4 * a * c >= 0;
}
bool RaySphereIntersect(in float3 orig, in float3 dir, in float radius, out float t0, out float t1)
{
// They ray dir is normalized so A = 1
float A = dot(dir, dir);
float B = 2 * dot(dir, orig);
float C = dot(orig, orig) - radius * radius;
if (!solveQuadratic(A, B, C, t0, t1))
return false;
if (t0 > t1)
swap(t0, t1);
return true;
}
bool RaySphereIntersect(in float3 orig, in float3 dir, in float radius)
{
// They ray dir is normalized so A = 1
float A = dot(dir, dir);
float B = 2 * dot(dir, orig);
float C = dot(orig, orig) - radius * radius;
return solveQuadratic(A, B, C);
}
float2x3 ComputeIncidentLight(in Atmosphere atmosphere, in float3 orig, in float3 dir, float tmin, float tmax)
{
float t0, t1;
if (!RaySphereIntersect(orig, dir, atmosphere.atmosphereRadius, t0, t1) || t1 < 0)
return 0;
if (t0 > tmin && t0 > 0)
tmin = t0;
if (t1 < tmax)
tmax = t1;
uint numSamples = 16;
uint numSamplesLight = 8;
float segmentLength = (tmax - tmin) / numSamples;
float tCurrent = tmin;
float3 sumR = 0.;
float3 sumM = 0.; //mie and rayleigh contribution
float opticalDepthR = 0, opticalDepthM = 0;
float mu = dot(dir, atmosphere.sunDirection); //mu in the paper which is the cosine of the angle between the sun direction and the ray direction
float phaseR = 3.f / (16.f * PI) * (1 + mu * mu);
float g = 0.76f;
float phaseM = 3.f / (8.f * PI) * ((1.f - g * g) * (1.f + mu * mu)) / ((2.f + g * g) * pow(1.f + g * g - 2.f * g * mu, 1.5f));
for (uint i = 0; i < numSamples; ++i)
{
float3 samplePosition = orig + (tCurrent + segmentLength * 0.5f) * dir;
float height = length(samplePosition) - atmosphere.earthRadius;
// compute optical depth for light
float hr = exp(-height / atmosphere.Hr) * segmentLength;
float hm = exp(-height / atmosphere.Hm) * segmentLength;
opticalDepthR += hr;
opticalDepthM += hm;
// light optical depth
float t0Light, t1Light;
RaySphereIntersect(samplePosition, atmosphere.sunDirection, atmosphere.atmosphereRadius, t0Light, t1Light);
float segmentLengthLight = t1Light / numSamplesLight, tCurrentLight = 0;
float opticalDepthLightR = 0, opticalDepthLightM = 0;
uint j;
for (j = 0; j < numSamplesLight; ++j)
{
float3 samplePositionLight = samplePosition + (tCurrentLight + segmentLengthLight * 0.5f) * atmosphere.sunDirection;
float heightLight = length(samplePositionLight) - atmosphere.earthRadius;
if (heightLight < 0)
break;
opticalDepthLightR += exp(-heightLight / atmosphere.Hr) * segmentLengthLight;
opticalDepthLightM += exp(-heightLight / atmosphere.Hm) * segmentLengthLight;
tCurrentLight += segmentLengthLight;
}
if (j == numSamplesLight)
{
float3 tau = atmosphere.betaR * (opticalDepthR + opticalDepthLightR) + atmosphere.betaM * 1.1f * (opticalDepthM + opticalDepthLightM);
float3 attenuation = exp(-tau);
sumR += attenuation * hr;
sumM += attenuation * hm;
}
tCurrent += segmentLength;
}
// We use a magic number here for the intensity of the sun (20). We will make it more
// scientific in a future revision of this lesson/code
float3 tau = atmosphere.betaR * opticalDepthR + atmosphere.betaM * 1.1f * opticalDepthM;
return float2x3((sumR * atmosphere.betaR * phaseR + sumM * atmosphere.betaM * phaseM) * 20, exp(-tau));
}
//#define MIE_G 0.76
//#define SQR_G (MIE_G * MIE_G)
//#define RAYLEIGH_SCALE 8e3 /* Rayleigh scale height (m). */
//#define MIE_SCALE 1.2e3 /* Mie scale height (m). */
//float phase_rayleigh(float mu)
//{
//return 3.0f / (16.0f * PI) * (1.0f + mu * mu);
//}
//float phase_mie(float mu)
//{
//return (3.0f * (1.0f - SQR_G) * (1.0f + mu * mu)) /
//(8.0f * PI * (2.0f + SQR_G) * pow((1.0f + SQR_G - 2.0f * MIE_G * mu), 1.5));
//}
//float density_rayleigh(float height)
//{
//return exp(-height / RAYLEIGH_SCALE);
//}
//float density_mie(float height)
//{
//return exp(-height / MIE_SCALE);
//}
//float density_ozone(float height)
//{
//float den = 0.0f;
//if (height >= 10000.0f && height < 25000.0f)
//{
//den = 1.0f / 15000.0f * height - 2.0f / 3.0f;
//}
//else if (height >= 25000 && height < 40000)
//{
//den = -(1.0f / 15000.0f * height - 8.0f / 3.0f);
//}
//return den;
//}
///* Parameters for optical depth quadrature.
//* See the comment in ray_optical_depth for more detail.
//* Computed using sympy and following Python code:
//* # from sympy.integrals.quadrature import gauss_laguerre
//* # from sympy import exp
//* # x, w = gauss_laguerre(8, 50)
//* # xend = 25
//* # print([(xi / xend).evalf(10) for xi in x])
//* # print([(wi * exp(xi) / xend).evalf(10) for xi, wi in zip(x, w)])
//*/
//static const int quadrature_steps = 8;
//static const float quadrature_nodes[] =
//{
//0.006811185292f,
//0.03614807107f,
//0.09004346519f,
//0.1706680068f,
//0.2818362161f,
//0.4303406404f,
//0.6296271457f,
//0.9145252695f
//};
//static const float quadrature_weights[] =
//{
//0.01750893642f,
//0.04135477391f,
//0.06678839063f,
//0.09507698807f,
//0.1283416365f,
//0.1707430204f,
//0.2327233347f,
//0.3562490486f
//};
//float3 ray_optical_depth(in Atmosphere atmosphere, float3 ray_origin)
//{
///* This function computes the optical depth along a ray.
//* Instead of using classic ray marching, the code is based on Gauss-Laguerre quadrature,
//* which is designed to compute the integral of f(x)*exp(-x) from 0 to infinity.
//* This works well here, since the optical depth along the ray tends to decrease exponentially.
//* By setting f(x) = g(x) exp(x), the exponentials cancel out and we get the integral of g(x).
//* The nodes and weights used here are the standard n=6 Gauss-Laguerre values, except that
//* the exp(x) scaling factor is already included in the weights.
//* The parametrization along the ray is scaled so that the last quadrature node is still within
//* the atmosphere. */
//float t0, t1;
//float3 ray_dir = atmosphere.sunDirection;
//RaySphereIntersect(ray_origin, ray_dir, atmosphere.atmosphereRadius, t0, t1);
//float3 ray_end = ray_dir * t1;
//float ray_length = distance(ray_origin, ray_end);
//float3 segment = ray_length * ray_dir;
///* instead of tracking the transmission spectrum across all wavelengths directly,
//* we use the fact that the density always has the same spectrum for each type of
//* scattering, so we split the density into a constant spectrum and a factor and
//* only track the factors */
//float3 optical_depth = float3(0.0f, 0.0f, 0.0f);
//for (int i = 0; i < quadrature_steps; i++)
//{
//float3 P = ray_origin + quadrature_nodes[i] * segment;
///* height above sea level */
//float height = length(P) - atmosphere.earthRadius;
//float3 density =float3(
//density_rayleigh(height), density_mie(height), density_ozone(height));
//optical_depth += density * quadrature_weights[i];
//}
//return optical_depth * ray_length;
//}
//float3 ComputeIncidentLight(in Atmosphere atmosphere, in float3 orig, in float3 dir, float tmin, float tmax)
//{
//float t0, t1;
//if (!RaySphereIntersect(orig, dir, atmosphere.atmosphereRadius, t0, t1) || t1 < 0)
//return 0;
//if (t0 > tmin && t0 > 0)
//tmin = t0;
//if (t1 < tmax)
//tmax = t1;
//uint numSamples = 16;
//uint numSamplesLight = 8;
//float segmentLength = (tmax - tmin) / numSamples;
//float tCurrent = tmin;
//float mu = dot(dir, atmosphere.sunDirection); //mu in the paper which is the cosine of the angle between the sun direction and the ray direction
//float3 opticalDepth = float3(0., 0., 0.);
//float3 densityScale = float3(1., 1., 1.);
//float3 phaseFunction = float3(phase_rayleigh(mu), phase_mie(mu), 0.0f);
//float3 spectrum = float3(0., 0., 0.);
//for (uint i = 0; i < numSamples; ++i)
//{
//float3 samplePosition = orig + (tCurrent + segmentLength * 0.5f) * dir;
//float height = length(samplePosition) - atmosphere.earthRadius;
//float3 density = densityScale * float3(density_rayleigh(height),
//density_mie(height),
//density_ozone(height));
//opticalDepth += segmentLength * density;
////if (!RaySphereIntersect(samplePosition, atmosphere.sunDirection, atmosphere.earthRadius))
////{
//float3 light_optical_depth = densityScale * ray_optical_depth(atmosphere, samplePosition);
//float3 total_optical_depth = opticalDepth + light_optical_depth;
//for (int wl = 0; wl < 3; ++wl)
//{
//float3 extinction_density = total_optical_depth * float3(atmosphere.betaR[wl],
//1.11f * atmosphere.betaM[wl],
//0. /*ozone_coeff[wl]*/);
//float attenuation = exp(-dot(extinction_density, float3(1., 1., 1.)));
//float3 scattering_density = density * float3(atmosphere.betaR[wl], atmosphere.betaM[wl], 0.0f);
///* the total inscattered radiance from one segment is:
//* Tr(A<->B) * Tr(B<->C) * sigma_s * phase * L * segment_length
//*
//* These terms are:
//* Tr(A<->B): Transmission from start to scattering position (tracked in optical_depth)
//* Tr(B<->C): Transmission from scattering position to light (computed in
//* ray_optical_depth) sigma_s: Scattering density phase: Phase function of the scattering
//* type (Rayleigh or Mie) L: Radiance coming from the light source segment_length: The
//* length of the segment
//*
//* The code here is just that, with a bit of additional optimization to not store full
//* spectra for the optical depth
//*/
//spectrum[wl] += attenuation * dot(phaseFunction * scattering_density, float3(1., 1., 1.)) *
//atmosphere.irradiance[wl] * 7. * segmentLength;
//}
////}
//// compute optical depth for light
////float hr = exp(-height / atmosphere.Hr) * segmentLength;
////float hm = exp(-height / atmosphere.Hm) * segmentLength;
////opticalDepthR += hr;
////opticalDepthM += hm;
////// light optical depth
////float t0Light, t1Light;
////RaySphereIntersect(samplePosition, atmosphere.sunDirection, atmosphere.atmosphereRadius, t0Light, t1Light);
////float segmentLengthLight = t1Light / numSamplesLight, tCurrentLight = 0;
////float opticalDepthLightR = 0, opticalDepthLightM = 0;
////uint j;
////for (j = 0; j < numSamplesLight; ++j)
////{
//// float3 samplePositionLight = samplePosition + (tCurrentLight + segmentLengthLight * 0.5f) * atmosphere.sunDirection;
//// float heightLight = length(samplePositionLight) - atmosphere.earthRadius;
//// if (heightLight < 0)
//// break;
//// opticalDepthLightR += exp(-heightLight / atmosphere.Hr) * segmentLengthLight;
//// opticalDepthLightM += exp(-heightLight / atmosphere.Hm) * segmentLengthLight;
//// tCurrentLight += segmentLengthLight;
////}
////if (j == numSamplesLight)
////{
//// float3 tau = atmosphere.betaR * (opticalDepthR + opticalDepthLightR) + atmosphere.betaM * 1.1f * (opticalDepthM + opticalDepthLightM);
//// float3 attenuation = float3(exp(-tau.x), exp(-tau.y), exp(-tau.z));
//// sumR += attenuation * hr;
//// sumM += attenuation * hm;
////}
//tCurrent += segmentLength;
//}
//return spectrum;
//// We use a magic number here for the intensity of the sun (20). We will make it more
//// scientific in a future revision of this lesson/code
////return (sumR * atmosphere.betaR * phaseR + sumM * atmosphere.betaM * phaseM) * 20;
//}
void CalcAtmosphere(inout float3 color, in float3 viewDir, in float3 sunDir, in float altitude, in float far, in float3 sunColor)
{
//float3x3 RGB2CAS = float3x3(1., 0., 0., 0., 1., 0., 0., 0., 1.);
//float3x3 CAS2RGB = float3x3(1., 0., 0., 0., 1., 0., 0., 0., 1.);
Atmosphere atmosphere;
atmosphere.earthRadius = 6360.e3;
atmosphere.atmosphereRadius = 6420.e3;
atmosphere.Hr = 7994.;
atmosphere.Hm = 1200.;
atmosphere.betaR = float3(5.8e-6, 13.5e-6, 33.1e-6);
atmosphere.betaR = float3(7.2865e-6, 1.2863e-5, 2.7408e-5);
//atmosphere.betaR = float3(8.658882115850087e-6, 1.5119838839850446e-5, 3.15880085046994e-5);
atmosphere.betaM = (float3) 21.e-6;
//atmosphere.irradiance = float3(1.4662835232339013, 1.7557753193321481, 1.7149614362178376);
if (isinf(far))
{
color = 1.e-7;
}
else
{
color = mul(RGB2CAS, color);
}
atmosphere.sunDirection = sunDir;
float2x3 atm = ComputeIncidentLight(atmosphere, float3(0., atmosphere.earthRadius + altitude, 0.), viewDir, 0., far);
color = atm[0] + color * atm[1];
color = mul(CAS2RGB, color);
//light = light.x * float3(0.933600, 0.358317, 0.000000) + light.y * float3(0.413697, 0.910408, 0.003489) + light.z * float3(0.178639, 0.025820, 0.983576);
//light = light.x * float3(0.899398, 0.437131, 0.000009) + light.y * float3(0.256233, 0.966471, 0.016686) + light.z * float3(0.175352, 0.029461, 0.984065);
//light = max(0., light.x * float3(3.240970, -0.969244, 0.055630) + light.y * float3(-1.537383, 1.875968, -0.203977) + light.z * float3(-0.498611, 0.041555, 1.056972));
//return max(0., mul(float3x3(
// 1.6218, -0.4493, 0.0325,
//-0.0374, 1.0598, -0.0742,
//-0.0283, -0.1119, 1.0491), light));
//float3x3 transform = float3x3(1.6218, -.4493, .0325, -.0374, 1.0598, -.0742, -.0283, -.0119, 1.0491);
//return mul(light, transform);
//return max(0., light.x * float3(2.512071, -0.057897, -0.043816) + light.y * float3(-0.622965, 1.469301, -0.155101) + light.z * float3(0.063069, -0.143869, 2.035035));
//return light;
//return (light.x * float3(1., 0., 0.) + light.y * float3(0.64, 1., 0.) + light.z * float3(0., 0., 1.)) / 1.64;
//return (light.x * float3(0.12544465, 0., 0.) + light.y * float3(0., 1.4467391, 0.)
//+ light.z * float3(0.22230228, 0., 1.86127603));
//return float3(dot(light, float3(0.64,
//0.33,
//0.03)), dot(light, float3(0.30,
//0.60,
//0.10)),
//dot(light, float3(0.15,
//0.06,
//0.79)));
}
void Sky(inout float3 color, in float3 viewDir, in float3 sunDir, in float altitude)
{
CalcAtmosphere(color, viewDir, sunDir, altitude, 1.#INF, 1.);
}
#endif

View File

@@ -0,0 +1,107 @@
#ifndef SKY_HLSLI
#define SKY_HLSLI
#include "math.hlsli"
#include "sky_common.hlsli"
Texture2D g_Sky : register(t13);
Texture2DArray<float4> g_AerialPerspectiveLut : register(t14);
SamplerState g_SkySampler : register(s13);
#define SKY_INF 1. #INF
float3 CalcSphereNormal(in float3 viewDir, in float3 sphereDir, in float cosAngularSize);
void CalcSun(inout float3 color, in float3 viewDir, in float3 sunDir, in float altitude)
{
float3 ray_origin = float3(0.0, EARTH_RADIUS + max(altitude, 1.), 0.0);
float ground_dist = ray_sphere_intersection(ray_origin, viewDir, EARTH_RADIUS);
if (ground_dist < 0.0)
{
if (dot(viewDir, sunDir) > 0.99998869014)
{
color = linear_srgb_from_spectral_samples(sun_spectral_irradiance);
}
}
}
void CalcMoon(inout float3 color, in float3 viewDir, in float3 moonDir, in float3 sunDir, in float altitude)
{
float3 ray_origin = float3(0.0, EARTH_RADIUS + max(altitude, 1.), 0.0);
float ground_dist = ray_sphere_intersection(ray_origin, viewDir, EARTH_RADIUS);
if (ground_dist < 0.0)
{
float3 normal = CalcSphereNormal(viewDir, moonDir, 0.99998869014);
if (dot(normal, normal) > 0.)
{
color = .07 * max(dot(normal, sunDir), 0.) * linear_srgb_from_spectral_samples(sun_spectral_irradiance);
}
}
}
float3 CalcSphereNormal(in float3 viewDir, in float3 sphereDir, in float cosAngularSize)
{
float sqrSinAngularSize = 1. - cosAngularSize * cosAngularSize;
float tca = dot(sphereDir, viewDir);
float d2 = dot(sphereDir, sphereDir) - tca * tca;
if (d2 > sqrSinAngularSize)
return (float3)(0.);
float thc = sqrt(sqrSinAngularSize - d2);
float t = tca - thc;
if (t < 0.)
return (float3)(0.);
float3 p = viewDir * t;
return normalize(p - sphereDir);
}
float2 CalcEquirectangularCoords(in float3 viewDir, in float3 sunDir, in float2 size)
{
float2 uv_scale = float2(1., 1. - 1. / size.y);
// float lon = acos(clamp(dot(normalize(sunDir.xz), normalize(viewDir.xz)), -1., 1.));
float lon = atan2(-viewDir.z, -viewDir.x);
float lat = asin(viewDir.y);
return float2(InvLerp(-PI, PI, lon), InvLerp(PI_OVER_TWO, -PI_OVER_TWO, lat)) * uv_scale;
}
float2 CalcEquirectangularCoordsBottomClipped(in float3 viewDir, in float3 sunDir, in float2 size)
{
float2 uv_scale = float2(1., 1. - 1. / size.y);
// float lon = acos(clamp(dot(normalize(sunDir.xz), normalize(viewDir.xz)), -1., 1.));
float lon = atan2(-viewDir.z, -viewDir.x);
float lat = asin(viewDir.y);
return float2(InvLerp(-PI, PI, lon), InvLerp(PI_OVER_TWO, -.125 * PI_OVER_TWO, lat)) * uv_scale;
}
void CalcAtmosphere(inout float3 color, in float3 viewDir, in float3 sunDir)
{
uint2 size;
g_Sky.GetDimensions(size.x, size.y);
float4 sky = g_Sky.SampleLevel(g_SkySampler, CalcEquirectangularCoordsBottomClipped(viewDir, sunDir, size), 0.);
color = color * sky.a + sky.rgb;
}
void ApplyAerialPerspective(inout float3 color, in float alpha, in float3 viewDir, in float3 sunDir, in float depth)
{
uint3 size;
g_AerialPerspectiveLut.GetDimensions(size.x, size.y, size.z);
float slice = sqrt(depth) * size.z;
float slice_factor = frac(slice);
slice -= 1.;
float2 uv = CalcEquirectangularCoords(viewDir, sunDir, size.xy);
float4 sample_near = slice > 0. ? g_AerialPerspectiveLut.SampleLevel(g_SkySampler, float3(uv, floor(slice)), 0.) : float4(0., 0., 0., 1.);
float4 sample_far = g_AerialPerspectiveLut.SampleLevel(g_SkySampler, float3(uv, ceil(slice)), 0.);
float4 aerial = lerp(sample_near, sample_far, slice_factor);
color = color * aerial.a + aerial.rgb * alpha;
}
#endif

View File

@@ -0,0 +1,630 @@
#ifndef SKY_HLSLI
#define SKY_HLSLI
#include "math.hlsli"
// Public Domain under http://unlicense.org, see link for details.
// A Nishita93-style single scattering atmosphere approximation.
// This can serve as a (not quite drop-in) replacement for e.g.
// https://github.com/wwwtyro/glsl-atmosphere
// Inherently provides both the sky model and the aerial perspective
// in the same code.
// Works for both surface and space views.
// Consider a spherical planet of radius R, centered at origin.
// Assuming, for the purposes of scattering computations, a 2-component
// (Rayleigh and Mie particles) atmosphere with position-dependent
// relative (i.e. dimensionless) densities Wr(r), Wm(r), the optical
// depth on the (segment of) ray ro+t*rd is:
// D(ro,rd,l,h) = ∫ Ber*Wr(ro+t*rd)+Bem*Wm(ro+t*rd) dt on [l;h]
// where Ber and Bem are volume-extinction coefficients (in units of m^-1):
// Ber=Bsr+Bar
// Bem=Bsm+Bam
// Bar, Bam - Rayleigh and Mie volume-absorption coefficients, respectively
// Bsr, Bsm - Rayleigh and Mie volume-scattering coefficients, respectively
// at W=1.
// NOTE: the volume-scattering, etc., coefficients are typically
// denoted β in literature, and are written as B here to keep the
// identifiers in ASCII. Similarly, the densities are typically denoted ρ.
// The integral is typically taken either on [0;+∞), or until the ray
// hits an obstacle (e.g. the planet itself).
// From this, transmittance is calculated as usual via Beer-Lambert law:
// T=exp(-D)
// Consider a light source (effectively infinitely far away) in the
// direction ld (so that the light direction is -ld). In the
// single-scattering approximation we only consider the following scenario:
// photons from the light source arrive (attenuated by the atmosphere) from
// the direction ld at some point on view ray ro+t*rd, experience scattering
// and arrive (again attenuated by the atmosphere) at viewer. The radiance,
// in W*sr^-1*m^-2, that the viewer observes in direction rd, then is:
// L = E * ∫ S(ro+t*rd)*(Pr*Bsr*Wr(ro+t*rd)+Pm*Bsm*Wm(ro+t*rd))*exp(-(D(ro,rd,0,t)+D(ro+t*rd,ld,0,+∞))) dt
// where
// E - irradiance from the light source, in W*m^-2
// Pr, Pm - Rayleigh and Mie phase functions, dimensionless
// S(r) - shadow/obstruction (0 if a ray from r towards ld intersects an obstacle, 1 otherwise)
// The irradiance is assumed constant in vicinity of the planet,
// and the light is only coming from ld (otherwise we would need to
// integrate the incoming radiance from all possible ld; here
// we basically assume that Lin(dir)=E*delta(dir-ld)).
// The phase functions depend on rd and ld, and describe the directional
// dependence of scattering.
// If the light is not monochromatic, the quantities L, E, and B* are actually
// spectral: they depend on the wavelength (or frequency, depending on how
// you style the spectrum). The units of L and E also change accordingly (e.g.
// E in W*sr^-1*m^-3 and E in W*m^-3 for spectrum-as-a-function-of-wavelength),
// but B* is still in m^-1.
// The entire double integral (since D(...) inside are themselves integrals)
// for the single-scattering approximation, in TeX:
// L_\lambda=E_\lambda \int_0^\infty S(\vec{r}_o+t \vec{r}_d) \left( P_r(\vec{r}_d,\vec{l}_d) B_{s r \lambda} W_r(\vec{r}_o+t \vec{r}_d) + P_m(\vec{r}_d,\vec{l}_d) B_{s m \lambda} W_m(\vec{r}_o+t \vec{r}_d) \right) \exp \left( -\left( \int_0^t \left( (B_{s r \lambda}+B_{a r \lambda}) W_r(\vec{r}_o+s \vec{r}_d)+ (B_{s m \lambda}+B_{a m \lambda}) W_m(\vec{r}_o+s \vec{r}_d) \right) \, ds + \int_0^\infty \left((B_{s r \lambda}+B_{a r \lambda}) W_r(\vec{r}_o+t \vec{r}_d+s \vec{l}_d)+ (B_{s m \lambda}+B_{a m \lambda}) W_m(\vec{r}_o+t \vec{r}_d+s \vec{l}_d) \right) \, ds \right) \right)
// NOTE: this can be also expressed as a differential equation. In this
// formulation, the computation is similar to alpha-blending infinitely many
// [t;dt] elements. This blending corresponds to exp(-D(ro,rd,0,t)) part in the
// integral, which is, therefore, absent from the differential equation.
//
// Computing this integral analytically is challenging, even for simple setups,
// so it is typically tackled numerically.
// The inner integral, D(...), however, admits fairly reasonable analytical
// approximations (whether they are faster than numerical methods is a separate
// question - we can get away with surprisingly few samples for
// numerical integration). A single-component version is sufficient to explore,
// since the multi-component version (Rayleigh and Mie) is a simple sum (this
// is NOT the case for the entire integral L, which is part of the reason
// for difficulty of computing it).
// Consider density exponentially decreasing with the altitude (a simple
// version of barometric formula)
// W=exp(-h/H)=exp(-(|r|-R)/H)
// where H is the scale height (generally, Rayleigh and Mie particles have
// different scale heights).
// This, basically, reduces the integral D(...) to the Chapman function, for
// which several approximations exist.
// We can approximate the density as a quadratic (gaussian) exponent:
// W≈exp(-(r^2-R^2)/(2*R*H))=exp(-h/H-h^2/(2*R*H))
// In this case, the single-component integral is:
// D = ∫ B*W(ro+t*rd) dt = B*sqrt(pi*R*H/2)*exp((R^2+dot(ro,rd)^2-dot(ro,ro))/(2*R*H))*erf((t+dot(ro,rd))/sqrt(2*R*H))
// This, actually, corresponds an existing approximation of the Chapman
// function. It has the advantage of being relatively simple and fully analytic.
// For Earth (R/H≈700) the error is reported ~0.1%.
// The expression for D above is useful both as a building block for the full
// radiance integral L, and on its own: for computing the optical depth, e.g.
// for aerial perspective and irradiance (both magnitude and color!) on a surface.
// Once we have this, the rest is just numerical intergration
// to compute L. We use an effective atmosphere radius (tuning parameter, currently
// const*(Hr+Hm), which balances errors from large integration step vs discarding
// larger part of atmosphere; more samples incentivize larger effective height) to only
// integrate on a finite segment where the atmospheric density is considered high
// enough. We then potentially split the integration segment by the planet's shadow,
// to only integrate on the lit parts of the ray (an alternative, to integrate on the
// entire segment, but skip the obstructed samples, seems less attractive, especially
// if the number of samples is low), since that are the only parts that contribute
// the light in our model. The integration itself is done using rectangle rule with
// a specifically chosen offset. This was found to provide better results than
// the commonly used midpoint rule (i.e. offset=1/2), while being similarly as simple.
// Notably, same as midpoint rule, there are no samples on segment endpoints,
// so the (sub-)segments don't share samples, which simplifies the computation.
// The results may look passable with as few as 4 samples total. The default in this
// shader is 8.
// Literature:
// https://en.wikipedia.org/wiki/Optical_depth
// https://en.wikipedia.org/wiki/Rayleigh_scattering
// https://en.wikipedia.org/wiki/Mie_scattering
// https://en.wikipedia.org/wiki/Beer%E2%80%93Lambert_law
// https://en.wikipedia.org/wiki/Sunlight
// https://en.wikipedia.org/wiki/Barometric_formula
// https://en.wikipedia.org/wiki/Scale_height
// https://en.wikipedia.org/wiki/Chapman_function
// https://pbr-book.org/4ed/Volume_Scattering
// https://www.scratchapixel.com/lessons/procedural-generation-virtual-worlds/simulating-sky/simulating-colors-of-the-sky.html
// http://www.thetenthplanet.de/archives/4519
// https://zero-radiance.github.io/post/analytic-media/
// Tomoyuki Nishita, Takao Sirai, Katsumi Tadamura, and Eihachiro Nakamae. 1993. Display of the earth taking into account atmospheric scattering. In Proceedings of the 20th annual conference on Computer graphics and interactive techniques (SIGGRAPH '93). Association for Computing Machinery, New York, NY, USA, 175182. https://doi.org/10.1145/166117.166140
// http://nishitalab.org/user/nis/cdrom/sig93_nis.pdf
// Eric Bruneton. 2016. A qualitative and quantitative evaluation of 8 clear sky models. IEEE transactions on visualization and computer graphics 23, 12 (2016), 2641--2655.
// https://arxiv.org/pdf/1612.04336.pdf
// Joseph T. Kider, Daniel Knowlton, Jeremy Newlin, Yining Karl Li, and Donald P. Greenberg. 2014. A framework for the experimental comparison of solar and skydome illumination. ACM Trans. Graph. 33, 6, Article 180 (November 2014), 12 pages. https://doi.org/10.1145/2661229.2661259
// https://spectralskylight.github.io/Kider2014.pdf
// Lukas Hosek and Alexander Wilkie. 2012. An analytic model for full spectral sky-dome radiance. ACM Trans. Graph. 31, 4, Article 95 (July 2012), 9 pages. https://doi.org/10.1145/2185520.2185591
// https://cgg.mff.cuni.cz/projects/SkylightModelling/HosekWilkie_SkylightModel_SIGGRAPH2012_Preprint.pdf
// A. J. Preetham, Peter Shirley, and Brian Smits. 1999. A practical analytic model for daylight. In Proceedings of the 26th annual conference on Computer graphics and interactive techniques (SIGGRAPH '99). ACM Press/Addison-Wesley Publishing Co., USA, 91100. https://doi.org/10.1145/311535.311545
// https://courses.cs.duke.edu/fall01/cps124/resources/p91-preetham.pdf
// Naty hoffman & Arcot J. Preetham, Photorealistic Real-Time Outdoor Light Scattering, GDC 2002
// https://renderwonk.com/publications/gdm-2002/GDM_August_2002.pdf
// https://drivers.amd.com/developer/gdc/GDC02_HoffmanPreetham.pdf
// https://github.com/ebruneton/clear-sky-models/tree/master
// https://github.com/spectralskylight
// https://github.com/wwwtyro/glsl-atmosphere
//==============================================================================
// Mathematical constants.
static const float INF = 1e17;
static const float pi = 3.14159265358979;
//==============================================================================
// Physical constants.
// Astrophysical Constants
// Source: https://pdg.lbl.gov/2023/reviews/rpp2022-rev-astrophysical-constants.pdf
static const float au = 149597870700.0; // Astronomical unit, in m.
static const float Rs = 6.957e8; // Nominal Solar equatorial radius, in m.
static const float Sa = 6.79431e-5; // Solar solid angle from 1 a.u., in sr (=2*pi*(1-sqrt(au^2-rs^2)/au))
static const float Re = 6.3781e6; // Nominal Earth equatorial radius, in m.
static const float Esc = 128e3; // Solar illuminance constant, in lux (https://en.wikipedia.org/wiki/Sunlight).
// NOTE: you want to compute such constants yourself (perhaps you are dealing with
// other planet and/or star), you would need the light source's spectrum (spectral
// radiance, L(λ) in W*sr^-1*m^-2*nm^-1, e.g. AM0 for Sun, see
// https://www.shadertoy.com/view/4XjGDK), its solid angle S as seen from viewer's
// position, and, if you want to work with photometric (perceptual) quantities,
// the lumious efficiency function (see http://www.cvrl.org/cie.htm and
// e.g. https://www.shadertoy.com/view/msXyDH).
// Then
// Ee = S*∫ L(λ)*dλ - irradiance (radiometric quantity, in W/m^2)
// Ev = 863*S*∫V(λ)*L(λ)*dλ - illuminance (photometric quantity, in lux)
// and the full color (using color matching functions,
// http://www.cvrl.org/cie.htm), also in lux, is:
// X = 863*S*∫x(λ)*dλ
// Y = 863*S*∫y(λ)*dλ
// Z = 863*S*∫z(λ)*dλ
// NOTE: V(λ)=y(λ). The relative, dimensionless color would be (X,Y,Z)/Ev.
// If you don't need photometric quantities (e.g. you are doing spectral
// render), then you just integrate L(λ) as is. E.g. for the AM0 solar spectrum
// this yields Ee equal to the solar constant Gsc=1361 W/m^2. The spectral
// radiance can be approximated as blackbody spectrum (see e.g.
// https://www.shadertoy.com/view/4XjGDK).
static const float Lsun = Esc / Sa; // Solar luminance, in cd/m^2.
//==============================================================================
// Colorspace.
// Following http://www.thetenthplanet.de/archives/4519, we adopt
// a chromatic adapation space (which we will unimaginatively call
// CAS) with the aim to make the scattering computations *independent* per
// channel, see the post for details.
// The channels in this space ~correspond to point sampling
// the spectrum at 615 nm, 535 nm, and 445 nm, respectively.
static const float3x3 CAS2RGB = float3x3(
1.6218, -0.4493, 0.0325,
-0.0374, 1.0598, -0.0742,
-0.0283, -0.1119, 1.0491);
static const float3x3 RGB2CAS = float3x3(
6.2267e-1, 2.6392e-1, -6.2375e-4,
2.3324e-2, 9.6056e-1, 6.7215e-2,
1.9285e-2, 1.0958e-1, 9.6035e-1);
static const float3 CASwavelength = float3(615, 535, 445);
// Sun color in CAS, normalized to unit luminance.
static const float3 CASsun = float3(0.9420, 1.0269, 1.0241);
// Rayleigh volume-scattering coefficient for
// CAS wavelengths, in m^-1, as per
// https://www.shadertoy.com/view/43j3zm.
static const float3 CASrayleigh = float3(7.2865e-6, 1.2863e-5, 2.7408e-5);
//==============================================================================
float2 quadratic_solve(float a, float b, float c)
{
float d = b * b - a * c;
#if 1
return d > 0.0 ? (-b + sqrt(d) * float2(-1, +1)) / a : float2(+INF, -INF);
#else
// Expected to be more accurate.
if(!(d>0.0)) return float2(+INF,-INF);
float q=-b+(b<0.0?sqrt(d):-sqrt(d)),l=c/q,h=q/a; // NOT sign(b), in case b=0.
return float2(min(l,h),max(l,h));
#endif
}
//==============================================================================
// Approximations for erf and erfcx.
#if 0
// From https://www.shadertoy.com/view/ml3yWj
// Eabs ~ 2.8e-8
// Erel ~ 2.8e-8
float erf(float x)
{
x=clamp(x,-4.0,+4.0);
float x2=x*x;
return tanh(x*(1.12837919+x2*(0.275732946+x2*(0.0408672727+x2*0.00200393011)))/(1.0+x2*(0.153282651+x2*(0.0224402472+x2*0.000285807058))));
}
// Erel=2.89760318e-08
float erfcx(float x)
{
float q=1.0+sqrt(pi)*abs(x);
float t=abs(x)>1e9?1.0:abs(x)/(1.0+abs(x));
float p=(1.0+t*(-1.83923738+t*(0.716286914+t*(0.798943258+t*(-0.777020726+t*0.196551435)))))/(1.0+t*(-2.48331366+t*(2.67169669+t*(-1.47041333+t*(0.417875249+t*-0.0403214433)))));
float y=p/q;
return x<0.0?2.0*exp(x*x)-y:y;
}
#else
// See https://www.shadertoy.com/view/4XSGzW
// Eabs<0.0037.
float erf(float x)
{
x = clamp(x, -3.0, +3.0);
return (1.13072 * x) / (1.0 + (x * x) * (0.357055 + (x * x) * -0.01014));
}
// Eabs<0.00150
// Erel<0.01100
float erfcx(float x)
{
const float a = 0.4956;
float y = (1.0 + a * abs(x)) / (1.0 + abs(x) * ((2.0 / sqrt(pi) + a) + sqrt(pi) * a * abs(x)));
return x >= 0.0 ? y : 2.0 * exp(x * x) - y;
}
#endif
//==============================================================================
// Gaussian integrals.
// Source: https://www.shadertoy.com/view/4XSGzW.
// Compute exp(x)*(erf(z)-erf(y)).
float gauss_segment(float x, float y, float z)
{
return y * z < 0.0 ?
exp(x) * (erf(z) - erf(y)) :
sign(y + z) * (exp(x - y * y) * erfcx(abs(y)) - exp(x - z * z) * erfcx(abs(z)));
}
// ∫ ρ(R,H,k,r(t)) dt on [l;h]
// where
// r(t) = ro+t*rd
// ρ = k*exp(-(r^2-R^2)/(2*R*H)) = k*exp(-h/H-h^2/(2*H*R))
// h = |r|-R
float density_integral(float R, float H, float k, float3 ro, float3 rd, float l, float h)
{
float A = 0.5 / (R * H);
float B = dot(ro, rd) / (R * H);
float C = 0.5 * (dot(ro, ro) - R * R) / (R * H);
float W = 0.25 * B * B / A - C;
return 0.5 * k * sqrt(pi / A) * gauss_segment(W, sqrt(A) * l + 0.5 * B / sqrt(A), sqrt(A) * h + 0.5 * B / sqrt(A));
}
//==============================================================================
// Phase functions.
// We follow the approach from
// Johannes Jendersie and Eugene d'Eon. 2023. An Approximate Mie Scattering Function for Fog and Cloud Rendering. In ACM SIGGRAPH 2023 Talks (SIGGRAPH '23). Association for Computing Machinery, New York, NY, USA, Article 47, 12. https://doi.org/10.1145/3587421.3595409
// https://research.nvidia.com/labs/rtr/approximate-mie/publications/approximate-mie.pdf
// and approximate the Mie phase function as a weighted sum of Henyey-Greenstein and Draine
// phase functions.
// NOTE: this reduces to pure Henyey-Greenstein for phase_mie(float4(g,0,0,0),x).
// Draine phase function.
// B.T. Draine. 2003. Scattering by interstellar dust grains. I. Optical and ultraviolet. The Astrophysical Journal 598, 2 (2003), 1017. https://doi.org/10.1086/379118
// NOTE: this reduces to Henyey-Greenstein for a=0,
// to Rayleigh for g=0, a=1 and to Cornette-Shanks for a=1.
float phase_draine(float a, float g, float x)
{
float d = 1.0 + g * g - 2.0 * g * x;
return 1.0 / (4.0 * pi) * (1.0 - g * g) / (1.0 + a * (1.0 + 2.0 * g * g) / 3.0) * // <-- constant factor.
(1.0 + a * x * x) / (d * sqrt(d));
}
float phase_rayleigh(float x)
{
return phase_draine(1.0, 0.0, x);
}
// Parametrization by droplet size.
// Input: droplet size, in m.
// Valid range: 5e-6<d<50e-6.
float4 phase_params_mie(float d)
{
d *= 1e6; // Convert to micrometres.
return float4(
exp(-0.0990567 / (d - 1.67154)), // gHG
exp(-2.20679 / (d + 3.91029) - 0.428934), // gD
exp(3.62489 - 8.29288 / (d + 5.52825)), // alpha
exp(-0.599085 / (d - 0.641583) - 0.665888) // wD
);
}
float phase_mie(float4 M, float x)
{
return lerp(phase_draine(0.0, M.x, x), phase_draine(M.z, M.y, x), M.w);
}
//==============================================================================
// Convert turbidity (used e.g. by Preetham and Hosek sky models) to Mie
// volume-scattering coefficient.
// The turbidity is defined as the ratio of vertical optical depths:
// T=(t_mie+t_rayleigh)/t_rayleigh
// at wavelength 550 nm (since, at least, t_rayleigh is wavelength-dependent).
// From this we get t_mie=(T-1)*t_rayleigh.
// Assuming volume-absorption coefficients are negligible, we get an
// expression for Mie volume-scattering coefficient in terms of
// Rayleigh volume-scattering coefficient (at 550 nm), scale heights
// and turbidity.
// To quote "An analytic model for full spectral sky-dome radiance"
// by Hosek and Wilkie:
// "This allows the user to easily define sky appearance without worrying
// too much about the intricacies of meteorology: T = 2 yields a very clear,
// Arctic-like sky, T = 3 a clear sky in a temperate climate, T = 6 a sky
// on a warm, moist day, T = 10 a slightly hazy day, and values of T above
// 50 represent dense fog."
// NOTE: this description doesn't appear to match the visuals: the
// sky at T>=2 is considerably hazier than the description suggests.
// The reason is not clear, though it may be partially explained by
// Hosek model using a simple Henyey-Greenstein phase function for Mie
// scattering.
// Returns Mie volume-scattering coefficient (scalar, assumed
// wavelength-independent), in m^-1.
float turbidity2mie(float B550, float Hr, float Hm, float T)
{
// NOTE: the vertical optical depth is just B*H.
return (T - 1.0) * B550 * Hr / Hm;
}
// Using the reference value of Rayleigh volume-scattering
// coefficient at 550 nm.
// See https://www.shadertoy.com/view/43j3zm.
float turbidity2mie(float Hr, float Hm, float T)
{
return turbidity2mie(1.149e-5, Hr, Hm, T);
}
//==============================================================================
// Atmospheric scattering.
// The number of samples is 2*NUM_STEPS.
static const int NUM_STEPS = 4;
// Returns color (C, in [0]), and alpha (A, in [1]), of the computed atmospheric scattering,
// which then can be used in the typical premultiplied alpha blending:
// Cresult=C+(1-A)*Cbackground
// Aresult=A+(1-A)*Abackground
// except alpha is per-component.
// NOTE: the units are for the reference. You can use other
// set of units, as long as it is consistent (everything in
// parsecs, etc.) - the internals of this function are not
// tied to SI.
float2x3 atmosphere(
float3 ro, // Ray origin.
float3 rd, // Ray direction.
float3 ld, // Direction toward light source.
float E, // Illuminance from light source, in lux.
float R, // Planet radius, in m.
float3 Bsr, float3 Bar, // Rayleigh volume-scattering and volume-absorption coefficients, in m^-1.
float3 Bsm, float3 Bam, // Mie volume-scattering and volume-absorption coefficients, in m^-1.
float Hr, float Hm, // Rayleigh and Mie scale height, in m.
float4 M, // Mie phase function parameters.
float2 s // Initial integration segment.
)
{
ld = normalize(ld);
rd = normalize(rd);
float Ra = R + 7.5 * (Hr + Hm); // Effective radius of atmosphere.
float2 sp = quadratic_solve(1.0, dot(ro, rd), dot(ro, ro) - R * R); // Segment inside planet.
float2 sa = quadratic_solve(1.0, dot(ro, rd), dot(ro, ro) - Ra * Ra); // Segment inside atmosphere.
if (sp.x < sp.y && sp.y > s.x)
s.y = min(s.y, sp.x); // Clamp segment by planet.
float3 Ber = Bsr + Bar;
float3 Bem = Bsm + Bam;
// Compute alpha. Needs to happen while the segment is still only
// clamped by planet.
float3 alpha = 1.0 - exp(-(
Ber * density_integral(R, Hr, 1.0, ro, rd, s.x, s.y) +
Bem * density_integral(R, Hm, 1.0, ro, rd, s.x, s.y)));
s = float2(max(s.x, sa.x), min(s.y, sa.y)); // Clamp segment by atmosphere.
float3 po = ro - dot(ro, ld) * ld;
float3 pd = rd - dot(rd, ld) * ld;
float2 ss = quadratic_solve(dot(pd, pd), dot(po, pd), dot(po, po) - R * R); // Segment inside planet's shadow.
if (dot(rd, ld) > 0.0)
ss.y = min(ss.y, -dot(ro, ld) / dot(rd, ld)); // Shadow is only half
else
ss.x = max(ss.x, -dot(ro, ld) / dot(rd, ld)); // of the cylinder.
float2 sl = s, sh = s.y;
if (ss.x < ss.y)
sl = float2(s.x, min(s.y, ss.x)); // Possibly split segment
if (ss.x < ss.y)
sh = float2(max(s.x, ss.y), s.y); // into two by shadow.
// If there's only one real segment, split it in half, so there
// are always 2.
if (!(sl.x < sl.y))
{
sl = float2(sh.x, 0.5 * (sh.x + sh.y));
sh = float2(0.5 * (sh.x + sh.y), sh.y);
}
else if (!(sh.x < sh.y))
{
sh = float2(0.5 * (sl.x + sl.y), sl.y);
sl = float2(sl.x, 0.5 * (sl.x + sl.y));
}
float3 Cr = 0;
float3 Cm = 0;
float Pr = phase_rayleigh(dot(rd, ld));
float Pm = phase_mie(M, dot(rd, ld));
// Integrate both segments.
//s = sl;
for (int k = 0; k < 1; ++k)
{
if (!(s.x < s.y))
continue;
// Find the step size, and estimate the optimized offset.
float dt = (s.y - s.x) / float(NUM_STEPS);
float3 B = Ber + Bem;
float X = max(max(B.x, B.y), B.z) * dt;
// The offset is chosen to provide exact answer
// for integrating exp(-(X/dt)*t) on [0;dt].
float offset = (X < 0.25 ? 0.5 - X / 24.0 : log(X / (1.0 - exp(-X))) / X);
if (k == 1)
offset = 1.0 - offset;
for (int i = 0; i < NUM_STEPS; ++i)
{
float t = s.x + (float(i) + offset) * dt;
float3 r = ro + rd * t;
float h = (dot(r, r) - R * R) / (2.0 * R);
float Dr = density_integral(R, Hr, 1.0, ro, rd, 0.0, t); // + density_integral(R, Hr, 1.0, r, ld, 0.0, INF);
float Dm = density_integral(R, Hm, 1.0, ro, rd, 0.0, t);// + density_integral(R, Hm, 1.0, r, ld, 0.0, INF);
float3 a = exp(-(Ber * Dr + Bem * Dm));
Cr += exp(-h / Hr) * a * dt;
Cm += exp(-h / Hm) * a * dt;
}
s = sh;
}
return float2x3(E * (Pr * Bsr * Cr + Pm * Bsm * Cm), alpha);
}
//==============================================================================
// Main image.
//void mainImage(out float4 fragColor, in float2 fragCoord)
//{
// float F = 1.5;
// float R = Re;
// float3 ro = float3(0, R + 80.0, 5.0);
// float3 rd = normalize(float3((2.0 * fragCoord - iResolution.xy) / iResolution.y, -F));
// float3 md = normalize(float3((2.0 * iMouse.xy - iResolution.xy) / iResolution.y, -F));
// float3 ld = normalize(float3(1, sin(0.5 * iTime) + 0.875, -2));
// if (length(iMouse.xy) > 16.0)
// ld = md;
// float Hr = 7.994e3, Hm = 1.2e3;
// if (texture(iChannel0, float2(49.5, 0.5) / float2(256, 3)).x > 0.0)
// ro.z = sqrt(R * Hr);
// if (texture(iChannel0, float2(50.5, 0.5) / float2(256, 3)).x > 0.0)
// ro.z = R / 16.0;
// if (texture(iChannel0, float2(51.5, 0.5) / float2(256, 3)).x > 0.0)
// ro.z = R / 8.0;
// if (texture(iChannel0, float2(52.5, 0.5) / float2(256, 3)).x > 0.0)
// ro.z = 4.0 * R;
// float T = 1.375; // Turbidity.
// float D = 17e-6; // Droplet size, in m.
// float3 Bsr = CASrayleigh, Bar = float3(0);
// float3 Bsm = float3(turbidity2mie(Hr, Hm, T)), Bam = float3(0); // NOTE: glsl-atmosphere uses Bsm=21e-6
// float4 M = phase_params_mie(D);
// // From https://en.wikipedia.org/wiki/Orders_of_magnitude_(illuminance)#Luminance
// // 5e3 - Typical photographic scene in full sunlight
// // 7e3 - Average clear sky
// float Lref = 7e3; // Scale luminance, in cd/m^2.
// // NOTE: we work entirely in CAS, and only convert to RGB right
// // before the tonemapping.
// float3 col = float3(1e-7); // Airglow.
// float2 s = quadratic_solve(1.0, dot(ro, rd), dot(ro, ro) - R * R);
// // Render Sun.
// if (true)
// col += Lsun * CASsun * smoothstep(1.0 - (0.5 * Rs * Rs / (au * au)), 1.0, dot(rd, ld));
// else
// col += Lsun * CASsun * exp(-(1.0 - dot(rd, ld)) / (2.0 * Rs * Rs / (au * au)));
// // Render planet.
// if (s.x < s.y && s.y > 0.0)
// {
// float3 a = exp(-(Bsr + Bar) * density_integral(R, Hr, 1.0, ro, rd, 0.0, s.x) + (Bsm + Bam) * density_integral(R, Hm, 1.0, ro, rd, 0.0, s.x));
// col = Esc / pi * CASsun * a * (1.0 / 128.0 + max(dot(ld, normalize(ro + s.x * rd)), 0.0)) * (RGB2CAS * float3(0.125, 0.25, 0.0625));
// }
// float2x3 m = atmosphere(ro, normalize(rd), normalize(ld), Esc, R, Bsr, Bar, Bsm, Bam, Hr, Hm, M, float2(0.0, INF));
// col = CASsun * m[0] + col * (1.0 - m[1]);
// col = CAS2RGB * col;
// // Apply exposure.
// col = 1.0 - exp(-col / Lref);
// col = lerp(12.92 * col, 1.055 * pow(col, float3(1.0 / 2.4)) - 0.055, step(0.0031308, col)); // sRGB
//
// float3x3 t = float3x3(0., 0., 0., 0., 0., 0., 0.1, 0.5, 0.3);
//
// fragColor = float4(col, 1.0);
// //fragColor.rgb = t * float3(0., 0., 1.);
//}
void CalcAtmosphere(inout float3 color, in float3 viewDir, in float3 sunDir, in float altitude, in float far, in float3 sunColor)
{
float R = Re;
float3 ro = float3(0, R + 80.0, 5.0);
float3 rd = normalize(viewDir);
float3 ld = normalize(sunDir);
float Hr = 7.994e3, Hm = 1.2e3;
float T = 1.375; // Turbidity.
float D = 17e-6; // Droplet size, in m.
float3 Bsr = CASrayleigh, Bar = 0;
float3 Bsm = turbidity2mie(Hr, Hm, T), Bam = 0; // NOTE: glsl-atmosphere uses Bsm=21e-6
float4 M = phase_params_mie(D);
// From https://en.wikipedia.org/wiki/Orders_of_magnitude_(illuminance)#Luminance
// 5e3 - Typical photographic scene in full sunlight
// 7e3 - Average clear sky
float Lref = 7e3; // Scale luminance, in cd/m^2.
// NOTE: we work entirely in CAS, and only convert to RGB right
// before the tonemapping.
float3 col = 1e-7; // Airglow.
float2 s = quadratic_solve(1.0, dot(ro, rd), dot(ro, ro) - R * R);
// Render Sun.
//if (true)
// col += Lsun * CASsun * smoothstep(1.0 - (0.5 * Rs * Rs / (au * au)), 1.0, dot(rd, ld));
//else
// col += Lsun * CASsun * exp(-(1.0 - dot(rd, ld)) / (2.0 * Rs * Rs / (au * au)));
// Render planet.
if (far < INF)
{
col = mul(RGB2CAS, color);
}
else
{
far = INF;
}
//if (s.x < s.y && s.y > 0.0)
//{
// float3 a = exp(-(Bsr + Bar) * density_integral(R, Hr, 1.0, ro, rd, 0.0, s.x) + (Bsm + Bam) * density_integral(R, Hm, 1.0, ro, rd, 0.0, s.x));
// col = Esc / pi * CASsun * a * (1.0 / 128.0 + max(dot(ld, normalize(ro + s.x * rd)), 0.0)) * (mul(RGB2CAS, float3(0.125, 0.25, 0.0625)));
//}
float2x3 m = atmosphere(ro, normalize(rd), normalize(ld), 4., R, Bsr, Bar, Bsm, Bam, Hr, Hm, M, float2(0.0, far));
col = sunColor * m[0] + col * (1.0 - m[1]);
col = mul(CAS2RGB, col);
// Apply exposure.
//col = 1.0 - exp(-col / Lref);
color = col;
//col = lerp(12.92 * col, 1.055 * pow(col, float3(1.0 / 2.4)) - 0.055, step(0.0031308, col)); // sRGB
//
//float3x3 t = float3x3(0., 0., 0., 0., 0., 0., 0.1, 0.5, 0.3);
//Atmosphere atmosphere;
//atmosphere.earthRadius = 6360.e3;
//atmosphere.atmosphereRadius = 6420.e3;
//atmosphere.Hr = 7994.;
//atmosphere.Hm = 1200.;
//atmosphere.betaR = float3(5.8e-6, 13.5e-6, 33.1e-6);
//atmosphere.betaR = float3(7.2865e-6, 1.2863e-5, 2.7408e-5);
////atmosphere.betaR = float3(8.658882115850087e-6, 1.5119838839850446e-5, 3.15880085046994e-5);
//atmosphere.betaM = (float3) 21.e-6;
////atmosphere.irradiance = float3(1.4662835232339013, 1.7557753193321481, 1.7149614362178376);
//
//atmosphere.sunDirection = sunDir;
//float3 light = ComputeIncidentLight(atmosphere, float3(0., atmosphere.earthRadius + 2500. + altitude, 0.), viewDir, 0., far);
////light = light.x * float3(0.933600, 0.358317, 0.000000) + light.y * float3(0.413697, 0.910408, 0.003489) + light.z * float3(0.178639, 0.025820, 0.983576);
////light = light.x * float3(0.899398, 0.437131, 0.000009) + light.y * float3(0.256233, 0.966471, 0.016686) + light.z * float3(0.175352, 0.029461, 0.984065);
////light = max(0., light.x * float3(3.240970, -0.969244, 0.055630) + light.y * float3(-1.537383, 1.875968, -0.203977) + light.z * float3(-0.498611, 0.041555, 1.056972));
//return max(0., mul(float3x3(
// 1.6218, -0.4493, 0.0325,
//-0.0374, 1.0598, -0.0742,
//-0.0283, -0.1119, 1.0491), light));
//float3x3 transform = float3x3(1.6218, -.4493, .0325, -.0374, 1.0598, -.0742, -.0283, -.0119, 1.0491);
//return mul(light, transform);
//return max(0., light.x * float3(2.512071, -0.057897, -0.043816) + light.y * float3(-0.622965, 1.469301, -0.155101) + light.z * float3(0.063069, -0.143869, 2.035035));
//return light;
//return (light.x * float3(1., 0., 0.) + light.y * float3(0.64, 1., 0.) + light.z * float3(0., 0., 1.)) / 1.64;
//return (light.x * float3(0.12544465, 0., 0.) + light.y * float3(0., 1.4467391, 0.)
//+ light.z * float3(0.22230228, 0., 1.86127603));
//return float3(dot(light, float3(0.64,
//0.33,
//0.03)), dot(light, float3(0.30,
//0.60,
//0.10)),
//dot(light, float3(0.15,
//0.06,
//0.79)));
}
void Sky(inout float3 color, in float3 viewDir, in float3 sunDir, in float altitude)
{
CalcAtmosphere(color, viewDir, sunDir, altitude, INF, 1.);
}
#endif

View File

@@ -0,0 +1,116 @@
#include "sky_common.hlsli"
#include "sky_inscattering.hlsli"
cbuffer DispatchConstants : register(b0)
{
float4x4 g_InverseView;
float4x4 g_InverseProjection;
float3 g_SunDir;
float g_Altitude;
float3 g_MoonDir;
float g_MaxDepth;
}
RWTexture2DArray<float4> g_AerialLut : register(u0);
RWTexture2D<float4> g_Sky : register(u1);
Texture2D<float4> g_TransmittanceLut : register(t13);
SamplerState g_TransmittanceLutSampler : register(s13);
[numthreads(8, 8, 1)] void CS_AerialLUT(uint3 PixCoord : SV_DispatchThreadID)
{
uint3 texture_size;
g_AerialLut.GetDimensions(texture_size.x, texture_size.y, texture_size.z);
float moon_phase = acos(dot(g_SunDir, g_MoonDir)) * ONE_OVER_PI;
uint2 pix_coord = PixCoord.xy;
float start_depth = 0.;
float4 L = 0.;
float4 M = 0.;
float4 transmittance = 1.;
float4 transmittance_m = 1.;
float4 position_ndc = float4(((pix_coord + 0.5) / texture_size.xy - .5) * float2(2., -2.), 1., 1.);
position_ndc = mul(g_InverseProjection, position_ndc);
position_ndc /= position_ndc.w;
float lat = lerp(PI_OVER_TWO, -PI_OVER_TWO, pix_coord.y / (texture_size.y - 1.));
float lon = lerp(0., TWO_PI, pix_coord.x / float(texture_size.x));
float3 ray_dir;
ray_dir.y = sin(lat);
ray_dir.x = cos(lon) * cos(lat);
ray_dir.z = sin(lon) * cos(lat);
float3 sun_dir;
sun_dir.y = g_SunDir.y;
sun_dir.x = sqrt(1 - sun_dir.y * sun_dir.y);
sun_dir.z = 0.;
// float3 ray_dir = mul((float3x3)g_InverseView, normalize(position_ndc.xyz));
float3 ray_origin = float3(0.0, EARTH_RADIUS + max(g_Altitude, 1.), 0.0);
float t_d = g_MaxDepth;
float atmos_dist = ray_sphere_intersection(ray_origin, ray_dir, ATMOSPHERE_RADIUS);
float ground_dist = ray_sphere_intersection(ray_origin, ray_dir, EARTH_RADIUS);
// We are inside the atmosphere
if (ground_dist < 0.0)
{
// No ground collision, use the distance to the outer atmosphere
t_d = min(t_d, atmos_dist);
}
else
{
// We have a collision with the ground, use the distance to it
t_d = min(t_d, ground_dist);
}
float cos_theta = dot(-ray_dir, g_SunDir);
float cos_theta_moon = dot(-ray_dir, g_MoonDir);
float molecular_phase = molecular_phase_function(cos_theta);
float aerosol_phase = aerosol_phase_function(cos_theta);
float molecular_phase_moon = molecular_phase_function(cos_theta_moon);
float aerosol_phase_moon = aerosol_phase_function(cos_theta_moon);
for (uint i = 0; i < texture_size.z; ++i)
{
float t = (i + 1.) / (float)texture_size.z;
float end_depth = min(t_d, t * t * g_MaxDepth);
compute_inscattering(g_TransmittanceLut, g_TransmittanceLutSampler, molecular_phase, aerosol_phase, 5, ray_origin, ray_dir, start_depth, end_depth, g_SunDir, L, transmittance);
compute_inscattering(g_TransmittanceLut, g_TransmittanceLutSampler, molecular_phase_moon, aerosol_phase_moon, 5, ray_origin, ray_dir, start_depth, end_depth, g_MoonDir, M, transmittance_m);
g_AerialLut[uint3(pix_coord, i)] = float4(linear_srgb_from_spectral_samples(L + .07 * moon_phase * M) * exp2(EXPOSURE), dot(transmittance, .25));
start_depth = end_depth;
}
}
[numthreads(8, 8, 1)] void CS_Sky(uint3 PixCoord : SV_DispatchThreadID)
{
uint2 texture_size;
g_Sky.GetDimensions(texture_size.x, texture_size.y);
float moon_phase = acos(dot(g_SunDir, g_MoonDir)) * ONE_OVER_PI;
uint2 pix_coord = PixCoord.xy;
float lat = lerp(PI_OVER_TWO, -.125 * PI_OVER_TWO, pix_coord.y / (texture_size.y - 1.));
float lon = lerp(0., TWO_PI, pix_coord.x / float(texture_size.x-1.));
float3 ray_dir;
ray_dir.y = sin(lat);
ray_dir.x = cos(lon) * cos(lat);
ray_dir.z = sin(lon) * cos(lat);
float3 sun_dir;
sun_dir.y = g_SunDir.y;
sun_dir.x = sqrt(1 - sun_dir.y * sun_dir.y);
sun_dir.z = 0.;
float4 position_ndc = float4(((pix_coord + 0.5) / texture_size.xy - .5) * float2(2., -2.), 1., 1.);
position_ndc = mul(g_InverseProjection, position_ndc);
position_ndc /= position_ndc.w;
// float3 ray_dir = mul((float3x3)g_InverseView, normalize(position_ndc.xyz));
g_Sky[pix_coord] = get_inscattering(g_TransmittanceLut, g_TransmittanceLutSampler, 32, g_Altitude, ray_dir, g_MaxDepth, 1.#INF, g_SunDir);
g_Sky[pix_coord] += .07 * moon_phase * get_inscattering(g_TransmittanceLut, g_TransmittanceLutSampler, 32, g_Altitude, ray_dir, g_MaxDepth, 1.#INF, g_MoonDir);
}

View File

@@ -0,0 +1,328 @@
#ifndef SKY_COMMON_HLSLI
#define SKY_COMMON_HLSLI
#include "color_transform.hlsli"
#include "math.hlsli"
cbuffer SkyConstants : register(b13){
float4 g_AerosolAbsorptionCrossSection;
float4 g_AerosolScatteringCrossSection;
float4 g_GroundAlbedo;
float g_AerosolHeightScale;
float g_AerosolTurbidity;
float g_AerosolBaseDensity;
float g_AerosolBackgroundDividedByBaseDensity;
float g_OzoneMean;
float g_FogDensity;
float g_FogHeightOffset;
float g_FogHeightScale;
}
// Configurable parameters
// #define ANIMATE_SUN 1
// // 0=equirectangular, 1=fisheye, 2=projection
// #define CAMERA_TYPE 2
// // 0=Background, 1=Desert Dust, 2=Maritime Clean, 3=Maritime Mineral,
// // 4=Polar Antarctic, 5=Polar Artic, 6=Remote Continental, 7=Rural, 8=Urban
//#define AEROSOL_TYPE 8
// static const float SUN_ELEVATION_DEGREES = 0.0; // 0=horizon, 90=zenith
// static const float EYE_ALTITUDE = 0.5; // km
//static const int MONTH = 0; // 0-11, January to December
//static const float AEROSOL_TURBIDITY = 1.0;
//static const float4 GROUND_ALBEDO = 0.3;
// // Ray marching steps. More steps mean better accuracy but worse performance
// static const int TRANSMITTANCE_STEPS = 32;
// static const int IN_SCATTERING_STEPS = 32;
// // Camera settings
static const float EXPOSURE = -4.0;
// // For the "projection" type camera
// static const float CAMERA_FOV = 90.0;
// static const float CAMERA_YAW = 15.0;
// static const float CAMERA_PITCH = -12.0;
// static const float CAMERA_ROLL = 0.0;
// // Debug
// #define ENABLE_SPECTRAL 1
// #define ENABLE_MULTIPLE_SCATTERING 1
// #define ENABLE_AEROSOLS 1
// #define SHOW_RELATIVE_LUMINANCE 0
// #define TONEMAPPING_TECHNIQUE 0 // 0=ACES, 1=simple
//-----------------------------------------------------------------------------
// Constants
// All parameters that depend on wavelength (float4) are sampled at
// 630, 560, 490, 430 nanometers
//static const float PI = 3.14159265358979323846;
static const float INV_PI = 0.31830988618379067154;
static const float INV_4PI = 0.25 * INV_PI;
static const float PHASE_ISOTROPIC = INV_4PI;
static const float RAYLEIGH_PHASE_SCALE = (3.0 / 16.0) * INV_PI;
static const float g = 0.8;
static const float gg = g*g;
static const float EARTH_RADIUS = 6371.0e3; // km
static const float ATMOSPHERE_THICKNESS = 100.0e3; // km
static const float ATMOSPHERE_RADIUS = EARTH_RADIUS + ATMOSPHERE_THICKNESS;
// static const float EYE_DISTANCE_TO_EARTH_CENTER = EARTH_RADIUS + EYE_ALTITUDE;
// static const float SUN_ZENITH_COS_ANGLE = cos(radians(90.0 - SUN_ELEVATION_DEGREES));
// static const float3 SUN_DIR = float3(-sqrt(1.0 - SUN_ZENITH_COS_ANGLE*SUN_ZENITH_COS_ANGLE), 0.0, SUN_ZENITH_COS_ANGLE);
// Extraterrestial Solar Irradiance Spectra, units W * m^-2 * nm^-1
// https://www.nrel.gov/grid/solar-resource/spectra.html
static const float4 sun_spectral_irradiance = float4(1.679, 1.828, 1.986, 1.307);
// Rayleigh scattering coefficient at sea level, units km^-1
// "Rayleigh-scattering calculations for the terrestrial atmosphere"
// by Anthony Bucholtz (1995).
static const float4 molecular_scattering_coefficient_base = float4(6.605e-6, 1.067e-5, 1.842e-5, 3.156e-5);
// Fog scattering/extinction cross section, units m^2 / molecules
// Mie theory results for IOR of 1.333. Particle size is a log normal
// distribution of mean diameter=15 and std deviation=0.4
static const float4 fog_scattering_cross_section = float4(5.015e-10, 4.987e-10, 4.966e-10, 4.949e-10);
// Ozone absorption cross section, units m^2 / molecules
// "High spectral resolution ozone absorption cross-sections"
// by V. Gorshelev et al. (2014).
static const float4 ozone_absorption_cross_section = float4(3.472e-25, 3.914e-25, 1.349e-25, 11.03e-27);
// Mean ozone concentration in Dobson for each month of the year.
// static const float ozone_mean_monthly_dobson[] = {
// 347.0, // January
// 370.0, // February
// 381.0, // March
// 384.0, // April
// 372.0, // May
// 352.0, // June
// 333.0, // July
// 317.0, // August
// 298.0, // September
// 285.0, // October
// 290.0, // November
// 315.0 // December
// };
/*
* Every aerosol type expects 5 parameters:
* - Scattering cross section
* - Absorption cross section
* - Base density (km^-3)
* - Background density (km^-3)
* - Height scaling parameter
* These parameters can be sent as uniforms.
*
* This model for aerosols and their corresponding parameters come from
* "A Physically-Based Spatio-Temporal Sky Model"
* by Guimera et al. (2018).
*/
// #if AEROSOL_TYPE == 0 // Background
// static const float4 aerosol_absorption_cross_section = float4(4.5517e-19, 5.9269e-19, 6.9143e-19, 8.5228e-19);
// static const float4 aerosol_scattering_cross_section = float4(1.8921e-26, 1.6951e-26, 1.7436e-26, 2.1158e-26);
// static const float aerosol_base_density = 2.584e14;
// static const float aerosol_background_density = 2e3;
// #elif AEROSOL_TYPE == 1 // Desert Dust
// static const float4 aerosol_absorption_cross_section = float4(4.6758e-16, 4.4654e-16, 4.1989e-16, 4.1493e-16);
// static const float4 aerosol_scattering_cross_section = float4(2.9144e-16, 3.1463e-16, 3.3902e-16, 3.4298e-16);
// static const float aerosol_base_density = 1.8662e15;
// static const float aerosol_background_density = 2e3;
// static const float aerosol_height_scale = 2.0;
// #elif AEROSOL_TYPE == 2 // Maritime Clean
// static const float4 aerosol_absorption_cross_section = float4(6.3312e-19, 7.5567e-19, 9.2627e-19, 1.0391e-18);
// static const float4 aerosol_scattering_cross_section = float4(4.6539e-26, 2.721e-26, 4.1104e-26, 5.6249e-26);
// static const float aerosol_base_density = 2.0266e14;
// static const float aerosol_background_density = 2e6;
// static const float aerosol_height_scale = 0.9;
// #elif AEROSOL_TYPE == 3 // Maritime Mineral
// static const float4 aerosol_absorption_cross_section = float4(6.9365e-19, 7.5951e-19, 8.2423e-19, 8.9101e-19);
// static const float4 aerosol_scattering_cross_section = float4(2.3699e-19, 2.2439e-19, 2.2126e-19, 2.021e-19);
// static const float aerosol_base_density = 2.0266e14;
// static const float aerosol_background_density = 2e3;
// static const float aerosol_height_scale = 2.0;
// #elif AEROSOL_TYPE == 4 // Polar Antarctic
// static const float4 aerosol_absorption_cross_section = float4(1.3399e-16, 1.3178e-16, 1.2909e-16, 1.3006e-16);
// static const float4 aerosol_scattering_cross_section = float4(1.5506e-19, 1.809e-19, 2.3069e-19, 2.5804e-19);
// static const float aerosol_base_density = 2.3864e13;
// static const float aerosol_background_density = 2e3;
// static const float aerosol_height_scale = 30.0;
// #elif AEROSOL_TYPE == 5 // Polar Arctic
// static const float4 aerosol_absorption_cross_section = float4(1.0364e-16, 1.0609e-16, 1.0193e-16, 1.0092e-16);
// static const float4 aerosol_scattering_cross_section = float4(2.1609e-17, 2.2759e-17, 2.5089e-17, 2.6323e-17);
// static const float aerosol_base_density = 2.3864e13;
// static const float aerosol_background_density = 2e3;
// static const float aerosol_height_scale = 30.0;
// #elif AEROSOL_TYPE == 6 // Remote Continental
// static const float4 aerosol_absorption_cross_section = float4(4.5307e-18, 5.0662e-18, 4.4877e-18, 3.7917e-18);
// static const float4 aerosol_scattering_cross_section = float4(1.8764e-18, 1.746e-18, 1.6902e-18, 1.479e-18);
// static const float aerosol_base_density = 6.103e15;
// static const float aerosol_background_density = 2e3;
// static const float aerosol_height_scale = 0.73;
// #elif AEROSOL_TYPE == 7 // Rural
// static const float4 aerosol_absorption_cross_section = float4(5.0393e-23, 8.0765e-23, 1.3823e-22, 2.3383e-22);
// static const float4 aerosol_scattering_cross_section = float4(2.6004e-22, 2.4844e-22, 2.8362e-22, 2.7494e-22);
// static const float aerosol_base_density = 8.544e15;
// static const float aerosol_background_density = 2e3;
// static const float aerosol_height_scale = 0.73;
// #elif AEROSOL_TYPE == 8 // Urban
// static const float4 aerosol_absorption_cross_section = float4(2.8722e-24, 4.6168e-24, 7.9706e-24, 1.3578e-23);
// static const float4 aerosol_scattering_cross_section = float4(1.5908e-22, 1.7711e-22, 2.0942e-22, 2.4033e-22);
// static const float aerosol_base_density = 1.3681e17;
// static const float aerosol_background_density = 2e3;
// static const float aerosol_height_scale = 0.73;
// #endif
//static const float aerosol_background_divided_by_base_density = aerosol_background_density / aerosol_base_density;
// static const float fog_density = 29811535.0;
// static const float fog_scale_height = .01;
// static const float fog_height_offset = 0.325;
//-----------------------------------------------------------------------------
/*
* Helper function to obtain the transmittance to the top of the atmosphere
* from Buffer A.
*/
float4 transmittance_from_lut(Texture2D lut, SamplerState lut_sampler, float cos_theta, float normalized_altitude)
{
float u = saturate(cos_theta * 0.5 + 0.5);
float v = saturate(normalized_altitude);
return lut.SampleLevel(lut_sampler, float2(u, v), 0.);
}
/*
* Returns the distance between ro and the first intersection with the sphere
* or -1.0 if there is no intersection. The sphere's origin is (0,0,0).
* -1.0 is also returned if the ray is pointing away from the sphere.
*/
float ray_sphere_intersection(float3 ro, float3 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));
}
/*
* Rayleigh phase function.
*/
float molecular_phase_function(float cos_theta)
{
return RAYLEIGH_PHASE_SCALE * (1.0 + cos_theta*cos_theta);
}
/*
* Henyey-Greenstrein phase function.
*/
float aerosol_phase_function(float cos_theta)
{
float den = 1.0 + gg + 2.0 * g * cos_theta;
return INV_4PI * (1.0 - gg) / (den * sqrt(den));
}
float4 get_multiple_scattering(Texture2D transmittance_lut, SamplerState lut_sampler, float cos_theta, float normalized_height, float d)
{
// Solid angle subtended by the planet from a point at d distance
// from the planet center.
float omega = 2.0 * PI * (1.0 - sqrt(max(0., d*d - EARTH_RADIUS*EARTH_RADIUS)) / d);
float4 T_to_ground = transmittance_from_lut(transmittance_lut, lut_sampler, cos_theta, 0.0);
float4 T_ground_to_sample =
transmittance_from_lut(transmittance_lut, lut_sampler, 1.0, 0.0) /
transmittance_from_lut(transmittance_lut, lut_sampler, 1.0, normalized_height);
// 2nd order scattering from the ground
float4 L_ground = PHASE_ISOTROPIC * omega * (g_GroundAlbedo / PI) * T_to_ground * T_ground_to_sample * cos_theta;
// Fit of Earth's multiple scattering coming from other points in the atmosphere
float4 L_ms = 0.02 * float4(0.217, 0.347, 0.594, 1.0) * (1.0 / (1.0 + 5.0 * exp(-17.92 * cos_theta)));
return L_ms + L_ground;
}
/*
* Return the molecular volume scattering coefficient (km^-1) for a given altitude
* in kilometers.
*/
float4 get_molecular_scattering_coefficient(float h)
{
return molecular_scattering_coefficient_base * exp(-0.07771971 * pow(h, 1.16364243));
}
/*
* Return the molecular volume absorption coefficient (km^-1) for a given altitude
* in kilometers.
*/
float4 get_molecular_absorption_coefficient(float h)
{
h += 1e-4; // Avoid division by 0
float t = log(h) - 3.22261;
float density = 3.78547397e17 * (1.0 / h) * exp(-t * t * 5.55555555);
return ozone_absorption_cross_section * g_OzoneMean * density;
}
/*
* Return the fog volume scattering coefficient (m^-1) for a given altitude in
* kilometers.
*
* Fog (or mist, depending on density) is a special kind of aerosol consisting
* of water droplets or ice crystals. Visibility is mostly dependent on fog.
*/
float4 get_fog_scattering_coefficient(float h)
{
if (g_FogDensity > 0.0) {
return fog_scattering_cross_section * g_FogDensity
* min(1.0, exp((-h + g_FogHeightOffset) / g_FogHeightScale));
} else {
return 0.0;
}
}
float get_aerosol_density(float h)
{
return g_AerosolBaseDensity * (exp(-h / g_AerosolHeightScale)
+ g_AerosolBackgroundDividedByBaseDensity);
}
/*
* Get the collision coefficients (scattering and absorption) of the
* atmospheric medium for a given point at an altitude h.
*/
void get_atmosphere_collision_coefficients(in float h,
out float4 aerosol_absorption,
out float4 aerosol_scattering,
out float4 molecular_absorption,
out float4 molecular_scattering,
out float4 fog_scattering,
out float4 extinction)
{
h = max(h, 1.e-3); // In case height is negative
h *= 1.e-3;
float aerosol_density = get_aerosol_density(h);
aerosol_absorption = g_AerosolAbsorptionCrossSection * aerosol_density * g_AerosolTurbidity;
aerosol_scattering = g_AerosolScatteringCrossSection * aerosol_density * g_AerosolTurbidity;
molecular_absorption = get_molecular_absorption_coefficient(h);
molecular_scattering = get_molecular_scattering_coefficient(h);
fog_scattering = get_fog_scattering_coefficient(h);
extinction = aerosol_absorption + aerosol_scattering + molecular_absorption + molecular_scattering + fog_scattering;
}
//-----------------------------------------------------------------------------
// Spectral rendering stuff
static const float4x3 M = float4x3(
137.672389239975, -8.632904716299537, -1.7181567391931372,
32.549094028629234, 91.29801417199785, -12.005406444382531,
-38.91428392614275, 34.31665471469816, 29.89044807197628,
8.572844237945445, -11.103384660054624, 117.47585277566478
);
float3 linear_srgb_from_spectral_samples(float4 L)
{
return REC709_to_XYZ(mul(L, M));
}
#endif

View File

@@ -0,0 +1,95 @@
#ifndef SKY_INSCATTERING_HLSLI
#define SKY_INSCATTERING_HLSLI
/*
* Buffer B: Sky texture
*
* "A Scalable and Production Ready Sky and Atmosphere Rendering Technique"
* by Sébastien Hillaire (2020).
*
* We render the sky to a texture instead of raymarching on the entire screen.
* This is not very useful in Shadertoy, but very useful for someone looking
* to implement this on a real application.
*
* It is important to note that quality decreases significantly when rendering
* space views. To avoid this, the compute_inscattering() function can be used
* directly when rendering to a fullscreen quad.
*/
#include "sky_common.hlsli"
void compute_inscattering(Texture2D transmittance_lut, SamplerState lut_sampler, in float molecular_phase, in float aerosol_phase, int steps, float3 ray_origin, float3 ray_dir, float t_min, float t_max, float3 sun_dir, inout float4 L_inscattering, inout float4 transmittance)
{
float dt = (t_max - t_min) / float(steps);
for (int i = 0; i < steps; ++i) {
float t = t_min + (float(i) + 0.5) * dt;
float3 x_t = ray_origin + ray_dir * t;
float distance_to_earth_center = length(x_t);
float3 zenith_dir = x_t / distance_to_earth_center;
float altitude = distance_to_earth_center - EARTH_RADIUS;
float normalized_altitude = altitude / ATMOSPHERE_THICKNESS;
float sample_cos_theta = dot(zenith_dir, sun_dir);
float4 aerosol_absorption, aerosol_scattering;
float4 molecular_absorption, molecular_scattering;
float4 fog_scattering;
float4 extinction;
get_atmosphere_collision_coefficients(
altitude,
aerosol_absorption, aerosol_scattering,
molecular_absorption, molecular_scattering,
fog_scattering,
extinction);
float4 transmittance_to_sun = transmittance_from_lut(
transmittance_lut, lut_sampler, sample_cos_theta, normalized_altitude);
float4 ms = get_multiple_scattering(
transmittance_lut, lut_sampler, sample_cos_theta, normalized_altitude,
distance_to_earth_center);
float4 S = sun_spectral_irradiance *
(molecular_scattering * (molecular_phase * transmittance_to_sun + ms) +
(aerosol_scattering + fog_scattering) * (aerosol_phase * transmittance_to_sun + ms));
float4 step_transmittance = exp(-dt * extinction);
// Energy-conserving analytical integration
// "Physically Based Sky, Atmosphere and Cloud Rendering in Frostbite"
// by Sébastien Hillaire
float4 S_factor = (1. - step_transmittance) / max(extinction, 1e-7);
float4 S_int = S * S_factor;
L_inscattering += transmittance * S_int;
transmittance *= step_transmittance;
}
}
float4 get_inscattering(Texture2D transmittance_lut, SamplerState lut_sampler, int steps, float altitude, float3 ray_dir, float t_min, float t_max, float3 sun_dir) {
float cos_theta = dot(-ray_dir, sun_dir);
float molecular_phase = molecular_phase_function(cos_theta);
float aerosol_phase = aerosol_phase_function(cos_theta);
float3 ray_origin = float3(0.0, EARTH_RADIUS + max(altitude, 1.), 0.0);
float atmos_dist = ray_sphere_intersection(ray_origin, ray_dir, ATMOSPHERE_RADIUS);
float ground_dist = ray_sphere_intersection(ray_origin, ray_dir, EARTH_RADIUS);
// We are inside the atmosphere
if (ground_dist < 0.0) {
// No ground collision, use the distance to the outer atmosphere
t_max = min(t_max, atmos_dist);
} else {
// We have a collision with the ground, use the distance to it
t_max = min(t_max, ground_dist);
}
if(t_min >= t_max) return float4(0., 0., 0., 1.);
float4 L = 0.;
float4 transmittance = 1.;
compute_inscattering(transmittance_lut, lut_sampler, molecular_phase, aerosol_phase, 32, ray_origin, ray_dir, t_min, t_max, sun_dir, L, transmittance);
return float4(linear_srgb_from_spectral_samples(L) * exp2(EXPOSURE), dot(transmittance, .25));
}
#endif

View File

@@ -0,0 +1,94 @@
#ifndef SKY_INSCATTERING_HLSLI
#define SKY_INSCATTERING_HLSLI
/*
* Buffer B: Sky texture
*
* "A Scalable and Production Ready Sky and Atmosphere Rendering Technique"
* by Sébastien Hillaire (2020).
*
* We render the sky to a texture instead of raymarching on the entire screen.
* This is not very useful in Shadertoy, but very useful for someone looking
* to implement this on a real application.
*
* It is important to note that quality decreases significantly when rendering
* space views. To avoid this, the compute_inscattering() function can be used
* directly when rendering to a fullscreen quad.
*/
#include "sky_common.hlsli"
void compute_inscattering(Texture2D transmittance_lut, SamplerState lut_sampler, in float molecular_phase, in float aerosol_phase, int steps, float3 ray_origin, float3 ray_dir, float t_min, float t_max, float3 sun_dir, inout float4 L_inscattering, inout float4 transmittance)
{
float dt = (t_max - t_min) / float(steps);
for (int i = 0; i < steps; ++i) {
float t = t_min + (float(i) + 0.5) * dt;
float3 x_t = ray_origin + ray_dir * t;
float distance_to_earth_center = length(x_t);
float3 zenith_dir = x_t / distance_to_earth_center;
float altitude = distance_to_earth_center - EARTH_RADIUS;
float normalized_altitude = altitude / ATMOSPHERE_THICKNESS;
float sample_cos_theta = dot(zenith_dir, sun_dir);
float4 aerosol_absorption, aerosol_scattering;
float4 molecular_absorption, molecular_scattering;
float4 fog_scattering;
float4 extinction;
get_atmosphere_collision_coefficients(
altitude,
aerosol_absorption, aerosol_scattering,
molecular_absorption, molecular_scattering,
fog_scattering,
extinction);
float4 transmittance_to_sun = transmittance_from_lut(
transmittance_lut, lut_sampler, sample_cos_theta, normalized_altitude);
float4 ms = get_multiple_scattering(
transmittance_lut, lut_sampler, sample_cos_theta, normalized_altitude,
distance_to_earth_center);
float4 S = sun_spectral_irradiance *
(molecular_scattering * (molecular_phase * transmittance_to_sun + ms) +
(aerosol_scattering + fog_scattering) * (aerosol_phase * transmittance_to_sun + ms));
float4 step_transmittance = exp(-dt * extinction);
// Energy-conserving analytical integration
// "Physically Based Sky, Atmosphere and Cloud Rendering in Frostbite"
// by Sébastien Hillaire
float4 S_int = (S - S * step_transmittance) / max(extinction, 1e-7);
L_inscattering += transmittance * S_int;
transmittance *= step_transmittance;
}
}
float4 get_inscattering(Texture2D transmittance_lut, SamplerState lut_sampler, int steps, float altitude, float3 ray_dir, float t_min, float t_max, float3 sun_dir) {
float cos_theta = dot(-ray_dir, sun_dir);
float molecular_phase = molecular_phase_function(cos_theta);
float aerosol_phase = aerosol_phase_function(cos_theta);
float3 ray_origin = float3(0.0, EARTH_RADIUS + max(altitude, 1.), 0.0);
float atmos_dist = ray_sphere_intersection(ray_origin, ray_dir, ATMOSPHERE_RADIUS);
float ground_dist = ray_sphere_intersection(ray_origin, ray_dir, EARTH_RADIUS);
// We are inside the atmosphere
if (ground_dist < 0.0) {
// No ground collision, use the distance to the outer atmosphere
t_max = min(t_max, atmos_dist);
} else {
// We have a collision with the ground, use the distance to it
t_max = min(t_max, ground_dist);
}
if(t_min >= t_max) return float4(0., 0., 0., 1.);
float4 L = 0.;
float4 transmittance = 1.;
compute_inscattering(transmittance_lut, lut_sampler, molecular_phase, aerosol_phase, 32, ray_origin, ray_dir, t_min, t_max, sun_dir, L, transmittance);
return float4(linear_srgb_from_spectral_samples(L) * exp2(EXPOSURE), dot(transmittance, .25));
}
#endif

View File

@@ -0,0 +1,40 @@
#include "sky_common.hlsli"
#define TRANSMITTANCE_STEPS 32
float4 main(in float2 uv : Position) : SV_Target0
{
float sun_cos_theta = uv.x * 2.0 - 1.0;
float3 sun_dir = float3(-sqrt(1.0 - sun_cos_theta*sun_cos_theta), 0.0, sun_cos_theta);
float distance_to_earth_center = lerp(EARTH_RADIUS, ATMOSPHERE_RADIUS, uv.y);
float3 ray_origin = float3(0.0, 0.0, distance_to_earth_center);
float t_d = ray_sphere_intersection(ray_origin, sun_dir, ATMOSPHERE_RADIUS);
float dt = t_d / float(TRANSMITTANCE_STEPS);
float4 result = 0.0;
for (int i = 0; i < TRANSMITTANCE_STEPS; ++i) {
float t = (float(i) + 0.5) * dt;
float3 x_t = ray_origin + sun_dir * t;
float altitude = length(x_t) - EARTH_RADIUS;
float4 aerosol_absorption, aerosol_scattering;
float4 molecular_absorption, molecular_scattering;
float4 fog_scattering;
float4 extinction;
get_atmosphere_collision_coefficients(
altitude,
aerosol_absorption, aerosol_scattering,
molecular_absorption, molecular_scattering,
fog_scattering,
extinction);
result += extinction * dt;
}
float4 transmittance = exp(-result);
return transmittance;
}

View File

@@ -0,0 +1,23 @@
#ifndef VIEW_DATA_HLSLI
#define VIEW_DATA_HLSLI
cbuffer DrawConstants : register(b2) {
float4x4 g_InverseModelView;
float4x4 g_InverseProjection;
float3 g_LightDir;
float g_Altitude;
float3 g_LightColor;
float g_Time;
}
float2 PixelToCS(in float2 pixel, in float2 size) {
return ((pixel + .5) / size - .5) * float2(2., -2.);
}
float3 ReconstructPos(in float2 cs, in float depth) {
float4 ndc = float4(cs, depth, 1.);
ndc = mul(g_InverseProjection, ndc);
return ndc.xyz / ndc.w;
}
#endif