Cycles: improve sampling of ellipse area light with spread

**Problem**:
Area lights in Cycles have spread angle, in which case some part of the area light might be invisible to a shading point. The current implementation samples the whole area light, resulting some samples invisible and thus simply discarded. A technique is applied on rectangular light to sample a subset of the area light that is potentially visible (rB3f24cfb9582e1c826406301d37808df7ca6aa64c), however, ellipse (including disk) area lights remained untreated. The purpose of this patch is to apply a techniques to ellipse area light.
**Related Task**:
T87053
**Results**:
These are renderings before and after the patch:
|16spp|Disk light|Ellipse light|Square light (for reference, no changes)
|Before|{F13996789}|{F13996788}|{F13996822}
|After|{F13996759}|{F13996787}|{F13996852}
**Explanation**:
The visible region on an area light is found by drawing a cone from the shading point to the plane where the area light lies, with the aperture of the cone being the light spread.
{F13990078,height=200}
Ideally, we would like to draw samples only from the intersection of the area light and the projection of the cone onto the plane (forming a circle). However, the shape of the intersection is often irregular and thus hard to sample from directly.
{F13990104,height=200}
Instead, the current implementation draws samples from the bounding rectangle of the intersection. In this case, we still end up with some invalid samples outside of the circle, but already much less than sampling the original area light, and the bounding rectangle is easy to sample from.
{F13990125}
The above technique is only applied to rectangle area lights, ellipse area light still suffers from poor sampling. We could apply a similar technique to ellipse area lights, that is, find the
smallest regular shape (rectangle, circle, or ellipse) that covers the intersection (or maybe not the smallest but easy to compute).
For disk area light, we consider the relative position of both circles. Denoting `dist` as the distance between the centre of two circles, and `r1`, `r2` their radii. If `dist > r1 + r2`, the area light is completely invisible, we directly return `false`. If `dist < abs(r1 - r2)`, the smaller circle lies inside the larger one, and we sample whichever circle is smaller. Otherwise, the two circles intersect, we compute the bounding rectangle of the intersection, in which case `axis_u`, `len_u`, `axis_v`, `len_v` needs to be computed anew. Depending on the distance between the two circles, `len_v` is either the diameter of the smaller circle or the length of the common chord.
|{F13990211,height=195}|{F13990225,height=195}|{F13990274,height=195}|{F13990210,height=195}
|`dist > r1 + r2`|`dist < abs(r1 - r2)`|`dist^2 < abs(r1^2 - r2^2)`|`dist^2 > abs(r1^2 - r2^2)`
For ellipse area light, it's hard to find the smallest bounding shape of the intersection, therefore, we compute the bounding rectangle of the ellipse itself, then treat it as a rectangle light.
|{F13990386,height=195}|{F13990385,height=195}|{F13990387,height=195}
We also check the areas of the bounding rectangle of the intersection, the ellipse (disk) light, and the spread circle, then draw samples from the smallest shape of the three. For ellipse light, this also detects where one shape lies inside the other. I am not sure if we should add this measure to rectangle area light and sample from the spread circle when it has smaller area, as we seem to have a better sampling technique for rectangular (uniformly sample the solid angle). Maybe we could add [area-preserving parameterization for spherical
ellipse](https://arxiv.org/pdf/1805.09048.pdf) in the future.
**Limitation**:
At some point we switch from sampling the ellipse to sampling the rectangle, depending on the area of the both, and there seems to be a visible line (with |slope| =1) on the final rendering
which demonstrate at which point we switch between the two methods. We could see that the new sampling method clearly has lower variance near the boundaries, but close to that visible line,
the rectangle sampling method seems to have larger variance. I could not spot any bug in the implementation, and I am not sure if this happens because different sampling patterns for ellipse and rectangle are used.
|Before (256spp)|After (256spp)
|{F13996995}|{F13996998}

Differential Revision: https://developer.blender.org/D16694
This commit is contained in:
Weizhen Huang 2022-12-05 14:35:25 +01:00
parent 96b6ea9ee6
commit 53ef52f165
5 changed files with 180 additions and 86 deletions

View File

@ -75,19 +75,19 @@ void BlenderSync::sync_light(BL::Object &b_parent,
switch (b_area_light.shape()) {
case BL::AreaLight::shape_SQUARE:
light->set_sizev(light->get_sizeu());
light->set_round(false);
light->set_ellipse(false);
break;
case BL::AreaLight::shape_RECTANGLE:
light->set_sizev(b_area_light.size_y());
light->set_round(false);
light->set_ellipse(false);
break;
case BL::AreaLight::shape_DISK:
light->set_sizev(light->get_sizeu());
light->set_round(true);
light->set_ellipse(true);
break;
case BL::AreaLight::shape_ELLIPSE:
light->set_sizev(b_area_light.size_y());
light->set_round(true);
light->set_ellipse(true);
break;
}
light->set_light_type(LIGHT_AREA);

View File

@ -102,49 +102,121 @@ ccl_device float area_light_spread_attenuation(const float3 D,
return max((1.0f - (cot_half_spread * tan_a)) * normalize_spread, 0.0f);
}
/* Compute subset of area light that actually has an influence on the shading point, to
* reduce noise with low spread. */
ccl_device bool area_light_spread_clamp_area_light(const float3 P,
const float3 lightNg,
ccl_private float3 *lightP,
const float3 axis_u,
ccl_private float *len_u,
const float3 axis_v,
ccl_private float *len_v,
const float cot_half_spread)
/* Compute the minimal rectangle, circle or ellipse that covers the valid sample region, to reduce
* noise with low spread. */
ccl_device bool area_light_spread_clamp_light(const float3 P,
const float3 lightNg,
ccl_private float3 *lightP,
ccl_private float3 *axis_u,
ccl_private float *len_u,
ccl_private float3 *axis_v,
ccl_private float *len_v,
const float cot_half_spread,
ccl_private bool *sample_rectangle)
{
/* Closest point in area light plane and distance to that plane. */
const float3 closest_P = P - dot(lightNg, P - *lightP) * lightNg;
const float t = len(closest_P - P);
/* Radius of circle on area light that actually affects the shading point. */
const float radius = t / cot_half_spread;
const float r_spread = t / cot_half_spread;
/* Local uv coordinates of closest point. */
const float closest_u = dot(axis_u, closest_P - *lightP);
const float closest_v = dot(axis_v, closest_P - *lightP);
const float spread_u = dot(*axis_u, closest_P - *lightP);
const float spread_v = dot(*axis_v, closest_P - *lightP);
/* Compute rectangle encompassing the circle that affects the shading point,
* clamped to the bounds of the area light. */
const float min_u = max(closest_u - radius, -*len_u * 0.5f);
const float max_u = min(closest_u + radius, *len_u * 0.5f);
const float min_v = max(closest_v - radius, -*len_v * 0.5f);
const float max_v = min(closest_v + radius, *len_v * 0.5f);
const bool is_round = !(*sample_rectangle) && (*len_u == *len_v);
/* Skip if rectangle is empty. */
if (min_u >= max_u || min_v >= max_v) {
return false;
/* Whether we should sample the spread circle. */
bool sample_spread;
if (is_round) {
/* Distance between the centers of the disk light and the valid region circle. */
const float dist = len(make_float2(spread_u, spread_v));
/* Radius of the disk light. */
const float r = *len_u * 0.5f;
if (dist >= r + r_spread) {
/* Two circles are outside each other or touch externally. */
return false;
}
sample_spread = (dist <= fabsf(r - r_spread)) && (r_spread < r);
if (dist > fabsf(r - r_spread)) {
/* Two circles intersect. Find the smallest rectangle that covers the intersection */
const float len_u_ = r + r_spread - dist;
const float len_v_ = (fabsf(sqr(r) - sqr(r_spread)) >= sqr(dist)) ?
2.0f * fminf(r, r_spread) :
sqrtf(sqr(2.0f * r_spread) -
sqr(dist + (sqr(r_spread) - sqr(r)) / dist));
const float rect_area = len_u_ * len_v_;
const float circle_area = M_PI_F * sqr(r);
const float spread_area = M_PI_F * sqr(r_spread);
/* Sample the shape with minimal area. */
if (rect_area < fminf(circle_area, spread_area)) {
*sample_rectangle = true;
*axis_u = normalize(*lightP - closest_P);
*axis_v = rotate_around_axis(*axis_u, lightNg, M_PI_2_F);
*len_u = len_u_;
*len_v = len_v_;
*lightP = 0.5f * (*lightP + closest_P + *axis_u * (r_spread - r));
return true;
}
sample_spread = (spread_area < circle_area);
}
}
else {
/* Compute rectangle encompassing the circle that affects the shading point,
* clamped to the bounds of the area light. */
const float min_u = max(spread_u - r_spread, -*len_u * 0.5f);
const float max_u = min(spread_u + r_spread, *len_u * 0.5f);
const float min_v = max(spread_v - r_spread, -*len_v * 0.5f);
const float max_v = min(spread_v + r_spread, *len_v * 0.5f);
/* Skip if rectangle is empty. */
if (min_u >= max_u || min_v >= max_v) {
return false;
}
const float rect_len_u = max_u - min_u;
const float rect_len_v = max_v - min_v;
const float rect_area = rect_len_u * rect_len_v;
const float ellipse_area = (*sample_rectangle) ? FLT_MAX : M_PI_4_F * (*len_u) * (*len_v);
const float spread_area = M_PI_F * sqr(r_spread);
/* Sample the shape with minimal area. */
/* NOTE: we don't switch to spread circle sampling for rectangle light because rectangle light
* supports solid angle sampling, which has less variance than sampling the area. If ellipse
* area light also supports solid angle sampling, `*sample_rectangle ||` could be deleted. */
if (*sample_rectangle || rect_area < fminf(ellipse_area, spread_area)) {
*sample_rectangle = true;
/* Compute new area light center position and axes from rectangle in local
* uv coordinates. */
const float new_center_u = 0.5f * (min_u + max_u);
const float new_center_v = 0.5f * (min_v + max_v);
*len_u = rect_len_u;
*len_v = rect_len_v;
*lightP = *lightP + *axis_u * new_center_u + *axis_v * new_center_v;
return true;
}
*sample_rectangle = false;
sample_spread = (spread_area < ellipse_area);
}
/* Compute new area light center position and axes from rectangle in local
* uv coordinates. */
const float new_center_u = 0.5f * (min_u + max_u);
const float new_center_v = 0.5f * (min_v + max_v);
*len_u = max_u - min_u;
*len_v = max_v - min_v;
*lightP = *lightP + new_center_u * axis_u + new_center_v * axis_v;
if (sample_spread) {
*lightP = *lightP + *axis_u * spread_u + *axis_v * spread_v;
*len_u = r_spread * 2.0f;
*len_v = r_spread * 2.0f;
return true;
}
/* Don't clamp. */
return true;
}
@ -159,13 +231,7 @@ ccl_device_inline bool area_light_sample(const ccl_global KernelLight *klight,
{
ls->P = klight->co;
const float3 axis_u = klight->area.axis_u;
const float3 axis_v = klight->area.axis_v;
const float len_u = klight->area.len_u;
const float len_v = klight->area.len_v;
float3 Ng = klight->area.dir;
float invarea = fabsf(klight->area.invarea);
bool is_round = (klight->area.invarea < 0.0f);
if (!in_volume_segment) {
if (dot(ls->P - P, Ng) > 0.0f) {
@ -173,40 +239,65 @@ ccl_device_inline bool area_light_sample(const ccl_global KernelLight *klight,
}
}
const float3 axis_u = klight->area.axis_u;
const float3 axis_v = klight->area.axis_v;
const float len_u = klight->area.len_u;
const float len_v = klight->area.len_v;
float invarea = fabsf(klight->area.invarea);
bool is_ellipse = (klight->area.invarea < 0.0f);
bool sample_rectangle = !is_ellipse;
float3 inplane;
if (is_round || in_volume_segment) {
if (in_volume_segment) {
/* FIXME: handle rectangular light. */
inplane = ellipse_sample(axis_u * len_u * 0.5f, axis_v * len_v * 0.5f, randu, randv);
ls->P += inplane;
ls->pdf = invarea;
}
else {
inplane = ls->P;
float3 old_P = ls->P;
float3 sample_axis_u = axis_u;
float3 sample_axis_v = axis_v;
float sample_len_u = len_u;
float sample_len_v = len_v;
if (!in_volume_segment && klight->area.cot_half_spread > 0.0f) {
if (!area_light_spread_clamp_area_light(P,
Ng,
&ls->P,
axis_u,
&sample_len_u,
axis_v,
&sample_len_v,
klight->area.cot_half_spread)) {
if (klight->area.cot_half_spread > 0.0f) {
if (!area_light_spread_clamp_light(P,
Ng,
&ls->P,
&sample_axis_u,
&sample_len_u,
&sample_axis_v,
&sample_len_v,
klight->area.cot_half_spread,
&sample_rectangle)) {
return false;
}
}
ls->pdf = area_light_rect_sample(
P, &ls->P, axis_u, sample_len_u, axis_v, sample_len_v, randu, randv, true);
inplane = ls->P - inplane;
if (sample_rectangle) {
ls->pdf = area_light_rect_sample(
P, &ls->P, sample_axis_u, sample_len_u, sample_axis_v, sample_len_v, randu, randv, true);
}
else {
ls->P += ellipse_sample(
sample_axis_u * sample_len_u * 0.5f, sample_axis_v * sample_len_v * 0.5f, randu, randv);
ls->pdf = 4.0f * M_1_PI_F / (sample_len_u * sample_len_v);
}
inplane = ls->P - old_P;
}
const float light_u = dot(inplane, axis_u) / len_u;
const float light_v = dot(inplane, axis_v) / len_v;
/* Sampled point lies outside of the area light. */
if (is_ellipse && (sqr(light_u) + sqr(light_v) > 0.25f)) {
return false;
}
if (!is_ellipse && (fabsf(light_u) > 0.5f || fabsf(light_v) > 0.5f)) {
return false;
}
/* NOTE: Return barycentric coordinates in the same notation as Embree and OptiX. */
ls->u = light_v + 0.5f;
ls->v = -light_u - light_v;
@ -222,7 +313,7 @@ ccl_device_inline bool area_light_sample(const ccl_global KernelLight *klight,
ls->D, ls->Ng, klight->area.cot_half_spread, klight->area.normalize_spread);
}
if (is_round) {
if (!sample_rectangle) {
ls->pdf *= lamp_light_pdf(Ng, -ls->D, ls->t);
}
@ -252,7 +343,7 @@ ccl_device_inline bool area_light_intersect(const ccl_global KernelLight *klight
{
/* Area light. */
const float invarea = fabsf(klight->area.invarea);
const bool is_round = (klight->area.invarea < 0.0f);
const bool is_ellipse = (klight->area.invarea < 0.0f);
if (invarea == 0.0f) {
return false;
}
@ -281,7 +372,7 @@ ccl_device_inline bool area_light_intersect(const ccl_global KernelLight *klight
t,
u,
v,
is_round);
is_ellipse);
}
ccl_device_inline bool area_light_sample_from_intersection(
@ -303,32 +394,35 @@ ccl_device_inline bool area_light_sample_from_intersection(
ls->D = ray_D;
ls->Ng = Ng;
const bool is_round = (klight->area.invarea < 0.0f);
if (is_round) {
ls->pdf = invarea * lamp_light_pdf(Ng, -ray_D, ls->t);
float3 sample_axis_u = klight->area.axis_u;
float3 sample_axis_v = klight->area.axis_v;
float sample_len_u = klight->area.len_u;
float sample_len_v = klight->area.len_v;
bool is_ellipse = (klight->area.invarea < 0.0f);
bool sample_rectangle = !is_ellipse;
if (klight->area.cot_half_spread > 0.0f) {
if (!area_light_spread_clamp_light(ray_P,
Ng,
&light_P,
&sample_axis_u,
&sample_len_u,
&sample_axis_v,
&sample_len_v,
klight->area.cot_half_spread,
&sample_rectangle)) {
return false;
}
}
if (sample_rectangle) {
ls->pdf = area_light_rect_sample(
ray_P, &light_P, sample_axis_u, sample_len_u, sample_axis_v, sample_len_v, 0, 0, false);
}
else {
const float3 axis_u = klight->area.axis_u;
const float3 axis_v = klight->area.axis_v;
float sample_len_u = klight->area.len_u;
float sample_len_v = klight->area.len_v;
if (klight->area.cot_half_spread > 0.0f) {
if (!area_light_spread_clamp_area_light(ray_P,
Ng,
&light_P,
axis_u,
&sample_len_u,
axis_v,
&sample_len_v,
klight->area.cot_half_spread)) {
return false;
}
}
ls->pdf = area_light_rect_sample(
ray_P, &light_P, axis_u, sample_len_u, axis_v, sample_len_v, 0, 0, false);
ls->pdf = 4.0f * M_1_PI_F / (sample_len_u * sample_len_v) * lamp_light_pdf(Ng, -ray_D, ls->t);
}
ls->eval_fac = 0.25f * invarea;
if (klight->area.cot_half_spread > 0.0f) {

View File

@ -106,7 +106,7 @@ NODE_DEFINE(Light)
SOCKET_FLOAT(sizeu, "Size U", 1.0f);
SOCKET_VECTOR(axisv, "Axis V", zero_float3());
SOCKET_FLOAT(sizev, "Size V", 1.0f);
SOCKET_BOOLEAN(round, "Round", false);
SOCKET_BOOLEAN(ellipse, "Ellipse", false);
SOCKET_FLOAT(spread, "Spread", M_PI_F);
SOCKET_INT(map_resolution, "Map Resolution", 0);
@ -911,7 +911,7 @@ void LightManager::device_update_lights(Device *device, DeviceScene *dscene, Sce
float3 axis_u = normalize_len(extentu, &len_u);
float3 axis_v = normalize_len(extentv, &len_v);
float area = len_u * len_v;
if (light->round) {
if (light->ellipse) {
area *= -M_PI_4_F;
}
float invarea = (area != 0.0f) ? 1.0f / area : 1.0f;
@ -1032,7 +1032,7 @@ void LightManager::device_update_lights(Device *device, DeviceScene *dscene, Sce
float3 axis_u = normalize_len(extentu, &len_u);
float3 axis_v = normalize_len(extentv, &len_v);
float area = len_u * len_v;
if (light->round) {
if (light->ellipse) {
area *= -M_PI_4_F;
}
float invarea = (area != 0.0f) ? 1.0f / area : 1.0f;

View File

@ -45,7 +45,7 @@ class Light : public Node {
NODE_SOCKET_API(float, sizeu)
NODE_SOCKET_API(float3, axisv)
NODE_SOCKET_API(float, sizev)
NODE_SOCKET_API(bool, round)
NODE_SOCKET_API(bool, ellipse)
NODE_SOCKET_API(float, spread)
NODE_SOCKET_API(Transform, tfm)

@ -1 +1 @@
Subproject commit 0b0052bd53ad8249ed07dfb87705c338af698bde
Subproject commit fdfd24de034d4bba4fb67731d0aae81dc4940239