🔭 astrophysics

cTreeBalls: a fast 3-point correlation function code for clustering measurements

cTreeBalls is a high-performance Python/C package that leverages octree and kd-tree algorithms to efficiently compute 3-point correlation functions for massive datasets (over 200 million pixels) in under 10 minutes, specifically enabling feasible clustering analyses for upcoming LSST data.

Original authors: Mario A. Rodriguez-Meza, Eladio Moreno, Alejandro Aviles, Gustavo Niz

Published 2026-04-13

📖 4 min read☕ Coffee break read

CC BY 4.0

Original authors: Mario A. Rodriguez-Meza, Eladio Moreno, Alejandro Aviles, Gustavo Niz

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine you are a detective trying to solve a cosmic mystery. You have a map of the universe filled with billions of stars and galaxies. Your goal is to understand how these objects are arranged. Do they cluster together in groups? Do they form specific shapes?

For a long time, astronomers have looked at pairs of objects (like two stars) to understand the universe. This is like looking at how often two people stand next to each other at a party. It's useful, but it only tells you half the story.

Sometimes, the real secrets are hidden in triangles. If you look at three galaxies at once, you can see if they form a tight little group, a long line, or a wide-open triangle. This "three-point" information reveals the hidden, non-random forces that shaped our universe.

However, there's a massive problem: The math is too hard.

If you have 100 galaxies, checking every possible pair is easy. But checking every possible triangle among 100 galaxies is hard. Now, imagine you have 200 million galaxies (which is what the next generation of telescopes will see). Trying to check every possible triangle with a standard computer is like trying to count every grain of sand on a beach by picking them up one by one. It would take longer than the age of the universe to finish the job.

Enter cTreeBalls (or "cBalls")

The authors of this paper have built a new super-tool called cTreeBalls. Think of it as a "smart organizer" for the universe. Here is how it works, using some everyday analogies:

1. The "Smart Filing Cabinet" (Tree Structures)

Instead of looking at every single galaxy one by one, cTreeBalls organizes the universe into a giant, nested filing cabinet (called an Octree).

Imagine a big box containing the whole universe.
Inside that box, you have 8 smaller boxes.
Inside those, 8 even smaller boxes, and so on.
If a small box is empty (a void in space), the computer skips it immediately. It only opens the boxes that actually contain galaxies. This saves a massive amount of time.

2. The "Shared Neighborhood" Trick

The code is also incredibly clever about how it searches.

Imagine you are at a party. If you want to find everyone standing within 5 feet of you, you scan the crowd.
Now, imagine your friend moves just one step to the left. Do you need to scan the entire crowd again? No! You already know who is nearby. You just need to check the tiny difference.
cTreeBalls uses this logic. If two galaxies are close to each other, they likely share the same neighbors. The code calculates the neighbors once and shares that list, rather than recalculating it for every single galaxy.

3. The "Magic Translator" (Harmonic Decomposition)

This is the secret sauce. Usually, calculating triangles involves complex math with sine and cosine waves (trigonometry), which is slow for computers.

cTreeBalls translates these complex waves into a different language called Chebyshev polynomials.
Think of it like translating a difficult poem into a simple song. The computer can "sing" the answer (calculate the result) much faster than it can "read" the complex poem.

Why Does This Matter?

The Rubin Observatory (a giant new telescope) is about to start taking pictures of the entire sky. It will find 20 billion objects.

Without cTreeBalls: Analyzing the "triangle" patterns of this data would take centuries.
With cTreeBalls: The authors tested it on a simulation with 200 million pixels (a full sky). It finished the job in less than 10 minutes on a single powerful computer.

The Bottom Line

cTreeBalls is a high-speed engine that allows astronomers to finally look at the "triangles" of the universe. By organizing data smartly, sharing work, and using a mathematical shortcut, it turns a task that was previously impossible into a routine calculation. This will help scientists unlock the secrets of Dark Energy and understand exactly how our universe grew from a smooth cloud of gas into the complex web of galaxies we see today.

Based on the provided paper, here is a detailed technical summary of cTreeBalls (cBalls):

1. Problem Statement

The analysis of large-scale structure in the universe relies heavily on correlation functions to extract cosmological information. While two-point statistics (2PCF) are standard, they fail to capture non-Gaussian information generated by non-linear cosmic evolution and primordial non-Gaussianities. The three-point correlation function (3PCF) is the necessary next step to capture this information (encoding spatial triangles in the distribution).

