Next Shader Development

[nani]

When the angles is pi/3 it can be considered a straight line so you can do an if (x > pi/3 ) then line interpolation etc so then you dont have that 0.02 error

[DanW58]

I'm COMPLETELY lost, Nani.  Please explain the whole thing, not just one part.

Edited March 2, 2021 by DanW58

[nani]

https://en.wikipedia.org/wiki/Trigonometric_functions#Power_series_expansion

put the cos x  expression from above in the f(x) = 1 - cos(x) and in our case using the first three terms of the expansion: 1 - x^2/2! + x^4/4!  but you can add more if you want more precision, then you do algebra simplification and stuff and gives you the result

[DanW58]

I still don't understand.  My input is in the range from 1 to 5.

IOW, I'm trying to compute % of diffuse light escaping as a function of refractive index.

What you're talking about is a trigonometry approximation;  no?

I think you're going from angle to solid angle;  I need a solution from refractive index.  There's a cosine function in glsl;  I don't need taylor expansion for that.

 

Edited March 2, 2021 by DanW58

[nani]

Ahhh, then how is it computed in the spreadsheet you posted? I only saw 1-cos(x) so i guessed you wanted that.

What is the complete function you want to approximate?

[DanW58]

Hahaha,  well,  that's the problem;  I don't know what the function is.  I made the huge table in page 1 of the spreadsheet to find out the angle of semi-reflection for each index, put the results into a table in page 2;  now the trick is inventing a function that fits the data.

y = 0.514 / x^7  +  0.486 / x^2

works nicely;  but maybe something not involving powers would be nicer, performance wise.

FresnelTable.ods

Edited March 2, 2021 by DanW58

[nani]

Uff, two dimensional data + many terms, I guess you could try using some toolbox from matlab (or octave) and see what it gives, havent tried that

[DanW58]

Well, no, it's not 2-dimensional;  it is just a linear function:  solid_angle = f( refractive_index ).

Inverting the Fresnel function would take me centuries;  that's why I did it this way.

I have some software that's much like Matlab, but never got around to learn it.

Many years ago I was using a beautiful app called DataFit;  but it is pay-for software, and I got no money.

This is the data:

refr.   hemispheres
index   solid angle

1.0       1.0000
1.5       0.2499
2.0       0.1314
2.5       0.0819
3.0       0.0562
3.5       0.0410
4.0       0.0312
4.5       0.0247
5.0       0.0200

y = 0.514 / x^7 + 0.486 / x^2   does it, but something cheaper to compute would be nice.  For example a trig function multiplied by a constant ... something along such lines ...

What's 1/RI?  1.0 ... 0.2

1-that = 1-(1/RI) = 0 ... 0.8

cos( (1-(1/RI)) * pi/2 ) = ...

Nah, doesn't work.

Edited March 2, 2021 by DanW58

[nani]

With matlab I found a good one:
f(x) = 0.6217/( x^2 + 1.207*x - 1.585 )

has an R^2 of 1  
 

[badosu]

as mentioned before, for low theta:

for theta ~ pi/2 should be similar (inverting the functions and adjusting approximations).

A great solution (from my limited experience with audio synthesis), if the function is periodic, is to use a wavetable and interpolate between points. i.e. use a sample step, store all function results to a hash, access the hash for the previous and next value for a particular input, interpolate linearly between.

PS: Ugh, nevermind it seems the function is more complicated than that. Still at least you can get the trig approximations that nani suggested.

Edited March 2, 2021 by badosu

[DanW58]

Thanks!  Back from a nap and you guys solved all my problems.  Well, a table is not an option with a shader.  Actually,  it is;  I have seen people commit complex things, such as the Fresnel formula, to a little texture, and read it for an answer.  But I don't like that;  it is slow in its own way, ocupies a texture unit, and has a precision problem;  but if many cases of complex functions crop up it might make sense to put them all into a single 16-bit texture and have lookups for several functions.

But definitely  f(x) = 0.6217/( x^2 + 1.207*x - 1.585 )  is a lot simpler to compute than my powerz abomination.

And it works like a charm and a half, too!

So, in glsl:
 

float DiffuseEscapeFractionFromDielectric( float RI )
{
  float temp = (2.207 * RI) + (RI * RI) - 1.585;
  return 0.6217 / temp;
}

