EEVEE: Implement the missing Sky texture

I'm not sure if the Sky was deliberately left out or was just waiting for a
better moment, but so many I was disappointed that Sky in EEVEE is
completely white.

There are already 2 implementations (osl and gpu) so this is the third one.
Looking at other cases it seems that we are not supposed to share sources
between cycles and the rest? So the new util_sky_model files are just
copies of what is already in cycles, except that the data file uses the RGB
variant of the Hosek/Wilkie model, because we output RGB anyway (but can be
easily changed to XYZ if desired - the results are nearly identical).
I am not sure if it is okay to pass 3*9 float values as 3 mat4 uniforms (I
wanted to use mat3 but it does not work).
Also, should I cache the sky model data between renders if the parameters
do not change?

Reviewed By: fclem, brecht

Differential Revision: https://developer.blender.org/D7108

5 files changed, 306 insertions, 6 deletions

diff --git a/source/blender/gpu/shaders/material/gpu_shader_material_tex_sky.glsl b/source/blender/gpu/shaders/material/gpu_shader_material_tex_sky.glsl
index 981d17b4283..b2f6147b629 100644
--- a/source/blender/gpu/shaders/material/gpu_shader_material_tex_sky.glsl
+++ b/source/blender/gpu/shaders/material/gpu_shader_material_tex_sky.glsl
@@ -1,4 +1,150 @@
-void node_tex_sky(vec3 co, out vec4 color)
+float sky_angle_between(float thetav, float phiv, float theta, float phi)
+{
+  float cospsi = sin(thetav) * sin(theta) * cos(phi - phiv) + cos(thetav) * cos(theta);
+
+  if (cospsi > 1.0) {
+    return 0.0;
+  }
+  if (cospsi < -1.0) {
+    return M_PI;
+  }
+
+  return acos(cospsi);
+}
+
+vec3 sky_spherical_coordinates(vec3 dir)
+{
+  return vec3(acos(dir.z), atan(dir.x, dir.y), 0);
+}
+
+/* Preetham */
+/* lam03+lam4: 5 floats passed as vec4+float */
+float sky_perez_function(vec4 lam03, float lam4, float theta, float gamma)
+{
+  float ctheta = cos(theta);
+  float cgamma = cos(gamma);
+
+  return (1.0 + lam03[0] * exp(lam03[1] / ctheta)) *
+         (1.0 + lam03[2] * exp(lam03[3] * gamma) + lam4 * cgamma * cgamma);
+}
+
+vec3 xyY_to_xyz(float x, float y, float Y)
+{
+  float X, Z;
+
+  if (y != 0.0) {
+    X = (x / y) * Y;
+  }
+  else {
+    X = 0.0;
+  }
+
+  if (y != 0.0 && Y != 0.0) {
+    Z = ((1.0 - x - y) / y) * Y;
+  }
+  else {
+    Z = 0.0;
+  }
+
+  return vec3(X, Y, Z);
+}
+
+void node_tex_sky_preetham(vec3 co,
+                           vec4 config_Y03,
+                           float config_Y4,
+                           vec4 config_x03,
+                           float config_x4,
+                           vec4 config_y03,
+                           float config_y4,
+                           vec2 sun_angles,
+                           vec3 radiance,
+                           vec3 xyz_to_r,
+                           vec3 xyz_to_g,
+                           vec3 xyz_to_b,
+                           out vec4 color)
+{
+  /* convert vector to spherical coordinates */
+  vec3 spherical = sky_spherical_coordinates(co);
+  float theta = spherical[0];
+  float phi = spherical[1];
+
+  float suntheta = sun_angles[0];
+  float sunphi = sun_angles[1];
+
+  /* angle between sun direction and dir */
+  float gamma = sky_angle_between(theta, phi, suntheta, sunphi);
+
+  /* clamp theta to horizon */
+  theta = min(theta, M_PI_2 - 0.001);
+
+  /* compute xyY color space values */
+  float Y = radiance[0] * sky_perez_function(config_Y03, config_Y4, theta, gamma);
+  float x = radiance[1] * sky_perez_function(config_x03, config_x4, theta, gamma);
+  float y = radiance[2] * sky_perez_function(config_y03, config_y4, theta, gamma);
+
+  /* convert to RGB */
+  vec3 xyz = xyY_to_xyz(x, y, Y);
+  color = vec4(dot(xyz_to_r, xyz), dot(xyz_to_g, xyz), dot(xyz_to_b, xyz), 1);
+}
+
+/* Hosek / Wilkie */
+float sky_radiance_hosekwilkie(
+    vec4 config03, vec4 config47, float config8, float theta, float gamma)
+{
+  float ctheta = cos(theta);
+  float cgamma = cos(gamma);
+
+  float expM = exp(config47[0] * gamma);
+  float rayM = cgamma * cgamma;
+  float mieM = (1.0 + rayM) / pow((1.0 + config8 * config8 - 2.0 * config8 * cgamma), 1.5);
+  float zenith = sqrt(ctheta);
+
+  return (1.0 + config03[0] * exp(config03[1] / (ctheta + 0.01))) *
+         (config03[2] + config03[3] * expM + config47[1] * rayM + config47[2] * mieM +
+          config47[3] * zenith);
+}
+
+void node_tex_sky_hosekwilkie(vec3 co,
+                              vec4 config_x03,
+                              vec4 config_x47,
+                              vec4 config_y03,
+                              vec4 config_y47,
+                              vec4 config_z03,
+                              vec4 config_z47,
+                              vec3 config_xyz8,
+                              vec2 sun_angles,
+                              vec3 radiance,
+                              vec3 xyz_to_r,
+                              vec3 xyz_to_g,
+                              vec3 xyz_to_b,
+                              out vec4 color)
+{
+  /* convert vector to spherical coordinates */
+  vec3 spherical = sky_spherical_coordinates(co);
+  float theta = spherical[0];
+  float phi = spherical[1];
+
+  float suntheta = sun_angles[0];
+  float sunphi = sun_angles[1];
+
+  /* angle between sun direction and dir */
+  float gamma = sky_angle_between(theta, phi, suntheta, sunphi);
+
+  /* clamp theta to horizon */
+  theta = min(theta, M_PI_2 - 0.001);
+
+  vec3 xyz;
+  xyz.x = sky_radiance_hosekwilkie(config_x03, config_x47, config_xyz8[0], theta, gamma) *
+          radiance.x;
+  xyz.y = sky_radiance_hosekwilkie(config_y03, config_y47, config_xyz8[1], theta, gamma) *
+          radiance.y;
+  xyz.z = sky_radiance_hosekwilkie(config_z03, config_z47, config_xyz8[2], theta, gamma) *
+          radiance.z;
+
+  color = vec4(dot(xyz_to_r, xyz), dot(xyz_to_g, xyz), dot(xyz_to_b, xyz), 1);
+}
+
+void node_tex_sky_nishita(vec3 co, out vec4 color)
 {
   color = vec4(1.0);
 }
