use geometry package octree; improve gsdfaux.RenderShader3D logging (#13)

Author: soypat

glrender/dual_contour.go

@@ -4,15 +4,17 @@ import (
 	"errors"
 
 	"github.com/chewxy/math32"
+	"github.com/soypat/geometry/i3"
 	"github.com/soypat/geometry/ms3"
 	"github.com/soypat/gsdf/gleval"
 )
 
 type DualContourRenderer struct {
 	sdf       gleval.SDF3
+	oct       ms3.Octree
 	contourer DualContourer
-	cubeMap   map[ivec]int
-	cubebuf   []icube
+	cubeMap   map[i3.Vec]int
+	cubebuf   []i3.Cube
 	posbuf    []ms3.Vec
 	distbuf   []float32
 	// cubeinfo stores both cube and edge information for dual contouring algorithm.
@@ -33,15 +35,16 @@ func (dcr *DualContourRenderer) Reset(sdf gleval.SDF3, res float32, vertexPlacer
 	if err != nil {
 		return err
 	}
-	nCubes := int(topCube.decomposesTo(1))
+	nCubes := int(topCube.DecomposesTo(1))
 	if cap(dcr.cubebuf) < nCubes*2 {
-		dcr.cubebuf = make([]icube, 0, nCubes*2) // Want to fully decompose cube.
+		dcr.cubebuf = make([]i3.Cube, 0, nCubes*2) // Want to fully decompose cube.
 	}
 	var ok bool
-	dcr.cubebuf, ok = octreeDecomposeBFS(dcr.cubebuf[:0], topCube, 1)
+	dcr.oct = ms3.Octree{Resolution: res, Origin: origin}
+	dcr.cubebuf, ok = dcr.oct.DecomposeBFS(dcr.cubebuf[:0], topCube, 1)
 	if !ok {
 		return errors.New("unable to decompose top level cube")
-	} else if dcr.cubebuf[0].lvl != 1 {
+	} else if dcr.cubebuf[0].Level != 1 {
 		return errors.New("short buffer decomposing all cubes")
 	} else if len(dcr.cubebuf) != nCubes {
 		panic("failed to decompose all edges?")
@@ -54,7 +57,7 @@ func (dcr *DualContourRenderer) Reset(sdf gleval.SDF3, res float32, vertexPlacer
 	posbuf := dcr.posbuf[:cap(dcr.posbuf)]
 
 	// Keep cubes that contain a surface and their neighbors.
-	dcr.cubebuf, _, err = octreePrune(sdf, dcr.cubebuf, origin, res, posbuf, dcr.distbuf[:len(posbuf)], userData, 2, true)
+	dcr.cubebuf, _, err = octreePrunea(sdf, dcr.cubebuf, origin, res, posbuf, dcr.distbuf[:len(posbuf)], userData, 2, true)
 	if err != nil {
 		return err
 	}
@@ -64,7 +67,7 @@ func (dcr *DualContourRenderer) Reset(sdf gleval.SDF3, res float32, vertexPlacer
 		dcr.cubeinfo = make([]DualCube, nCubes)
 	}
 	if dcr.cubeMap == nil {
-		dcr.cubeMap = make(map[ivec]int)
+		dcr.cubeMap = make(map[i3.Vec]int)
 	}
 	clear(dcr.cubeMap)
 	dcr.res = res
@@ -84,8 +87,9 @@ func (dcr *DualContourRenderer) RenderAll(dst []ms3.Triangle, userData any) ([]m
 	cubes := dcr.cubeinfo[:len(edges)]
 
 	for e, edge := range edges {
-		sz := edge.size(res)
-		edgeOrig := edge.origin(origin, sz)
+
+		sz := dcr.oct.CubeSize(edge)
+		edgeOrig := dcr.oct.CubeOrigin(edge, sz)
 		// posbuf has edge origin, edge extremes in x,y,z and the center of the voxel.
 		posbuf = append(posbuf,
 			edgeOrig,
@@ -94,7 +98,7 @@ func (dcr *DualContourRenderer) RenderAll(dst []ms3.Triangle, userData any) ([]m
 			ms3.Add(edgeOrig, ms3.Vec{Z: sz}),
 		)
 		cubes[e].FinalVertex = edgeOrig // By default set to center.
-		cubeMap[edge.ivec] = e
+		cubeMap[edge.Vec] = e
 	}
 
 	lenPos := len(posbuf)
@@ -108,7 +112,7 @@ func (dcr *DualContourRenderer) RenderAll(dst []ms3.Triangle, userData any) ([]m
 	posbuf = posbuf[:0]
 	for e, edge := range edges {
 		// First for loop accumulates edge biases into voxels/cubes.
-		cube := makeDualCube(edge.ivec, distbuf[e*4:])
+		cube := makeDualCube(edge.Vec, distbuf[e*4:])
 		cubes[e] = cube
 		if !cube.IsActive() {
 			continue
@@ -180,7 +184,7 @@ func (dcr *DualContourRenderer) RenderAll(dst []ms3.Triangle, userData any) ([]m
 	return dst, nil
 }
 
-func makeDualCube(ivec ivec, data []float32) DualCube {
+func makeDualCube(ivec i3.Vec, data []float32) DualCube {
 	if len(data) < 4 {
 		panic("short dual cube info buffer")
 	}
@@ -196,7 +200,7 @@ func makeDualCube(ivec ivec, data []float32) DualCube {
 // DualCube corresponds to a voxel anmd contains both cube and edge data.
 type DualCube struct {
 	// ivec stores the octree index of the cube, used to find neighboring cube ivec indices and the absolute position of the cube.
-	ivec ivec
+	ivec i3.Vec
 	// Neighbors contains neighboring index into dualCube buffer and contributing edge intersect axis.
 	//  - Neighbors[0]: Index into dualCube buffer to cube neighbor with edge.
 	//  - Neighbors[1]: Intersecting axis. 0 is x; 1 is y; 2 is z.
@@ -212,7 +216,9 @@ type DualCube struct {
 }
 
 func (dc *DualCube) SizeAndOrigin(res float32, octreeOrigin ms3.Vec) (float32, ms3.Vec) {
-	return res, icube{ivec: dc.ivec, lvl: 1}.origin(octreeOrigin, res)
+	cb := i3.Cube{Vec: dc.ivec, Level: 1}
+	oct := ms3.Octree{Resolution: res, Origin: octreeOrigin}
+	return res, oct.CubeOrigin(cb, oct.CubeSize(cb))
 }
 
 func (dc *DualCube) IsActive() bool {
@@ -235,22 +241,22 @@ func (dc *DualCube) FlipX() bool           { return dc.XDist-dc.OrigDist < 0 }
 func (dc *DualCube) FlipY() bool           { return dc.YDist-dc.OrigDist < 0 }
 func (dc *DualCube) FlipZ() bool           { return dc.ZDist-dc.OrigDist < 0 }
 
-func (dc *DualCube) EdgeNeighborsX() [4]ivec {
+func (dc *DualCube) EdgeNeighborsX() [4]i3.Vec {
 	const sub = -(1 << minIcubeLvl)
 	v := dc.ivec
-	return [4]ivec{v.Add(ivec{y: sub, z: sub}), v.Add(ivec{z: sub}), v, v.Add(ivec{y: sub})}
+	return [4]i3.Vec{v.Add(i3.Vec{Y: sub, Z: sub}), v.Add(i3.Vec{Z: sub}), v, v.Add(i3.Vec{Y: sub})}
 }
 
-func (dc *DualCube) EdgeNeighborsY() [4]ivec {
+func (dc *DualCube) EdgeNeighborsY() [4]i3.Vec {
 	const sub = -(1 << minIcubeLvl)
 	v := dc.ivec
-	return [4]ivec{v.Add(ivec{x: sub, z: sub}), v.Add(ivec{x: sub}), v, v.Add(ivec{z: sub})}
+	return [4]i3.Vec{v.Add(i3.Vec{X: sub, Z: sub}), v.Add(i3.Vec{X: sub}), v, v.Add(i3.Vec{Z: sub})}
 }
 
-func (dc *DualCube) EdgeNeighborsZ() [4]ivec {
+func (dc *DualCube) EdgeNeighborsZ() [4]i3.Vec {
 	const sub = -(1 << minIcubeLvl)
 	v := dc.ivec
-	return [4]ivec{v.Add(ivec{x: sub, y: sub}), v.Add(ivec{y: sub}), v, v.Add(ivec{x: sub})}
+	return [4]i3.Vec{v.Add(i3.Vec{X: sub, Y: sub}), v.Add(i3.Vec{Y: sub}), v, v.Add(i3.Vec{X: sub})}
 }
 
 // minecraftRender performs a minecraft-like render of the SDF using a dual contour method.
@@ -261,21 +267,27 @@ func minecraftRender(dst []ms3.Triangle, sdf gleval.SDF3, res float32) ([]ms3.Tr
 	if err != nil {
 		return dst, err
 	}
-	decomp := topCube.decomposesTo(1)
-	cubes := make([]icube, 0, decomp*2)
-	cubes, ok := octreeDecomposeBFS(cubes, topCube, 1)
+	oct := ms3.Octree{
+		Resolution: res,
+		Origin:     origin,
+	}
+
+	decomp := topCube.DecomposesTo(1)
+	cubes := make([]i3.Cube, 0, decomp*2)
+
+	cubes, ok := oct.DecomposeBFS(cubes, topCube, 1)
 	if !ok {
 		return dst, errors.New("unable to decompose top level cube")
-	} else if cubes[0].lvl != 1 {
+	} else if cubes[0].Level != 1 {
 		return dst, errors.New("short buffer decomposing all cubes")
 	}
 
 	var posbuf []ms3.Vec
 	iCubes := 0
 	for ; iCubes < len(cubes); iCubes++ {
 		cube := cubes[iCubes]
-		sz := cube.size(res)
-		cubeOrig := cube.origin(origin, sz)
+		sz := oct.CubeSize(cube)
+		cubeOrig := oct.CubeOrigin(cube, sz)
 		// Append origin and edge-end vertices.
 		posbuf = append(posbuf,
 			cubeOrig,
@@ -292,9 +304,9 @@ func minecraftRender(dst []ms3.Triangle, sdf gleval.SDF3, res float32) ([]ms3.Tr
 	}
 	for j := 0; j < iCubes; j++ {
 		cube := cubes[j]
-		sz := cube.size(res)
-		srcOrig := cube.origin(origin, sz)
-		dci := makeDualCube(cube.ivec, distbuf[j*4:])
+		sz := oct.CubeSize(cube)
+		srcOrig := oct.CubeOrigin(cube, sz)
+		dci := makeDualCube(cube.Vec, distbuf[j*4:])
 		origOff := sz
 		if dci.ActiveX() {
 			xOrig := ms3.Add(srcOrig, ms3.Vec{X: origOff})

glrender/glrender.go

@@ -35,110 +35,8 @@ func RenderAll(r Renderer, userData any) ([]ms3.Triangle, error) {
 	return result, err
 }
 
-type ivec struct {
-	x int
-	y int
-	z int
-}
-
-func (a ivec) Add(b ivec) ivec         { return ivec{x: a.x + b.x, y: a.y + b.y, z: a.z + b.z} }
-func (a ivec) AddScalar(f int) ivec    { return ivec{x: a.x + f, y: a.y + f, z: a.z + f} }
-func (a ivec) MulScalar(f int) ivec    { return ivec{x: a.x * f, y: a.y * f, z: a.z * f} }
-func (a ivec) DivScalar(f int) ivec    { return ivec{x: a.x / f, y: a.y / f, z: a.z / f} }
-func (a ivec) ShiftRight(lo int) ivec  { return ivec{x: a.x >> lo, y: a.y >> lo, z: a.z >> lo} }
-func (a ivec) ShiftLeft(hi int) ivec   { return ivec{x: a.x << hi, y: a.y << hi, z: a.z << hi} }
-func (a ivec) Sub(b ivec) ivec         { return ivec{x: a.x - b.x, y: a.y - b.y, z: a.z - b.z} }
-func (a ivec) Vec() ms3.Vec            { return ms3.Vec{X: float32(a.x), Y: float32(a.y), Z: float32(a.z)} }
-func (a ivec) AndScalar(f int) ivec    { return ivec{x: a.x & f, y: a.y & f, z: a.z & f} }
-func (a ivec) OrScalar(f int) ivec     { return ivec{x: a.x | f, y: a.y | f, z: a.z | f} }
-func (a ivec) XorScalar(f int) ivec    { return ivec{x: a.x ^ f, y: a.y ^ f, z: a.z ^ f} }
-func (a ivec) AndnotScalar(f int) ivec { return ivec{x: a.x &^ f, y: a.y &^ f, z: a.z &^ f} }
-
 const minIcubeLvl = 1
 
-type icube struct {
-	ivec
-	lvl int
-}
-
-func (c icube) isSmallest() bool       { return c.lvl == 1 }
-func (c icube) isSecondSmallest() bool { return c.lvl == 2 }
-
-// decomposesTo returns the amount of cubes generated from decomposing the cube down to cubes of the argument target level.
-func (c icube) decomposesTo(targetLvl int) uint64 {
-	if targetLvl > c.lvl {
-		panic("invalid targetLvl to icube.decomposesTo")
-	}
-	return pow8(c.lvl - targetLvl)
-}
-
-func (c icube) size(baseRes float32) float32 {
-	dim := 1 << (c.lvl - 1)
-	return float32(dim) * baseRes
-}
-
-// supercube returns the icube's parent octree icube.
-func (c icube) supercube() icube {
-	upLvl := c.lvl + 1
-	bitmask := (1 << upLvl) - 1
-	return icube{
-		ivec: c.ivec.AndnotScalar(bitmask),
-		lvl:  upLvl,
-	}
-}
-
-func (c icube) box(origin ms3.Vec, size float32) ms3.Box {
-	origin = c.origin(origin, size) // Replace origin with icube origin.
-	return ms3.Box{
-		Min: origin,
-		Max: ms3.AddScalar(size, origin),
-	}
-}
-
-func (c icube) origin(origin ms3.Vec, size float32) ms3.Vec {
-	idx := c.lvlIdx()
-	return ms3.Add(origin, ms3.Scale(size, idx.Vec()))
-}
-
-func (c icube) lvlIdx() ivec {
-	return c.ivec.ShiftRight(c.lvl) // icube indices per level in the octree.
-}
-
-func (c icube) center(origin ms3.Vec, size float32) ms3.Vec {
-	return c.box(origin, size).Center() // TODO(soypat): this can probably be optimized.
-}
-
-// corners returns the cube corners. Be aware size is NOT the minimum cube resolution but
-// can be calculated with the [icube.size] method using resolution. If [icube.lvl]==1 then size is resolution.
-func (c icube) corners(origin ms3.Vec, size float32) [8]ms3.Vec {
-	origin = c.origin(origin, size)
-	return [8]ms3.Vec{
-		ms3.Add(origin, ms3.Vec{X: 0, Y: 0, Z: 0}),
-		ms3.Add(origin, ms3.Vec{X: size, Y: 0, Z: 0}),
-		ms3.Add(origin, ms3.Vec{X: size, Y: size, Z: 0}),
-		ms3.Add(origin, ms3.Vec{X: 0, Y: size, Z: 0}),
-		ms3.Add(origin, ms3.Vec{X: 0, Y: 0, Z: size}),
-		ms3.Add(origin, ms3.Vec{X: size, Y: 0, Z: size}),
-		ms3.Add(origin, ms3.Vec{X: size, Y: size, Z: size}),
-		ms3.Add(origin, ms3.Vec{X: 0, Y: size, Z: size}),
-	}
-}
-
-func (c icube) octree() [8]icube {
-	lvl := c.lvl - 1
-	s := 1 << lvl
-	return [8]icube{
-		{ivec: c.Add(ivec{0, 0, 0}), lvl: lvl},
-		{ivec: c.Add(ivec{s, 0, 0}), lvl: lvl},
-		{ivec: c.Add(ivec{s, s, 0}), lvl: lvl},
-		{ivec: c.Add(ivec{0, s, 0}), lvl: lvl},
-		{ivec: c.Add(ivec{0, 0, s}), lvl: lvl},
-		{ivec: c.Add(ivec{s, 0, s}), lvl: lvl},
-		{ivec: c.Add(ivec{s, s, s}), lvl: lvl},
-		{ivec: c.Add(ivec{0, s, s}), lvl: lvl},
-	}
-}
-
 func min(a, b int) int {
 	if a < b {
 		return a