So, I'm going to make a drawing for diffuse case;  so I'll be back in 2 hours or so. 

Edited March 2, 2021 by DanW58

[DanW58]

Okay, so, the case for a dielectric, transparent, glossy layer over a diffuse base begins similarly to the case of dielectric over specular metal.  We have a sunray coming in, 'a',  that hits our dielectric and splits into a reflected ray, 'b' (the name of which got cut off in the pic), and a refracted ray 'c'.  We already have the math for all that,  except for one detail:  My calculation of reflected light by the diffuse base was based on ndotl (normal dot lightray vector, or, the cosine of the angle), but light refracts and changes angle, so I should consider the refracted lightray, versus the original lightray.  The good news is that the new Fresnel function gives me the cosine of the refracted angle already, as it needs to compute it internally anyways;  so I don't need a special calculation;  I just need to use that output of the function.

So, our concern here is what happens from the moment 'c' is about to hit the diffuse layer.

When 'c' hits the diffuse base, it explodes in all directions evenly.  At least that is the simplified model of diffuse optics.

Rays shoot out in a semi-sphere of even distribution.  Some of the rays, shown in green, make it out of the dielectric layer;  and some of them, shown in dark red, reflect off the dielectric interface, heading downwards again.

This is simplified because the rays escaping are actually split, with a portion reflecting;  however the portion is miniscule and can be ignored here.

I put gray lines between the green and red, forming a triangle.  This is the escape cone, whose half-angle (radius), which I marked in the drawing as 'x',  and therefore solid angle, are a function of the material's refractive index alone.  The higher the RI, the narrower the escape cone gets.  THIS is what the function in the previous posts was all about.

But note, as an aside, that the distribution of rays coming out of the dielectric and into the air is not limited to the cone of escape, because they refract, expanding again to a semi-spherical distribution.

With this we can calculate how much light will refract out, and how much light will reflect back in to produce yet another diffuse spherical distribution, and so on and so forth.  So, given a 'c' ray, let's follow the saga in pseudo-code:

void saga_of_a_photon( vec3 ray_c, vec3 normal, vec3 MatDiffuseRGB, vec3 RefractiveIndex
                           vec3 RefractingOutRGBlight, vec3 ReflectingBackInRGBlight )
{
  vec3 ExplodingRGBlight = dot(ray_c,normal) * MatDiffuseRGB; //diffuse explosion
  RefractingOutRGBlight = ExplodingRGBlight * RefractiveEscapeFunction( RefractiveIndex );
  ReflectingBackInRGBlight = ExplodingRGBlight - RefractingOutRGBlight;
}

Before you jump horrified that RefractiveIndex is a vec3 instead of a float, there's a reason for it:  most materials' refractive index changes with wavelength, so it is not the same for R, G or B.  In the latest water shader I have the index of refraction for water as a vec3, not that I'm going to boast of noticing the difference in the results, but every little bit helps.  And you might ask, how am I going to specify these "anisochromatic" refractive indices with only one texture channel, the definitive answer is I don't know yet.  But seriously, chromium oxide, whose properties are what causes iridescence in chrome-plated Harley exhaust pipes, has steeply changing refractive indexes.  To model iridescence I will need to model them correctly.

There is one problem that arises from the picture and the pseudocode.  I was going to grab the middle of the 'c' ray as my starting point to define the repeating cycle to then solve as standard geometric series.  But now I realize I cannot do that, because the angle that 'c' is coming at is input light vector-dependent, whereas the reflecting-back-in rays are not.

In fact, I would say that the average angle of light reflecting back in is (90 degrees - x) / 2.  Don't forget we need to compute dot(ray,normal) at each diffuse bounce.  So, my repeating cycle really begins AFTER the first diffuse explosion.

Time for a break;  I'll work on the other shaders for a bit.

Edited March 3, 2021 by DanW58

[DanW58]

I just had an idea for how to solve the problem of only having one channel for refractive index but needing separate RGB refractive indexes for many materials.  The idea is to organize all materials from the available databases online into a big spreadsheet.  Where two materials have the exact same refractive index, use the more commonly found in nature.  Sort the spreadsheet by refractive index on the green channel.  Then make a texture where the green channel increases linearly across and the red and blue channels represent deviation from the green x16, for better resolution.  So, start by putting the red and blue channels all at zero, then for each known material along the green line write the red and blue to encode its deviation from green. And there will be many empty spaces in the end, as there are, for example, no known materials with refractive index of 1.1, or 1.2;  they jump from.  1.05 to 1.25.  Similarly, there's not too many materials of refractivity around 2, except boron nitride.  So, for the empty spaces, might as well interpolate values between the nearest actual materials.