diff --git a/source/blender/imbuf/IMB_colormanagement.h b/source/blender/imbuf/IMB_colormanagement.h
index 3158e3419b0..4530f6c9fc0 100644
--- a/source/blender/imbuf/IMB_colormanagement.h
+++ b/source/blender/imbuf/IMB_colormanagement.h
@@ -72,6 +72,7 @@ BLI_INLINE float IMB_colormanagement_get_luminance(const float rgb[3]);
 BLI_INLINE unsigned char IMB_colormanagement_get_luminance_byte(const unsigned char[3]);
 BLI_INLINE void IMB_colormangement_xyz_to_rgb(float rgb[3], const float xyz[3]);
 BLI_INLINE void IMB_colormangement_rgb_to_xyz(float xyz[3], const float rgb[3]);
+const float *IMB_colormangement_get_xyz_to_rgb(void);
 
 /* ** Color space transformation functions ** */
 void IMB_colormanagement_transform(float *buffer,
diff --git a/source/blender/imbuf/intern/colormanagement.c b/source/blender/imbuf/intern/colormanagement.c
index 3f5a0f25cc5..6341706bede 100644
--- a/source/blender/imbuf/intern/colormanagement.c
+++ b/source/blender/imbuf/intern/colormanagement.c
@@ -1456,6 +1456,11 @@ bool IMB_colormanagement_space_name_is_data(const char *name)
   return (colorspace && colorspace->is_data);
 }
 
