Julia实战：从零实现3D Gaussian Splatting渲染与优化

3D Gaussian Splatting（3DGS）是当前三维视觉和实时渲染领域关注度很高的方向，它的核心思路是用一组带位置、形状、不透明度和颜色的三维高斯原语表达场景，通过可微光栅化把三维场景投到二维图像上，再根据真实照片反向优化这组高斯原语的参数。相比传统 NeRF，它的训练速度快、渲染速度快，因而在三维重建、新视角合成、动态场景编辑等场景里被广泛应用。Julia 是一门兼顾动态语言开发效率和静态性能的编程语言，它适合用来实现 3DGS 的原理验证、课程作业、小型实验，以及探索算法改进。这篇文章会从 3DGS 的核心原理讲起，再带你在 Julia 里搭建环境、设计数据结构、实现最小可运行的光栅化和训练流程，然后深入分析 Julia 性能优化和内存管理需要注意的关键点，最后给出调试排查路径和生产环境建议。

全文面向有一定 Julia 基础、正在学习 3DGS 原理，或者想把 3DGS 官方 CUDA 实现迁移到 Julia 进行算法实验的开发者。阅读完你可以掌握 3DGS 的渲染链路，理解 Julia 中“类型稳定”“减少堆分配”“StaticArrays 使用”这些优化手段在图形学代码里的实际作用，并具备继续扩展成一个完整可训练项目的能力。

1. 3D Gaussian Splatting 到底做了什么

1.1 一句话说清楚核心思路

3D Gaussian Splatting 要解决的核心问题是：如何用一组可微分的三维高斯分布来表示一个场景，并且让这组分布在给定相机视角下能够快速渲染成一张高质量图像。

如果把场景比作一幅拼图，那么传统点云用离散的点描述表面，NeRF 用神经网络隐式描述连续空间，而 3DGS 用的是无数个小椭圆体。每个椭圆体都有自己的中心位置、三轴缩放、旋转角度、不透明度和颜色。渲染时，从相机位置观察这些椭圆体，把它们投影到屏幕上形成半透明的二维椭圆，最后按深度从小到大做 alpha 混合，得到最终像素颜色。

这个过程和传统 Splatting 思路一脉相承，但有两个重要改进：所有参数都可以微分，因此可以直接用梯度下降优化；光栅化时做了基于 tile 的前向排序和 alpha 合成，让渲染速度达到实时。理解了这一点，后续看代码和参数就不会迷路。

1.2 渲染流程拆解：初始化、优化、密度控制

3DGS 的完整训练循环可以拆成四个阶段：

1. 初始化：从稀疏点云或随机位置生成初始高斯原语，每个高斯带一个协方差矩阵、不透明度、RGB 颜色和球谐系数。
2. 前向渲染：把三维高斯协方差投影到二维图像平面，结合相机内外参计算每个高斯在屏幕上的椭圆位置和大小。
3. 损失计算：比较渲染图像与真实照片的 L1、SSIM 或两者组合损失。
4. 反向更新：梯度回传到每个高斯的位置、协方差、不透明度、颜色；同时每隔一定迭代步执行密度控制：删除不透明度过低的高斯，克隆位置梯度大的高斯，分裂尺度较大的高斯。

这个流程和 NeRF 的端到端训练类似，但区别在于 3DGS 没有神经网络编码器，它直接用可优化参数作为场景表示。因此每一轮迭代都是在直接调整“椭圆体的摆放、大小、颜色和透明度”。

1.3 为什么选择 Julia 而不是 Python/CUDA 原型

官方 3DGS 实现基于 CUDA 和 PyTorch，性能很高，但对想理解算法细节的开发者来说，代码里大量并行逻辑会和核心数学表达式混在一起。Python 纯 NumPy 版本容易看懂，但逐像素光栅化太慢，实验迭代很不舒服。

Julia 在表现力上接近 Python，在数值计算的运行时性能上又接近 C。它的多重派发特性很适合为不同数据类型写同一套算法：CPU 实验用 Float32 数组，GPU 实验用 CuArray，算法主体可以保持一致。对于 3DGS 这种“数学表达简单、工程实现复杂”的算法，Julia 是很好的实验语言：先写 CPU 版本验证公式，再逐步替换为 GPU kernel 或 CUDA.jl 扩展。

