nx_isosurf.nx source
↩ module page · 144 lines · 7166 B
1// nx_isosurf.nx -- ★F0 UNIFY implicit + explicit (from the graphics foundation census): a SURFACE-NETS polygonizer that
2// turns a signed-distance FIELD into a real TRIANGLE MESH. This is the bridge that lets us MESH THE BEING (or any SDF):
3// sample the SDF on a grid, place one vertex per surface-straddling cell (averaged edge crossings), connect the 4 cells
4// around every sign-changing grid edge into a quad. Per-vertex normals come from the SDF GRADIENT (so smooth Gouraud
5// shading is winding-independent). Fills the nx_trimesh buffers -> render with the existing rasterizer. license_tier: ORIGINAL
6import "nx_syscalls.nx"
7import "nx_trimesh.nx"
8const K_MAGIC_1024: i64 = 1024
9const K_MAGIC_4096: i64 = 4096
10
11func iso_isqrt(v: i64) -> i64 { if v <= 0 { return 0 } var x: i64 = v; var y: i64 = (x+1)/2; while y < x { x=y; y=(x+v/x)/2 } return x }
12
13// ★polygonize a corner-SDF grid (dims (GX+1)x(GY+1)x(GZ+1)) at iso-level `iso`. Cells GX x GY x GZ. Origin ox,oy,oz;
14// corner (i,j,k) is at (ox+i*cell, oy+j*cell, oz+k*cell). grid[(i*(GY+1)+j)*(GZ+1)+k] = SDF value. Emits into trimesh.
15func surface_nets(grid: *i64, GX: i64, GY: i64, GZ: i64, ox: i64, oy: i64, oz: i64, cell: i64, iso: i64, col: i64) -> i64 {
16 let GY1: i64 = GY+1; let GZ1: i64 = GZ+1
17 let cv: *i64 = sys_mmap(GX*GY*GZ*8) as *i64
18 let cval: *i64 = sys_mmap(8*8) as *i64 // scratch, reused per cell (NOT allocated in the hot loop)
19 let ea: *i64 = sys_mmap(12*8) as *i64
20 let eb: *i64 = sys_mmap(12*8) as *i64
21 ea[0]=0;eb[0]=1; ea[1]=2;eb[1]=3; ea[2]=4;eb[2]=5; ea[3]=6;eb[3]=7
22 ea[4]=0;eb[4]=2; ea[5]=1;eb[5]=3; ea[6]=4;eb[6]=6; ea[7]=5;eb[7]=7
23 ea[8]=0;eb[8]=4; ea[9]=1;eb[9]=5; ea[10]=2;eb[10]=6; ea[11]=3;eb[11]=7
24 // ---- pass 1: one vertex per surface-straddling cell ----
25 var ci: i64 = 0
26 while ci < GX {
27 var cj: i64 = 0
28 while cj < GY {
29 var ck: i64 = 0
30 while ck < GZ {
31 cv[(ci*GY+cj)*GZ+ck] = 0-1
32 var bit: i64 = 0
33 while bit < 8 {
34 let dx: i64 = bit&1; let dy: i64 = (bit>>1)&1; let dz: i64 = (bit>>2)&1
35 cval[bit] = grid[((ci+dx)*GY1 + (cj+dy))*GZ1 + (ck+dz)]
36 bit = bit + 1
37 }
38 var sx: i64 = 0; var sy: i64 = 0; var sz: i64 = 0; var cnt: i64 = 0
39 var e: i64 = 0
40 while e < 12 {
41 let a: i64 = ea[e]; let b: i64 = eb[e]
42 var ina: i64 = 0; if cval[a] <= iso { ina = 1 } // inside = val<=iso (handles exact-0 corners)
43 var inb: i64 = 0; if cval[b] <= iso { inb = 1 }
44 if ina != inb {
45 let den: i64 = cval[b]-cval[a]
46 var t: i64 = 0
47 if den != 0 { t = (iso-cval[a])*K_MAGIC_1024/den }
48 let ax: i64 = ox+(ci+(a&1))*cell; let ay: i64 = oy+(cj+((a>>1)&1))*cell; let az: i64 = oz+(ck+((a>>2)&1))*cell
49 let bx: i64 = ox+(ci+(b&1))*cell; let by: i64 = oy+(cj+((b>>1)&1))*cell; let bz: i64 = oz+(ck+((b>>2)&1))*cell
50 sx = sx + ax + (bx-ax)*t/K_MAGIC_1024
51 sy = sy + ay + (by-ay)*t/K_MAGIC_1024
52 sz = sz + az + (bz-az)*t/K_MAGIC_1024
53 cnt = cnt + 1
54 }
55 e = e + 1
56 }
57 if cnt > 0 {
58 let vx: i64 = sx/cnt; let vy: i64 = sy/cnt; let vz: i64 = sz/cnt
59 // normal = SDF gradient at the base corner (central differences, clamped)
60 var ai: i64 = ci+1; if ai>GX {ai=GX}
61 var bi: i64 = ci-1; if bi<0 {bi=0}
62 var aj: i64 = cj+1; if aj>GY {aj=GY}
63 var bj: i64 = cj-1; if bj<0 {bj=0}
64 var ak: i64 = ck+1; if ak>GZ {ak=GZ}
65 var bk: i64 = ck-1; if bk<0 {bk=0}
66 let gx: i64 = grid[(ai*GY1+cj)*GZ1+ck] - grid[(bi*GY1+cj)*GZ1+ck]
67 let gy: i64 = grid[(ci*GY1+aj)*GZ1+ck] - grid[(ci*GY1+bj)*GZ1+ck]
68 let gz: i64 = grid[(ci*GY1+cj)*GZ1+ak] - grid[(ci*GY1+cj)*GZ1+bk]
69 let gl: i64 = iso_isqrt(gx*gx+gy*gy+gz*gz) + 1
70 let idx: i64 = tm_vert_n(vx, vy, vz, gx*K_MAGIC_4096/gl, gy*K_MAGIC_4096/gl, gz*K_MAGIC_4096/gl)
71 cv[(ci*GY+cj)*GZ+ck] = idx
72 }
73 ck = ck + 1
74 }
75 cj = cj + 1
76 }
77 ci = ci + 1
78 }
79 // ---- pass 2: quad per sign-changing INTERIOR corner-grid edge (4 surrounding cells) ----
80 // x-edges
81 var i: i64 = 0
82 while i < GX {
83 var j: i64 = 1
84 while j < GY {
85 var k: i64 = 1
86 while k < GZ {
87 var q0: i64 = 0; if grid[(i*GY1+j)*GZ1+k] <= iso { q0 = 1 }
88 var q1: i64 = 0; if grid[((i+1)*GY1+j)*GZ1+k] <= iso { q1 = 1 }
89 if q0 != q1 {
90 let m: i64 = cv[(i*GY+(j-1))*GZ+(k-1)]; let n: i64 = cv[(i*GY+j)*GZ+(k-1)]; let o: i64 = cv[(i*GY+j)*GZ+k]; let p: i64 = cv[(i*GY+(j-1))*GZ+k]
91 if m>=0 { if n>=0 { if o>=0 { if p>=0 {
92 if q0 == 1 { tm_quad(m,n,o,p,col) } else { tm_quad(m,p,o,n,col) } // sign-consistent outward winding
93 } } } }
94 }
95 k = k + 1
96 }
97 j = j + 1
98 }
99 i = i + 1
100 }
101 // y-edges
102 i = 1
103 while i < GX {
104 var j: i64 = 0
105 while j < GY {
106 var k: i64 = 1
107 while k < GZ {
108 var q0: i64 = 0; if grid[(i*GY1+j)*GZ1+k] <= iso { q0 = 1 }
109 var q1: i64 = 0; if grid[(i*GY1+(j+1))*GZ1+k] <= iso { q1 = 1 }
110 if q0 != q1 {
111 let m: i64 = cv[((i-1)*GY+j)*GZ+(k-1)]; let n: i64 = cv[(i*GY+j)*GZ+(k-1)]; let o: i64 = cv[(i*GY+j)*GZ+k]; let p: i64 = cv[((i-1)*GY+j)*GZ+k]
112 if m>=0 { if n>=0 { if o>=0 { if p>=0 {
113 if q0 == 1 { tm_quad(m,n,o,p,col) } else { tm_quad(m,p,o,n,col) } // sign-consistent outward winding
114 } } } }
115 }
116 k = k + 1
117 }
118 j = j + 1
119 }
120 i = i + 1
121 }
122 // z-edges
123 i = 1
124 while i < GX {
125 var j: i64 = 1
126 while j < GY {
127 var k: i64 = 0
128 while k < GZ {
129 var q0: i64 = 0; if grid[(i*GY1+j)*GZ1+k] <= iso { q0 = 1 }
130 var q1: i64 = 0; if grid[(i*GY1+j)*GZ1+(k+1)] <= iso { q1 = 1 }
131 if q0 != q1 {
132 let m: i64 = cv[((i-1)*GY+(j-1))*GZ+k]; let n: i64 = cv[(i*GY+(j-1))*GZ+k]; let o: i64 = cv[(i*GY+j)*GZ+k]; let p: i64 = cv[((i-1)*GY+j)*GZ+k]
133 if m>=0 { if n>=0 { if o>=0 { if p>=0 {
134 if q0 == 1 { tm_quad(m,n,o,p,col) } else { tm_quad(m,p,o,n,col) } // sign-consistent outward winding
135 } } } }
136 }
137 k = k + 1
138 }
139 j = j + 1
140 }
141 i = i + 1
142 }
143 return 0
144}