However, calculating the 3PCF for modern and future surveys (e.g., the Vera C. Rubin Observatory's LSST, which will contain billions of objects) presents a severe computational bottleneck:

Naïve Approach: A brute-force estimator scales as $O(N^3)$ , where $N$ is the number of data points. For $N > 10^9$ , this is computationally prohibitive.
Tree-Based Approach: Using spatial tree structures (like kd-trees or octrees) to find neighbors reduces scaling to $O(N^2 \log N)$ , which remains too slow for full-sky, high-resolution simulations (e.g., $N_{side}=4096$ HEALPix maps containing millions of pixels).

2. Methodology

The authors developed cTreeBalls, a hybrid Python/C package designed to efficiently compute 2PCF and 3PCF for projected scalar fields on a sphere (and periodic boxes for 2PCF). The code achieves high performance through several key algorithmic innovations:

Harmonic Decomposition & Basis Transformation:
Instead of calculating raw 3PCF counts, the code maps the 3PCF signal onto an orthogonal basis using harmonic decomposition of one of the triangle angles. This reduces the problem complexity from counting triplets to calculating pair counts, theoretically scaling as $O(N \log N)$ , similar to 2PCF codes.
Tree Structures (Octree & kd-tree):
The spatial search is optimized using octrees (Barnes & Hut) and kd-trees (Bentley).
- Pruning: The code implements "oc-ball-trees" by pruning empty cubic cells, ensuring the search only visits regions populated by galaxies.
- Shared Neighbor Lists: To further reduce redundancy, the algorithm identifies "close" pivots (galaxies) that share similar neighbor lists. It defines four length scales (root cube, smallest cell, bin width, mean separation) to automatically determine a pivot radius, allowing the reuse of neighbor lists for nearby pivots.
Chebyshev Polynomial Approximation:
The calculation of multipoles requires evaluating sine and cosine functions for triangle configurations. cBalls replaces these expensive trigonometric evaluations with Chebyshev polynomials, which are significantly faster to compute.
Parallelization & Architecture:
- Language: The core search engine is written in C for speed, wrapped in Python via Cython (cyballs) for user accessibility.
- Threading: The search loop is parallelized using OpenMP, distributing pivot galaxies across multiple CPU threads.
- I/O Flexibility: Supports various formats including CFITSIO, HEALPix, GADGET-2, and ASCII (e.g., ROCKSTAR). It can handle masks and cross-correlations between multiple maps.

3. Key Contributions

Scalability: The code can process >200 million HEALPix pixels (a full-sky simulation with $N_{side}=4096$ ) in under 10 minutes on a single high-performance computing node.
LSST Readiness: It is specifically engineered to handle the upcoming LSST data volume (planned >20 billion objects), making 3-point statistics feasible for future cosmological analyses.
Complementarity: Unlike existing tools like TreeCorr (ball-tree based) or Triumvirate (tri-polar spherical harmonic decomposition), cBalls uses a specific harmonic decomposition of the 3PCF angles. This offers a complementary approach with potentially different constraining power on cosmological parameters.
Versatility: It supports both projected fields on a sphere (using the Limber approximation) and cubic box simulations with periodic boundaries. It also includes edge corrections for partial-sky surveys.

4. Results & Performance Metrics

Scaling Law: Empirical testing on HEALPix maps ( $N_{side}$ from 256 to 4096) revealed a power-law scaling of $N^{1.1}$ , confirming the code behaves nearly linearly with the number of objects.
Multipole Scaling: The time required to compute multipoles scales linearly with the maximum multipole number ( $m_{max}$ $m_{ma x}$ ). The relationship is approximated as $t(m) = 1.85m + 63.7$ $t (m) = 1.85 m + 63.7$ seconds for $N_{side}=1024$ $N_{s i d e} = 1024$ .
- The slope (~2) corresponds to the operations needed for recursive Chebyshev relations.
- Computing up to $m=8$ is only ~20% slower than computing up to $m=1$ .
Large-Scale Benchmark: A full-sky analysis with $N_{side}=8192$ (~800 million pixels) and a search radius of 200 arcminutes took approximately 40 CPU minutes on 128 threads of a Perlmutter-NERSC node. The primary bottleneck identified was scanning the tree and updating histograms, not the mathematical decomposition.

5. Significance

cTreeBalls fills a critical gap in the analysis pipeline for next-generation cosmological surveys. By reducing the computational cost of 3-point statistics from prohibitive ( $O(N^3)$ ) to feasible ( $O(N^{1.1})$ ), it enables the extraction of non-Gaussian information from massive datasets like LSST. This capability is essential for:

Breaking degeneracies in cosmological parameters that 2-point statistics cannot resolve.
Testing models of primordial non-Gaussianity.
Integrating into massive workflows like TXPipe for the Rubin Observatory.

The software is open-source (CC BY 4.0), well-documented, and designed for integration into modern data analysis pipelines, marking a significant step forward in the practical application of higher-order clustering statistics.