Trace reflection paths

Trace specular reflection chains, solve exact paths for a known reflector sequence, evaluate reflected fields, and accumulate them on a receiver grid.

RayD has four reflection operations. They go from geometry only to a full field on a grid:

Operation	Input	Output
`trace_reflections`	rays	the chain of specular hits for each ray
`trace_refl_epc`	a source, a receiver, and the reflectors to use	the exact reflection points of that path (equivalent-path correction, EPC)
`trace_refl_epc_field`	the same, plus per-slot material	the complex electric field at the receiver
`accumulate_reflections`	rays from a transmitter, a grid, per-primitive material	field and power summed into grid cells

All examples use a single floor at z = 0:

import math
import drjit as dr
import rayd.drjit as rd
from drjit.cuda import Array3f, Array3i, Float, Int, Bool
 
floor = rd.Mesh(Array3f([-5.0, 5.0, 5.0, -5.0], [-5.0, -5.0, 5.0, 5.0], [0.0] * 4),
                Array3i([0, 0], [1, 2], [2, 3]))
scene = rd.Scene()
scene.add_mesh(floor)
scene.build()

Trace specular chains #

trace_reflections(ray, max_bounces) follows each ray through up to max_bounces mirror reflections.

ray = rd.Ray(Array3f([0.0], [0.0], [1.0]), Array3f([0.3], [0.0], [-1.0]))
 
trace = scene.trace_reflections(ray, 2)               # symbolic=True: ReflectionTrace
print(trace.bounce_count)                             # [1]
print(trace.bounce(0).hit_points)                     # [[0.3, 0, 0]]
 
chain = scene.trace_reflections(ray, 2, symbolic=False)   # flat ReflectionChain
print(chain.hit_points)       # [[0.3, 0, 0], [0, 0, 0]]: slot 1 is unused
print(chain.image_sources)    # [[0, 0, -1], [0, 0, 0]]
print(chain.trailing_dir)     # [[0.3, 0, 1]]: the ray after its last bounce

ReflectionTrace gives one ReflectionBounce per bounce through bounce(i) or bounces. ReflectionChain stores the same data flat, ray_count * max_bounces slots, ray-major. In both, a slot with prim_ids < 0 is unused; is_valid() returns that mask. The trailing_* fields describe the segment after the last reflection.

The image source of a bounce is the mirror image of the previous source across the reflector. It is what image-method RF solvers use to compute path length and delay.

Deduplicate paths #

Many rays from one transmitter find the same reflector sequence. In Dr.Jit, pass deduplicate=True (or a ReflectionTraceOptions) to keep one representative ray per sequence. canonical_prim_table maps primitives to a shared id, for example all triangles of one wall, and image_source_tolerance sets how close two image sources must be to merge. The result reports dedup_keep_mask, representative_ray_index, and discovery_count (how many rays found each path).

trace = scene.trace_reflections(ray, 2, deduplicate=True)

ReflectionTraceOptions.export_mode trims the result: rd.REFLECTION_EXPORT_FULL (default), REFLECTION_EXPORT_MINIMAL, or REFLECTION_EXPORT_COUNT_ONLY. return_trailing=False skips the segment after the last bounce.

Solve an exact path #

Ray tracing finds which reflectors a path uses, but its hit points are only as exact as the ray that found them. trace_refl_epc(source, receiver, max_bounces, options) solves the exact specular path between a source and a receiver that reflects from the given primitives, then checks that every segment is visible.

opt = rd.ReflEpcOptions()
opt.expected_prim_ids = Int([0])    # one slot per ray and bounce
seed = rd.Ray(Array3f([0.0], [0.0], [1.0]), Array3f([0.5], [0.0], [-1.0]))
rx = Array3f([1.0], [0.0], [1.0])
 
path = scene.trace_refl_epc(seed, rx, 1, opt)
print(path.valid, path.path_length)   # [True] [2.23607]
print(path.reflection_points)         # [[0.5, 0, 0]]

expected_prim_ids has ray_count * max_bounces entries. Surface groups (surface_group_id, surface_group_size, surface_group_members, surface_max_group_size) let a path move to a neighbouring triangle of the same flat surface; set visibility_ignore_mode = "surface_group" to ignore the whole group in the visibility checks.

Evaluate the reflected field #

