Better Gaussian Splatting in Julia: A Vibe-Coding Guide to Cutting-Edge 3D Rendering

Three years after the original 3D Gaussian Splatting paper transformed real-time radiance field rendering, the technology has moved from research labs into production pipelines, AR/VR experiences, and even smartphone apps. Yet the dominant implementations remain locked in Python prototypes and C++ engines. For developers who prefer a language that combines mathematical clarity with native performance, Julia offers a compelling alternative. And with the rise of "vibe coding" — the art of generating and refining code through AI assistants — recreating and improving Gaussian splatting has never been more accessible. In this guide, we'll explore how to leverage Julia's expressive type system, GPU programming capabilities, and the emerging ecosystem to build better Gaussian splatting pipelines, while using AI-based tools to iterate faster.

Understanding 3D Gaussian Splatting

Gaussian splatting, or 3D Gaussian Splatting (3DGS), is a point-based rendering technique introduced by Kerbl et al. in their SIGGRAPH 2023 paper "3D Gaussian Splatting for Real-Time Radiance Field Rendering" [1]. Unlike neural radiance fields (NeRFs) that query a neural network along each ray, 3DGS represents a scene implicitly as a collection of millions of anisotropic 3D Gaussian functions. Each Gaussian is characterized by:

  • Center μ: a 3D position.
  • Covariance matrix Σ: a symmetric positive-definite 3×3 matrix that encodes the ellipsoid's orientation and scaling.
  • Color c: usually per-point RGB values modulated by spherical harmonics (SH) for view-dependent effects.
  • Opacity α: a scalar controlling transparency.

To render an image, all Gaussians are projected to the 2D image plane, sorted by depth, and then alpha-blended. The projection of a 3D Gaussian under the pinhole model yields a 2D Gaussian whose covariance is Σ' = J Σ J^T, where J is the Jacobian of the perspective projection. The final pixel color C is the weighted sum of the ordered Gaussians that cover that pixel:

C = Σ c_i α_i ∏_{j<i} (1 - α_j)

The learning process optimizes the parameters (center, covariance, colors, opacity) to minimize photometric loss against ground truth images. This is done via differentiable rasterization and stochastic gradient descent. The original paper achieves real-time frame rates on modern GPUs, a significant improvement over NeRFs.

Why Julia for Gaussian Splatting?

Julia was designed to solve the "two-language problem": high-level languages are easy to develop in but slow, while low-level languages offer performance but are painful to write. Julia's LLVM-based JIT compiler allows dynamic typing and multiple dispatch while generating machine code comparable to C. The official Julia micro-benchmarks show that Julia often matches C and Fortran in numerical tasks [2].

For Gaussian splatting, several features make Julia particularly attractive:

  • Multiple dispatch: Define methods for different Gaussian types (e.g., Float32, Float64, Dual numbers for autodiff) without code duplication.
  • StaticArrays.jl: Fixed-size vectors and matrices are stack-allocated. When looping over 500,000 Gaussians, this eliminates heap allocations and boosts cache performance.
  • Native GPU programming: CUDA.jl allows writing CUDA kernels in Julia, with access to shared memory, atomics, and custom vectorization.
  • Automatic differentiation: Zygote.jl and ForwardDiff.jl provide source-to-source and forward-mode AD, so you can differentiate your rasterizer directly without a separate framework like PyTorch.
  • Metaprogramming: Macros generate specialized kernel code for different hardware, a significant advantage when targeting new GPU architectures.

These capabilities allow a single Julia implementation to match hand-tuned C++ performance while remaining as readable as Python.

Vibe Coding: A New Superpower

The term "vibe coding" was coined by Andrej Karpathy in early 2025. It describes a workflow where you rely on an AI assistant to generate and iterate on code, focusing on high-level intent rather than every syntax detail. Vibe coding exploded in popularity because AI assistants, trained on massive codebases, now produce correct and idiomatic code for most common tasks.

