Odd behavior for sphere intersections

[memamsaleh]

Summary

I get self intersections when I sample outgoing rays from the BSDF off the surface of a small sphere. In addition I noticed some odd behavior regarding the first bounce as well.

System configuration

System information:

OS: Windows-10 CPU: AMD64 Family 25 Model 33 Stepping 2, AuthenticAMD GPU: NVIDIA GeForce RTX 4090 Python: 3.11.5 | packaged by Anaconda, Inc. | (main, Sep 11 2023, 13:26:23) [MSC v.1916 64 bit (AMD64)] NVidia driver: 552.44 CUDA: 12.2.91 LLVM: -1.-1.-1

Dr.Jit: 0.4.4 Mitsuba: 3.5.0 Is custom build? False Compiled with: MSVC 19.34.31937.0 Variants: scalar_rgb cuda_ad_mono cuda_mono

Description

I tested with a small scene to reproduce the issue, the scene only contains one sphere however when I sample from the BSDF off the surface of the sphere I still get a few intersections.

one_sphere_scene = mi.load_dict({
    'type' : 'scene',
    'sphere1': {
        'type': 'sphere',
        'center': [0.0, 0.083, 0.4],
        'radius': 0.015,
        'bsdf': {
            'type': 'roughconductor',
            'alpha': 0.2,
            'distribution': 'ggx'
        }
    }
})

Here are the relevant parts to compute the intersections.

def intersect_scene(rays, scene):
    si = scene.ray_intersect(rays)

    n_intersections = dr.count(si.is_valid())
    if config.print_stats:
        print("Ray intersections: %3d" % n_intersections)
    return si, n_intersections

def get_bounce(si, tp, scene):
    bsdf_ctx = mi.BSDFContext()

    indices = dr.compress(si.is_valid())
    bsdf = si.bsdf()

    sp, tp_pdf = bsdf.sample(bsdf_ctx, si, sampler.next_1d(), sampler.next_2d())
    wo = si.to_world(sp.wo)

    wo = dr.gather(mi.Vector3f, wo, indices)
    tp_pdf = dr.gather(mi.Color1f, tp_pdf, indices)

    new_tp = dr.gather(mi.Float, tp, indices) * tp_pdf
    pts = dr.gather(mi.Vector3f, si.p, indices)
    rays = mi.Ray3f(o=pts, d=wo)

    new_si, n_intersections = intersect_scene(rays, scene)
    if n_intersections == 0:
        return None, tp
    return new_si, new_tp

I took 50 samples off the surface of the sphere to compute the incoming rays for the initial bounce leading to 50 intersections

However I get intersections when computing the first bounce off the surface, the number depends on the incoming directions and the setup. Here is a visualization of 27 rays leading to intersections for this particular example (shown in blue):

Additional tests:

I thought it might be related false intersections on the inside of the sphere, so I checked the dot product for the incoming ray with the normal at the point of intersections dr.count(dr.dot(si.n, -si.wi)>=0.0) and got only 35 when visualizing the 15 problematic rays they seem to all lie on one side of the sphere which I find odd.

I also noticed that the indices of those rays are always a subset of the ones leading to false intersections. However, when checking the normals with an estimate they seem to be correctly pointing outwards so it might not be the case.

n_est = dr.normalize(si.p - sphere_pos)
print(dr.count(dr.dot(si.n, n_est) > 0.99))

[njroussel]

Hi @memamsaleh

I appreciate the detailed investigation.

Self-intersections are very common, the only way to avoid the is to slightly offset your ray from the surface of the previous bounce. This is done by the spawn_ray method, take a look at this tutorial which briefly goes over simply rendering operations like ray intersections. You'll find the details of the offset in the source code.

[memamsaleh]

Hi @njroussel

Thank you, it did help. I still get exactly one intersection for some reason but it is easy to filter out for now as the outgoing direction seems to be facing away from the normal in this particular case.

[njroussel]

Hmm, that's still bizarre and shouldn't happen. Can you double-check that prior to the world coordinate transform, the outgoing direction is still in top upper hemisphere (i.e. sp.wo.z > 0) ?

[memamsaleh]