The way this works is, the standard texture set specifies an exact refractive index.  If it is 1.4603, we compute an offset into this texture and read the value, and it comes back (after scaling) as (1.4591,1.4603,1.4672).  If instead the main texture specifies 1.5464, a lookup in this texture returns (1.545,1.5464,1.5544).  The first is for fused quartz;  the second is for natural quartz.  A read of 1.77 will return (1.77,1.77,1.78), aluminium oxide.  A read of 2.15 returns (2.05,2.15,2.3), the wildly steep refractive indices of chromium's most common oxide.

And I know I was just ranting against the use of textures as lookup tables, but that is because I don't feel that a lookup is justified for functions that can be computed.  But in the case of materials there's a lot to be said for having a lookup table for them ...  Materials are as computable as phone-books.

Edited March 4, 2021 by DanW58

[DanW58]

Alright, taking a break from the other shaders.  Back to the world of Sanity! 

So as I finished my second-last post I was saying that the repeating cycle for diffuse reflections does NOT start with ray 'c' going towards the diffuse layer.  This is because the angle matters, and it will be a constant dependent on index of refraction alone after internal reflections begin, but at the stage of ray 'c' it is dependent on the incident light's angle.

So, our photon saga now looks like this (separating the repeating part):

//where..
void RealFresnel( vec3 NdotL, vec3 IOR_RGB, out vec3 Frefl, inout vec3 sin_t, inout vec3 cos_t )
{
  vec3 Z2 = WHITE / IOR_RGB;           // Assumes n1 = 1 thus Z1 = 1.
  vec3 cos_i = NdotL;                  // assignnment for name's sake
  vec3 sin_i = sqrt( WHITE - cos_i*cos_i );
  sin_t = min(WHITE, sin_i * Z2);      // Outputs sin(refraction angle).
  cos_t = sqrt( WHITE - sin_t*sin_t ); // Outputs cos(refraction angle).
  vec3 Rs = (Z2*cos_i-cos_t) / (Z2*cos_i+cos_t);
  vec3 Rp = (Z2*cos_t-cos_i) / (Z2*cos_t+cos_i);
  Frefl = mix( Rs*Rs, Rp*Rp, 0.5 );    // Outputs reflectivity.
}

//The saga, from the beginning (from sun-ray):
void saga_of_a_photon( vec3 ray_A, vec3 normal, vec3 MatDiffuseRGB, vec3 RefractiveIndex,
                                  out vec3 specularFactorRGB, out vec3 diffuseFactorRGB )
{
  //Getting ready ...
  float NdotL = max( 0.0, dot( ray_A, normal ) );
  vec3  ReflectivityRGB;
  vec3  sinRefrAngle;
  vec3  cosRefrAngle;
  RealFresnel( vec3(NdotL), RefractiveIndex, Reflectivity, sinRefrAngle, cosRefrAngle );
  float EscapeFraction = getDiffuseEscapeFractionFromDielectric( RefractiveIndex );
  float cosAvgReflAngle = getCosOfAvgReReflectionAngle( RefractiveIndex );
  vec3  diffuseFactorRGBacc = vec3(0.0);
  //repeating cycle vars:
  vec3  ExplodingRGBlight
  vec3  LightEscapingRGB;
  vec3  LightIncomingRGB;

  //First reflection (ray_B):
  specularFactorRGB = Reflectivity;

  //First refraction (ray_C):
  LightIncomingRGB = vec3(1.0) - Reflectivity;

  //First diffuse bounce
  vec3 ExplodingRGBlight = LightIncoming * cosRefrAngle * MatDiffuseRGB;

  //<Sound of trumpets here.>

  //First iteration of our repeating cycle:
  LightEscapingRGB = ExplodingRGBlight * EscapeFraction;
  diffuseFactorRGBacc += LightEscapingRGB;
  LightIncomingRGB = ExplodingRGBlight - LightEscapingRGB;
  vec3 ExplodingRGBlight = LightIncoming * cosAvgReflAngle * MatDiffuseRGB;

  //Second iteration of our repeating cycle:
  LightEscapingRGB = ExplodingRGBlight * EscapeFraction;
  diffuseFactorRGBacc += LightEscapingRGB;
  LightIncomingRGB = ExplodingRGBlight - LightEscapingRGB;
  vec3 ExplodingRGBlight = LightIncoming * cosAvgReflAngle * MatDiffuseRGB;

  //..............................................
  // repeat infinite times; then ...
  //..............................................

  diffuseFactorRGB = diffuseFactorRGBacc;
}

