Doppler spectrum (OHD, offline)¶

This tutorial renders the Doppler spectrum that a coherent lidar with optical heterodyne detection measures, with DopplerHistogramPathTracerInline. The tall box of the Cornell box comes towards the camera at 20 m/s and the short box moves away at 20 m/s. With a 1550 nm laser and the light next to the camera, light reflected once by a box is shifted by \(2v/\lambda \approx \pm 26\) MHz; light that bounces between the boxes and the walls spreads over other shifts. Nothing in the scene actually moves: each path’s shift follows from the velocities of its vertices.

It uses the scene described in Doppler rendering.

The Cornell box's steady image, mean Doppler shift and three pixels' spectra

1. Load the scene¶

The scene is loaded with DontMergeMaterials: Falcor otherwise merges identical materials, and every white surface, the two boxes included, would be called Floor. With it, each object keeps the name of its material, which the velocities refer to.

import falcorcomp as falcor
import numpy as np
import matplotlib
matplotlib.use("Agg")  # no window
import matplotlib.pyplot as plt

testbed = falcor.Testbed(create_window=False)
# Keep every material separate, so that each object keeps its name (identical materials are merged by default).
testbed.load_scene("cornell-box/scene-v4-nolight.pbrt", falcor.SceneBuilderFlags.DontMergeMaterials)
testbed.resize_frame_buffer(256, 256)
testbed.scene.camera.aspectRatio = 1.0
testbed.clock.pause()

2. Build the render graph¶

LaserLight puts a point light at the camera (isLightSourceLaser = false, laserCollocated = true). The tracer’s spectrum is set by:

wavelength: the laser wavelength, in nm, which turns path velocities into Doppler shifts.
frequencyMin, frequencyMax, frequencyBin: 256 bins from -100 to 100 MHz.
velocities: the motion of each object, by name; the camera looks along \(-z\), so \(+z\) is towards it.
accumulate: sum the frames in the spectrum, as the transient histogram tracer does.

computeDirect is on: the single reflection off each box is the strongest part of the signal.

F_MIN, F_MAX, BINS = -100.0, 100.0, 256  # Doppler shifts, MHz

graph = testbed.create_render_graph("DopplerSpectrum")
graph.create_pass("VBuffer", "VBufferRT", {"samplePattern": "Center", "sampleCount": 1})
# A point light at the camera, as in a lidar.
graph.create_pass("Light", "LaserLight", {
    "isLightSourceLaser": False, "laserCollocated": True, "laserPower": [10.0, 10.0, 10.0],
})
graph.create_pass("Tracer", "DopplerHistogramPathTracerInline", {
    "samplesPerPixel": 64, "maxBounces": 3, "computeDirect": True,
    "wavelength": 1550.0,  # nm
    "frequencyMin": F_MIN, "frequencyMax": F_MAX, "frequencyBin": BINS,
    # The camera looks along -z: the tall box comes towards it, the short box moves away. Nothing moves in the
    # scene; the velocities only decide the Doppler shift of each path.
    "velocities": {"TallBox": {"linear": [0.0, 0.0, 20.0]}, "ShortBox": {"linear": [0.0, 0.0, -20.0]}},
    "accumulate": True,
})
graph.add_edge("VBuffer.vbuffer", "Tracer.vbuffer")
graph.add_edge("VBuffer.viewW", "Tracer.viewW")
graph.add_edge("Light", "Tracer")
graph.mark_output("Tracer.spectrum")
testbed.render_graph = graph

3. Render¶

FRAMES = 64  # 64 frames x 64 spp
for _ in range(FRAMES):
    testbed.frame()

4. Read the spectrum¶

to_numpy() returns the spectrum as bins x height x width x RGBA, in radiance per MHz.

# With accumulate, the spectrum is the sum of the frames, per MHz: divide by the frame count for the mean.
spectrum = graph.get_output("Tracer.spectrum").to_numpy()[..., :3].mean(-1) / FRAMES  # (bins, height, width)
np.save("doppler_spectrum.npy", spectrum)
frequencies = F_MIN + (np.arange(BINS) + 0.5) * (F_MAX - F_MIN) / BINS  # bin centers, MHz