需要注意，Julia 不是 3DGS 官方生态的主流语言，很多现成工具如相机参数读取、COLMAP 点云转换、网格导出都需要自己封装。文章后面的示例定位是“学习用最小实现”，面向公式验证和算法研究，生产级大规模训练仍要考虑 CUDA 实现或调用官方工程链路。

2. 搭建 Julia 实验环境

2.1 Julia 版本与关键依赖

建议安装 Julia 1.9 以上版本。较早版本虽然也能运行，但泛型多重派发、任务并行和包解析行为与当前版本有差异，排查问题时会增加不必要的干扰。

核心依赖按使用场景划分如下：

在 Julia 中创建项目并添加依赖：

julia -e 'using Pkg; Pkg.activate("gaussian_splatting_julia"); Pkg.add(["StaticArrays", "Zygote", "CUDA", "Plots"])'

说明：这里使用Pkg.activate创建独立环境，避免污染全局环境。如果网络环境下载慢，可以先只添加必要依赖，CUDA 等扩展包放到 GPU 实验阶段再安装。

2.2 最小项目结构

推荐按模块拆分文件，每个文件只负责一个职责：

gaussian_splatting_julia/ ├── Project.toml ├── src/ │ ├── GaussianSplatting.jl │ ├── gaussian.jl │ ├── camera.jl │ └── rasterizer.jl ├── test/ │ └── run_basic.jl └── data/ └── 示例图片或点云

Project.toml由Pkg.generate创建，或者直接用文本编辑器维护。最小内容示例：

name = "GaussianSplattingJulia" uuid = "xxxxxxxx-xxxx-xxxx-xxxx-xxxxxxxxxxxx" authors = ["yourname"] version = "0.1.0" [deps] LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f"

注意：uuid必须在Pkg.generate时自动生成，不要手写固定值。如果手动维护Project.toml，请用Pkg.add保持格式一致。

2.3 环境检查清单

开始写算法前，先确认环境是否正常：

1. Julia 能正常启动：终端执行julia --version。
2. 包能正确解析：在项目目录执行julia --project -e 'using StaticArrays, LinearAlgebra'。
3. 自动微分包能加载：执行using Zygote，确认没有依赖冲突。
4. CUDA 包（如需要）能识别显卡：执行using CUDA; CUDA.functional()。

后面所有代码默认在--project=.激活环境下运行。如果加载某个包报错，最优先检查Project.toml里依赖列表和Manifest.toml是否一致，命令是：

julia --project -e 'using Pkg; Pkg.resolve(); Pkg.instantiate()'

3. 用 Julia 实现一个最小 Gaussian Splatting 流程

这一节会实现一个可运行的最小流程，重点展示 3DGS 的数学链路和 Julia 的数据结构设计。这里不追求与官方 CUDA 实现完全一致，而是用可读的 Julia 代码跑通“定义高斯 -> 投影 -> 光栅化 -> 计算损失”的闭环。

3.1 高斯原语的数据结构设计

每个三维高斯可以用以下参数描述：

• mean：三维中心点。
• scale：三轴缩放，对应协方差矩阵的特征值平方根。
• rotation：旋转四元数 (w, x, y, z)。
• opacity：不透明度，经 sigmoid 映射到 (0, 1) 区间。
• color：RGB 颜色。
• sh_coeffs：球谐系数，学习阶段可选。

在 Julia 中，用具体结构体而非抽象字段定义 Gaussian：

using StaticArrays struct Gaussian3D{T<:Real} mean::SVector{3,T} scale::SVector{3,T} rotation::SVector{4,T} opacity::T color::SVector{3,T} end

这里把字段类型写死为SVector{3,T}，关键目的是保证 Julia 编译器能推断出结构体内存的精确布局。如果把字段类型写成Vector或AbstractVector，每个字段都是一个需要动态分配内存的数组，性能和类型稳定性都会变差。

从四元数构造旋转矩阵：

function quaternion_to_matrix(q::SVector{4,T}) where T w, x, y, z = q return @SMatrix [ 1 - 2*(y*y + z*z) 2*(x*y - z*w) 2*(x*z + y*w) 2*(x*y + z*w) 1 - 2*(x*x + z*z) 2*(y*z - x*w) 2*(x*z - y*w) 2*(y*z + x*w) 1 - 2*(x*x + y*y) ] end