The function getCosOfAvgReReflectionAngle( RefractiveIndex ) is the next problem to solve.

What it represents is the cosine of the average angle to the normal for light reflecting back in.  It is really a function of material refractivity alone, same as the escape cone.  It is the cosine of half the complement of the cone of escape's radius.  Time for some spreadsheet work...

 

EDIT:  Don't worry about the "repeat infinite times";  that's just rethorical.  Once we get our per-cycle factor x, we have a case of

 f(x) = 1 + x + x^2 + x^3 + x^4 + ...  which is a geometric series and solves as

 = 1/(1-x)  for absolute values of x less than 1;  or...

f(x) = x + x^2 + x^3 + x^4 + ... =  1/(1-x) - 1;

Edited March 4, 2021 by DanW58

[DanW58]

@nani  Any chance you could do me one more matlab favor?  This is the angle of incidence after every diffuse bounce and internal reflection, as a function of refractive index:

refractive    cos of average
  index       reflect angle
==========     ===========
   1.0            0.0000
   1.5            0.4115
   2.0            0.5023
   2.5            0.5494
   3.0            0.5786
   3.5            0.5985
   4.0            0.6132
   4.5            0.6242
   5.0            0.6329

Thanks in advance.

Edited March 4, 2021 by DanW58

[DanW58]

In fact, it is not important that it include refractive index of 1.0, as that will never be used.  Instead, it would be more important that it include water's refractivity...

refractive    cos of average
  index       reflect angle
==========     ===========
   1.335 (H2O)    0.3576
   1.5            0.4115
   2.0            0.5023
   2.5            0.5494
   3.0            0.5786
   3.5            0.5985
   4.0            0.6132
   4.5            0.6242
   5.0            0.6329

You might ask why ...  Because this modeling of glossy paint is applicable to meters-thick dielectrics as much as it is to hair-thin top-coats.  It will help make water even more realistic looking ...

((In fact, the water shader will be my first test of this.))

EDIT:  Actually, the first test will be a math test:  With a refractive index of 1.0 I should get the same result as the diffuse base, namely zero specular and standard diffuse lighting.  This has to work because, to this shader, a simple material without a dielectric coat is simply a two-layer material whose coat's refractive index is 1.

EDIT2:

Here's a spreadsheet work update.  The values come from the rightmost column on sheet 2.

FresnelTable.ods

Edited March 4, 2021 by DanW58

[nani]

1) f(x) = (0.6821*x - 0.6818)/(x - 0.6562)
2) f(x) = (0.7012*x - 0.6062)/(x - 0.4146)

[DanW58]

Many thanks, friend!

EDIT:  I would never have found that in a million years!

EDIT2:  ABSOLUTELY PERFECT!

Edited March 4, 2021 by DanW58

[nani]

To be honest matlab curve fit toolbox is quite amazing, haven't found anything similar to it. Xd

[DanW58]

Seems to be as good as my old beloved DataFit, which never thought would be surpassed or reached.  It had whole families of forms of formulas to look for, like if you preferred a straight polynomial it tried that way, or if you didn't mind an inverse of a polynomial, or something with logs... you could pick and choose the type of formula you wanted, or just tell it to minimize computation costs, and it would go to work.  Seems from what you found that matlab has similar capabilities.

 

[DanW58]

Okay, this is the whole thing put together now:

const vec3 WHITE = vec3(1.0);
const vec3 BLACK = vec3(0.0);

float getDiffuseEscapeFractionFromDielectric( float RI )
{
  float temp = (2.207 * RI) + (RI * RI) - 1.585;
  return 0.6217 / temp;
}

float getCosOfAvgReReflectionAngle( float RI )
{
  return (0.7012 * RI - 0.6062) / (RI - 0.4146);
}