Julia is exceptionally well-suited to vibe coding. Its syntax is concise and mirrors mathematical notation, which reduces the chance of AI generating semantically wrong code. Interactive notebooks like Pluto.jl automatically execute cells when dependencies change, giving immediate visual feedback. For 3D rendering, Pluto is ideal: you can adjust a Gaussian's covariance and watch the rendered image change in real time.

A typical vibe-coding workflow for a Gaussian splatting task:

  1. Launch Pluto.
  2. Load a small dataset of camera images.
  3. Ask an AI: "Write Julia code to initialize a grid of 1000 Gaussians with random colors and covariances from point cloud data."
  4. Generate a rasterizer with AI: "Write a basic rasterizer with depth sorting and alpha blending for the Gaussians."
  5. Run the code, inspect the image, refine prompts like "add a scale regularization term" or "enable adaptive density control".
  6. Use @time and JET.jl to profile and ask the AI to optimize the bottleneck.

AI assistants like GitHub Copilot have become standard. ASI Biont supports integration with GitHub Copilot via its AI-assisted workflow — learn more at asibiont.com/courses.

Building a Minimal Gaussian Splatter in Julia

Let's create a simple, educational implementation. We'll assume Julia 1.10+ and the packages StaticArrays, LinearAlgebra, and CUDA. This example focuses on clarity rather than production-grade optimization.

Define the Gaussian Structure

using StaticArrays, LinearAlgebra

struct Gaussian{T}
    center::SVector{3,T}
    cov::SMatrix{3,3,T}
    color::SVector{3,T}
    opacity::T
end

We store the full covariance matrix for simplicity; a production version would store only the six unique parameters.

Projection to 2D

The projection of a 3D Gaussian under a pinhole camera involves transforming the Gaussian into camera coordinates and then linearly approximating perspective projection. The Jacobian J at point (x_c, y_c, z_c) is:

J = [ f_x/z_c 0 -f_x * x_c / z_c^2
0 f_y/z_c -f_y * y_c / z_c^2 ]

In Julia:

function project(g::Gaussian{T}, R::SMatrix{3,3,T}, t::SVector{3,T}, K::SMatrix{3,3,T}) where T
    cam_center = R * g.center + t
    x, y, z = cam_center
    J = SMatrix{2,3,T}(
        K[1,1]/z, 0, -K[1,1]*x/(z*z),
        0, K[2,2]/z, -K[2,2]*y/(z*z)
    )
    cov2d = J * R * g.cov * R' * J'
    return (u = x/z, v = y/z, cov2d = cov2d, color = g.color, opacity = g.opacity, depth = z)
end

This projected Gaussian can now be rasterized.

CPU Rasterization and Alpha Blending

For simplicity, we'll loop over all Gaussians and all pixels in the bounding box of each Gaussian. The blending equation is straightforward:

function rasterize_cpu(gaussians::Vector{Gaussian{T}}, R, t, K, width, height) where T
    img = zeros(RGB{Float32}, height, width)
    depth_buffer = fill(Inf32, height, width)
    for g in gaussians
        pr = project(g, R, t, K)
        # Compute pixel radius from the 2D covariance
        evals = LinearAlgebra.eigvals(Symmetric(pr.cov2d))
        radius = 3 * sqrt(maximum(evals))
        xmin = max(1, floor(Int, pr.u - radius))
        xmax = min(width, ceil(Int, pr.u + radius))
        ymin = max(1, floor(Int, pr.v - radius))
        ymax = min(height, ceil(Int, pr.v + radius))
        for y in ymin:ymax, x in xmin:xmax
            dx = x - pr.u; dy = y - pr.v
            # Compute Gaussian falloff
            inv_cov = inv(Symmetric(pr.cov2d))
            alpha = pr.opacity * exp(-0.5 * (dx^2 * inv_cov[1,1] + 2*dx*dy*inv_cov[1,2] + dy^2 * inv_cov[2,2]))
            if alpha > 0.01  # ignore negligible contributions
                idx = y + (x-1)*height
                if pr.depth < depth_buffer[idx]
                    img[idx] = (1 - alpha) * img[idx] + alpha * pr.color
                    depth_buffer[idx] = pr.depth
                end
            end
        end
    end
    return img