旋转矩阵作用于原有缩放矩阵后，得到三维协方差矩阵：

function cov_matrix(g::Gaussian3D{T}) where T R = quaternion_to_matrix(g.rotation) S = Diagonal(g.scale) Σ = R * S * S * R' return Σ end

这里把协方差矩阵拆成旋转和缩放两个部分，比直接优化协方差矩阵的六个独立元素更稳定。直接优化协方差矩阵容易产生非半正定矩阵，渲染时会出现椭圆面积异常甚至 NaN。

3.2 相机模型与投影：3D 协方差到二维椭圆

相机模型使用简化的针孔模型。相机参数包括位置、从世界坐标系到相机坐标系的旋转矩阵、焦距和主点：

struct PinholeCamera{T} position::SVector{3,T} orientation::SMatrix{3,3,T} fx::T fy::T cx::T cy::T end

把三维高斯投影到二维屏幕坐标的步骤：

1. 把世界坐标系中的高斯中心点转换到相机坐标系。
2. 用透视投影公式计算屏幕坐标。
3. 用仿射变换近似透视投影的雅可比矩阵，把三维协方差投影成二维协方差。

实现代码：

function project_gaussian(g::Gaussian3D{T}, cam::PinholeCamera{T}) where T cam_center = cam.orientation * (g.mean - cam.position) depth = cam_center[3] depth <= 0 && return nothing u = cam.fx * cam_center[1] / depth + cam.cx v = cam.fy * cam_center[2] / depth + cam.cy J = @SMatrix [ cam.fx/depth 0.0 -(cam.fx * cam_center[1]) / (depth^2) 0.0 cam.fy/depth -(cam.fy * cam_center[2]) / (depth^2) ] world_to_cam = cam.orientation cov3d = cov_matrix(g) cov2d = J * world_to_cam * cov3d * world_to_cam' * J' return (u=u, v=v, cov2d=cov2d, depth=depth, opacity=g.opacity, color=g.color) end

需要注意，depth <= 0的剔除是必要步骤。在针孔模型中，位于相机后方的高斯不能投影到成像平面，不剔除会导致屏幕上出现反向矩阵或错误叠加。

3.3 逐像素 alpha blending 光栅化

得到每个高斯的二维椭圆参数后，就可以把椭圆叠加到图像上。真实 3DGS 使用基于 tile 的并行排序光栅化，理解起来较复杂。学习用最小版本采用逐像素 alpha blending，逻辑更直观。

首先计算二维高斯在某个像素处的权重：

function gaussian_weight_2d(dx::T, dy::T, cov2d::SMatrix{2,2,T}) where T det = cov2d[1,1] * cov2d[2,2] - cov2d[1,2] * cov2d[2,1] det <= eps(T) && return zero(T) inv_det = inv(det) a = cov2d[2,2] * inv_det b = -cov2d[1,2] * inv_det c = cov2d[1,1] * inv_det exponent = -0.5 * (a * dx^2 + 2 * b * dx * dy + c * dy^2) return clamp(exponent, -30, 30) |> exp end

这里省略了归一化系数，因为后续 alpha 混合还会做归一化，保持相对权重即可。clamp(exponent, -30, 30)防止指数溢出，这是数值稳定性的关键。

光栅化主循环：

struct ProjectedGaussian{T} u::T v::T cov2d::SMatrix{2,2,T} depth::T opacity::T color::SVector{3,T} end function rasterize_frame(gaussians::Vector{Gaussian3D{T}}, cam::PinholeCamera{T}; width::Int, height::Int) where T projected = ProjectedGaussian{T}[] for g in gaussians proj = project_gaussian(g, cam) proj === nothing && continue push!(projected, ProjectedGaussian(proj.u, proj.v, proj.cov2d, proj.depth, proj.opacity, proj.color)) end sort!(projected, by=p -> p.depth, rev=true) image = zeros(RGB{Float64}, height, width) for py in 1:height for px in 1:width acc = RGB{Float64}(0, 0, 0) alpha_accum = 0.0 for p in projected alpha_accum >= 0.999 && break dx = px - 0.5 - p.u dy = py - 0.5 - p.v contrib = gaussian_weight_2d(dx, dy, p.cov2d) alpha = p.opacity * contrib acc = acc + (1 - alpha_accum) * alpha * p.color alpha_accum += (1 - alpha_accum) * alpha end image[py, px] = acc end end return image end