+const float *IMB_colormangement_get_xyz_to_rgb()
+{
+  return &imbuf_xyz_to_rgb[0][0];
+}
+
 /*********************** Threaded display buffer transform routines *************************/
 
 typedef struct DisplayBufferThread {
diff --git a/source/blender/nodes/CMakeLists.txt b/source/blender/nodes/CMakeLists.txt
index 55c5009ea56..2c4ca087efa 100644
--- a/source/blender/nodes/CMakeLists.txt
+++ b/source/blender/nodes/CMakeLists.txt
@@ -38,6 +38,7 @@ set(INC
   ../render/extern/include
   ../../../intern/glew-mx
   ../../../intern/guardedalloc
+  ../../../intern/sky/include
 )
 
 set(INC_SYS
@@ -298,6 +299,7 @@ set(SRC
 
 set(LIB
   bf_functions
+  bf_intern_sky
 )
 
 if(WITH_PYTHON)
diff --git a/source/blender/nodes/shader/nodes/node_shader_tex_sky.c b/source/blender/nodes/shader/nodes/node_shader_tex_sky.c
index 0daa948c139..70f3f80eb4b 100644
--- a/source/blender/nodes/shader/nodes/node_shader_tex_sky.c
+++ b/source/blender/nodes/shader/nodes/node_shader_tex_sky.c
@@ -18,6 +18,7 @@
  */
 
 #include "../node_shader_util.h"
+#include "sky_model.h"
 
 /* **************** OUTPUT ******************** */
 
@@ -53,19 +54,164 @@ static void node_shader_init_tex_sky(bNodeTree *UNUSED(ntree), bNode *node)
   node->storage = tex;
 }
 
+typedef struct SkyModelPreetham {
+  float config_Y[5], config_x[5], config_y[5]; /* named after xyY color space */
+  float radiance[3];
+} SkyModelPreetham;
+
+typedef struct XYZ_to_RGB /* transposed imbuf_xyz_to_rgb, passed as 3x vec3 */
+{
+  float r[3], g[3], b[3];
+} XYZ_to_RGB;
+
+static float sky_perez_function(float *lam, float theta, float gamma)
+{
+  float ctheta = cosf(theta);
+  float cgamma = cosf(gamma);
+
+  return (1.0 + lam[0] * expf(lam[1] / ctheta)) *
+         (1.0 + lam[2] * expf(lam[3] * gamma) + lam[4] * cgamma * cgamma);
+}
+
+static void sky_precompute_old(SkyModelPreetham *sunsky, float sun_angles[], float turbidity)
+{
+  float theta = sun_angles[0];
+  float theta2 = theta * theta;
+  float theta3 = theta2 * theta;
+  float T = turbidity;
+  float T2 = T * T;
+  float chi = (4.0f / 9.0f - T / 120.0f) * (M_PI - 2.0f * theta);
+
+  sunsky->radiance[0] = (4.0453f * T - 4.9710f) * tanf(chi) - 0.2155f * T + 2.4192f;
+  sunsky->radiance[0] *= 0.06f;
+
+  sunsky->radiance[1] = (0.00166f * theta3 - 0.00375f * theta2 + 0.00209f * theta) * T2 +
+                        (-0.02903f * theta3 + 0.06377f * theta2 - 0.03202f * theta + 0.00394f) *
+                            T +
+                        (0.11693f * theta3 - 0.21196f * theta2 + 0.06052f * theta + 0.25886f);
+
+  sunsky->radiance[2] = (0.00275f * theta3 - 0.00610f * theta2 + 0.00317f * theta) * T2 +
+                        (-0.04214f * theta3 + 0.08970f * theta2 - 0.04153f * theta + 0.00516f) *
+                            T +
+                        (0.15346f * theta3 - 0.26756f * theta2 + 0.06670f * theta + 0.26688f);
+
+  sunsky->config_Y[0] = (0.1787f * T - 1.4630f);
+  sunsky->config_Y[1] = (-0.3554f * T + 0.4275f);
+  sunsky->config_Y[2] = (-0.0227f * T + 5.3251f);
+  sunsky->config_Y[3] = (0.1206f * T - 2.5771f);
+  sunsky->config_Y[4] = (-0.0670f * T + 0.3703f);
+
+  sunsky->config_x[0] = (-0.0193f * T - 0.2592f);
+  sunsky->config_x[1] = (-0.0665f * T + 0.0008f);
+  sunsky->config_x[2] = (-0.0004f * T + 0.2125f);
+  sunsky->config_x[3] = (-0.0641f * T - 0.8989f);
+  sunsky->config_x[4] = (-0.0033f * T + 0.0452f);
+
+  sunsky->config_y[0] = (-0.0167f * T - 0.2608f);
+  sunsky->config_y[1] = (-0.0950f * T + 0.0092f);
+  sunsky->config_y[2] = (-0.0079f * T + 0.2102f);
+  sunsky->config_y[3] = (-0.0441f * T - 1.6537f);
+  sunsky->config_y[4] = (-0.0109f * T + 0.0529f);
+
+  sunsky->radiance[0] /= sky_perez_function(sunsky->config_Y, 0, theta);
+  sunsky->radiance[1] /= sky_perez_function(sunsky->config_x, 0, theta);
+  sunsky->radiance[2] /= sky_perez_function(sunsky->config_y, 0, theta);
+}
+
+static void get_XYZ_to_RGB_for_gpu(XYZ_to_RGB *data)
+{
+  const float *xyz_to_rgb = IMB_colormangement_get_xyz_to_rgb();
+  data->r[0] = xyz_to_rgb[0];
+  data->r[1] = xyz_to_rgb[3];
+  data->r[2] = xyz_to_rgb[6];
+  data->g[0] = xyz_to_rgb[1];
+  data->g[1] = xyz_to_rgb[4];
+  data->g[2] = xyz_to_rgb[7];
+  data->b[0] = xyz_to_rgb[2];
+  data->b[1] = xyz_to_rgb[5];
+  data->b[2] = xyz_to_rgb[8];
+}
+
 static int node_shader_gpu_tex_sky(GPUMaterial *mat,
                                    bNode *node,
                                    bNodeExecData *UNUSED(execdata),
                                    GPUNodeStack *in,
                                    GPUNodeStack *out)
 {
-  if (!in[0].link) {
-    in[0].link = GPU_attribute(mat, CD_ORCO, "");
-  }
-
+  node_shader_gpu_default_tex_coord(mat, node, &in[0].link);
   node_shader_gpu_tex_mapping(mat, node, in, out);
+  NodeTexSky *tex = (NodeTexSky *)node->storage;
+  float sun_angles[2]; /* [0]=theta=zenith angle  [1]=phi=azimuth */
+  sun_angles[0] = acosf(tex->sun_direction[2]);
+  sun_angles[1] = atan2f(tex->sun_direction[0], tex->sun_direction[1]);
 
-  return GPU_stack_link(mat, node, "node_tex_sky", in, out);
+  if (tex->sky_model == 0) {
+    /* Preetham */
+    SkyModelPreetham sunsky;
+    sky_precompute_old(&sunsky, sun_angles, tex->turbidity);
+    XYZ_to_RGB xyz_to_rgb;
+    get_XYZ_to_RGB_for_gpu(&xyz_to_rgb);
+    return GPU_stack_link(mat,
+                          node,
+                          "node_tex_sky_preetham",
+                          in,
+                          out,
+                          /* Pass config_Y/x/y as 3x(vec4+float) */
+                          GPU_uniform(&sunsky.config_Y[0]),
+                          GPU_uniform(&sunsky.config_Y[4]),
+                          GPU_uniform(&sunsky.config_x[0]),
+                          GPU_uniform(&sunsky.config_x[4]),
+                          GPU_uniform(&sunsky.config_y[0]),
+                          GPU_uniform(&sunsky.config_y[4]),
+                          GPU_uniform(sun_angles),
+                          GPU_uniform(sunsky.radiance),
+                          GPU_uniform(xyz_to_rgb.r),
+                          GPU_uniform(xyz_to_rgb.g),
+                          GPU_uniform(xyz_to_rgb.b));
+  }
+  else if (tex->sky_model == 1) {
+    /* Hosek / Wilkie */
+    sun_angles[0] = fmin(M_PI_2, sun_angles[0]); /* clamp to horizon */
+    SKY_ArHosekSkyModelState *sky_state = SKY_arhosek_xyz_skymodelstate_alloc_init(
+        tex->turbidity, tex->ground_albedo, fmax(0.0, M_PI_2 - sun_angles[0]));
+    /* Pass sky_state->configs[3][9] as 3*(vec4+vec4)+vec3 */
+    float config_x07[8], config_y07[8], config_z07[8], config_xyz8[3];
+    for (int i = 0; i < 8; ++i) {
+      config_x07[i] = (float)sky_state->configs[0][i];
+      config_y07[i] = (float)sky_state->configs[1][i];
+      config_z07[i] = (float)sky_state->configs[2][i];
+    }
+    for (int i = 0; i < 3; ++i) {
+      config_xyz8[i] = (float)sky_state->configs[i][8];
+    }
+    float radiance[3];
+    for (int i = 0; i < 3; i++) {
+      radiance[i] = sky_state->radiances[i] * (2 * M_PI / 683);
+    }
+    SKY_arhosekskymodelstate_free(sky_state);
+    XYZ_to_RGB xyz_to_rgb;
+    get_XYZ_to_RGB_for_gpu(&xyz_to_rgb);
+    return GPU_stack_link(mat,
+                          node,
+                          "node_tex_sky_hosekwilkie",
+                          in,
+                          out,
+                          GPU_uniform(&config_x07[0]),
+                          GPU_uniform(&config_x07[4]),
+                          GPU_uniform(&config_y07[0]),
+                          GPU_uniform(&config_y07[4]),
+                          GPU_uniform(&config_z07[0]),
+                          GPU_uniform(&config_z07[4]),
+                          GPU_uniform(config_xyz8),
+                          GPU_uniform(sun_angles),
+                          GPU_uniform(radiance),
+                          GPU_uniform(xyz_to_rgb.r),
+                          GPU_uniform(xyz_to_rgb.g),
+                          GPU_uniform(xyz_to_rgb.b));
+  }
+  else {
+    return GPU_stack_link(mat, node, "node_tex_sky_nishita", in, out);
+  }
 }
 
 /* node type definition */