Run GPU-Powered Black Hole Interferometry with NVIDIA cuPhoton

AI-Generated Summary NVIDIA cuPhoton.xBeam accelerates radio interferometry correlation for the Event Horizon Telescope by processing voltage recordings into visibilities on GPUs. The FX correlator ingests VDIF frames, applies delay and fringe corrections, channelizes signals via FFT, cross-multiplies station pairs, and writes visibilities in DiFX SWIN or xbeam-v1 format. On an NVIDIA Vera Rubin NVL72 system with 64 GPUs and GPUDirect Storage, xBeam correlated a 10 TB dataset in 70 seconds, achieving a 460x speedup over an x86 CPU baseline of 31,992 seconds. Scaling splits observation time across GPUs using --devices on a single node or --shard-count and --shard-index across nodes, with merge-shards combining partial SWIN outputs into a complete visibility set. AI agents can orchestrate end-to-end pipelines by parsing VDIF metadata, generating JSON job specifications, launching xBeam correlation across GPU clusters, and monitoring execution. Next Steps Explore the NVIDIA cuPhoton repository for imaging modules, xBeam, and sample recordings with job files. Read the Faster Scientific Image Analysis with NVIDIA cuPhoton walkthrough for an imaging-focused tutorial. Review the NVIDIA cuPhoton product page for additional resources and related examples. Imaging black holes is challenging because they are distant and appear tiny in the sky—for example, M87*’s observed ring spans just 42 microarcseconds. The black hole itself emits no light, so telescopes capture radio emission from superheated gas swirling outside the event horizon. The black hole’s intense gravity bends the paths of this light, shaping the bright ring that frames the central shadow. The Event Horizon Telescope (EHT) addresses this challenge by linking 13 radio observatory stations across locations, including Hawaii, Chile, Spain, Mexico, and the South Pole, through very-long-baseline interferometry (VLBI). These stations record millimeter and submillimeter radio signals as high-speed voltage streams, which are synchronized and combined to reconstruct detailed images. The next-generation EHT (ngEHT) aims to add observing sites to improve image quality and track changes around black holes over time. This post walks through GPU-based interferometry with cuPhoton.xBeam, from reading telescope recordings to generating visibilities and scaling across multiple GPUs. These measurements serve as inputs to downstream calibration and image reconstruction, which are outside the scope of this tutorial. NVIDIA cuPhoton.xBeam performs GPU-based correlation for interferometry and organizes recorded data for parallel reads. It correlates signals from pairs of telescopes to produce one complex visibility for each pair of stations, each frequency channel, and each short time interval. The next step can involve constructing a 3D gas model from these visibilities via a distinct downstream process, such as Gaussian splatting. This method models the emitting gas as numerous compact, glowing blobs moving along orbits near a spinning black hole, optimized to match the observed visibilities. cuPhoton does not currently incorporate this second component. The NVIDIA cuPhoton team collaborates with the Harvard-Smithsonian Center for Astrophysics and the MIT Haystack Observatory to accelerate data ingestion and correlation using xBeam. Running on an NVIDIA Vera Rubin NVL72 system with 64 GPUs and GPUDirect Storage (GDS) in compatibility mode, xBeam processed a 10 TB dataset in 70 seconds—about 460x faster than an x86 CPU baseline of 31,992 seconds (8.89 hours). The sky is an inverse problem The EHT records radio signals in two polarizations at up to 64 gigabits per second per telescope, generating petabytes of data. Correlating signals between telescope pairs produces visibilities, which sample Fourier components of the sky’s brightness. Earth’s rotation traces sparse measurement tracks, and reconstruction algorithms combine these samples into an image. Rapid changes in sources such as Sgr A* complicate this process because observations capture an evolving source over time. An Earth-sized array observing at 1.3 mm only collects (u, v) measurements along the elliptical tracks drawn by Earth’s rotation. Because no physical lens forms the image directly, reconstructive algorithms must solve an inverse problem to synthesize the picture from sparse samples. For targets like Sgr A*, the source changes during an observation session, complicating reconstruction methods that assume a static source (Focus on First Sgr A* Results). Figure 2 shows the imaging workflow using a synthetic source. The left panel shows the (u, v) tracks that an array laid out like the 2017 EHT samples over one night as Earth turns. The next panel shows visibility amplitudes, |V(u, v)| plotted against baseline length, followed by a dirty image formed from the sparse samples. The right panel shows a crescent model, a simple ring-shaped image that researchers can fit to those amplitudes. Light passing near a black hole can orbit it before escaping, forming a thin photon ring near the shadow’s edge. Its size and shape can help constrain the black hole’s mass, spin, and viewing angle, while its interferometric signal produces characteristic oscillations. The proposed Black Hole Explorer (BHEX) aims to measure this ring, and future ngEHT observations at 345 GHz could detect its signature through polarization measurements. What cuPhoton.xBeam computes cuPhoton.xBeam is the FX correlator in the cuPhoton toolkit. It restructures the recordings so they can be read in parallel, builds a .vdix index so later reads jump straight to each frame, and loads the voltages into GPU memory. Then it does the FX work, like DiFX (CPU-based, open-source software correlator): it Fourier-transforms each station’s voltage into frequency channels (F) and cross-multiplies the station and polarization pairs at each frequency channel requested by the job (X). Those steps run on the GPU. In order, xBeam does five things: 1. Ingest. Read the VDIF frames and place each one at its correct time. xBeam builds a small index file, the .vdix sidecar, that records where every frame is without copying the data, so later reads jump straight to the right place. xBeam indexes VDIF frames by timestamp, thread, and physical location. Its .vdix sidecars let subsequent reads locate the required data without rewriting the recording. Missing or invalid data remain explicitly tracked. 2. Delay and fringe. Shift each station’s samples by the delay in the job’s model, so the same wavefront lines up at every station, and remove the phase drift caused by Earth’s motion. Delay correction spans both sides of the fast Fourier transform (FFT). Integer-sample alignment and fringe rotation occur before channelization; a frequency-dependent phase correction handles the remaining fractional-sample delay afterward. 3. Channelize. Run an FFT on each station’s voltage stream to split it into frequency channels. 4. Cross-multiply and accumulate. For every pair of stations, multiply each channel’s spectrum from one station by the complex conjugate of the corresponding spectrum from the other station and sum the products over the integration time. Invalid input is masked, and the accumulated products are normalized using validity-derived exposure weights. The selected processing profile controls additional admission and reported-weight conventions. 5. Write visibilities. Save the result in DiFX’s SWIN format or in xBeam’s own xbeam-v1 format, whichever the job asks for. The delay xBeam applies comes from the job file. It is a sum of different delays: the geometric delay, which is parameterized as a polynomial in time around a reference time t0, plus a delay resulting from the offset and linear drift of each station’s clock, and a fixed clock offset for each frequency band: Within one FFT, xBeam either holds the delay constant or approximates it with a straight line or parabola, which is order 0, 1, or 2 in the setting fringe_approximation_order. The phase correction applied to the samples is where νref is the band’s reference frequency, and νLO is the local-oscillator frequency term. For stations a and b at frequency ν, the product that is summed is the spectrum of station a times the complex conjugate of the spectrum of station b: xBeam does not compute delays from station positions. It uses the delays given in the job. It also does not trace light paths around a black hole. DiFX is the CPU correlator used as the reference for validating xBeam results. Tools such as HOPS, rPICARD, and eht-imaging take the visibilities from there for calibration and imaging. Data moves from disk into GPU memory in one of two ways: standard POSIX reads through pinned host memory, or reads through cuFile, the GPUDirect Storage API. Whether cuFile uses a direct storage-to-GPU path depends on the configuration and runtime conditions. The phase correction, the FFT, and the cross-multiply all run on the GPU. Jobs that use several GPUs split the scan in time (--shard-count, --shard-index), and merge-shards joins the partial SWIN outputs. Run one job xBeam (part of cuPhoton v0.1.3) takes a normalized JSON job, schema version 1, or a supported DiFX .input file. A sample recording and matching job file can be found in the NVIDIA/cuPhoton repository. The job names stations, clocks, frequencies, VDIF datastreams, baselines, and scans. Point JOB.json at the recordings. The output directory must be new. The run writes visibilities in the job’s output format, swin or xbeam-v1. Readers can verify a successful run by checking that the output matches the expected SWIN visibility records and validating the completion headers. Inspect a recording first. inspect-vdif reads headers only: cuphoton xbeam inspect-vdif RECORDING.vdif --json VDIF files may contain interleaved threads, missing frames, or frames stored out of timestamp order. Build .vdix sidecars before correlation so xBeam can locate the requested timestamp ranges without rewriting the recordings. Reuse them while the source files remain unchanged. Every later read then jumps straight to the right frame. The builder uses cuFile by default; pass --backend posix if the machine does not have GPUDirect Storage: cuphoton xbeam build-vdif-sidecar RECORDING-*.vdif \ --index-dir SIDECARS \ --backend cufile \ --devices 0,1,2,3 Run the correlation on the GPU. –backend cupy selects the GPU path. By default, xBeam reads data through pinned host memory (POSIX). –vdif-io cufile selects the cuFile I/O path through GPUDirect Storage. That is the mode used in the Vera Rubin runs, and it needs no other flags. cuphoton xbeam correlate JOB.json \ --output RUN \ --backend cupy \ --vdif-io cufile \ --scan SCAN_ID \ --first-integration 0 \ --integration-count N \ --index-dir SIDECARS The same job can run from Python. load_job returns the same job object the command line uses. correlate returns visibility records one at a time; writing SWIN files is left to the command-line tool. from cuphoton.xbeam import CorrelationWindow, correlate, load_job job = load_job("JOB.json") window = CorrelationWindow("SCAN_ID", first_integration=0, integration_count=8) for record in correlate( job, backend="cupy", window=window, index_dir="SIDECARS", vdif_io="cufile", ): spectrum = record.spectrum # complex visibilities for one baseline product Gaps in a recording affect how the correlator assigns times and weights to subsequent samples. xBeam offers an explicit buffered-reader policy through --profile interop-real-usb-buffered-window-replay-v1. Within a retained reader buffer, this policy continues consuming recorded samples using the buffer’s existing time reference until the next allocation establishes a new anchor. It also determines which FFT blocks contribute and their reported visibility weights. It does not recover missing samples, and a full reported weight does not mean every timestamp was recorded. Without replay, xBeam aligns subsequent frames by their header timestamps and marks missing slots invalid, which can produce different visibilities and weights. Scale across GPUs and nodes cuPhoton.xBeam scales by splitting the scan in time across GPUs. Each GPU correlates its own stretch of time. On one node, xBeam can start several GPUs at once. Across many nodes, each GPU runs one piece, called a shard, and merge-shards joins the SWIN outputs. On one node, --devices starts one correlator process per GPU, gives each its share of the integrations, writes the results under RUN/shards/, and merges them: cuphoton xbeam correlate JOB.json \ --output RUN \ --backend cupy \ --vdif-io cufile \ --scan SCAN_ID \ --first-integration 0 \ --integration-count N \ --devices 0,1,2,3 \ --index-dir SIDECARS Across nodes, each process runs one time shard. Set --shard-count and --shard-index together, with the same scan ID, first integration, and integration count for every shard. For 64 GPUs, use --shard-count 64 with indices 0 through 63. In the following command, set SHARDS=64 and set SHARD to the unique index for each process: cuphoton xbeam correlate JOB.json \ --output RUN/shard-$SHARD \ --backend cupy \ --vdif-io cufile \ --scan SCAN_ID \ --first-integration 0 \ --integration-count N \ --shard-count $SHARDS \ --shard-index $SHARD \ --index-dir SIDECARS When every shard has finished, merge them into one set of visibilities. Before copying anything, the merge checks that all shards came from the same job, that together they cover the full time range, and that their SWIN headers agree: cuphoton xbeam merge-shards JOB.json RUN/shard-* \ --output RUN/complete --devices is the shortcut for one node. --shard-count and --shard-index with merge-shards is the path across nodes. Both run delay and phase corrections, FFTs, and cross-multiplications on the GPUs. Figure 4 shows the result of scaling xBeam on multi-node GPUs. On Vera Rubin NVL72, one 10 TB recording is split across 1 to 64 GPUs on several nodes, and the shards are merged into one set of visibilities. The reported speedup continues to increase over the tested range. It is measured on total processing time: reading the data, preparing it, and correlating it, end to end. On 64 GPUs that total is 70 s, against about 32,000 s for the x86 CPU baseline, which makes it about 460x faster. With GPUDirect Storage, read bandwidth on the same 10 TB recording grows with the number of GPUs, from 20 GB/s on 4 GPUs to 170 GB/s on 64 GPUs. Read bandwidth increases with GPU count over the tested range. AI agents and the rest of the toolkit AI agents offer a powerful framework for autonomous science by acting as automated assistants that orchestrate data ingestion, job configuration, parallel execution, and visibility extraction. In this workflow, an agent can autonomously parse incoming VDIF metadata, generate normalized JSON job specs (schema version 1), invoke cuphoton xbeam correlate across GPU clusters, and monitor execution. While current benchmarks focus on GPU acceleration rather than full agentic control, integrating AI agents with xBeam offers a compelling potential application to automate complex multi-telescope data pipelines end-to-end. The rest of the cuPhoton toolkit applies this same GPU-first pattern across domain data. For optical imaging, xDataReader loads image datasets, while xScan, xPois, and xFit analyze them. xRay processes detector frames for X-ray laser experiments directly in GPU memory. Researchers working with materials, quantum optics, and fusion can apply similar cuPhoton processing pipelines on multi-node, multi-GPU systems for data analysis and simulation. Get started The public cuPhoton repository provides the imaging modules and xBeam. See the public repository and related posts for available tools and examples. Repository: NVIDIA cuPhoton Imaging walkthrough: Faster Scientific Image Analysis with NVIDIA cuPhoton. Related examples: Accelerated X-Ray Analysis for Nanoscale Imaging (XANI) of Novel Materials and Using Accelerated Computing to Live-Steer Scientific Experiments at Massive Research Facilities. Learn more: NVIDIA cuPhoton About the Authors Comments