这个代码块有三个关键点：

1. sort!(projected, by=p -> p.depth, rev=true)表示从远到近排序，这样近处的高斯后绘制，能正确遮挡远处高斯。
2. alpha 混合公式里(1 - alpha_accum)表示当前像素还剩多少透明度没有被覆盖。
3. dx = px - 0.5 - p.u中的 0.5 是像素中心偏移，习惯用像素中心作为采样点，能避免半像素偏移误差。

3.4 最小训练循环与损失函数

在 3DGS 中，颜色参数和球谐系数可以直接优化，位置、缩放、旋转、透明度也可以用梯度更新。最小示例使用 MSE 损失，并借助 Zygote 做反向传播：

using Zygote function loss_function(gaussians::Vector{Gaussian3D{T}}, cam::PinholeCamera{T}, target::Matrix{RGB{Float64}}; width::Int, height::Int) where T rendered = rasterize_frame(gaussians, cam; width=width, height=height) diff = Float64.(red.(rendered)) .- Float64.(red.(target)) return sum(diff .^ 2) / length(diff) end

如果用 Zygote 计算梯度，需要注意rasterize_frame中所有操作必须可微。StaticArrays 和广播运算在 Zygote 下兼容性较好，但涉及到排序的判断分支可能产生梯度中断。学习阶段的替代方案是手写梯度或简化光栅化细节。

一个简化训练循环：

function train_step!(gaussians::Vector{Gaussian3D{Float32}}, cam::PinholeCamera{Float32}, target::Matrix{RGB{Float32}}; lr=0.01) grads = Zygote.gradient(g -> loss_function(g, cam, target, width=256, height=256), gaussians) for (g, grad) in zip(gaussians, grads[1]) g.mean .-= lr .* grad.mean g.scale .-= lr .* grad.scale g.rotation .-= lr .* grad.rotation g.opacity -= lr * grad.opacity g.color .-= lr .* grad.color end end

完整训练还需要加入熵正则、透明度截断、密度控制等策略，这里先用最简流程验证公式是否正确。把训练迭代和渲染结果保存下来，观察损失是否下降、图像是否逐渐接近目标即可。

4. 深入优化 Julia 性能与内存管理

3DGS 的核心计算是大量三维高斯的投影和像素合成。用 Julia 实现时，性能瓶颈往往不是语言本身，而是代码没有做到类型稳定、分配了过多临时数组、或者在没有必要的地方用了抽象类型。这一节是 Julia 侧性能优化的重点。

4.1 类型稳定是第一优先级

Julia 能接近 C 性能的前提是类型稳定：编译器在运行时能够推断每个变量的类型，从而生成高效的机器码。

类型不稳定在 3DGS 代码里最常见的表现是：

• 结构体字段是抽象类型，例如Vector{Real}。
• 函数返回类型在分支中不一致，例如一个分支返回SVector，另一个分支返回nothing。
• 全局变量的类型不稳定。

以ProjectedGaussian{T}为例，之所以把cov2d写成SMatrix{2,2,T}，而不是写Matrix，就是为了让光栅化循环能够精确知道协方差矩阵的大小和元素类型。一旦cov2d是动态大小的矩阵，每个gaussian_weight_2d调用都可能触发动态分派，性能会明显下降。

检查类型稳定性的命令：

using InteractiveUtils @code_warntype project_gaussian(example_gaussian, example_camera)

如果输出里有红色Union或Any，说明返回类型没有被完全推断，需要检查分支条件或字段类型。

4.2 用 StaticArrays 消除堆分配

Julia 的普通数组Vector{Float64}在堆上分配，字段访问和矩阵运算会产生指针间接引用。对于每帧数千到数百万个高斯的数据结构，堆分配会带来大量垃圾回收压力。

SVector和SMatrix是固定长度、固定类型的栈上数组。对于三维点、四元数、2x2 协方差矩阵这种固定大小结构，使用SVector是最合适的选择。

错误写法：

# 不推荐：每个字段都独立分配 struct BadGaussian mean::Vector{Float32} scale::Vector{Float32} rotation::Vector{Float32} color::Vector{Float32} end

