Cycles: more efficient and better-behaved sampling of spherical triangles

uses less inverse trigonometric functions and normalizations.
Also `q` is now guaranteed to be smaller than or equal to 1.

Pull Request: https://projects.blender.org/blender/blender/pulls/108380
This commit is contained in:
Weizhen Huang
2023-06-02 16:57:50 +02:00
committed by Weizhen Huang
parent c8bd998944
commit 97d9bbbc97

View File

@@ -72,18 +72,13 @@ ccl_device_forceinline float triangle_light_pdf(KernelGlobals kg,
/* sd contains the point on the light source /* sd contains the point on the light source
* calculate Px, the point that we're shading */ * calculate Px, the point that we're shading */
const float3 Px = sd->P + sd->wi * t; const float3 Px = sd->P + sd->wi * t;
const float3 v0_p = V[0] - Px;
const float3 v1_p = V[1] - Px;
const float3 v2_p = V[2] - Px;
const float3 u01 = safe_normalize(cross(v0_p, v1_p)); const float3 A = safe_normalize(V[0] - Px);
const float3 u02 = safe_normalize(cross(v0_p, v2_p)); const float3 B = safe_normalize(V[1] - Px);
const float3 u12 = safe_normalize(cross(v1_p, v2_p)); const float3 C = safe_normalize(V[2] - Px);
const float alpha = fast_acosf(dot(u02, u01)); const float solid_angle = 2.0f * fast_atanf(fabsf(dot(A, cross(B, C))) /
const float beta = fast_acosf(-dot(u01, u12)); (1.0f + dot(B, C) + dot(A, C) + dot(A, B)));
const float gamma = fast_acosf(dot(u02, u12));
const float solid_angle = alpha + beta + gamma - M_PI_F;
/* distribution_pdf_triangles is calculated over triangle area, but we're not sampling over /* distribution_pdf_triangles is calculated over triangle area, but we're not sampling over
* its area */ * its area */
@@ -160,59 +155,42 @@ ccl_device_forceinline bool triangle_light_sample(KernelGlobals kg,
float distance_to_plane = fabsf(dot(N0, V[0] - P) / dot(N0, N0)); float distance_to_plane = fabsf(dot(N0, V[0] - P) / dot(N0, N0));
if (!in_volume_segment && (longest_edge_squared > distance_to_plane * distance_to_plane)) { if (!in_volume_segment && (longest_edge_squared > distance_to_plane * distance_to_plane)) {
/* see James Arvo, "Stratified Sampling of Spherical Triangles" /* A modified version of James Arvo, "Stratified Sampling of Spherical Triangles"
* http://www.graphics.cornell.edu/pubs/1995/Arv95c.pdf */ * http://www.graphics.cornell.edu/pubs/1995/Arv95c.pdf */
/* project the triangle to the unit sphere /* Project the triangle to the unit sphere and calculate the three unit vector that spans the
* and calculate its edges and angles */ * spherical triangle. */
const float3 v0_p = V[0] - P; const float3 A = safe_normalize(V[0] - P);
const float3 v1_p = V[1] - P; const float3 B = safe_normalize(V[1] - P);
const float3 v2_p = V[2] - P; const float3 C = safe_normalize(V[2] - P);
const float3 u01 = safe_normalize(cross(v0_p, v1_p)); const float cos_a = dot(B, C);
const float3 u02 = safe_normalize(cross(v0_p, v2_p)); const float cos_b = dot(A, C);
const float3 u12 = safe_normalize(cross(v1_p, v2_p));
const float3 A = safe_normalize(v0_p);
const float3 B = safe_normalize(v1_p);
const float3 C = safe_normalize(v2_p);
const float cos_alpha = dot(u02, u01);
const float cos_beta = -dot(u01, u12);
const float cos_gamma = dot(u02, u12);
/* calculate dihedral angles */
const float alpha = fast_acosf(cos_alpha);
const float beta = fast_acosf(cos_beta);
const float gamma = fast_acosf(cos_gamma);
/* the area of the unit spherical triangle = solid angle */
const float solid_angle = alpha + beta + gamma - M_PI_F;
/* precompute a few things
* these could be re-used to take several samples
* as they are independent of `rand` */
const float cos_c = dot(A, B); const float cos_c = dot(A, B);
const float sin_alpha = fast_sinf(alpha); const float sin_b_sin_c_2 = (1.0f - sqr(cos_b)) * (1.0f - sqr(cos_c));
const float product = sin_alpha * cos_c;
/* Select a random sub-area of the spherical triangle const float mixed_product = fabsf(dot(A, cross(B, C)));
* and calculate the third vertex C_ of that new triangle */
const float phi = rand.x * solid_angle - alpha;
float s, t;
fast_sincosf(phi, &s, &t);
const float u = t - cos_alpha;
const float v = s + product;
const float3 U = safe_normalize(C - dot(C, A) * A); /* The area of the spherical triangle is equal to the subtended solid angle. */
const float solid_angle = 2.0f * fast_atanf(mixed_product / (1.0f + cos_a + cos_b + cos_c));
float q = 1.0f; /* Select a random sub-area of the spherical triangle and calculate the third vertex C_ of that
const float det = ((v * s + u * t) * sin_alpha); * new triangle. */
if (det != 0.0f) { const float A_hat = rand.x * solid_angle;
q = ((v * t - u * s) * cos_alpha - v) / det; float sin_A_hat, cos_A_hat;
} fast_sincosf(A_hat, &sin_A_hat, &cos_A_hat);
const float temp = max(1.0f - q * q, 0.0f);
const float3 C_ = safe_normalize(q * A + sqrtf(temp) * U); /* These values lack a `sin_b * sin_c` factor, will divide when computing `temp`. */
const float cos_alpha = cos_a - cos_b * cos_c;
const float sin_alpha = mixed_product;
const float t = cos_A_hat * cos_alpha + sin_A_hat * sin_alpha;
const float temp = (cos_c - 1.0f) * t * cos_alpha / sin_b_sin_c_2;
const float q = (cos_A_hat - cos_c + temp) / (1.0f - cos_A_hat * cos_c + temp);
const float3 U = safe_normalize(C - cos_b * A);
const float3 C_ = safe_normalize(q * A + sin_from_cos(q) * U);
/* Finally, select a random point along the edge of the new triangle /* Finally, select a random point along the edge of the new triangle
* That point on the spherical triangle is the sampled ray direction */ * That point on the spherical triangle is the sampled ray direction */