end

This naive implementation works for a few hundred Gaussians, but becomes slow for larger scenes. The real power comes from GPU acceleration.

GPU Kernel with CUDA.jl

On a GPU, we parallelize over tiles or pixels. A minimal CUDA kernel that splats all Gaussians to a single pixel (without sorting) is:

using CUDA

function splat_kernel!(img, depth, gaussians, R, t, K, width, height)
    idx = (blockIdx().x - 1) * blockDim().x + threadIdx().x
    if idx <= length(gaussians)
        g = gaussians[idx]
        pr = project(g, R, t, K)
        # Accumulate the Gaussian's contribution to the center pixel
        # For a real implementation, compute bounding box and loop
        # Use atomic operations to avoid race conditions when blending
    end
    return
end

The actual 3DGS implementation uses a tile-based approach with depth sorting per tile. Writing that in Julia is entirely feasible, and the readability of the kernel code is far superior to raw CUDA C++.

Techniques for "Better" Gaussian Splatting

"Better" can mean higher visual quality, fewer artifacts, or faster rendering. Here are three areas where Julia gives you an edge.

1. Adaptive Density Control

The original 3DGS paper uses a heuristic to densify regions with large positional gradients and prune low-opacity Gaussians. This adaptive control is crucial for capturing fine details. In Julia, you can implement it with simple array operations and views. The entire gradient computation can be written in a differentiable way, and you can use Zygote to derive the gradients of the loss with respect to the Gaussian parameters.

2. Regularization to Prevent Floaters

A common problem is the appearance of "floaters" — semi-transparent Gaussians that hover in the air. Researchers have proposed adding regularization terms to the loss, such as penalizing the trace of the covariance matrix (encouraging smaller Gaussians) or adding a prior on opacity. Julia's native AD makes it trivial to add such terms:

using Zygote

function regularization_loss(g::Gaussian)
    # Encourage small, isotropic Gaussians
    return tr(g.cov) + 0.5 * log(g.opacity)
end

# In the training loop
grads = Zygote.gradient((g) -> photometric_loss(g) + λ * regularization_loss(g), g)[1]

3. Anti-Aliasing with Mip-Splatting

Mip-Splatting [3] introduces a low-pass filter to reduce aliasing at different scales. The technique convolves the Gaussian's spectral content with a footprint function. Implementing this in Julia involves Fourier transforms and custom kernels, but the high-level syntax makes the math directly translatable.

4. Spherical Harmonics for View-Dependent Effects

Using SH coefficients allows the model to capture specular highlights. Julia's StaticArrays and fast basis function implementations make SH evaluation and gradient computation very efficient. You can even use macros to generate optimized code for different SH degrees.

Performance: Julia vs. Python vs. C++

It's tempting to compare Julia against PyTorch in a general sense, but Gaussian splatting is a custom computation. The baseline implementations in the original repository are written in CUDA C++, while many research forks use PyTorch for the optimization loop and CUDA for rasterization. Julia offers a unified language that can handle both.

Feature Julia (CUDA.jl) Python (PyTorch) C++ (CUDA)
Kernel authoring effort Low (native Julia) Moderate (C++ extensions) High
Autodiff integration Native Zygote PyTorch autograd Manual/gradients
Memory safety Automatic GC Automatic GC Manual
Type stability Fully controllable Dynamic Static
Metaprogramming Excellent Limited Preprocessor only
Learning curve Moderate Low High