推荐写法见 3.1 节。使用SVector后，结构体本身可以通过malloc或栈分配，并且整个数组可以作为连续内存块存储。

检查垃圾回收压力的方法：

@time rasterize_frame(gaussians, cam; width=256, height=256)

输出中allocations如果非常大，说明有大量临时数组被创建。重点检查SMatrix乘法、Diagonal构造和sort是否创建了不必要的大对象。

4.3 避免临时数组与 in-place 优化

光栅化主循环里最常见的问题是坐标和颜色计算产生临时数组。例如：

delta = SVector(dx, dy) # 固定的 SVector，可以 weight = gaussian_weight_2d(delta[1], delta[2], cov2d) tmp_color = (1 - alpha_accum) * alpha * p.color # 仍是 SVector，可以

但如果写成：

tmp = Vector{Float64}(undef, 2) tmp[1] = dx tmp[2] = dy

就会在每次像素循环里分配数组，光栅化图像越大，性能越差。

对图像缓冲区，可以先预分配一次，而不是在每次光栅化时重新创建：

function rasterize_frame!(buffer::Matrix{RGB{Float64}}, gaussians::Vector{Gaussian3D{T}}, cam::PinholeCamera{T}) where T # 直接在 buffer 中写入像素 end

这样训练循环中同一块缓冲区反复使用，减少分配。

多线程化是另一个方向。逐像素循环天然可并行，可以使用Threads.@threads把像素行拆分：

Threads.@threads for py in 1:height for px in 1:width # ... end end

使用多线程前确保 Julia 以JULIA_NUM_THREADS=4或julia -t 4启动。线程数不是越多越好，要根据 CPU 核心数和缓存情况调整。

4.4 广播融合与 SIMD

Julia 的广播语法可以将多个数组操作融合成单个循环，而不是为每个中间结果生成临时数组。例如计算协方差矩阵投影时：

cov2d .= J * world_to_cam * cov3d * world_to_cam' * J'

如果cov2d是固定大小的SMatrix，编译器更容易做优化。广播语法.=避免在嵌套调用中产生中间临时矩阵。

Simd 指令由编译器自动生成，但前提是数据在内存中连续并且循环边界明确。也就是说，应当尽量让所有高斯数据存储在Vector{Gaussian3D{Float32}}这种同构连续数组里，而不是Vector{Any}或Vector{Gaussian3D}的非具体参数化版本。

4.5 GPU 路径：CUDA.jl 的接入思路

如果实验规模变大，CPU 光栅化会成为瓶颈。CUDA.jl 提供了 Julia 直接编写 CUDA kernel 的能力。接入思路通常分两步：

1. 把高斯参数转换为适合 GPU 并行处理的结构，例如CuVector{NTuple{12,Float32}}或拆成多个CuVector{Float32}。
2. 在 GPU kernel 中实现透视投影和 tile 级并行光栅化，用CuArray存渲染结果。

基本示例：

using CUDA function gpu_rasterize_wrapper(gaussians_cu::CuDeviceArray, params..., out_cu::CuDeviceArray) i = (blockIdx().x - 1) * blockDim().x + threadIdx().x if i <= length(out_cu) # 读取高斯、投影、写像素 end return nothing end

这一部分比 CPU 版本复杂很多。初学者建议先完成 CPU 版本并验证公式正确，再迁移 GPU。不要一开始就写 CUDA，否则调试成本和测试成本都会成倍增加。

5. 调试、验证与常见问题排查

5.1 可视化三个关键中间结果

调试 3DGS 时，最怕看到一张全黑或全白的图却不知道问题出在哪一层。推荐按以下顺序检查中间结果：

1. 屏幕坐标投影：画出每个高斯的(u, v)分布，确认点云是否正确投影到画面范围内。
2. 二维协方差椭圆：在图上叠加二维椭圆，确认长轴、短轴和旋转方向是否符合相机视角。
3. alpha 累积图：单独输出每个像素的累计 alpha 值，正常情况下应该是一张平滑的深度遮罩。

用 Plots 画投影点：

using Plots function visualize_projection(gaussians, cam; width=800, height=600) xs = Float64[] ys = Float64[] for g in gaussians proj = project_gaussian(g, cam) proj === nothing && continue push!(xs, proj.u) push!(ys, proj.v) end scatter(xs, ys, markersize=1, label="projected centers") end