void RealFresnel( vec3 NdotL, vec3 IOR_RGB, out vec3 Frefl, inout vec3 sin_t, inout vec3 cos_t )
{
  vec3 Z2 = WHITE / IOR_RGB;           // Assumes n1 = 1 thus Z1 = 1.
  vec3 cos_i = NdotL;                  // assignnment for name's sake
  vec3 sin_i = sqrt( WHITE - cos_i*cos_i );
  sin_t = min(WHITE, sin_i * Z2);      // Outputs sin(refraction angle).
  cos_t = sqrt( WHITE - sin_t*sin_t ); // Outputs cos(refraction angle).
  vec3 Rs = (Z2*cos_i-cos_t) / (Z2*cos_i+cos_t);
  vec3 Rp = (Z2*cos_t-cos_i) / (Z2*cos_t+cos_i);
  Frefl = mix( Rs*Rs, Rp*Rp, 0.5 );    // Outputs reflectivity.
}

//The saga, from the beginning (from sun-ray):
void saga_of_a_photon( vec3 ray_A, vec3 normal, vec3 MatDiffuseRGB, vec3 RefractiveIndex,
            vec3 AmbientRGBlight, out vec3 specularFactorRGB, out vec3 diffuseFactorRGB )
{
  //Getting ready ...
  float NdotL = max( 0.0, dot( ray_A, normal ) );
  vec3  ReflectivityRGB;
  vec3  sinRefrAngle;
  vec3  cosRefrAngle;
  RealFresnel( vec3(NdotL), RefractiveIndex, Reflectivity, sinRefrAngle, cosRefrAngle );
  float EscapeFraction = getDiffuseEscapeFractionFromDielectric( RefractiveIndex );
  float cosAvgReflAngle = getCosOfAvgReReflectionAngle( RefractiveIndex );
  vec3  diffuseFactorRGBacc = BLACK;
  //repeating cycle vars:
  vec3  ExplodingRGBlight
  vec3  LightEscapingRGB;
  vec3  LightIncomingRGB;

  //First reflection (ray_B):
  specularFactorRGB = Reflectivity;

  //First refraction (ray_C):
  LightIncomingRGB = WHITE - Reflectivity;
  
  //First diffuse bounce
  vec3 ExplodingRGBlight = LightIncoming * cosRefrAngle;
  //Adding ambient light here:
  ExplodingRGBlight += (vec3(0.5(EscapeFraction+1.0))*AmbientRGBlight);
  ExplodingRGBlight *= MatDiffuseRGB;

  //<Sound of trumpets here.>

  //First iteration of our repeating cycle:
  LightEscapingRGB = ExplodingRGBlight * EscapeFraction;
  diffuseFactorRGBacc += LightEscapingRGB;
  LightIncomingRGB = ExplodingRGBlight - LightEscapingRGB;
  ExplodingRGBlight = LightIncoming * cosAvgReflAngle * MatDiffuseRGB;

  //Second iteration of our repeating cycle:
  LightEscapingRGB = ExplodingRGBlight * EscapeFraction;
  diffuseFactorRGBacc += LightEscapingRGB;
  LightIncomingRGB = ExplodingRGBlight - LightEscapingRGB;
  ExplodingRGBlight = LightIncoming * cosAvgReflAngle * MatDiffuseRGB;

  //..............................................
  // repeat infinite times; then ...
  //..............................................

  diffuseFactorRGB = diffuseFactorRGBacc;
}

Note that I added ambient light there;  that's because ambient light follows the same multiple bounce rule as the multi-bounce cycle of diffuse lighting.  And the choice of 0.5*(EscapeFraction+1.0) as a factor for the ambient light, it is a bit of a gut feeling...  You would think the cone of escape is the cone of ambient light that gets in there.  By that reasoning you would use EscapeFraction as the factor.  However, light refracts from 180 degrees all around INTO that cone;  so by this reasoning you'd say the factor should be 1.0.  However, not all of the light refracts in;  some of it reflects out.  So I figured somewhere between EscapeFraction and 1.0 lies the true answer.

Next I'm going to try to simplify the repeating cycle to a single factor.

Edited March 4, 2021 by DanW58

[DanW58]

In the repeating cycle,

  LightEscapingRGB = ExplodingRGBlight * EscapeFraction;
  diffuseFactorRGBacc += LightEscapingRGB;
  LightIncomingRGB = ExplodingRGBlight - LightEscapingRGB;
  ExplodingRGBlight = LightIncoming * cosAvgReflAngle * MatDiffuseRGB;

  LightEscapingRGB = ExplodingRGBlight * EscapeFraction;
  //.....