I checked and prior to the world coordinate transform the z coordinate is always negative for the specific index causing the issue when it occurs.

[njroussel]

Something is wrong.

Can you share how you set the incident direction si.wi ? Can you also check that prior to calling bsdf.sample its z coordinate is positive ?

[memamsaleh]

For setting the incoming direction in the current example I was sampling points on the sphere surface and computing the direction from the eye, however I tested now and I still get the same problems when using different sampling methods like sampling from a cone around the eye.

def sample_target_shape(rng, scene, idx=0):
    sample_1, sample_2 = rng.next_float32(), rng.next_float32()
    sample = scene.shapes()[idx].sample_position(0, [sample_1, sample_2])
    return sample.p

def sample_cone(rng):
    sample_1, sample_2 = rng.next_float32(), rng.next_float32()
    return mi.warp.square_to_uniform_cone([sample_1, sample_2], 0.75)

def sample_dir(rng, eye, scene, mode="shape"):
        fixed = (mode == "shape")
        if fixed:
            pos = md[mode](rng, scene)
            dir = dr.normalize(pos-eye)
        else:
            dir = md[mode](rng)
        return dir

The rays are then set using the eye as the origin and the resulting directions

origins = mi.Vector3f(np.repeat(Ti.x, N_Paths), np.repeat(Ti.y, N_Paths), np.repeat(Ti.z, N_Paths))
ray_d = smp.sample_dir(rng, Ti, scene, mode="cone")
rays = mi.Ray3f(o=origins, d=ray_d)
si, n_intersections = intersect_scene(rays, scene)

As for the z coordinate of the incident ray prior to the bsdf.sample it is always positive.

One additional note on the sampling I use rng = mi.PCG32(size=N_Paths) for the initial sampling shown above but for the bsdf sampling I use a different independent sampler

sampler = mi.load_dict({
    'type': 'independent'
})

sampler.seed(0xBADCAFE, N_Paths)

The issue still holds if I use the same random number generator (rng) for both.

I reseed both at the beginning of the whole simulation to help with debugging and the issue usually occur with 1-5 rays as I move the origin for incident rays through the simulation but I couldn't detect a pattern for the incoming directions that causes the issue.

I also quickly tested with a large rectangle object using the same bsdf and I noticed I do still get a few cases 1-5 where the local outgoing direction lies in the bottom hemisphere but somehow they don't lead to self intersections which I also find odd (as far as I understand spawn_ray seems to shift the point of the intersection along the normal which should lead to intersections in this case).

[njroussel]

I don't see what's wrong with your code, if you could give us a reproducer it would greatly help.

Here's something I quickly put together and that works as expected:

import mitsuba as mi
import drjit as dr

mi.set_variant('cuda_ad_rgb')

N = 10_000_000

# Empty scene with a small off-centered sphere
scene = {
    'type': 'scene',
    'sphere': {
        'type': 'sphere',
        'radius': 0.01, # small
        'center': mi.ScalarPoint3f(0.01, 0.02, -0.01), # not exactly centered
        'bsdf': {
            'type': 'diffuse'
        }
    }
}
scene = mi.load_dict(scene)

# Sample positions on the sphere
sphere = scene.shapes()[0]
sampler = mi.load_dict({'type': 'independent'})
sampler.seed(0xBADCAFE, N)
ps = sphere.sample_position(0, sampler.next_1d())

# Build rays with origin (1, 1, 1) pointing to samples positions on the sphere
ray_o = mi.Point3f(1)
ray_d = dr.normalize(ps.p - ray_o)
ray = mi.Ray3f(o=ray_o, d=ray_d)

# Intersect the scene, some might miss the sphere
si = scene.ray_intersect(ray)
print(f"Valid primary rays: {dr.count(si.is_valid())}") # Should be close to N

# Sample the BSDF
ctx = mi.BSDFContext()
bs, _ = si.bsdf().sample(ctx, si, sampler.next_1d(), sampler.next_2d())
wo_world = si.to_world(bs.wo)

# Use the sampled outgoing direction to intersect the scene again
ray = si.spawn_ray(wo_world)
si = scene.ray_intersect(ray, si.is_valid())
print(f"Valid secondary rays: {dr.count(si.is_valid())}") # Should be 0

