nx_gsplat_aniso_gate.nx source

1// nx_gsplat_aniso_gate.nx -- ★ANISOTROPIC surface splatting (the defining 3DGS feature). Oriented ellipse splats 2// (surfels aligned to the surface normal) instead of isotropic blobs -> smooth coverage, no dotty grid. 3// T1 ★ANISOTROPY: an EDGE-ON surfel (normal in the screen plane) renders a THIN, elongated ellipse; a FACE-ON 4// surfel (normal toward camera) renders a round splat -- the aspect ratios differ sharply (orientation is real) 5// T2 the anatomical mesh, splatted as surface-aligned surfels, covers the face densely (no isotropic-grid holes) 6// + eyeball PNG knowledge/nx_gsplat_face_aniso.png 7// T3 determinism 8// license_tier: ORIGINAL expect_exit: 0 9import "nx_syscalls.nx" 10import "nx_faceanat.nx" 11import "nx_meshgen.nx" 12import "nx_gsplat.nx" 13import "nx_png.nx" 14 15func hw(s: *u8) -> i64 { var n: i64=0; while s[n]!=(0 as u8){n=n+1} sys_write(1,s,n); return 0 } 16func pn(v: i64) -> i64 { let b: *u8=sys_mmap(32) as *u8; var x: i64=v; var ng: i64=0; if x<0{ng=1;x=0-x} var i: i64=31; if x==0{b[i]=48 as u8;i=i-1} while x>0{b[i]=(48+x%10) as u8;x=x/10;i=i-1} if ng==1{b[i]=45 as u8;i=i-1} sys_write(1,(b as i64+i+1) as *u8,31-i); return 0 } 17func lum(px: i64) -> i64 { return (px & 255) + ((px >> 8) & 255) + ((px >> 16) & 255) } 18 19func main() -> i64 { 20 hw("=== nx_gsplat_aniso_gate -- anisotropic surface splatting (oriented ellipse surfels) ===\n" as *u8) 21 var fails: i64 = 0 22 let W: i64 = gs_w() 23 let H: i64 = gs_h() 24 let npx: i64 = W*H 25 let NGX: i64 = 6000 26 let gauss: *i64 = sys_mmap(NGX*12*8) as *i64 27 let fb: *i64 = sys_mmap(npx*8) as *i64 28 let acc: *i64 = sys_mmap(npx*3*8) as *i64 29 let trans: *i64 = sys_mmap(npx*8) as *i64 30 let depth: *i64 = sys_mmap(NGX*8) as *i64 31 let sxb: *i64 = sys_mmap(NGX*8) as *i64 32 let syb: *i64 = sys_mmap(NGX*8) as *i64 33 let pa: *i64 = sys_mmap(NGX*8) as *i64 34 let pb: *i64 = sys_mmap(NGX*8) as *i64 35 let pc: *i64 = sys_mmap(NGX*8) as *i64 36 let pdet: *i64 = sys_mmap(NGX*8) as *i64 37 let order: *i64 = sys_mmap(NGX*8) as *i64 38 let count: *i64 = sys_mmap((gs_nb()+2)*8) as *i64 39 let explut: *i64 = sys_mmap(gs_expn()*8) as *i64 40 gs_build_explut(explut) 41 42 // measure the non-bg bounding box of the splat 43 // T1a EDGE-ON: normal = +x (in the screen plane) -> tangent plane holds z(depth)+y -> thin vertical ellipse 44 gs_set_aniso(gauss, 0, 0, 0, 0, 256, 0, 0, 140, 240, 240, 240, 256) 45 gs_render_aniso(gauss, 1, 0, 4, fb, acc, trans, depth, sxb, syb, pa, pb, pc, pdet, order, count, explut, 26, 28, 44) 46 var exminx: i64 = W; var exmaxx: i64 = 0; var exminy: i64 = H; var exmaxy: i64 = 0 47 var pi: i64 = 0 48 while pi < npx { if lum(fb[pi]) > 200 { let px: i64=pi%W; let py: i64=pi/W; if px<exminx{exminx=px} if px>exmaxx{exmaxx=px} if py<exminy{exminy=py} if py>exmaxy{exmaxy=py} } pi=pi+1 } 49 let ew: i64 = exmaxx - exminx 50 let eh: i64 = exmaxy - exminy 51 // T1b FACE-ON: normal = -z (toward camera at yaw0) -> tangent plane = x-y -> round splat 52 gs_set_aniso(gauss, 0, 0, 0, 0, 0, 0, 0-256, 140, 240, 240, 240, 256) 53 gs_render_aniso(gauss, 1, 0, 4, fb, acc, trans, depth, sxb, syb, pa, pb, pc, pdet, order, count, explut, 26, 28, 44) 54 var fminx: i64 = W; var fmaxx: i64 = 0; var fminy: i64 = H; var fmaxy: i64 = 0 55 pi = 0 56 while pi < npx { if lum(fb[pi]) > 200 { let px: i64=pi%W; let py: i64=pi/W; if px<fminx{fminx=px} if px>fmaxx{fmaxx=px} if py<fminy{fminy=py} if py>fmaxy{fmaxy=py} } pi=pi+1 } 57 let fw: i64 = fmaxx - fminx 58 let fh: i64 = fmaxy - fminy 59 hw(" edge-on splat wxh="); pn(ew); hw("x"); pn(eh); hw(" face-on wxh="); pn(fw); hw("x"); pn(fh); hw("\n" as *u8) 60 var t1: i64 = 0 61 // edge-on: strongly elongated (height >> width); face-on: roughly round (w ~ h) 62 if eh > ew * 2 { if fw*10 > fh*7 { if fh*10 > fw*7 { t1 = 1 } } } 63 if t1 == 1 { hw("T1 PASS anisotropy: edge-on surfel = thin ellipse, face-on = round -> orientation is real\n" as *u8) } 64 else { fails=fails+1; hw("T1 FAIL aspect ratios\n" as *u8) } 65 66 // T2 the anatomical mesh as surface-aligned surfels 67 let base: i64 = sys_mmap(sdf_bytes()) as i64 68 let Fg: *i64 = sys_mmap((MG_N+1)*(MG_N+1)*(MG_N+1)*8) as *i64 69 let cubevi: *i64 = sys_mmap(MG_N*MG_N*MG_N*8) as *i64 70 let vbuf: *i64 = sys_mmap(MG_MAXV*3*8) as *i64 71 let fbuf: *i64 = sys_mmap(MG_MAXF*3*8) as *i64 72 let mo: *i64 = sys_mmap(16) as *i64 73 faceanat_build(base) 74 mg_build(base, Fg, cubevi, vbuf, fbuf, mo) 75 let nv: i64 = mo[0] 76 var vi: i64 = 0 77 while vi < nv { 78 let x: i64 = vbuf[vi*3] 79 let y: i64 = vbuf[vi*3+1] 80 let z: i64 = vbuf[vi*3+2] 81 let e: i64 = 12 82 var gx: i64 = sdf_eval(base, x+e, y, z) - sdf_eval(base, x-e, y, z) 83 var gy: i64 = sdf_eval(base, x, y+e, z) - sdf_eval(base, x, y-e, z) 84 var gz: i64 = sdf_eval(base, x, y, z+e) - sdf_eval(base, x, y, z-e) 85 let gl: i64 = gs_isqrt(gx*gx + gy*gy + gz*gz) 86 var nx: i64 = 0; var ny: i64 = 256; var nz: i64 = 0 87 if gl > 0 { nx = gx*256/gl; ny = gy*256/gl; nz = gz*256/gl } 88 var sh: i64 = 130 + (nx*(0-121) + ny*191 + nz*(0-121))/256 * 90 / 256 89 if sh < 70 { sh = 70 } 90 if sh > 255 { sh = 255 } 91 gs_set_aniso(gauss, vi, x, y, z, nx, ny, nz, 34, sh*236/255, sh*180/255, sh*156/255, 235) 92 vi = vi + 1 93 } 94 gs_render_aniso(gauss, nv, 200, 4, fb, acc, trans, depth, sxb, syb, pa, pb, pc, pdet, order, count, explut, 26, 28, 44) 95 var face: i64 = 0 96 pi = 0 97 while pi < npx { if lum(fb[pi]) > 130 { face = face + 1 } pi = pi + 1 } 98 hw(" aniso face: "); pn(nv); hw(" surfels -> "); pn(face); hw(" covered px\n" as *u8) 99 write_png(fb, W, H, "knowledge/nx_gsplat_face_aniso.png" as *u8) 100 var t2: i64 = 0 101 if face > 12000 { t2 = 1 } 102 if t2 == 1 { hw("T2 PASS surface-aligned surfels densely cover the face (PNG knowledge/nx_gsplat_face_aniso.png)\n" as *u8) } 103 else { fails=fails+1; hw("T2 FAIL coverage="); pn(face); hw("\n" as *u8) } 104 105 // T3 determinism 106 let fb2: *i64 = sys_mmap(npx*8) as *i64 107 gs_render_aniso(gauss, nv, 200, 4, fb2, acc, trans, depth, sxb, syb, pa, pb, pc, pdet, order, count, explut, 26, 28, 44) 108 var diff: i64 = 0 109 pi = 0 110 while pi < npx { if fb[pi] != fb2[pi] { diff = diff + 1 } pi = pi + 1 } 111 var t3: i64 = 0 112 if diff == 0 { t3 = 1 } 113 if t3 == 1 { hw("T3 PASS deterministic\n" as *u8) } 114 else { fails=fails+1; hw("T3 FAIL diff="); pn(diff); hw("\n" as *u8) } 115 116 if fails == 0 { hw("GSPLAT-ANISO-GATE 3/3 GREEN -- anisotropic surface splatting: oriented ellipse surfels, dense surface coverage, integer + deterministic (the defining 3DGS feature, sovereign)\n" as *u8); sys_exit(0); return 0 } 117 hw("GSPLAT-ANISO-GATE RED fails="); pn(fails); hw("\n" as *u8) 118 sys_exit(1) 119 return 1 120}