all that matters to us is light escaping.  The ratio of light escaping's for subsequent cycles is constant;  we just have to work out what it is, now.

The second line is just output to an accumulator.

The third line is a bit of a problem, as it relies on the previous cycle's ExplodingLight and LightEscaping to obtain light coming back in for another diffuse bounce.  However, we can also compute incoming light from LightEscaping.

How?

LightEscaping is previous Exploding times a constant for this call, EscapeFraction.

LightIncoming is previous Exploding times (1.0-EscapeFraction).

So LightIncoming = ((1.0 - EscapeFraction)/EscapeFraction)*LightEscaping

Yes?

Then we have the next Exploding bounce = LightIncoming * cosAvgReflAngle * MatDiffuseRGB.

And finally the next escape:  LightEscaping = Exploding * EscapeFraction.

Putting it all together:

NextEscape = ((1.0-EscapeFraction)/EscapeFraction)*ThisEscape*cosAvgReflAngle*MatDiffuseRGB*EscapeFraction.

But now, if EscapeFraction is a constant, then (1-EscapeFraction)/EscapeFraction is a constant too.

Furthermore, far to the right we are multiplying by EscapeFraction, so divide and multiply simplify out.

So, let DontEscapeFraction = 1.0-EscapeFraction, and we can write,

NextEscape = DontEscapeFraction*ThisEscape*cosAvgReflAngle*MatDiffuseRGB

So, NextEscape/ThisEscape = DontEscapeFraction*cosAvgReflAngle*MatDiffuseRGB

And so we have that after we calculate light escaping from the first diffuse bounce, call it FirstDiffuse, the subsequent escapes will be,

float r = (1.0-EscapeFraction)*cosAvgReflAngle*MatDiffuseRGB;

TotalDiffuse = FirstDiffuse * (1 + r + r^2 + r^3 + ... ) = FirstDiffuse * (1/(1-r))

Next post I write the simplified code.

 

EDIT:  Some astute reader might ask "Haven't you forgot specular ambient light?

This is one HUGE mistake often made by shader hackers...

Specularly reflected ambient light ***IS*** the environment mapping.  Environment mapping accounts for 100% of specularly reflected ambient light.  In fact, ambient light should not exist in a shader that has a sky-box, as readings of ambient light can be taken from the sky-box itself, along the object's normal (or bent normal).  Ambient light is an average of the environment;  not something that in any way "adds to" the environment.

Edited March 4, 2021 by DanW58

[DanW58]

Here it is:

float getDiffuseEscapeFractionFromDielectric( float RI )
{
  float temp = (2.207 * RI) + (RI * RI) - 1.585;
  return 0.6217 / temp;
}

float getCosOfAvgReReflectionAngle( float RI )
{
  return (0.7012 * RI - 0.6062) / (RI - 0.4146);
}

void RealFresnel
(
  vec3 NdotL,
  vec3 IOR_RGB,
  out vec3 Frefl,
  inout vec3 sin_t,
  inout vec3 cos_t
)
{
  vec3 Z2 = WHITE / IOR_RGB;           // Assumes n1 = 1 thus Z1 = 1.
  vec3 cos_i = NdotL;                  // assignnment for name's sake
  vec3 sin_i = sqrt( WHITE - cos_i*cos_i );
  sin_t = min(WHITE, sin_i * Z2);      // Outputs sin(refraction angle).
  cos_t = sqrt( WHITE - sin_t*sin_t ); // Outputs cos(refraction angle).
  vec3 Rs = (Z2*cos_i-cos_t) / (Z2*cos_i+cos_t);
  vec3 Rp = (Z2*cos_t-cos_i) / (Z2*cos_t+cos_i);
  Frefl = mix( Rs*Rs, Rp*Rp, 0.5 );    // Outputs reflectivity.
}