[memamsaleh]

It seems the issue is in the bsdf, I just changed the bsdf in your code to a roughconductor and already got lots of valid secondary rays.

import mitsuba as mi
import drjit as dr

mi.set_variant('cuda_mono')

N = 10_000_000

# Empty scene with a small off-centered sphere
scene = {
    'type': 'scene',
    'sphere': {
        'type': 'sphere',
        'radius': 0.01, # small
        'center': mi.ScalarPoint3f(0.01, 0.02, -0.01), # not exactly centered
        'bsdf': {
            'type': 'roughconductor',
            'alpha': 0.2,
            'distribution': 'ggx'
        }
    }
}
scene = mi.load_dict(scene)

# Sample positions on the sphere
sphere = scene.shapes()[0]
sampler = mi.load_dict({'type': 'independent'})
sampler.seed(0xBADCAFE, N)
ps = sphere.sample_position(0, sampler.next_1d())

# Build rays with origin (1, 1, 1) pointing to samples positions on the sphere
ray_o = mi.Point3f(1)
ray_d = dr.normalize(ps.p - ray_o)
ray = mi.Ray3f(o=ray_o, d=ray_d)

# Intersect the scene, some might miss the sphere
si = scene.ray_intersect(ray)
print(f"Valid primary rays: {dr.count(si.is_valid())}") # Should be close to N

# Sample the BSDF
ctx = mi.BSDFContext()
bs, _ = si.bsdf().sample(ctx, si, sampler.next_1d(), sampler.next_2d())
wo_world = si.to_world(bs.wo)

# Use the sampled outgoing direction to intersect the scene again
ray = si.spawn_ray(wo_world)
si = scene.ray_intersect(ray, si.is_valid())
print(f"Valid secondary rays: {dr.count(si.is_valid())}") # Should be 0

Edit: It seems also to be tied with the roughness, the rougher the material (larger alpha) the more self intersections we get. When using a conductor material instead of rough conductor there are still some self intersections but they are much less in comparison.

[njroussel]

Oh my bad, I thought that you were exclusively using the diffuse BSDF. This is expected behavior! The sample() method does not fail, if you're feeding it wrong/impossible configurations it will still produce some BSDFSample3f object, however from the backside (the weight it returns should be 0. For example, if you hit the backside of your surface the ray should terminate. Any BSDF will still return a valid bs.wo direciton but the weight associated to that interaction will be exactly 0.

Why does this matter here? Well, first of all, the ray intersection accelerators we use are not perfect. Sometimes, they will return an intersection on the backside (i.e si.wi.z < 0). See here:

import mitsuba as mi
import drjit as dr

mi.set_variant('cuda_ad_rgb')

N = 10_000_000

# Empty scene with a small off-centered sphere
scene = {
    'type': 'scene',
    'sphere': {
        'type': 'sphere',
        'radius': 0.01, # small
        'center': mi.ScalarPoint3f(0.01, 0.02, -0.01), # not exactly centered
        'bsdf': {
            'type': 'roughconductor',
            'alpha': 0.2,
            'distribution': 'ggx'
        }
    }
}
scene = mi.load_dict(scene)

# Sample positions on the sphere
sphere = scene.shapes()[0]
sampler = mi.load_dict({'type': 'independent'})
sampler.seed(0xBADCAFE, N)
ps = sphere.sample_position(0, sampler.next_1d())

# Build rays with origin (1, 1, 1) pointing to samples positions on the sphere
ray_o = mi.Point3f(1)
ray_d = dr.normalize(ps.p - ray_o)
ray = mi.Ray3f(o=ray_o, d=ray_d)

# Intersect the scene, some might miss the sphere
si = scene.ray_intersect(ray)
print(f"Valid primary rays: {dr.count(si.is_valid())}") # Should be close to N
print(f"{dr.count(si.wi.z < 0)=}") # Greater than zero due to precision issues in OptiX

Secondly, for something like a roughconductor, the sampling scheme might still produce "invalid" BSDF samples. For example, the sampling of the microfact normal might lead to a outgoing direction that is below the surface. Once again, this is signalled implicitly by having output value be 0.