In benchmarks comparing custom CUDA kernels, the original 3DGS paper reports rendering at 100+ FPS on an RTX 3090 [1]. A Julia implementation, once properly optimized, can approach that rate because the generated PTX instructions are nearly identical to hand-written CUDA. The key advantage is development speed: you iterate in a high-level language while the compiler handles low-level optimizations.

Ecosystem and Resources for Julia Graphics

The Julia ecosystem has a solid set of tools for 3D rendering:

  • Makie.jl: High-performance interactive 3D plotting for visualizing point clouds and Gaussian distributions.
  • GeometryBasics.jl: Defines meshes and geometric primitives.
  • JET.jl: Static analysis to catch type instabilities, which is essential for performance.
  • LoopVectorization.jl: SIMD optimizations for CPU loops.
  • CUDA.jl: First-class GPU development, including support for writing custom kernels and using shared memory.

For datasets, you can use the original 3DGS scenes from the official project page [1]. The Julia community on Discourse [4] is active, and there are emerging GitHub repositories exploring Gaussian splatting in Julia.

A Step-by-Step Vibe Coding Workflow

Let's put everything together into a concrete workflow that you can replicate today.

  1. Set Up Pluto: Pkg.add("Pluto"), then load a notebook.
  2. Load Images: Use Images.jl to read a series of photos of a static scene.
  3. Initialize Gaussians: Ask an AI for code that converts a depth map or sparse point cloud into Gaussian objects.
  4. Render an Initial View: Generate a rasterizer and display the output in Pluto.
  5. Define the Loss: Use Zygote to compute gradients of a photometric loss (e.g., mean squared error against a reference image).
  6. Optimize: Write a simple SGD or Adam loop. In Pluto, re-rendering happens automatically each iteration.
  7. Refine: Add features like adaptive densification, SH coefficients, or regularization. Use prompts to the AI to generate code snippets.
  8. Profile: Use @time, JET.@profile, or CUDA.@profile to find bottlenecks. Ask the AI to optimize the inner loops.

This iterative vibe-coding approach lets you focus on the science rather than boilerplate code. The tight feedback loop of Pluto combined with Julia's dynamic nature makes exploration a joy.

Challenges and Limitations

While Julia is a great choice for exploring Gaussian splatting, there are some hurdles:

  • Ecosystem maturity: Julia's graphics ecosystem is smaller than Python's. You may need to build some tools yourself.
  • Compilation latency: The first run of a Julia function can be slow due to JIT compilation. This can disrupt the vibe flow, but PackageCompiler can precompile fixed pipelines.
  • Community: The number of researchers using Julia for 3DGS is still small, so you'll often be breaking new ground.

Despite these, the speed of development and execution makes Julia a strong candidate for the next generation of 3D rendering tools.

Conclusion

Gaussian splatting in Julia is not a pipe dream. With its powerful type system, GPU integration, and mathematical clarity, Julia offers a clean way to build high-performance, custom 3D reconstruction pipelines. The vibe-coding trend lowers the barrier to entry, enabling developers to generate intermediate implementations and explore novel ideas quickly. Whether you're a graphics researcher looking to prototype new optimizations or a hobbyist who wants to understand photorealistic 3D, Julia is a compelling choice.

The future of 3D rendering is undoubtedly being splatted, and those who combine the expressiveness of Julia with the velocity of AI-assisted development will be the ones to make it better.

References

[1] Kerbl, B., Kopanas, G., Leimkühler, T., Drettakis, G. "3D Gaussian Splatting for Real-Time Radiance Field Rendering." ACM Transactions on Graphics, SIGGRAPH 2023. https://repo-sam.inria.fr/fungraph/3d-gaussian-splatting/

[2] Julia Language Benchmarks. https://julialang.org/benchmarks/

[3] Barron, J. T., Mildenhall, B., Hedman, P., et al. "Mip-Splatting: Aliasing-Free 3D Gaussian Splatting." 2023. https://arxiv.org/abs/2311.16477

[4] Julia Discourse. https://discourse.julialang.org/

← All posts

Comments