void DielectricOverDiffuse
(
  vec3 ray_A,
  vec3 normal,
  vec3 MatDiffuseRGB,
  vec3 RefractiveIndex,
  vec3 AmbientRGBlight,
  out vec3 specularFactorRGB,
  out vec3 diffuseFactorRGB
)
{
  //Getting ready ...
  float NdotL = max( 0.0, dot( ray_A, normal ) );
  vec3  ReflectivityRGB;
  vec3  sinRefrAngle;
  vec3  cosRefrAngle;
  RealFresnel( vec3(NdotL), RefractiveIndex, ReflectivityRGB, sinRefrAngle, cosRefrAngle );
  float EscapeFraction = getDiffuseEscapeFractionFromDielectric( RefractiveIndex );
  float cosAvgReflAngle = getCosOfAvgReReflectionAngle( RefractiveIndex );

  //First reflection (ray_B):
  specularFactorRGB = ReflectivityRGB; //straight output

  //First refraction and diffuse bounce
  vec3 temp3 = (WHITE-ReflectivityRGB) * cosRefrAngle;
  //Adding ambient light here:
  temp3 += (vec3(0.5(EscapeFraction+1.0))*AmbientRGBlight);
  //Multiplying by the diffuse color
  temp3 *= MatDiffuseRGB;
  //And computing the first diffuse escape:
  temp3 *= EscapeFraction;
  
  //And now the great simplification for infinite bounces:
  vec3 r = vec3((1.0-EscapeFraction) * cosAvgReflAngle) * MatDiffuseRGB;
  diffuseFactorRGB = temp3 / (WHITE-r);
}

Next post I will try to combine specular and diffuse into one function.  I usually keep them separate, but the problem is that in the case of a two-layer material many of the calculations are common to specular and diffuse.  And as you saw, we computed outgoing Fresnel spec, even though this was supposed to be the diffuse routine.

There's another advantage in combining diffuse and specular in one function:  I'm not sure how good the compiler is at optimizing a shader on load, but it can't hurt to pre-optimize by interspersing specular and diffuse calculations to reduce dependencies.

EDIT:  I'm having second thoughts about adding ambient in this routine, though.  In a multiple-light situation, where this routine would be called for each light, doing so would add ambient light multiple times.  Better leave it out.

EDIT: OTOH, the fact that ambient light bounces multiple times too, using the same instructions instead of repeating them, is a good thing.  In the case of multiple lights, a user of this shader could add ambient to the first call and pass BLACK as argument for ambient light in subsequent calls.

I'll postpone this decision for now.

Edited March 4, 2021 by DanW58

[DanW58]

Quoting myself,

On 27/02/2021 at 10:05 PM, DanW58 said:

So, at each iteration, re-reflected light is RRratioRGB = refract_factor * specular_color.  So the total light re-reflected after infinite bounces is RRratioRGB^1 + RRratioRGB^2 + RRratioRGB^3 + ...

This is solving an infinite geometric series.

Where r meets the criteria −1 < r < 1,  a + ar + ar² + ar³ + ... = a/(1−r)

So  total_re_reflected_light = c_RGB / ( RGB(1,1,1) - RRratioRGB )

and the total emerging light (rays f, i, etc.) is that times ((1/fresnel_reflection)-1)

but since we don't care about the light that dies bouncing back and forth, we ONLY care about the light that we see, we can simply put the formulas together and say that, for the specular base case of specColRGB, we proceed as follows, in pseudo-glsl now:
 

void FresnelSplit( vec3 ray, vec3 normal, float n1, float n2, vec3 specColRGB,
                   out vec3 bounce1, out vec3 bounceN )
{
  float raydotnormal = max( 0.0, dot( normalize(ray), normalize(normal)) );
  float R0 = (n2-n1)/(n2+n1);  R0 = R0*R0;
  float angle_part = pow( 1.0 - raydotnormal, 5 );
  reflectivity = R0 + (1.0-R0) * angle_part;
  bounce1 = vec3(reflectivity);
  float refractivity = 1.0 - reflectivity;
  float RRratio = vec3(reflectivity) * specColRGB;
  vec3 dead_light = vec3(refractivity) / ( vec3(1.0) - RRratio );
  bounceN = vec3((1.0/reflectivity)-1.0) * dead_light;
}

That's where we left off the specular case.

I don't remember why I called this "FresnelSplit";  it should be called "DielectricOverMSpec()".

Anyways, one thing to say is, I am not sure yet, in the case of dielectric over specular, when a specular power (smoothness) is specified, whether it applies to the metallic specularity underneath or to the dielectric surface above  (or both).  And necessity may force the choice...

