nx_gsplat_density.nx source
↩ module page · 215 lines · 10426 B
1// nx_gsplat_density.nx -- ADAPTIVE SPLAT DENSITY (split / clone / prune), rung 4 of the fitting ladder.
2//
3// The fitter (gs_pos_grad_step_aniso) already measures, per primitive, everything this decision needs:
4// gradbuf[g*3], gradbuf[g*3+1] -- the screen-space residual gradient
5// gradbuf[g*3+2], wnorm[g] -- the Gauss-Newton curvature denominators
6// pa/pc[g] -- the conic, hence the screen footprint via gs_splat_rx/ry
7// depth[g] >= 0 -- visible in this view
8// So densification is a CLASSIFIER over state that already exists. It derives no new quantity and
9// re-renders nothing; a second measurement here would be a second ruler.
10//
11// EVERY THRESHOLD IS DERIVED FROM THE POPULATION OR THE FORMAT. NONE IS PICKED.
12// The published 3DGS rule densifies on a view-space gradient above a FIXED tau (0.0002 in the paper).
13// That constant is calibrated to their units, their scene scale and their float pipeline; copying it
14// into an integer fixed-point renderer would be a picked number wearing a citation. So the bar here is
15// the population's OWN MEAN, recomputed every step:
16// DENSIFY iff |gx| + |gy| > mean(|gx| + |gy|) over VISIBLE primitives
17// A mean needs no parameter, adapts to scene scale and to fitting progress by construction, and cannot
18// go stale. L1 is used rather than L2 because this is an ORDERING against a mean of the same measure --
19// no square root, exact in integers, and monotone in the same direction as L2.
20//
21// SPLIT vs CLONE, also against the population's own mean, and the reasoning is the paper's:
22// footprint > mean footprint -> SPLIT (over-reconstruction: one large splat spanning detail it
23// cannot represent; two smaller ones can)
24// footprint <= mean footprint -> CLONE (under-reconstruction: the region is under-covered, so add
25// a second primitive of the same size)
26// Footprint is the screen half-extent gs_splat_rx(ca) + gs_splat_ry(cc) -- the renderer's OWN support
27// radius, so a splat is judged by exactly the extent it actually paints.
28//
29// PRUNE, derived from the FIXED-POINT FORMAT rather than from the population:
30// a splat contributes c*T*alpha/GFXA to a channel; peak over T<=GFXA and alpha<=op gives c*op/GFXA.
31// With the framebuffer quantised to 8 bits (c <= 255), that peak is below one quantisation step when
32// op * 255 < GFXA
33// Such a primitive CANNOT change any pixel at any depth or colour. It is invisible by the arithmetic
34// of the format, not by a judgement about smallness.
35// A view-dependent prune was CONSIDERED AND REJECTED, and the reason is multi-view correctness.
36// Pruning on zero accumulated footprint in this view, or on depth[g] < 0, would delete primitives that
37// are culled HERE and visible from another camera -- and the fitter is now multi-view (it cycles yaw
38// across steps). A prune rule that is correct per-view and wrong per-scene is the worst kind, because
39// every single-view test passes. Opacity is view-independent, so the shipped rule is.
40//
41// OFFSET DIRECTION IS THE FITTER'S OWN, NOT A SECOND CONVENTION. gs_pos_grad_step_aniso moves a splat
42// by x -= gradbuf[0]*... and y += gradbuf[1]*... so descent is -sign(gx) in world X and +sign(gy)
43// in world Y. New primitives are offset along that same direction, so the densifier and the fitter can
44// never disagree about which way downhill is. Offset MAGNITUDE is rtan/2 -- the primitive's own
45// radius, so split children land inside the parent's support and the clone lands adjacent to it.
46//
47// CAPACITY IS A CALLER PARAMETER AND OVERFLOW REFUSES BY NAME. ngmax is the caller's allocation, not a
48// module constant -- there is no MAX_PRIMITIVES here and there must never be. If the classified output
49// would exceed it the organ REFUSES (-1) and leaves gauss untouched, because a silently truncated
50// densification is a fit that thinks it densified and did not.
51
52import "nx_syscalls.nx"
53import "nx_gsplat.nx"
54
55const GD_S_IN: i64 = 0
56const GD_S_KEPT: i64 = 1
57const GD_S_SPLIT: i64 = 2
58const GD_S_CLONED: i64 = 3
59const GD_S_PRUNED: i64 = 4
60const GD_S_OUT: i64 = 5
61const GD_S_GMEAN: i64 = 6
62const GD_S_RMEAN: i64 = 7
63const GD_S_VIS: i64 = 8
64const GD_S_REFUSED: i64 = 9
65// GD_S_ABOVEBAR -- how many primitives EXCEEDED the gradient bar. Added because a MUTATION BITE proved
66// the stats could not distinguish two very different states: "nothing needed densifying" and "the
67// criterion is broken and fires on nothing" BOTH report out == in, and the mutant's arms produced
68// byte-identical losses. With this counter the two are separable at a glance:
69// out == in AND above_bar == 0 -> a legitimate no-op; the population genuinely had no outlier
70// out == in AND above_bar > 0 -> a DEFECT SIGNATURE: the bar fired and the emit dropped it
71const GD_S_ABOVEBAR: i64 = 10
72const GD_STATN: i64 = 11
73
74const GD_C_PRUNE: i64 = 0
75const GD_C_KEEP: i64 = 1
76const GD_C_SPLIT: i64 = 2
77const GD_C_CLONE: i64 = 3
78
79func gd_stat_count() -> i64 { return GD_STATN }
80func gd_abs(v: i64) -> i64 { if v < 0 { return 0 - v } return v }
81
82// THE VISIBILITY QUANTUM OF THE FORMAT: op*255 < GFXA means the peak contribution to any channel at
83// any depth is below one 8-bit step. Derived, not chosen.
84func gd_invisible(op: i64) -> i64 { if op * 255 < gs_fxa() { return 1 } return 0 }
85
86// population means over VISIBLE primitives. Returns the visible count so the caller can bind every
87// downstream claim to a non-empty denominator -- a mean over zero primitives must never be used.
88func gd_means(ng: i64, depth: *i64, pa: *i64, pc: *i64, gradbuf: *i64, out2: *i64) -> i64 {
89 var vis: i64 = 0
90 var gsum: i64 = 0
91 var rsum: i64 = 0
92 var g: i64 = 0
93 while g < ng {
94 if depth[g] >= 0 {
95 vis = vis + 1
96 gsum = gsum + gd_abs(gradbuf[g*3]) + gd_abs(gradbuf[g*3+1])
97 rsum = rsum + gs_splat_rx(pa[g]) + gs_splat_ry(pc[g])
98 }
99 g = g + 1
100 }
101 out2[0] = 0
102 out2[1] = 0
103 if vis > 0 { out2[0] = gsum / vis; out2[1] = rsum / vis }
104 return vis
105}
106
107// classify every primitive into exactly one class, so the partition sums by construction.
108// The class is parked in gradbuf[g*3+2], the Gauss-Newton curvature denominator: it is DEAD the moment
109// the fitter's step returns (the step consumed it to scale its own move) and it is the only free
110// per-primitive slot here. Declared at the site rather than discovered later.
111func gd_classify(gauss: *i64, ng: i64, depth: *i64, pa: *i64, pc: *i64, gradbuf: *i64, gmean: i64, rmean: i64, stats: *i64) -> i64 {
112 let st: i64 = gs_stride_aniso()
113 var g: i64 = 0
114 while g < ng {
115 let op: i64 = gauss[g*st+10]
116 var cls: i64 = GD_C_KEEP
117 if gd_invisible(op) == 1 {
118 cls = GD_C_PRUNE
119 } else {
120 if depth[g] >= 0 {
121 let gm: i64 = gd_abs(gradbuf[g*3]) + gd_abs(gradbuf[g*3+1])
122 if gm > gmean {
123 stats[GD_S_ABOVEBAR] = stats[GD_S_ABOVEBAR] + 1
124 let fp: i64 = gs_splat_rx(pa[g]) + gs_splat_ry(pc[g])
125 if fp > rmean { cls = GD_C_SPLIT } else { cls = GD_C_CLONE }
126 }
127 }
128 }
129 gradbuf[g*3+2] = cls
130 if cls == GD_C_PRUNE { stats[GD_S_PRUNED] = stats[GD_S_PRUNED] + 1 }
131 if cls == GD_C_KEEP { stats[GD_S_KEPT] = stats[GD_S_KEPT] + 1 }
132 if cls == GD_C_SPLIT { stats[GD_S_SPLIT] = stats[GD_S_SPLIT] + 1 }
133 if cls == GD_C_CLONE { stats[GD_S_CLONED] = stats[GD_S_CLONED] + 1 }
134 g = g + 1
135 }
136 return 0
137}
138
139func gd_copy(src: *i64, si: i64, dst: *i64, di: i64, st: i64) -> i64 {
140 var k: i64 = 0
141 while k < st { dst[di*st+k] = src[si*st+k]; k = k + 1 }
142 return 0
143}
144
145// gs_densify -- classify, then emit. Returns the new primitive count, or -1 REFUSED on overflow.
146// outbuf is the caller's staging allocation (ngmax * gs_stride_aniso() i64). Building into it and
147// copying back keeps the emit trivially correct: an in-place expansion would overwrite unread inputs.
148func gs_densify(gauss: *i64, ng: i64, ngmax: i64, depth: *i64, pa: *i64, pc: *i64, gradbuf: *i64, outbuf: *i64, m2: *i64, stats: *i64) -> i64 {
149 let st: i64 = gs_stride_aniso()
150 var i: i64 = 0
151 while i < GD_STATN { stats[i] = 0; i = i + 1 }
152 stats[GD_S_IN] = ng
153 if ng <= 0 { stats[GD_S_OUT] = 0; return 0 }
154
155 let vis: i64 = gd_means(ng, depth, pa, pc, gradbuf, m2)
156 stats[GD_S_VIS] = vis
157 stats[GD_S_GMEAN] = m2[0]
158 stats[GD_S_RMEAN] = m2[1]
159
160 gd_classify(gauss, ng, depth, pa, pc, gradbuf, m2[0], m2[1], stats)
161
162 // a split or a clone each yield TWO primitives; a keep yields one; a prune yields none.
163 let out: i64 = stats[GD_S_KEPT] + 2*stats[GD_S_SPLIT] + 2*stats[GD_S_CLONED]
164 if out > ngmax {
165 stats[GD_S_REFUSED] = 1
166 stats[GD_S_OUT] = ng
167 return 0 - 1
168 }
169
170 var w: i64 = 0
171 var g: i64 = 0
172 while g < ng {
173 let cls: i64 = gradbuf[g*3+2]
174 if cls != GD_C_PRUNE {
175 let rtan: i64 = gauss[g*st+6]
176 var dx: i64 = 0
177 var dy: i64 = 0
178 if gradbuf[g*3] > 0 { dx = 0 - 1 }
179 if gradbuf[g*3] < 0 { dx = 1 }
180 if gradbuf[g*3+1] > 0 { dy = 1 }
181 if gradbuf[g*3+1] < 0 { dy = 0 - 1 }
182 let off: i64 = rtan / 2
183 if cls == GD_C_KEEP {
184 gd_copy(gauss, g, outbuf, w, st)
185 w = w + 1
186 }
187 if cls == GD_C_SPLIT {
188 gd_copy(gauss, g, outbuf, w, st)
189 outbuf[w*st+6] = rtan / 2
190 outbuf[w*st] = gauss[g*st] + dx*off
191 outbuf[w*st+1] = gauss[g*st+1] + dy*off
192 w = w + 1
193 gd_copy(gauss, g, outbuf, w, st)
194 outbuf[w*st+6] = rtan / 2
195 outbuf[w*st] = gauss[g*st] - dx*off
196 outbuf[w*st+1] = gauss[g*st+1] - dy*off
197 w = w + 1
198 }
199 if cls == GD_C_CLONE {
200 gd_copy(gauss, g, outbuf, w, st)
201 w = w + 1
202 gd_copy(gauss, g, outbuf, w, st)
203 outbuf[w*st] = gauss[g*st] + dx*off
204 outbuf[w*st+1] = gauss[g*st+1] + dy*off
205 w = w + 1
206 }
207 }
208 g = g + 1
209 }
210
211 i = 0
212 while i < w*st { gauss[i] = outbuf[i]; i = i + 1 }
213 stats[GD_S_OUT] = w
214 return w
215}