trace_refl_epc_field(source, receiver, max_bounces, options) solves the same path and evaluates the complex electric field at the receiver. The source is a Ray (whose first hit picks the reflector) or a transmitter position (then expected_prim_ids is required). Give the material and plane of each slot in a ReflEpcFieldOptions:

fo = rd.ReflEpcFieldOptions()
fo.slot_plane_normal = Array3f([0.0], [0.0], [1.0])   # required: the default is a zero normal
fo.slot_eta_r = Float([4.0])
fo.slot_mu_r = Float([1.0])
fo.slot_sigma = Float([0.0])
fo.slot_gain = Float([1.0])
fo.wavelength = 0.125
fo.omega = 2 * math.pi * 299792458.0 / 0.125
 
field = scene.trace_refl_epc_field(seed, rx, 1, fo)
print(field.valid, field.path_length)     # [True] [2.23607]
print(field.field_x_re, field.field_x_im) # [-0.000896339] [-0.000755447]
 
# Same path from a transmitter position: name the reflector explicitly.
fo.expected_prim_ids = Int([0])
field = scene.trace_refl_epc_field(Array3f([0.0], [0.0], [1.0]), rx, 1, fo)

tx_polarization sets the transmitter polarization (default (1, 0, 0)). The return_* flags choose which geometry the result carries. The field is differentiable in the Dr.Jit backend through ReflEpcFieldOptionsAD and AD inputs; in Torch this operation is forward-only.

Warning

If slot_plane_normal is left at its default zero vector, every path comes back invalid. Set the plane normal of each slot.

Accumulate onto a grid #

accumulate_reflections(ray, tx_position, grid, material, max_bounces, options) launches the given rays from a transmitter, follows their reflections, and sums field and power into the cells of an axis-aligned grid. It is the building block of a coverage map.

n = 4096
i = dr.arange(Float, n)
theta = 0.45 * dr.pi * (i + 0.5) / n
phi = 2.0 * dr.pi * (i * 0.618034 - dr.floor(i * 0.618034))
d = Array3f(dr.sin(theta) * dr.cos(phi), dr.sin(theta) * dr.sin(phi), -dr.cos(theta))
tx = Array3f(0.0, 0.0, 1.0)
rays = rd.Ray(dr.zeros(Array3f, n) + tx, d)
 
grid = rd.AccumGrid()
grid.axis = 2                 # plane normal along z
grid.position = 1.0           # the plane z = 1
grid.coord0_min, grid.coord0_max = -2.0, 2.0
grid.coord1_min, grid.coord1_max = -2.0, 2.0
grid.resolution0 = grid.resolution1 = 8
 
mat = rd.Material()           # one entry per global primitive
mat.eta_r = dr.full(Float, 5.0, 2)
mat.sigma = dr.full(Float, 0.01, 2)
mat.mu_r = dr.full(Float, 1.0, 2)
mat.gain = dr.full(Float, 1.0, 2)
mat.valid = dr.full(Bool, True, 2)
 
opts = rd.AccumOptions()
opts.wavelength = 0.125
opts.k = 2 * math.pi / 0.125          # wavenumber
opts.solid_angle_per_ray = 2 * math.pi / n
opts.cell_area = (4.0 / 8) ** 2
 
res = scene.accumulate_reflections(rays, tx, grid, mat, 1, opts)
print(res.grid_cell_count)            # 64
print(res.reflection_count)           # [2434]

reflection_power has one value per cell, and reflection_field_x, _y, _z hold the complex field per cell. Material fields are relative permittivity eta_r, conductivity sigma, relative permeability mu_r, an extra gain, and valid. AccumOptions also controls Russian roulette (rr_depth, rr_prob, stop_threshold) and can record diffraction-wedge events (collect_wedges, wedge_capacity) for a later diffraction pass.

In Torch, the grid is rayd.torch.AccumGrid, the material is ReflMaterial (use ReflMaterial.default(count, device=...) for a start), and the operation is forward-only.

Gradients #

trace_reflections is differentiable with respect to geometry and rays in both backends. trace_refl_epc_field and accumulate_reflections are differentiable in Dr.Jit and forward-only in Torch. trace_refl_epc has no derivative. Check Differentiation for the full table.

Trace specular chains#

Deduplicate paths#

Solve an exact path#

Evaluate the reflected field#

Accumulate onto a grid#

Gradients#

Related#