Thing is, roughness in the dielectric surface would affect the appearance of an ultra-smooth metallic reflector underneath by refraction.  And the simplest way to represent that is by giving the metallic specularity the same specular power as to the dielectric.  And this is really the simplest course of action:  take the two specularities together and Phong them as one.

So, next I will try to put specular and diffuse together, into a Reflectivity() function.

Just to remind us all, myself included, of where we are at in the shader:  This function I'm about to put together does not have sunlight color or any such thing.

Here sunlight is assumed at a generic 1.0 and output may be 0.6 for specular and 0.2 for diffuse.  RGB, of course;  I'm just simplifying.  Another thing to note is that environment mapping is not part of this, nor is ambient light (finally arrived to a decision).  This is only from light source.  It also does not compute phong, and therefore doesn't take the eye-vector as an input.  This simply returns specular reflectivity (combined metallic and dielectric for the given light vector) and the diffuse reflectivity;  both of them in RGB less than WHITE.  In fact, the sum of the two is guaranteed to be less than 1.0, per channel. 

So, to summarize, this is just a per-light source reflectivity calculator.

I'm also thinking now that I might as well dispense with the ray and normal inputs, and take a single NdotL float as input.  Reason being that other functions will need NdotL;  no need to recalculate if it is external.

Edited March 5, 2021 by DanW58

[DanW58]

Done!    (But some secretive cabal out there has blown its cover for accusing me of "not doing anything productive" ...  Figures!)

void Reflectivity( float raydotnormal, vec3 RefractiveIndex, vec3 specColRGB,
  vec3 MatDiffuseRGB, out vec3 specularFactorRGB, out vec3 diffuseFactorRGB )
{
  vec3  FReflectivityRGB;
  vec3  sinRefrAngle;
  vec3  cosRefrAngle;
  RealFresnel( vec3(raydotnormal), RefractiveIndex, ReflectivityRGB, sinRefrAngle, cosRefrAngle );
  vec3  FRefractivityRGB = WHITE-FReflectivityRGB;
  vec3  EscapeFraction = getDiffuseEscapeFractionFromDielectric( RefractiveIndex );
  float cosAvgReflAngle = getCosOfAvgReReflectionAngle( RefractiveIndex );

  //specular:
  vec3  RRratio = FReflectivityRGB * specColRGB;
  vec3  dead_light = FRefractivityRGB / ( WHITE - RRratio );
  //simplification for infinite specular bounces:
  vec3  nbouncesRGB = vec3((1.0/reflectivity)-1.0) * dead_light;
  specularFactorRGB = FReflectivityRGB + nbouncesRGB;
  
  //diffuse:
  //First refraction and diffuse bounce
  vec3 temp3 = (WHITE-ReflectivityRGB) * cosRefrAngle * MatDiffuseRGB;
  //And computing the first diffuse escape:
  temp3 *= EscapeFraction;
  //simplification for infinite diffuse bounces:
  vec3 r = (WHITE-EscapeFraction) * vec3(cosAvgReflAngle) * MatDiffuseRGB;
  diffuseFactorRGB = temp3 / (WHITE-r);
}

Next step is to test this in the water shader.

EDIT:  Wait a minute;  this is not going to be easy to test in the water shader...  It is doable, I'm sure;  but not easy.  First of all, I'm already computing Fresnel for sunlight going INTO the water, and I'm doing that in the vertex shader, to relieve the fragment shader of some burden.  Secondly, the refractive and reflective parts in that shader are separate functions.  The refractive function is HUGE, and I don't understand half of it.

I think I can do it,  but it's going to be a hell of a lot of work.  Won't be ready tonight, by any stretch of the imagination.

Or, how else can I test this?  In case you just stumbled here, this is (at this moment) a shader for glossy paint.  Why a shader for glossy paint for a first century game?  Well, a shader for glossy paint can easily be tuned to represent cheaper paints, plastics and most biological things.  A glossy paint shader is really the starting point for representing non-metals other than matte materials.  If a shader can represent a glossy paint well AND a metal well, it can represent just about anything that exists.  Well, transparent materials are a separate category.

Problem with testing is, the current texture set doesn't have a channel for gloss.

I suppose I could do a color selection, like make all green things look glossy...

Give all plants a refractive index of 2.5;  make them look like Diamond Age plants

Edited March 5, 2021 by DanW58