5. Show the image, the mean shift and three spectra¶

Summed over all shifts, the spectrum is the ordinary image. Its mean shift per pixel shows the tall box at about +26 MHz and the short box at about -26 MHz, and the walls next to each box tinted by light that reached them through the box. The three spectra show the single-reflection peaks of the two boxes and of the static back wall (at 0), and the spread of light that bounced between moving and static surfaces; the tall box’s spectrum reaches 50 to 65 MHz through paths that meet it twice.

image = spectrum.sum(0) * (F_MAX - F_MIN) / BINS  # the steady image: the spectrum summed over frequency
mean_shift = (spectrum * frequencies[:, None, None]).sum(0) / np.maximum(spectrum.sum(0), 1e-12)
pixels = {"tall box": (140, 90), "short box": (185, 175), "back wall": (60, 200)}

figure, axes = plt.subplots(1, 3, figsize=(15, 4.4))
axes[0].imshow(np.clip(image / np.percentile(image, 99), 0, 1) ** (1 / 2.2), cmap="gray")
axes[0].set_title("steady image")
shown = axes[1].imshow(mean_shift, cmap="bwr", vmin=-30, vmax=30)
axes[1].set_title("mean Doppler shift (MHz)")
figure.colorbar(shown, ax=axes[1])
for (name, (y, x)), color in zip(pixels.items(), ["C0", "C1", "C2"]):
    for axis in axes[:2]:
        axis.plot(x, y, "o", mfc="none", mec=color, mew=2)
    axes[2].semilogy(frequencies, spectrum[:, y, x] + 1e-6, color=color, label=name)
axes[2].set_xlabel("Doppler shift (MHz)")
axes[2].set_title("spectra (per MHz)")
axes[2].legend()
for axis in axes[:2]:
    axis.set_xticks([])
    axis.set_yticks([])
figure.tight_layout()
figure.savefig("doppler_spectrum.png", dpi=90)

The full script: doppler_spectrum_offline.py. It takes about 10 seconds.

A Doppler-gated image¶

A spectrum keeps every shift. To render only the light in one band of shifts, a single image instead of frequencyBin bins, replace the tracer with DopplerGatedPathTracerInline and give it a gate: a center and a width in MHz. The scene, the light and the velocities stay the same; its color output is averaged with AccumulatePass. A gate 2 MHz wide at +25.8 MHz keeps the light that the tall box reflected once (the Doppler-gated image of the overview):

graph.create_pass("Tracer", "DopplerGatedPathTracerInline", {
    "samplesPerPixel": 64, "maxBounces": 3, "computeDirect": True,
    "wavelength": 1550.0,  # nm
    "frequencyCenter": 25.8, "frequencyGateWindow": 2.0, "frequencyGateMode": "box",  # MHz
    "velocities": {"TallBox": {"linear": [0.0, 0.0, 20.0]}, "ShortBox": {"linear": [0.0, 0.0, -20.0]}},
})
graph.create_pass("Accumulate", "AccumulatePass", {"precisionMode": "SingleCompensated"})
graph.add_edge("VBuffer.vbuffer", "Tracer.vbuffer")
graph.add_edge("VBuffer.viewW", "Tracer.viewW")
graph.add_edge("Light", "Tracer")
graph.add_edge("Tracer.color", "Accumulate.input")
graph.mark_output("Accumulate.output")

Connections whose shift falls outside the gate are skipped before they are shaded, so a narrow gate renders faster than the spectrum. See DopplerGatedPathTracerInline for the gate kernels and the gate sweep.

An FMCW lidar¶

With chirpBandwidth, the same pass renders the up- and down-chirp spectra of an FMCW lidar, whose beat frequencies also depend on the path length; see the FMCW lidar tutorial.