Gaia DR3 Variable Star Finder

InterSystems Employee Programming Challenge #1

Identifies every astronomical source in the Gaia DR3 epoch-photometry archive whose Blue Photometer (BP) or Red Photometer (RP) flux changed by more than 100% across all valid observations.

Result: 57,099 variable sources across 20 input files, processed in ~1.3 seconds.

How It Works

Pipeline Overview

20x EpochPhotometry_*.csv.gz
        |
        | ThreadPoolExecutor (one thread per file)
        v
  _process_file() — Polars: decompress gzip, parse CSV,
                    extract bp_flux / rp_flux columns,
                    compute per-source min/max
        |
        | pl.concat() — merge all 20 DataFrames
        v
  filter: pct_change > 100%  →  write data/out/result.csv
        |
        v
  [benchmark timer stops here]
ai_summary() — bonus step, runs after timing, requires GEMINI_API_KEY
|
v
data/out/ai_summary.txt — Gemini interpretation of top 5 sources

Formula

For each source, percentage change is computed per band:

pct_change = (max_flux - min_flux) / min_flux × 100

The result is the larger of the two bands (BP and RP). Sources where this exceeds 100% are written to output.

Polars applies this as a vectorized column expression across all rows at once — no Python loop, no row-by-row iteration:

result = (
    pl.concat(parts)
    .with_columns([
        ((pl.col("bp_max") - pl.col("bp_min")) / pl.col("bp_min") * 100).alias("bp_pct"),
        ((pl.col("rp_max") - pl.col("rp_min")) / pl.col("rp_min") * 100).alias("rp_pct"),
    ])
    .with_columns(pl.max_horizontal("bp_pct", "rp_pct").fill_null(0).alias("pct"))
    .filter(pl.col("pct") > 100)
)

NaN, null, and empty flux values are ignored during parsing — a source only appears if a real band swung past 100%.

Parallel Ingestion with Polars

with ThreadPoolExecutor() as ex:
    parts = list(ex.map(_process_file, files))

Files are sorted largest-first so the biggest gzip stream starts immediately. Each thread handles one file independently: decompress → CSV parse → min/max extraction. Polars reads only the three needed columns (source_id, bp_flux, rp_flux) out of 48 — the other 45 columns are never decoded.

The per-row flux arrays arrive as strings ("[1820.8, 2013.8, NaN, ...]"). These are parsed entirely inside Polars’ Rust engine in a single expression chain:

pl.col("bp_flux")
  .str.replace_all("NaN", "")
  .str.strip_chars("[]")
  .str.split(",")
  .list.eval(pl.element().cast(pl.Float64, strict=False))

list.min() and list.max() are then computed in one scan per column, with no separate Python calls.

AI Summary (bonus)

After the benchmark timer stops, the top 5 most variable sources are sent to Google Gemini for astronomical interpretation. This runs outside the timed section and is skipped gracefully if GEMINI_API_KEY is not set.

The code first attempts the InterSystems AI Hub path (langchain_intersystems) — AI Hub manages LLM configurations (API keys, model names) inside IRIS via %ConfigStore.Configuration, providing a unified init_chat_model() interface. If AI Hub is unavailable, it falls back to calling Gemini directly via langchain-google-genai.

Performance Journey

We profiled the baseline to find the real bottleneck before optimizing:

Bottleneck Profile (3-file sample, gzip + csv.reader baseline)

This told us that json.loads — the obvious target — was only 6% of runtime. The real cost was I/O and row parsing.

What We Tried (20 files, 8 CPUs)

Why Polars

• Columnar parsing: only source_id, bp_flux, rp_flux are decoded — the other 45 CSV columns are skipped entirely.
• Rust gzip: decompression happens inside Polars’ Rust engine, not Python’s gzip module.
• No intermediate list: the flux string is split and cast to Float64 in one expression; list.min() and list.max() scan it once.

Why Numpy Vectorization Lost

Gaia flux arrays have variable length — from a handful to hundreds of values per row. To use numpy’s matrix operations, every row must be padded to the same length. The padding itself requires a Python loop, and np.genfromtxt on a large padded block costs more than the vectorization saves. Numpy vectorization is optimal for fixed-shape data (e.g. image pixels), not variable-length arrays.

Correctness Guarantee

Every optimization was verified against a fixed oracle:

python -m pytest test_runchallenge.py    # 27 unit tests — all pass

Installation

Prerequisites

• Docker Desktop
• git
• Google Gemini API key — free at aistudio.google.com (AI summary step only)

Environment Variables

Set on Windows:

setx GEMINI_API_KEY "AIzaSy..."

Run

git clone https://github.com/angela81ku/isc-gaia-challenge.git
cd isc-gaia-challenge
docker-compose up --build -d
docker-compose exec iris iris session iris -U USER
USER> do ^RunScript

Output

Project Structure

├── src/
│   ├── gaia.py              # Core logic — parallel parse, Polars, AI summary
│   ├── register_llm.py      # Registers Gemini config in IRIS AI Hub at startup
│   └── RunScript.mac        # IRIS entry point  (do ^RunScript)
├── data/
│   ├── in/                  # 20 EpochPhotometry_*.csv.gz input files
│   └── out/result.csv       # Output — 57,099 variable sources
├── Dockerfile               # IRIS Community 2026.3.0AI + Polars + langchain
├── docker-compose.yml
└── startup.sh               # Compiles RunScript.mac + registers Gemini config

Tech Stack

Data Attribution

Input: ESA Gaia DR3 Epoch Photometry.
Column definitions: Gaia Archive Documentation
How to cite Gaia: Credits