如果投影点大量堆在图像边界或者出现镜像分布，优先检查相机朝向矩阵和景深方向。

5.2 常见错误与排查顺序

排查顺序建议：

1. 检查输入高斯初值是否合理。
2. 检查投影结果是否在合理范围内。
3. 检查二维协方差矩阵是否半正定。
4. 检查排序方向是否正确。
5. 检查 alpha 混合公式是否写对。

如果这些都正确，再往训练和梯度方向排查。不要在损失不降时就怀疑自动微分，先用可视化确认前向链路没问题。

5.3 学习环境与生产环境的差异

学习环境只需要验证公式和流程，跑一个小场景即可。通常可以用 256x256 图像和几百个高斯做实验，即使 CPU 版本非常慢也不会影响学习节奏。

生产环境如果要处理真实照片、几十万高斯和实时渲染，还需要额外考虑：

学习环境里用一行println足够，但生产环境必须考虑异常恢复。建议训练过程中每 N 步保存一次高斯参数，至少记录损失、梯度范数和有效高斯数量。

6. 优化扩展与实践建议

6.1 工程上值得做的改进

当前最小实现仍有大量可以改进的性能点：

• 使用 tile 光栅化：把图像分成固定大小的 tile，每个 tile 只处理与其覆盖区域有交集的高斯，避免全图逐像素扫描全部高斯。
• 高斯排序优化：按 tile 分组排序，而不是全屏排序；全屏排序在大规模高斯下会浪费大量计算。
• SoA 布局：把结构体数组 AoS 改成数组结构 SoA，即所有高斯的 mean 单独一个矩阵、所有 opacity 单独一个向量。GPU 环境下 SoA 内存访问更友好。
• 混合精度训练：使用 Float32 训练，必要时用 Float16 存储中间量，降低显存占用。
• 球谐系数：增加球谐阶数后，颜色表达从固定 RGB 变成视角相关，可用于复杂反射和新视角一致性优化。

这些点和 Julia 本身的性能优化不冲突：CPU 版本先把类型稳定和预分配做好，GPU 版本再考虑并行度和内存布局。

6.2 学习路径建议：从原理速通到源码阅读

如果刚开始接触 3DGS，推荐以下顺序：

1. 先跑通本文最小实现，确认前向链路和 alpha 混合理解无误。
2. 在平面上验证单个高斯的投影：改变相机位置，观察椭圆大小、方向和叠加关系。
3. 阅读官方 3DGS 论文原图，理解 tile 排序、密度控制和球谐系数的损失函数。
4. 阅读开源实现中关于光栅化和 backward kernel 的部分，重点看协方差反向传播公式。
5. 尝试在 Julia 里用几万高斯做小规模训练，看看密度控制、学习率衰减和损失函数的配合效果。

这个顺序比直接读大型代码库更平滑。了解原理后，“迭代参数与渲染”就不再是黑盒，而是能根据现象定位问题。

6.3 可复用清单

发布前或进入下一阶段实验前，建议对照以下清单检查：

• [ ] 高斯结构体字段是否全部具体类型，无抽象类型字段。
• [ ] 投影函数是否处理了depth <= 0的剔除。
• [ ] 协方差矩阵是否经过半正定检查。
• [ ] 排序方向是否正确，远到近排序。
• [ ] 光栅化是否预分配了图像缓冲区。
• [ ]@code_warntype是否显示无红色Union输出。
• [ ]@time是否观察到异常大量 allocations。
• [ ] 训练时是否记录损失、梯度范数、有效高斯数量。
• [ ] 是否定期保存 checkpoint。
• [ ] 是否对 NaN 和高斯透明度爆炸做了保护。

这份清单同样适合以后迁移到 CUDA 或加入球谐扩展时复用。技术方案会变，但类型稳定、数值稳定、异常可排查这三个原则不会变。

3D Gaussian Splatting 把三维场景表示从离散点云推进到了可微高斯原语，Julia 的高性能数值计算特性让研究者可以快速验证算法公式并优化实现。后续往真实场景扩展时，最值得优先投入的方向是 tile 光栅化和 CUDA kernel，其次是密度控制策略和球谐颜色表达。建议初学者先把本文的最小 CPU 流程跑通，理解每个中间矩阵的来源，再决定是否向 GPU 和完整训练流程推进。